Open Access

Optimization of magnetic flux density for fast MREIT conductivity imaging using multi-echo interleaved partial fourier acquisitions

  • Munish Chauhan1,
  • Woo Chul Jeong1,
  • Hyung Joong Kim1Email author,
  • Oh In Kwon2 and
  • Eung Je Woo1
BioMedical Engineering OnLine201312:82

https://doi.org/10.1186/1475-925X-12-82

Received: 19 June 2013

Accepted: 21 August 2013

Published: 27 August 2013

Abstract

Background

Magnetic resonance electrical impedance tomography (MREIT) has been introduced as a non-invasive method for visualizing the internal conductivity and/or current density of an electrically conductive object by externally injected currents. The injected current through a pair of surface electrodes induces a magnetic flux density distribution inside the imaging object, which results in additional magnetic flux density. To measure the magnetic flux density signal in MREIT, the phase difference approach in an interleaved encoding scheme cancels out the systematic artifacts accumulated in phase signals and also reduces the random noise effect by doubling the measured magnetic flux density signal. For practical applications of in vivo MREIT, it is essential to reduce the scan duration maintaining spatial-resolution and sufficient contrast. In this paper, we optimize the magnetic flux density by using a fast gradient multi-echo MR pulse sequence. To recover the one component of magnetic flux density B z , we use a coupled partial Fourier acquisitions in the interleaved sense.

Methods

To prove the proposed algorithm, we performed numerical simulations using a two-dimensional finite-element model. For a real experiment, we designed a phantom filled with a calibrated saline solution and located a rubber balloon inside the phantom. The rubber balloon was inflated by injecting the same saline solution during the MREIT imaging. We used the multi-echo fast low angle shot (FLASH) MR pulse sequence for MRI scan, which allows the reduction of measuring time without a substantial loss in image quality.

Results

Under the assumption of a priori phase artifact map from a reference scan, we rigorously investigated the convergence ratio of the proposed method, which was closely related with the number of measured phase encode set and the frequency range of the background field inhomogeneity. In the phantom experiment with a partial Fourier acquisition, the total scan time was less than 6 seconds to measure the magnetic flux density B z data with 128×128 spacial matrix size, where it required 10.24 seconds to fill the complete k-space region.

Conclusion

Numerical simulation and experimental results demonstrated that the proposed method reduces the scanning time and provides the recovered B z data comparable to what we obtained by measuring complete k-space data.

Keywords

MREIT MRI Interleaved partial fourier acquisition Magnetic flux density Current density

Background

Magnetic resonance electrical impedance tomography (MREIT) utilizes a magnetic resonance imaging (MRI) scanner to measure magnetic flux density B z data inside an imaging object induced by the externally injected current. The internal current density distribution has been studied in magnetic resonance current density imaging (MRCDI) by measuring the whole magnetic flux density data B=(B x ,B y ,B z ) [1, 2]. Combining MRCDI and electrical impedance tomography (EIT) technique, MREIT provides the cross-sectional conductivity images of the object with high spatial resolution [310]. Since an MRI scanner measures only one component B z of B without rotating the imaging object, most MREIT algorithms assumed that the internal conductivity is isotropic and focused on visualizing its distribution by using one component of the magnetic flux density data B z of B[1118].

Recent MREIT imaging techniques have been developed with respect to both the capacity of measurement techniques and the numerical reconstruction algorithms. Experimental results from in vivo animal and human have been reported [19, 20] in MREIT. As an innovation of current MREIT, a fast MREIT imaging technique referring to the continuous monitoring of objects includes various wide application areas [21]. Recently, one of challenging problem in MREIT is to implement a new imaging technique with a very short acquisition time for the imaging of neural activities of brain related to the conductivity change.

Since current MREIT experiments suffer from poor SNR of the measured B z data under the typical data acquisition durations and a small amount of injected current, it is important to reduce the scan time, while maintaining the spatial-resolution and sufficient contrast, for practical implementations of in vivo MREIT. Recently, to reduce the scan time in MREIT, Hamamura et al[22] reconstructed the interior conductivity using a single-shot spin-echo echo planar imaging (SS-SEPI) pulse sequence and Muftuler et al[23] used a SENSE-accelerated imaging technique to acquire phase signal by the injected current. One of basic approaches for maintaining the spatial resolution is to reduce the number of phase encoding steps because each phase encoding step requires a certain amount of time for execution. Since the MREIT techniques use an interleaved phase encoding acquisition scheme to double the B z signal, Park et al[24] reconstructed the phase signal B z by filling the skipped k-space region using the interleaved measurement property.

To obtain the static conductivity image in MREIT [19, 20], a spin-echo based MREIT pulse sequence has been predominantly used to reduce the background artifact and to increase the imaging quality. In real situations, it is difficult to employ the fast conventional MR pulse sequences because the noise standard deviation of B z is inversely proportional to the width of injection current and the intensity of MR magnitude, simultaneously. An MREIT pulse sequence should be devised to enhance changes in MR phase images for given current amplitudes.

In this paper, we used the multi-echo fast low angle shot (FLASH) MR pulse sequence which allows the reduction of imaging time without any substantial loss in image quality. In addition, the multi-echo FLASH sequence maximizes the width of injection current extending the duration of injection current until the end of a readout gradient in MREIT. To reconstruct the internal conductivity distribution, most algorithms require at least two independent injection currents in an interleaved sense, which require relatively a long scanning duration [7, 11]. To reduce the scanning time, we adopt a partial phase encoding acquisition scheme using the multi-echo FLASH MREIT pulse sequence and rigorously investigate the relationship between the convergence ratio of the algorithm and the background field inhomogeneity [24]. We consider the discrete 2-norm to evaluate the convergence ratio.

