Newton Method to Recover the Phase Accumulated during MRI Data Acquisition

For an internal conductivity image, magnetic resonance electrical impedance tomography MREIT injects an electric current into an object and measures the induced magnetic flux density, which appears in the phase part of the acquired MR image data. To maximize signal intensity, the injected current nonlinear encoding ICNE method extends the duration of the current injection until the end of the MR data reading. It disturbs the usual linear encoding of the MR k-space data used in the inverse Fourier transform. In this study, we estimate the magnetic flux density, which is recoverable from nonlinearly encoded MR k-space data by applying a Newton method.


Introduction
An electric current injected into an electrically conducting object, such as the human body, induces an internal distribution of the magnetic flux density B B x , B y , B z .Magnetic resonance electrical impedance tomography MREIT visualizes the internal conductivity distribution from the z-component B z of B which can be measured in practice using an MRI scanner.This technique was originally proposed by Joy et al. in 1989 1 ; since then, several researchers 2-10 have investigated and further developed MREIT as well as magnetic resonance current density image MRCDI , which has similar modalities 11-13 .
The magnetic flux density B z induced by injecting current through the electrodes attached on the surface of a conducting object Ω accumulates its signals in the phase parts of acquired MR image data.The conventional current-injection method 1, 8 injects the current during time T c , between the end of the first RF pulse and the beginning of the reading gradient, in order to ensure gradient linearity.
Since the signal-to-noise ratio SNR of the MR magnitude depends on the echo time T E , it is impossible to increase both T c and the SNR of the MR magnitude simultaneously in order to reduce noise effects.As an attempt to reduce the noise level, the injected current nonlinear encoding ICNE method was developed in 2007 14, 15 ; it extends the duration of the injection current until the end of a reading gradient in order to maximize the signal intensity of B z .Then, it disturbs the usual linear encoding of MR k-space data used in the inverse Fourier transform.
For example, the one-dimensional inverse problem in the conventional acquisition method is to find the unknown discrete magnetic flux density b l , l 0, . . ., N − 1 from the N × N matrix A satisfying the following: where T c is a constant and ρ, S are known quantities that can be measured.By an inverse Fourier transform F −1 , the first equation in 1.1 becomes The unknown data b l is simply recovered from 1.2 .In the ICNE method, however, the matrix A is perturbed to A kl e −i 2π/N k l b l .Then, it becomes a system of nonlinear equations for unknown b l data, where the conventional inverse Fourier transform is no longer applicable.
In this paper, we prove a unique determination of the magnetic flux density from measured MR signal obtained by the ICNE acquisition method.Secondly, applying a Newton method, we suggest a bound of l 2 -norm for recoverable magnetic flux density from nonlinearly encoded MR k-space data.Numerical experiments show the feasibility of the proposed method.

ICNE Method and Invertibility
For a standard spin echo pulse sequence in MR imaging, the k-space MR signal In MREIT, we inject the current I through the electrodes attached on the threedimensional conducting object Ω, having conductivity distribution σ.The injection current I produces the internal current density J and the magnetic flux density B B x , B y , B z in Ω, satisfying the Ampère and Biot-Savart laws.Since an MRI scanner measures only the main magnetic field direction component of B, the z-component B z , we focus on the problem of measuring B z x, y B z x, y, z 0 , where z 0 is the center of the selected imaging slice.Since MREIT is a methodology for reconstructing the internal conductivity σ from B z data, it is important to measure B z more precisely.

Conventional B z Acquisition
For a conventional B z acquisition, current is not injected during T s of the MR data acquisition, ADC as shown in Figure 1.In this case, the induced magnetic flux density B z provides additional dephasing of spins, and, consequently, extra phase is accumulated during the total injection time T c .Then, the measured k-space data for the injection current I can be represented as follows: ρ x, y e iδ x,y e −iγT c B z x,y e −i2π k x x k y y dx dy, 2.3 where γ 26.75 × 10 7 rad/T•s is the gyromagnetic ratio of hydrogen.From the measured S and S I in 2.1 and 2.3 , by applying an inverse Fourier transform, we obtain ρ 0 in 2.2 and ρ I x, y ρ x, y e iδ x,y e −iγT c B z x,y .