To show the feasibility of the proposed algorithm, we performed numerical simulations and compared the performance to the simulated true B z data. We designed a cylindrical acrylic phantom filled with a calibrated saline solution and located a rubber balloon inside the phantom. The rubber balloon was inflated by injecting the same saline solution during the scan. The phantom was designed to provide a homogeneous magnitude image, but distinguishable signals of measured B z between the inside and outside the balloon. The phantom experiment demonstrated that the proposed method reduces the scanning time and recovers the reasonable resolution of B z , which is comparable to the recovered B z using the complete k-space data.

Method

k-space signal and B z data

In a conventional spin echo MREIT pulse sequence, both positive and negative currents of the same amplitude and duration are injected with reverse polarity. These injection currents with the pulse width of T c accumulate extra phases. Corresponding k-space MR signals can be described as
S ± ( k x , k y ) = Ω ρ ( x , y ) e ( x , y ) e ± B z ( x , y ) T c e i 2 π ( k x x + k y y ) dxdy
(1)
where ρ is the T2 weighted spin density, δ is any systematic phase artifact, and Ω is a field-of-view (FOV). Here, the superscript of S±(k x ,k y ) denotes a brief notation for S+(k x ,k y ) and S(k x ,k y ). For the standard coverage of k-space, we set
k x = γ 2 π G x ( nΔt T E ) for n = N x / 2 , , N x / 2 1 k y = γ 2 π G y T pe for m = N y / 2 , , N y / 2 1
(2)
where γ=26.75×107rad/T·s is the gyromagnetic ratio of hydrogen, Δ t is the time between samplings, G x is the frequency encoding gradient strength, T E is the echo time, Δ G y is the phase encoding step, and T pe is the phase encoding time. The induced magnetic flux density ±B z is generated by the positive and negative injection currents I±. Applying the inverse Fourier transform to the measured k-space data sets in (1), we can compute the magnetic flux density B z as
B z ( r ) = 1 2 γ T c tan 1 α ( r ) β ( r )
(3)

where α and β are the imaginary and real part of ρ e e B z T c / ρ e e B z T c , respectively [2].

The current MREIT method is based on the electromagnetic information embedded in the measured B z data in order to visualize the conductivity(or current density) based on Biot-Savart law:
B z ( r ) = μ 0 4 π Ω ( y y ) J x ( r ) ( x x ) J y ( r ) | r r | 3 d r , r = ( x , y , z ) , r = ( x , y , z )
(4)
where μ0=4π 10−7Tm/A is the magnetic permeability of the free space. The current density J in Ω is given for the isotropic conductivity σ in Ω
J ( r ) = σ u ( r )
(5)
and satisfies the following elliptic equation
· σ u = 0 in Ω σ u · ν = g on ∂Ω and ∂Ω u ds = 0
(6)

where ν is the outward unit normal vector on Ω and g is an applied current density on the surface.

Recovery of complex T 2 weighted spin density using interleaved partial Fourier acquisition

To simplify, we develop a theory using a conventional cartesian k-space that has three zones within the k y (phase-encode) domain; the central region ( P 0 ), the positive region ( P + ) and the negative region ( P ):
P 0 = { ( k x , k y ) k y = γ 2 π G y T pe , N c m N c } P + = { ( k x , k y ) k y = γ 2 π G y T pe , m > N c } P = { ( k x , k y ) k y = γ 2 π G y T pe , m < N c }
(7)
Here, N c denotes the number of partial Fourier over-sampling phase-encodes. The interleaved k-space data S p + ( k x , k y ) and S p - ( k x , k y ) as partially acquired for the phase k y = γ 2 π G y T pe can be expressed
S p ± ( k x , k y ) = S ± ( k x , k y ) , ( k x , k y ) P + P 0 0 , ( k x , k y ) P
(8)
where S p ± represents S p + and S p simultaneously. We propose an algorithm to determine the T 2 weighted spin density ρ :
ρ ± ( r ) : = ρ ( r ) e ( r ) e ± B z ( r ) T c = F T 1 ( S ± ( k x , k y ) ) ( r )
(9)

using the partially scanned k-space data S p ± .

The systematic phase artifact ei δ(x,y), unavoidable artifacts due to the main field inhomogeneity and the mismatch between the center of data acquisition interval and echo formation, arises from a low frequency field, which mainly belongs to the central region P 0 . Including a small perturbed phase artifact ei δ(x,y), we start with the initial guess S p ± and design an alternating procedure by updating the skipped k-space regions.

The recovered ρ p ± = F T 1 ( S p ± ) can be formally expressed by the equation
ρ p ± ( r ) = F T 1 ( S ± ( k x , k y ) ) ( r ) + F T P 1 ( S p ± ( k x , k y ) ) ( r ) F T P 1 ( S ± ( k x , k y ) ) ( r ) = ρ ± ( r ) + E p ± ( r )
(10)
where F T P 1 ( S p ± ( k x , k y ) ) denotes the inverse Fourier transform by zero-filling the k-space except in the region P and the remainder term E p ± is
E p ± ( r ) : = F T P 1 S p ± ( k x , k y ) S ± ( k x , k y ) ( r ) .
(11)
Let us define a support region D δ for a low spatially varying magnetic field due to background field inhomogeneities:
D δ : = { ( k x , k y ) | FT ( e ) ( k x , k y ) = 0 , k y = γ 2 π G y T pe , N δ m N δ }
(12)

Observation 1. If D δ P 0 , FT ( ρ p e 2 ) ( k x , k y ) = FT ( ρ e 2 ) ( k x , k y ) for ( k x , k y ) P + .

The proof of observation 1 is provided in the Appendix A. By using the observation 1, we fill the skipped k-space regions in S p ±
S u ± ( k x , k y ) = S ± ( k x , k y ) , ( k x , k y ) P + P 0 FT ( ρ p e 2 ) ¯ ( k x , k y ) , ( k x , k y ) P
(13)
Observation 2. If D δ P 0 , the k-space S u ± in recovers the T2 (or T 2 ) weighted spin density ρ± without loss of information.
F T 1 ( S u ± ) = ρ ± = ρ e e ± B z ( x , y ) T c

The proof of observation 2 is provided in the Appendix 2. The observation 2 shows that the skipped region in the measured k-space S+ can be recovered by using the interleaved acquired S and the estimated background sensitivity map.

Convergence characteristics

When the common measured P 0 does not cover the support region D δ , N δ >N c , for the phase-encode ( k x , k y ) P + , the Fourier transform of ρ p e 2 can be written by following the observation 1:
FT ( ρ p e 2 ) ( k x , k y ) = ( k x , k m ) P + P 0 FT ( ρ ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m )
(14)
From the relation (14), for the phase-encode ( k x , k y ) P + , we have
FT ( ρ e 2 ) ( k x , k y ) FT ( ρ p e 2 ) ( k x , k y ) = ( k x , k m ) P FT ( ρ ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m )
(15)
We set
[ ψ ] R : = F T 1 FT ( ψ ) | R
(16)

where R is a subregion of the k-space and FT ( ψ ) | R denotes the restriction of F T(ψ) to the region R . The discrete 2-norm of [ ψ ] R is equivalent to that of FT ( ψ ) | R , i.e., [ ψ ] R 2 = C FT ( ψ ) | R 2 , where the constant C is independent to the function ψ and the region R .

When the skipped k-space region, P , are filled the previously updated as in (13), we have
( ρ ρ n ) 2 = C FT ( ρ ρ n ) 2 = C FT ( ρ ρ n ) | P 2 = C FT ( ρ + e 2 ) ¯ FT ( ρ n 1 + e 2 ) ¯ | P ( k x , k y ) 2 = C FT ( ρ + e 2 ) FT ( ρ n 1 + e 2 ) | P + ( k x , k y ) 2 = C FT ( ( ρ + ρ n 1 + ) e 2 ) | P + ( k x , k y ) 2
(17)
For the k-space region ( k x , k y ) P + , we have the following identity
FT ( ρ + ρ n 1 + ) e 2 ) | P + ( k x , k y ) = ( k x , k m ) P FT ( ρ + ρ n 1 + ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) + Θ ( ρ + ρ n 1 + , P 0 P + )
(18)
where
Θ ( ρ + ρ n 1 + , P 0 P + ) : = ( k x , k m ) P 0 P + FT ( ρ + ρ n 1 + ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) .
Since the updated complex density FT ( ρ n ± ) | P + P 0 = FT ( ρ ± ) | P + P 0 , the remainder term Θ ( ρ + ρ n 1 + , P 0 P + ) = 0 . Thus, the discrete 2-norm of the difference between the true and the iteratively updated T 2 weighted spin density can be estimated
FT ( ( ρ + ρ n 1 + ) e 2 ) | P + 2 2 FT ( ( ρ + ρ n 1 + ) | P 2 2 k y P + ( k x , k m ) P | FT ( e 2 ) ( k x , k y k m ) | 2
(19)

Detailed estimates of 2-norm calculation are presented in the Appendix C Estimation of 2 -norm.

Define an estimator for the convergence of the proposed algorithm
Z 2 δ , P ± : = k y P + ( k x , k m ) P FT ( e 2 ) ( k x , k y k m ) 2
(20)
Using the same procedures of (17) and (19), we have
FT ( ( ρ + ρ n 1 + ) | P 2 2 = FT ( ( ρ ρ n 2 ) e 2 ) | P + ( k x , k y ) 2 Z 2 δ , P ± FT ( ( ρ ρ n 2 ) | P 2 2
(21)
The relations (17), (19) and (21) show that the convergence of the proposed method depends on the estimator Z 2 δ , P ± :
( ρ ρ n ) 2 Z 2 δ , P ± ( ρ ρ n 2 ) 2
(22)

Optimization of B z using gradient multi-echo data

Since the noise standard deviation s B z of the measured B z is inversely proportional to injection current duration T c and the SNR of MR magnitude image Υ M [25, 26] as
s B z ( r ) = 1 2 T c Υ M ( r )
(23)

To reduce the noise level of B z and the imaging time, we applied the proposed method to the gradient multi-echo pulse sequence as a fast MR imaging technique. Subsequently, by using the gradient multi-echo, it is possible to inject the current for a long duration to maximize the off-resonance phase.

When T E 1 denotes the first echo time and Δ T E is the echo spacing, the m-th echo time is T E m = T E 1 + ( m 1 ) Δ T E , m = 1 , , N E , where N E is the echo number. The phase artifact δ T E m depends on the echo time T E m . Using a priori estimation of the phase artifact map δ T E m from a reference scan, for a time varying functional MREIT technique, the measured k-space region including P 0 and P + can be determined by taking account of the imaging multi-echo times T E m , the echo number N E and the repetition time T R .

The recovered multiple T 2 -weighted complex densities ρ T E m ± , m = 1 , , N E , using the proposed algorithm in each echo time T E m can be optimized to generate a representative measured B z m data
B z m ( r ) = 1 2 γ T c tan 1 α m ( r ) β m ( r )
(24)
where α m and β m are the imaginary and real parts of ρ T E m + / ρ T E m , respectively. The recovered multiple B z m data include pixel-by-pixel different noise level depending on the different imaging time T E m and the width of injection current. The multiple measured B z m data are optimally combined to reduce the noise level of B z [27]:
B z ( r ) = m = 1 N E ω m ( r ) B z m ( r )
(25)
where the point-wise weighting factor ω m (r) is given as
ω m ( r ) = 1 / ( s B z m ( r ) ) 2 k = 1 N E 1 / ( s B z k ( r ) ) 2
(26)
Now, we setup an algorithm to reconstruct the magnetic flux density B z data using the partially measured k-space data S p ± and involving the following steps:
  1. 1.

    Take an initial guess ρ 0 ± = F T 1 ( S ± | P 0 P + ) .

     
  2. 2.

    Transform the n-th updated ρ n ± e 2 to k-space by taking the Fourier transform.

     
  3. 3.

    Update the measured k-space data S n + 1 ± by filling the skipped region in S p ± with the transformed FT ( ρ n ± e 2 ) data.

     
  4. 4.

    Update ρ n + 1 ± by taking the two-dimensional inverse Fourier transform.

     
  5. 5.

    Stop if ρ n + 1 ± ρ n ± Ω ρ n ± Ω ε where ε>0 is a given tolerance and · Ω is a standard L 2-norm in Ω. Otherwise, repeat the process.

     
  6. 6.

    Reconstruct B z m using (24) for m=1,2,,N E .

     
  7. 7.

    Determine a weighting factor map ω m , m=1,2,,N E using (26) and reconstruct an optimally weighted B z in (25).

     

Experimental setup

Numerical simulation setup

To validate the proposed algorithm, we performed numerical simulations with the two-dimensional finite-element model of a object 20×20 cm2 with 256×256 rectangular elements and with the origin at its bottom-left, shown in Figure 1. We added the different complex field inhomogeneity artifacts to the simulated spin density image in Figure 1(a). The target magnetic flux density B z in Figure 1(b) was generated by solving the elliptic equation (6) and by using the Biot-Savart law given by (4). Since the magnetic flux density B z in Figure 1(b) is continuous and has no abrupt changes, we used |B z | image displayed in Figure 1(c) to enhance image for the magnetic flux density. Figure 1(d) shows a partially measured k-space data corresponding to ρ+.
Figure 1

Simulation setup. a) T2-weighted spin density, b) simulated magnetic flux density B z image, c) intensity of |B z image, d) magnitude image of the k-space data on the measured region.