2.4
Then, the magnetic flux density B z is precisely computed as where α and β are the imaginary and real parts of ρ I /ρ 0 , respectively.

ICNE B z Acquisition
In the ICNE B z acquisition, in order to improve the SNR of B z , we prolong the current injection time T c until the end of the MR data acquisition, as shown in Figure 1.Then, since the induced B z disturbs the linearity of the reading gradient, the measured k-space data has lost the linear encoding characteristic as where G is a constant that denotes the strength of the magnetic reading gradient.The inverse problem arising in the ICNE method is to recover B z x, y from ρ 0 x, y obtained in 2.2 and the measured signal S C k x , k y in 2.6 .Although the inversion is not uniquely solvable for B z x, y in general, we can uniquely determine B z x, y by assuming that ϕ x : x B z x, y /G is monotone increasing in the following theorem.
Theorem 2.1.Let ρ x, y have a finite support Ω.If B z x, y is sufficiently small to guarantee that ϕ x x B z x, y /G is monotone increasing so that ϕ x > 0 for each y, then B z x, y is uniquely recovered in Ω from ρ 0 x, y in 2.2 and S C k x , k y in 2.6 .
Proof.We note that the linear encoding characteristic in the k y -variable remains unperturbed in 2.6 .Thus, by one-dimensional inverse Fourier transform, S C in 2.6 is reduced to S C in the k x , y -hybrid space as the following: , y e iδ x,y e −iγT c B z x,y e −i2π k x x B z x,y /G dx.

2.7
Then, the ICNE inverse problem suffices to consider the x-directional inversion of B z x, y from ρ 0 x, y in 2.2 and S C k x , y in 2.7 for each fixed y.By change of variables with ϕ x , 2.7 is changed into ,y e −iγT c B z x,y e −i2πk x ϕ dϕ.

2.8
From 2.8 , by inverse Fourier transform for the ϕ-variable, ϕ satisfies Φ ϕ ρ x, y ϕ x e iδ x,y e −iγT c B z x,y , 2.9 where Φ is a function defined with the given S C by The relation 2.9 gives us the simple ordinary differential equation as follows: Φ ϕ x ϕ x ρ x, y .

2.13
By the same argument, the reverse inequality is not possible.Thus, we have φ x 0 β.

2.14
By separation of variables, 2.11 and 2.14 lead us to For any given x, ϕ x is uniquely determined from 2.15 .It completes the proof.
Remark 2.2.In Theorem 2.1, we assume that ϕ x 1 1/G ∂B z x, y /∂x > 0. The magnetic flux density B z is smooth and its intensity is 10 −7 ∼ 10 −8 T in practical experimental environments.Furthermore, the usual range of the reading gradient G is 10 −3 ∼ 10 −4 T/m.Thus, the assumption of ϕ x > 0 is not severe in Theorem 2.1.

Discrete ICNE Inverse Problem
In a practical MRI scanner, the MR k-space data in 2.1 , 2.
In the rest of the paper, we assume that N is even and k, n 0, 1, . . ., N − 1 denote the row and column numbers, respectively.A matrix whose k, n entry is M kn is represented by M kn . 3.4 For a vector x x 0 , x 1 , . . ., x N−1 t , V x e −i 2π/N n x n k denotes a Vandermonde matrix as

Newton Iterations
Define a function F F 0 , F 1 , . . ., F N−1 t by The discrete ICNE inverse problem is to find the zero of F for s k given in 3.3 .
The Jacobian matrix DF x is composed of four parts as where D K , D x j , D ρ are diagonal matrices such that

3.8
Newton iterations to find the zero of F are as the following: x j 1 x j − real DF x j −1 F x j , 3.9 with an initial x 0 and the iterates The previous method in 14, 15 was based on the Taylor approximation, but as a coincidental result, it can be interpreted as the first Newton iterate x 1 in 3.9 with x 0 0. The following theorem suggests a condition for b in which the Newton iterations in 3.9 converge to b with the trivial initial guess 0 0, 0, . . ., 0 t .The proof is based on the Theorem 6.14 in 17 , which states a sufficient condition for the convergence of Newton iterations that

Convergence of Newton Iterations
where α DF 0 −1 F 0 ∞ and β, γ are the respective bounds of