The target conductivity distribution σ had different anomalies with different conductivity values and the amount of injection current was 10 mA. Set an applied current density g on the surface as
g ( x , y ) = 10 , if | y 5 | 2 , and x = 0 10 , if | y 5 | 2 , and x = 10 0 , otherwise.
(27)

Phantom imaging experimental setup

For the practical application of the proposed method as a fast MREIT imaging, we designed a cylindrical phantom filled with the saline solution of conductivity 1 S/m (shown in Figure 2(a)), including a rubber balloon for the visualization of isotropic conductivity excluding other artifacts by any concentration gradient in the phantom. The inside of balloon was filled with the same saline solution and the volume of balloon was controlled by injecting saline solution during the imaging experiment. After positioning the phantom inside the bore of 3T MR scanner (Achieva TX, Philips Medical Systems, Best, The Netherlands) with 8 channel RF coil, we collected k-space data using the gradient multi-echo injection current nonlinear encoding (ICNE) pulse sequence which was originated from FALSH (Figure 2(b)). To obtain the MR magnitude and magnetic flux density (B z ) images, it extends throughout the duration of injection current until the end of a readout gradient [28]. Since the multi-echo ICNE pulse sequence was synchronized to the injection currents with alternating polarity, it enabled to maximize the width of the injection currents and minimize the noise standard deviation of the measured B z data. The maximum amplitude of injection current was 5 mA and the total imaging time was 10.24 second to fill the k-space in the interleaved sense. Since the total imaging time was corresponding to the whole k-space scan, the actual imaging time would be reduced to 5.52 and 5.92 seconds for the number of partial region N ( P 0 P + ) = 69 and 74, respectively. The imaging parameters were followings: slice thickness 5 mm, number of imaging slices one, repetition time T R =40 ms, echo spacing Δ T E =6 ms, flip angle 40 degree, and multi-echo time T E m = 6 + ( m 1 ) × 6 ms for N E =4. The FOV was 160 ×160 mm2 with a matrix size of 128×128. The duration of current injection T c m was almost same to the multi-echo time T E m = 6 + ( m 1 ) × 6 , m = 1 , 2 , 3 , 4 , because the current was continuously injected until the end of the readout gradient.
Figure 2

Experimental setup. a) saline phantom with balloon, b) diagram of the ICNE-multi-echo MR pulse sequence based on a gradient echo.

Figure 3(a) shows the multiply acquired magnitude images | ρ T E m + | , m = 1 , , 4 , where ρ T E m + was the m-th measured T 2 weighted complex spin density, Figure 3(c) and (e) show the measured magnetic flux densities B z m and the absolute of B z m images at each echo m=1,,4, respectively. The slope of B z m reflecting the width of injected current linearly increased as the multi-echo time T E m was increasing. Figure 3(b), (d), and (f) are the averaged images corresponding to Figure 3(a), (c), and (e), respectively.
Figure 3

Phantom imaging. a) T 2 -weighted magnitude images | ρ T E m + | , m = 1 , , 4 , b) averaged T 2 weighted magnitude image | ρ ̄ + | : = 1 4 m = 1 4 | ρ T E m + | , c) measured B z m images at each echo m=1,,4, d) averaged B z image, B z : = 1 4 m = 1 4 B z m , e) | B z m | images at each echo m=1,,4, f) averaged |B z | image, | B z | : = 1 4 m = 1 4 | B z m | .

Results

Simulation results

We compared the reconstructed magnetic flux density B z using the complete k-space data to the B z achieved using the partial k-space data. To evaluate the convergence characteristics of the proposed algorithm, we define the relative 2-errors:
E ( ψ n ) : = ψ ψ n ψ
(28)

where ψ and ψ n are the recovered images with the complete k-space data and the n-th updated data using the partially measured k-space data, respectively, and · denotes the 2-norm.

To investigate the estimator Z 2 δ , P ± for the convergence of the proposed iterative algorithm to fill the sipped k-space region P , we fixed the number of partially measured k-space region as N ( P 0 P + ) = 138 , i.e., N ( P 0 ) = 20 and N ( P + ) = 118 , and changed the frequency range of background field inhomogeneity. We generated several background field inhomogeneities changing the phase frequency range in k-space region by taking N δ = 5 , 10 , 20 , 30 where the background field inhomogeneity δ satisfies F T(e i δ )(k x ,k y )=0 for | k y | > N δ .

Figure 4(a)-(d) show the background field inhomogeneities used in the reconstruction procedure and Table 1 shows the estimated Z 2 δ , P ± in each background field inhomogeneity for N δ = 5 , 10 , 20 , and 30, respectively.
Figure 4

Simulated background field inhomogeneity distributions. a-d) real part images of e i δ for N δ = 5 , 10 , 20 , 30 .

Table 1

Calculated estimator Z 2 δ , P ± in (20) for N δ = 5 , 10 , 20 , 30 and fixed measuredk-space data N ( P 0 P + ) = 138

 

N δ = 5

N δ = 10

N δ = 20

N δ = 30

Z 2 δ , P ±

0.0013

0.0749

0.4927

0.7738

Figure 5(a) shows the reconstructed |B z | images using the fixed background field inhomogeneity with N δ = 5 . From the top-left to the bottom-right, each image was corresponding to the j-th iterative updated |B z |, j=0,1,,10. The recovered B z using the partially measured k-space data included a large amount of artifacts (the top-left image in Figure 5(a)). However, since the value of Z 2 δ , P ± was small, the first updated magnetic flux density almost recovered the true B z .
Figure 5

Reconstructed|B z | images using the fixed background field inhomogeneity with N δ = 5 and N δ = 30 .a) reconstructed |B z | images using the fixed background field inhomogeneity with N δ = 5 . b) reconstructed |B z | images using the fixed background field inhomogeneity with N δ = 30 . From top-left to bottom-right, each image corresponds to the j-th iterative updated |B z |, j=0,1,,9.

Figure 5(b) shows reconstructed |B z | images using the fixed background field inhomogeneity with N δ = 30 corresponding to Figure 5(a). Since the value of Z 2 δ , P ± was 0.7738, the convergence ratio of ρ± was relatively slow comparing to the field inhomogeneity with N δ = 5 . Table 2 shows the relative discrete 2-errors of updated complex spin density for each iteration number and the convergence ratio of ρ± was depending on the number of N δ .
Table 2

Relative 2 -errors of updated complex spin density for each iteration number()

 

0

1

2

4

6

8

10

E ( ρ N δ = 5 )

0.06291

0.00098

0.00044

0.00022

0.00021

0.00021

0.00021

E ( ρ N δ = 10 )

0.08992

0.00282

0.00104

0.00039

0.00024

0.00020

0.00020

E ( ρ N δ = 20 )

0.20801

0.03346

0.00706

0.00088

0.00043

0.00029

0.00023

E ( ρ N δ = 30 )

0.22640

0.04623

0.01451

0.00309

0.00090

0.00041

0.00028

Table 3 shows the relative 2-errors of reconstructed B z for each iteration number depending on the number of N δ . The decay rates of the relative 2-error were very fast as the number of N δ was small, but we needed relatively many iterations to approach the required accuracy as the number of N δ was increase, even though the update procedure was rapidly computed by use of the fast Fourier transform.
Table 3

Relative discrete 2 -errors ofB z for each iteration number()

 

0

1

2

4

6

8

10

E ( B z N δ = 5 )

0.8755

0.0753

0.0462

0.0265

0.0214

0.0202

0.0199

E ( B z N δ = 10 )

1.0945

0.2189

0.0949

0.0552

0.0373

0.0288

0.0252

E ( B z N δ = 20 )

1.3132

0.6737

0.3039

0.1126

0.0730

0.0548

0.0435

E ( B z N δ = 30 )

1.3276

0.8300

0.4764

0.2242

0.1245

0.0823

0.0626

Phantom experimental results

For the phantom experiment, we changed N c =5,,10 for the set P 0 to investigate the convergence behavior with respect to a given background field inhomogeneity. Using the collected k-space data with 8 channel RF coil and the gradient multi-echo by alternating readout gradient, we measured the T 2 weighted complex densities ρ m ± n , n = 1 , , N CH , m = 1 , , N E , where N C H =8 denotes the coil number and N E =4 is the echo number. Figure 6 shows the measured background field inhomogeneities by displaying the real part of e 2 i δ mn corresponding to the n-th coil and the m-th echo image. According to the increase of echo number, the accumulated background field inhomogeneity also increased.
Figure 6

Measured background field inhomogeneity distributions. Real part of e 2 i δ mn , n = 1 , , N C H , m=1,,N E , where N C H =8 and N E =4 denote the coil and echo numbers, respectively.

Figure 7(a) and (c) show the measured T 2 weighted magnitude and magnetic flux density B z images at the time T E m , m = 1 , 2 , 3 , 4 , using partially acquired k-space region P 0 P + with N c =5 by a transversally injected current. Although the amount of accumulated phase signal by the injected current increased as the echo time varied from T E 1 to T E 4 , the magnitude image at the 4-th echo was more deteriorated comparing to the 1-st echo case. Figure 7(b) shows an averaged MR magnitude image at each echo time and Figure 7(d) is a weighted B z image depending on the width of injected current using the phase signal in Figure 7(c). Figure 7(e)-(h) shows the measured magnitude and magnetic flux density B z images corresponding to Figure 7(a)-(d) using partially acquired k-space region P 0 P + with N c =10.
Figure 7

Measured T 2 weighted magnitude and magnetic flux densityB z images using P 0 P + withN c =5 and N c =10.a) and e) T 2 weighted magnitude image at each echo time T E m , m = 1 , 2 , 3 , 4 , with N c =5 and N c =10, respectively. b) and f) combined T 2 weighted magnitude image with N c =5 and N c =10, respectively. c) and g) recovered B z image at each echo time T E m , m = 1 , 2 , 3 , 4 with N c =5 and N c =10, respectively. d) and h) weighted B z image using multiple B z image at each echo time with N c =5 and N c =10, respectively.