3.13
Then, starting with x 0 0, the Newton iterations in 3.9 are well defined and converge to b.One also has the following quadratic error estimate: , j 0, 1, 2, . . . .

3.14
As a consequence, the zero of F satisfying 3.13 is unique.

3.16
By 3.7 , 3.16 implies that For two constants in 3.17 and 3.21 , let From 3.15 and the condition 3.13 , we have q ≤ 1 2 .

3.19
Then, with the aid of the Theorem 6.14 in 17 , F has a unique zero x in the ball and the Newton iterations in 3.9 converge to x with x 0 0. Since b is a zero of F contained in B from 3.15 , we have x b.The quadratic error estimate in 3.14 also comes from the same theorem in 17 .
Proof.From 3.7 , we expand
Lemma 3.3.For the trivial initial x 0 0, one has

3.26
Thus, we expand where

3.28
From Lemma 3.4, we estimate that

3.29
Since V 0 is the matrix of the discrete Fourier transform, we have

3.32
Proof.If |θ| ≤ 3, we have 3.32 since the series in the following expansion is alternating.

Norm of Inverse of Vandermonde Matrix
The norms of inverses of Vandermonde matrices were estimated by Gautschi 18,19 .In some estimations there, the equality holds if all base points are on the same ray through the origin.Compared to 3.31 , for a small perturbation x, a bound of V x −1 ∞ is investigated in the following lemma.Since the norm estimation of inverse of Vandermonde matrix must be interesting, we separate the result in this subsection from other ingredients for Theorem 3.1. 3.35 Proof.Let I δ kn be the identity matrix and W V 0 −1 V x − I, whose entries are

3.37
The diagonal entries are estimated by the Cauchy mean value theorem as follows:

3.38
Regarding summations of |W kn |, the maximum occurs when k N/2 ± 1 from the symmetry of the sine function in 3.37 .In both cases, we have

Numerical Results
From the Biot-Savart law, we simulate the magnetic flux density B z n, m induced by a horizontal current through the Logan shape ρ as in Figure 2, where N 60.We obtain the simulated MR signals S C d k, l as depicted in Figure 3, through 3.1 with M 32, δ 0.  The log of error max n,m |B z n, m − B z n, m j | is given in Figure 4 a , which means that the error decay is quadratic.In the Newton iterations in 3.9 , we have to solve a Vandermonde system for V x , which may be consuming time or unstable.Instead of V x in 3.7 , we can fix V 0 and simplify the Newton iterations in 3.9 .Then, the error decay is reduced to be linear as in Figure 4 b .

Figure 1 :
Figure 1: Conventional and ICNE current injections in a spin echo pulse sequence.

2 For
3 , and 2.6 are acquired by finite sampling with a dwell time dt.If N is the reading time T s divided by dt, we have the following N × N discrete signals with dimensionless variables instead of those in 2.6 : m e iδ n,m e −i 2π/N MB z n,m e −i 2π/N k n B z n,m lm , 3.1 for a constant M > 0. The discrete ICNE inverse problem is to recover B z n, m from 3.1 with known a priori ρ n, m e iδ n,m and the measured signal S C d k, l , where n, m, k, l 0, 1, . . ., N − 1.By discrete inverse Fourier transform for l, 3.1 can be suppressed into m e iδ n,m e −i 2π/N M k B z n,m e −i 2π/N kn .3.each fixed m, let s k , ρ n , and b n denote S C d k, m , ρ n, m e iδ n,m , and B z n, m , respectively.Then, the discrete ICNE inversion problem is a system of N nonlinear equations for N unknowns, b 0 , b 1 , . . ., b N−1 such that

Let b b 0 the
, b 1 , . . ., b N−1 t , ρ ρ 0 , ρ 1 , . . ., ρ N−1 t , and ρ max max n |ρ n |, ρ min min n |ρ n |.If Jacobian DF b in 3.7 is invertible, since the Vandermonde matrix V b in 3.5 is based on the N distinct points.Thus, the Newton iterations in 3.9 converge to b for an initial x 0 which is sufficiently close to b 17 .

b
Logan shape ρ

Figure 4 :
Figure 4: log max n,m |B z n, m − B z n, m j | error decay.