Comparing to the measured images in Figure 7(a)-(d), in contrast to the recovery of low phase frequency information corresponding to P 0 , the increased background field inhomogeneity caused relatively high frequency artifacts.

Figure 8(a)-(d) shows iteratively updated T 2 weighted magnitude and magnetic flux density B z images using P 0 P + with N c =5. We fixed the update iteration number as 20 for all experiments. When we fixed N c =5, the 1-st and 2-nd recovered T 2 weighted complex densities in Figure 8(a) and (c) were relatively close to the recovered ones using the complete k-space data. However, as the phase artifact increased, the 3-rd and 4-th recovered T 2 weighted complex densities were deficient in reflecting full information of B z signal. Especially, the 4-th updated magnetic flux density B z image shows some defective region due to the insufficient recovery of T 2 weighted complex density.
Figure 8

Iteratively updated T 2 weighted magnitude and magnetic flux densityB z images using P 0 P + withN c =5 andN c =10.a) and e) recovered T 2 weighted magnitude image at each echo time T E m , m = 1 , 2 , 3 , 4 , with N c =5 and N c =10, respectively. b) and f) combined T 2 weighted magnitude image using the recovered magnitude image at each echo time with N c =5 and N c =10, respectively. c) and g) recovered B z image at each echo time T E m , m = 1 , 2 , 3 , 4 ,with N c =5 and N c =10, respectively. d) and h) weighted B z image using multiple B z image at each echo time with N c =5 and N c =10, respectively.

Figure 8(e)-(f) shows iteratively updated T 2 weighted magnitude and magnetic flux density B z images corresponding to Figure 8(a)-(d) using P 0 P + with N c =10. When we used N c =10, the updated T 2 weighted complex densities almost recovered the magnetic flux density B z data comparing to those using the complete k-space data.

Discussion

We used the gradient multi-echo MREIT pulse sequence to reduce the imaging time and to maximize injection current duration. Since the MREIT techniques utilize accumulated phase signal by the injected current, it requires enough repetition time T R to accumulate the phase signal. In this sense, the gradient multi-echo MREIT pulse sequence seems practical approach for the improvement of B z quality as well as reducing the imaging time. In this paper, we used a partially acquired k-space data in the phantom experiment by filling the k-space as much as 74 line by line, results in 5.92 second to image the resolution of 128×128. Experimental results show that the proposed interleaved partial Fourier strategy for MREIT has a potential to reduce scan times and maintain the information of B z data comparable to what is obtained with complete k-space data.

The convergence ratio of the iteratively updated phase signal heavily depends on the frequency of the background filed inhomogeneity and the number of half-Fourier over-sampling phase-encodes P 0 . Instead of the gradient multi-echo, if we use the spin multi-echo pulse sequence, the proposed iterative algorithm would rapidly recover T2-weighted complex spin density due to a small amount of background field inhomogeneity. However, in spite of some advantages of the spin multi-echo MREIT pulse sequence, for a real-time MREIT imaging, MR pulse sequence should be carefully investigate by taking into account of the width of injection current, the scan duration and the low SNR of measured Bz signal.

In this paper, we assumed a priori background field inhomogeneity which is typically used in the sensitivity encoding (SENSE) as a fast MRI measurement technique. Since the MREIT techniques typically used interleaved acquisition by injecting alternative currents, it may be possible to extract background field inhomogeneity information under a low frequency range assumption and by cancelation of Bz information:
ρ + ( x , y ) ρ ( x , y ) = ρ 2 ( x , y ) e 2 ( x , y ) ρ + ( x , y ) ρ ( x , y ) = e 2 B z ( x , y ) T c
(29)

Several studies reported for the feasibility of MREIT to detect neural activities in the brain, directly [29, 30]. Functional MREIT technique is suggested to image brain activity via conductivity change related to neural activity through the fast MREIT pulse sequence. Our future study will focus on applying the proposed method to produce functional conductivity images of animal and/or human brain to pursue rapidly changing conductivity associated with neural activities.

Conclusion

In MREIT, the inherent challenges are to reduce the scan time and maintain current injection duration to make it feasible for the clinical applications. We developed an iterative method to optimize the measured magnetic flux density Bz using the multi-echo interleaved partial Fourier acquisitions for fast imaging in MREIT. The proposed method used a fast gradient multi-echo MR pulse sequence to reduce the scan time and to maximize the phase signal by injection current. Under the assumption of a priori background field inhomogeneity map, we rigorously investigated the convergence ratio of the proposed method using the discrete 2-norm, which was closely related with the number of measured phase encode set and the frequency range of the background field inhomogeneity. To evaluate the proposed method, a specially designed conductivity phantom was used to provide a homogeneous magnitude, but it yielded distinguishable Bz signal inside and outside the anomaly. For the phantom experiment, total imaging time was 10.24 seconds to fill the complete k-space region in the interleaved sense and it was less than 6 seconds to fill the partial k-space region to implement the proposed method. The proposed interleaved partial Fourier strategy for the fast MREIT has a potential to reduce scan times and maintain the information of Bz data comparable to what is obtained with the complete k-space data.

Appendix

A Proof of Observation 1

For the phase-encode ( k x , k y ) P + , the Fourier transform of ρ p e 2 can be separated as
FT ( ρ p e 2 ) ( k x , k y ) = ( FT ( ρ p ) FT ( e 2 ) ) ( k x , k y ) = ( k x , k m ) P + P 0 P FT ( ρ p ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) = ( k x , k m ) P + P 0 FT ( ρ p ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) + ( k x , k m ) P FT ( ρ p ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m )
(30)
where denotes the convolution with respect to k y . Since the updated S p data conserve the measured data in P + P 0 , FT ( ρ p ) ( k x , k y ) = FT ( ρ ) ( k x , k y ) for ( k x , k y ) P + P 0 . Thus, we have
FT ( ρ p e 2 ) ( k x , k y ) = ( k x , k m ) P + P 0 FT ( ρ ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) + ( k x , k m ) P FT ( ρ p ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) = ( k x , k m ) P + P 0 FT ( ρ ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m )
(31)

Since the central phase-encode set P 0 includes all phase frequencies of the systematic phase artifact ei δ, the range of the phase frequency k y k m for ( k x , k y ) P + and ( k x , k m ) P is over 2Nδ. This means that F T(e−2i δ)(k x ,k y k m )=0. Thus, we have FT ( ρ p e 2 ) ( k x , k y ) = FT ( ρ e 2 ) ( k x , k y ) for ( k x , k y ) P + . The case for ρ p + e 2 is similar.

B Proof of Observation 2

Since FT ( ρ p e 2 ) ( k x , k y ) = FT ( ρ e 2 ) ( k x , k y ) for ( k x , k y ) P + due to the observation 1, we have
FT ( ρ + e 2 ) ( k x , k y ) = FT ( ρ e e B z T c ) ( k x , k y ) for ( k x , k y ) P +
(32)
From the relation (32), by taking the complex conjugate, we recover the skipped k-space region P
FT ( ρ ) ( k x , k y ) = FT ( ρ e e B z T c ) ( k x , k y ) = FT ( ρ e e B z T c ) ¯ ( k x , k y ) = FT ( ρ p + e 2 ) ¯ ( k x , k y )
(33)

C Estimation of 2 -norm

The discrete ℓ2-norm of the difference between the true and the iteratively updated T 2 weighted spin density can be estimated as following:
FT ( ( ρ + ρ n 1 + ) e 2 ) | P + 2 2 = k y P + | FT ( ρ + ρ n 1 + ) e 2 ( k x , k y ) | 2 = k y P + ( k x , k m ) P FT ( ρ + ρ n 1 + ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) 2 k y P + ( k x , k m ) P | FT ( ρ + ρ n 1 + ) ( k x , k m ) FT ( e 2 ) ( k x , k y k m ) | 2 k y P + FT ( ( ρ + ρ n 1 + ) | P 2 2 ( k x , k m ) P | FT ( e 2 ) ( k x , k y k m ) | 2 FT ( ( ρ + ρ n 1 + ) | P 2 2 k y P + ( k x , k m ) P | FT ( e 2 ) ( k x , k y k m ) | 2
(34)

Abbreviations

MRI: 

Magnetic resonance imaging

MREIT: 

Magnetic resonance electrical Impedance Tomography

MRCDI: 

Magnetic resonance current density imaging

EIT: 

Electrical impedance tomography

SNR: 

Signal-to-noise ratio

SS-SEPI: 

Single-shot Spin-echo Echo Planar Imaging

FOV: 

Field of view

ICNE: 

Injection current nonlinear encoding.

Declarations

Acknowledgements

E J Woo was supported by the National Research Foundation of Korea (NRF) grant funded by the Korea government (MSIP)(No. 2010-0018275). O I Kwon and H J Kim were supported by Basic Science Research Program through the National Research Foundation of Korea(NRF) funded by the Ministry of Education, Science and Technology (No. 2013R1A2A2A04016066, 2012R1A1A2008477).

Authors’ Affiliations

(1)
Department of Biomedical Engineering, Kyung Hee University
(2)
Department of Mathematics, Konkuk University

References

  1. Joy MLG, Scott GC, Henkelman RM: In vivo detection of applied electric currents by magnetic resonance imaging. Magn Reson Imag 1989, 7: 89–94. 10.1016/0730-725X(89)90328-7View ArticleGoogle Scholar
  2. Scott GC, Joy MLG, Armstrong RL, Henkelman RM: Measurement of nonuniform current density by magnetic resonance. IEEE Trans Med Imag 1991, 10: 362–374. 10.1109/42.97586View ArticleGoogle Scholar
  3. Ider YZ, Birgul O: Use of the magnetic field generated by the internal distribution of injected currents for Electrical Impedance Tomography (MR-EIT). Elektrik 1998, 6: 215–225.Google Scholar
  4. Eyuboglu M, Birgul O, IY Z: A dual modality system for high resolution-true conductivity imaging. Proc. XI Int. Conf. Electrical Bioimpedance (ICEBI) 2001, 409–413.Google Scholar
  5. Birgul O, Eyuboglu BM, Ider YZ: Current constrained voltage scaled reconstruction (CCVSR) algorithm for MR-EIT and its performance with different probing current patterns. Phys Med Biol 2003, 48: 653–671. 10.1088/0031-9155/48/5/307View ArticleGoogle Scholar
  6. Kwon O, Woo EJ, Yoon JR, Seo JK: Magnetic resonance electrical impedance tomography (MREIT): simulation study of J-substitution algorithm. IEEE Trans Biomed Eng 2002, 48: 160–167.View ArticleGoogle Scholar
  7. Kim YJ, Kwon O, Seo JK, Woo EJ: Uniqueness and convergence of conductivity image reconstruction in magnetic resonance electrical impedance tomography. Inverse Probl 2003, 19: 1213–1225. 10.1088/0266-5611/19/5/312MathSciNetView ArticleGoogle Scholar
  8. Ider YZ, Onart S, Lionheart WRB: Uniqueness and reconstruction in magnetic resonance-electrical impedance tomography(MR-EIT). Physiol Meas 2003, 24: 591–604. 10.1088/0967-3334/24/2/368View ArticleGoogle Scholar
  9. Muftuler L, Hamamura M, Birgul O, Nalcioglu O: Resolution and contrast in magnetic resonance electrical impedance tomography (MREIT) and its application to cancer imaging. Tech Cancer Res Treat 2004, 3: 599–609.View ArticleGoogle Scholar
  10. Ozdemir M, Eyuboglu BM, Ozbek O: Equipotential projection-based magnetic resonance electrical impedance tomography and experimental realization. Phys Med Biol 2004, 49: 4765–4783. 10.1088/0031-9155/49/20/008View ArticleGoogle Scholar
  11. Seo JK, Yoon JR, Woo EJ, Kwon O: Reconstruction of conductivity and current density images using only one component of magnetic field measurements. IEEE Trans Biomed Eng 2003, 50: 1121–1124. 10.1109/TBME.2003.816080View ArticleGoogle Scholar
  12. Oh SH, Lee BI, Woo EJ, Lee SY, Cho MH, Kwon O, Seo JK: Conductivity and current density image reconstruction using harmonic Bz algorithm in magnetic resonance electrical impedance tomography. Phys Med Biol 2003, 48: 3101–3116. 10.1088/0031-9155/48/19/001View ArticleGoogle Scholar
  13. Park C, Kwon O, Woo EJ, Seo JK: Electrical conductivity imaging using gradient Bz decomposition algorithm in magnetic resonance electrical impedance tomography (MREIT). IEEE Trans Med Imag 2004, 23: 388–394. 10.1109/TMI.2004.824228View ArticleGoogle Scholar
  14. Joy MLG: MR current density and conductivity imaging: the state of the art. In Proc. 26th Ann. Int. Conf. IEEE EMBS. San Francisco; 2004:5315–5319.Google Scholar
  15. Lee BI, Lee SH, Kim TS, Kwon O, Woo EJ, Seo JK: Harmonic decomposition in PDE-based denoising technique for magnetic resonance electrical impedance tomography. IEEE Trans Biomed Eng 2005, 52: 1912–1920. 10.1109/TBME.2005.856258View ArticleGoogle Scholar
  16. Oh SH, Lee BI, Woo EJ, Lee SY, Kim TS, Kwon O, Seo JK: Electrical conductivity images of biological tissue phantoms in MREIT. Physiol Meas 2005, 26: S279-S288. 10.1088/0967-3334/26/2/026View ArticleGoogle Scholar
  17. Gao N, Zhu SA, He BA: New magnetic resonance electrical impedance tomography (MREIT) algorithm: the RSM-MREIT algorithm with applications to estimation of human head conductivity. Phys Med Biol 2006, 51: 3067–3083. 10.1088/0031-9155/51/12/005View ArticleGoogle Scholar
  18. Birgul O, Hamamura M, Muftuler L, Nalcioglu O: Contrast and spatial resolution in MREIT using low amplitude current. Phys Med Biol 2006, 51: 5035–5049. 10.1088/0031-9155/51/19/020View ArticleGoogle Scholar
  19. Kim HJ, Oh TI, Kim YT, Lee BI, Woo EJ, Seo JK, Lee SY, Kwon O, Park C, Kang BT, Park HM: In vivo electrical conductivity imaging of a canine brain using a 3 T MREIT system. Physiol Meas 2008, 29: 1145–1155. 10.1088/0967-3334/29/10/001View ArticleGoogle Scholar
  20. Kim HJ, Kim YT, Minhas AS, Jeong WC, Woo EJ, Seo JK, Kwon OJ: In Vivo high-resolution conductivity imaging of the human leg using MREIT: the first human experiment. IEEE Trans Med Imag 2009, 99: 160–167.Google Scholar
  21. Woo EJ, Seo JK: Magnetic resonance electrical impedance tomography (MREIT) for high-resolution conductivity imaging. Physiol Meas 2008, 29: R1-R26. 10.1088/0967-3334/29/10/R01View ArticleGoogle Scholar
  22. Hamamura M, Muftuler L: Fast imaging for magnetic resonance electrical impedance tomography. Magn Reson Imaging 2008, 26: 739–745. 10.1016/j.mri.2008.01.031View ArticleGoogle Scholar
  23. Muftuler L, Chen G, Hamamura M, Ha SH: MREIT with SENSE acceleration using a dedicated RF coil design. Physiol Meas 2009, 30: 913–929. 10.1088/0967-3334/30/9/004View ArticleGoogle Scholar
  24. Park HM, Nam HS, Kwon O: Magnetic flux density reconstruction using interleaved partial Fourier acquisitions in MREIT. Phys Med Biol 2011, 56: 2059–2073. 10.1088/0031-9155/56/7/010View ArticleGoogle Scholar
  25. Scott GC, Joy MLG, Armstrong RL, Henkelman RM: Sensitivity of magnetic resonance current density imaging. J Magn Reson 1992, 97: 235–254.Google Scholar
  26. Sadleir R, Grant S, Zhang SU, Lee BI, Pyo HC, Oh SH, Park C, Woo EJ, Lee SY, Kwon O, Seo JK: Noise analysis in MREIT at 3 and 11 Tesla field strength. Physiol Meas 2005, 26: 875–884. 10.1088/0967-3334/26/5/023View ArticleGoogle Scholar
  27. Nam H, Kwon O: Optimization of multiply acquired magnetic flux density Bz using ICNE-Multiecho train in MREIT. Phys Med Biol 2010, 55: 2743–2759. 10.1088/0031-9155/55/9/021View ArticleGoogle Scholar
  28. Park C, Lee BI, Kwon O, Woo EJ: Measurement of induced magnetic flux density using injection current nonlinear encoding (ICNE) in MREIT. Physiol Meas 2007, 28: 117–127. 10.1088/0967-3334/28/2/001View ArticleGoogle Scholar
  29. Woo EJ: Functional brain imaging using MREIT and EIT: Requirements and feasibility. 8th Int. Conf. on Bioelectromagnetism 2011, 131–134.Google Scholar
  30. Sadleir RJ, Grant SC, Woo EJ: Can high-field MREIT be used to directly detect neural activity? Theoretical considerations. Neuroimage 2010, 52: 205–216. 10.1016/j.neuroimage.2010.04.005View ArticleGoogle Scholar

Copyright

© Chauhan et al.; licensee BioMed Central Ltd. 2013

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.