Home | About | Journals | Submit | Contact Us | Français |

**|**HHS Author Manuscripts**|**PMC2962531

Formats

Article sections

Authors

Related links

Phys Med Biol. Author manuscript; available in PMC 2010 November 1.

Published in final edited form as:

Published online 2009 November 4. doi: 10.1088/0031-9155/54/22/018

PMCID: PMC2962531

NIHMSID: NIHMS239940

Corresponding author: Y.-C. N. Cheng, MRI Center/Concourse Research, Harper University Hospital, 3990 John R Street, Detroit MI, 48201, Tel: 313-966-8220, Fax: 313-745-9182, Email: ude.enyaw@61cxy

The publisher's final edited version of this article is available at Phys Med Biol

See other articles in PMC that cite the published article.

A new procedure for accurately measuring effective magnetic moments of long cylinders is presented. Partial volume, dephasing, and phase aliasing effects are naturally included and overcome in our approach. Images from a typical gradient echo sequence at one single echo time are usually sufficient to quantify the effective magnetic moment of a cylindrical-like object. Only pixels in the neighborhood of the object are needed. Our approach can accurately quantify the magnetic moments and distinguish subpixel changes of cross sections between cylindrical objects. Uncertainties of our procedure are studied through the error propagation method. Images acquired with different parameters are used to test the robustness of our method. Alternate approaches and their limitations to extract magnet moments of objects with different orientations are also discussed. Our method has the potential to be applied to any long object whose cross section is close to a disk.

Changes of blood oxygenation cause a change in the magnetic susceptibility of venous blood. This is the basis of functional magnetic resonance imaging (fMRI) (e.g., see [1]). For this reason, there has been a strong interest to quantify the magnetic susceptibility of venous blood. As a forward problem, researchers in MRI usually model a blood vessel as an infinitely long cylinder and calculate its MR signal with a known susceptibility of the vessel obtained from *ex vivo* studies [2]. In order to avoid the partial volume effect in MRI but to obtain accurate results, an *ex vivo* study usually analyzes the phase information both inside and outside a cylindrical test tube whose diameter occupies at least 10 image pixels [2]. If the wall thickness of the test tube can be neglected, this approach can accurately recover the susceptibility difference between the materials inside and outside the cylinder. On the other hand, the *in vivo* measurement of the blood susceptibility has also been performed on vessels with large enough diameters (e.g., [3]). Nonetheless, the quantification of susceptibility of a cylindrical object with a small diameter has been a challenging problem and has not yet been tackled. Until recently, we [4] and Sedlacik et al. [5] have attempted to solve this inverse problem with two independent approaches. The former requires *a priori* knowledge of the cylinder radius while the latter fits several variables to MR images from multiple long echo times, which can become impractical in *in vivo* studies. In this paper, we will demonstrate a general approach in the quantification of the effective magnetic moments of narrow cylindrical objects from only one echo time of a gradient echo image. This new approach relies only on the known imaging parameters but not on the object size, which is usually not known in advance. From the magnetostatic theory [6], in the first order approximation, the magnetic moment per length of a long but narrow object is equal to that of an infinitely long cylinder. This concept makes our approach described in this paper more applicable to long objects whose cross sections may not be exactly a disk. Practically, our approach may be applied to clinical MR images and industrial problems [7].

In this paper, we will also show our improved method to identify the center of the cylindrical object. This improved method is robust as we will prove mathematically. Uncertainties of our methods including both the thermal noise from an object and the systematic uncertainty due to discrete pixels will be studied from the error propagation method [8]. By minimizing the uncertainty due to the thermal noise, we will derive criteria of how to effectively evaluate the magnetic moment of an object from images. Certain technical procedures that have been discussed in detail in [4] will be briefly mentioned here but will not be repeated. In the following sections, we will provide the general formulas but focus on null-signal objects perpendicular to the MRI main field in our phantom studies. Results from different imaging parameters will be shown. Objects with different orientations and limitations of our approach will be discussed toward the end of this paper.

The basic concept of our method is to add the complex MR signal of each voxel around a long cylindrical object of interest. As we model each object as an infinitely long cylinder, we only need to examine signals on a plane whose normal is parallel to the axis of the cylindrical object. The cross section of the cylinder on the plane is a disk of which the circumference is a circle. Thus, these terms will also be used in this paper. Theoretically, the MR signal around the infinitely long cylinder can be easily modeled. The total signal within a co-axial cylinder can be formulated and can be compared to the actual signal obtained from images, such that the effective magnetic moment can be solved. However, as different types of noise exist in images and the actual center of the object requires a subpixel determination, the challenge is to identify a set of robust procedures to fully solve this problem.

Applying an external magnetic field on an infinitely long cylindrical object with an absolute susceptibility *χ _{i}* embedded in an environment of susceptibility

$${\phi}_{\mathit{out}}=-\gamma \mathrm{\Delta}B{T}_{E}=-\gamma \frac{\mathrm{\Delta}\chi}{2}\frac{{a}^{2}}{{\rho}^{2}}{B}_{0}{T}_{E}cos2\psi {sin}^{2}\theta $$

(1)

where *γ* is the proton gyromagnetic ratio (2*π* · 42.58 MHz/T), *T _{E}* is the echo time, Δ

In the context of this paper, it is Δ*χ* that we wish to measure from experiments. Thus, the term susceptibility in this paper refers to Δ*χ* rather than the absolute susceptibility of the object. In addition, we have adopted SI units. From Eq. 1, we define the maximal (or minimal, depending on the sign of Δ*χ*) phase value *g* as

$$g\equiv 0.5\gamma \mathrm{\Delta}\chi {B}_{0}{T}_{E}$$

(2)

and the effective magnetic moment (or “magnetic moment, hereafter) as

$$p\equiv {ga}^{2}$$

(3)

Under the same consideration, the induced phase value inside the cylinder is [1]

$${\phi}_{in}=-\gamma \frac{\mathrm{\Delta}\chi}{6}(3{cos}^{2}\theta -1){B}_{0}{T}_{E}=\frac{g}{3}(1-3{cos}^{2}\theta )$$

(4)

If we consider a pseudo cylinder with radius *R* that is co-axial with the cylindrical object (see Fig. 2), for *g* ≠ 0, the overall MR signal of the pseudo cylinder from a gradient echo sequence is [4]

A schematic drawing shows the cross section of a cylindrical object with radius *a*, enclosed by a co-axial pseudo cylinder whose radius is *R*. The total MR signal within the pseudo cylinder can be formulated.

$$S=\ell {\rho}_{0}{\int}_{a}^{R}d\rho \phantom{\rule{0.16667em}{0ex}}\rho {\int}_{0}^{2\pi}d\psi \phantom{\rule{0.16667em}{0ex}}{e}^{i{\phi}_{\mathit{out}}}+\pi \ell {a}^{2}{\rho}_{0,c}{e}^{i{\phi}_{in}}=\pi \ell {\rho}_{0}p{\int}_{p/{R}^{2}}^{g}\frac{dx}{{x}^{2}}{J}_{0}(x{sin}^{2}\theta )+\pi \ell {a}^{2}{\rho}_{0,c}{e}^{i{\phi}_{in}}$$

(5)

where is an arbitrary length of the cylindrical object and can be the slice thickness of the image, *J*_{0} is the zeroth order Bessel function, *ρ*_{0,}* _{c}* is the effective spin density of the cylindrical object, and

Equation 5 is valid only when the center of the disk (i.e., the cross section of the cylindrical object) is identical to the center of the pseudo disk. Thus the first step is to identify the center of the disk. When a cylindrical object itself has no MR signal, i.e., *ρ*_{0,}* _{c}* = 0, we can identify the center of the object by minimizing the imaginary part of the signal (Eq. 5) within the pseudo disk [4].

When a cylindrical object has a non-zero *ρ*_{0,}* _{c}*, even if the pseudo circle is replaced by two pseudo concentric circles such that the imaginary part from the annular ring is theoretically zero, the above approach still fails. This is because when the object has an imaginary part of the MR signal, its point spread function leads to an additional imaginary signal in each pixel outside the object. On the other hand, although the “subvoxel shift approach by [5] is valid, it requires many Fourier transformations of an image and thus longer computing time. We offer two alternate approaches below.

The first alternate approach is to maximize the real part of the signal from an annular ring region where the disk is completely inside the smaller pseudo circle. Assume that the annular ring region is formed by two pseudo circles with radii *R*_{1} and *R*_{2} and *a* < *R*_{1} < *R*_{2}. The center of the object is located at coordinates (*x*_{0}, *y*_{0}) with
$\sqrt{{x}_{0}^{2}+{y}_{0}^{2}}+a\le {R}_{1}$. The MR signal from the annular ring is

$$S=\ell {\rho}_{0}{\int}_{{R}_{1}}^{{R}_{2}}d\rho \phantom{\rule{0.16667em}{0ex}}\rho {\int}_{0}^{2\pi}d\psi \phantom{\rule{0.16667em}{0ex}}exp\left\{-i\phantom{\rule{0.16667em}{0ex}}g\frac{{a}^{2}}{{\rho}^{2}}cos2\psi {sin}^{2}\theta \right\}$$

(6)

where the relations between *ρ*, *ρ*, *ψ*, and *ψ* in the integrand are given by *ρ*cos*ψ* = *x*_{0} + *ρ*cos*ψ* and *ρ*sin*ψ* = *y*_{0} + *ρ*sin*ψ*. One can analytically prove that at *x*_{0} = *y*_{0} = 0,
${\scriptstyle \frac{\partial S}{\partial {x}_{0}}}={\scriptstyle \frac{\partial S}{\partial {y}_{0}}}=0$, and
${\scriptstyle \frac{{\partial}^{2}S}{\partial {x}_{0}\partial {y}_{0}}}=0$. Appendix A proves that the real part of the signal is maximum at *x*_{0} = *y*_{0} = 0 under certain conditions.

The second approach is to minimize the real part of the signal within one pseudo circle with radius *R*. With the parameters defined above and with
$\sqrt{{x}_{0}^{2}+{y}_{0}^{2}}<a$, the MR signal of the pseudo disk is

$$S=\ell {\rho}_{0}{\int}_{0}^{2\pi}d\psi {\int}_{r(\psi )}^{R}d\rho \phantom{\rule{0.16667em}{0ex}}\rho \phantom{\rule{0.16667em}{0ex}}exp\left\{-i\phantom{\rule{0.16667em}{0ex}}g\frac{{a}^{2}}{{\rho}^{2}}cos2\psi {sin}^{2}\theta \right\}+\pi \ell {a}^{2}{\rho}_{0,c}{e}^{i{\phi}_{in}}$$

(7)

where the relations *ρ*cos*ψ* = *x*_{0}+ *ρ*cos*ψ* and *ρ*sin*ψ* = *y*_{0}+ *ρ*sin*ψ* still hold. In addition, *r*(*ψ*) is the minimal value of *ρ* and satisfies *r*(*ψ*) cos*ψ* = *x*_{0} + *a* cos*ψ* and *r*(*ψ*) sin*ψ* = *y*_{0} + *a* sin*ψ*. This leads to

$${(r(\psi )cos\psi {x}_{0})}^{2}+{(r(\psi )sin\psi {y}_{0})}^{2}={a}^{2}$$

(8)

and

$$\begin{array}{l}r(\psi )=({x}_{0}cos\psi +{y}_{0}sin\psi )+\sqrt{{({x}_{0}cos\psi +{y}_{0}sin\psi )}^{2}+({a}^{2}{x}_{0}^{2}{y}_{0}^{2})}\\ =\sqrt{{x}_{0}^{2}+{y}_{0}^{2}}cos(\psi -{\psi}_{0})+\sqrt{{a}^{2}({x}_{0}^{2}+{y}_{0}^{2}){sin}^{2}(\psi -{\psi}_{0})}\end{array}$$

(9)

where $cos{\psi}_{0}\equiv {x}_{0}/\sqrt{{x}_{0}^{2}+{y}_{0}^{2}}$ and $sin{\psi}_{0}\equiv {y}_{0}/\sqrt{{x}_{0}^{2}+{y}_{0}^{2}}$.

Again, one can analytically prove that at *x*_{0} = *y*_{0} = 0,
${\scriptstyle \frac{\partial S}{\partial {x}_{0}}}={\scriptstyle \frac{\partial S}{\partial {y}_{0}}}=0$, and
${\scriptstyle \frac{{\partial}^{2}S}{\partial {x}_{0}\partial {y}_{0}}}=0$. Appendix B proves that the real part of the signal is minimum at *x*_{0} = *y*_{0} = 0.

All these approaches will fail when *θ* = 0. In this special case, no phase value is outside the object (see Eq. 1). The MR signal becomes:

$$S=\ell {\rho}_{0}(A-\pi {a}^{2})+\pi \ell {a}^{2}{\rho}_{0,c}{e}^{-2i\phantom{\rule{0.16667em}{0ex}}g/3}$$

(10)

where *A* is an arbitrarily shaped area that encloses the disk. Only the size of the object can be determined in this case.

As shown in Fig. 3, if three arbitrary concentric circles with radii *R*_{1}, *R*_{2}, and *R*_{3} are chosen, then three complex signals *S*_{1}, *S*_{2}, and *S*_{3} can be calculated from MR images as described in [4]. In order to improve the accuracy of our approach, each image pixel is further divided into 100 subpixels in the calculations of the complex signals *S _{i}* [4]. From Eq. 5, we obtain

(a) Magnitude and (b) its associated phase image at an echo time of 5 ms show an air cylinder in the gel phantom. (c) Magnitude and (d) phase image at an echo time of 20 ms show the same air cylinder in the gel phantom. Although the actual radius of the **...**

$$({S}_{1}-{S}_{2}){\int}_{p/{R}_{2}^{2}}^{p/{R}_{3}^{2}}\frac{dx}{{x}^{2}}{J}_{0}(x{sin}^{2}\theta )=({S}_{2}-{S}_{3}){\int}_{p/{R}_{1}^{2}}^{p/{R}_{2}^{2}}\frac{dx}{{x}^{2}}{J}_{0}(x{sin}^{2}\theta )$$

(11)

where the effective magnetic moment, *p*, becomes the only unknown in the equation. Note that both ±*p* satisfy Eq. 11 but the correct sign of *p* may be determined from MR phase images with the help of Eq. 1. From Eq. 1, it is clear that
$p{sin}^{2}\theta /{R}_{i}^{2}$ is the maximal (or minimal) phase value at the circumference of the *i*-th circle (where *i* = 1, 2, 3). Therefore, if any *R _{i}* is chosen larger than the phase aliasing area, then
$\mid p{sin}^{2}\theta /{R}_{i}^{2}\mid $ will be always less than

The above discussion is generally valid except when the angle *θ* is close to zero. At *θ* = 0, no magnetic moment appears in the signal, which is given by Eq. 10.

It is important to study the uncertainty of the method when it is applied to actual images. In addition, in order to determine the optimal choices of three radii for the quantification of the magnetic moment, the study of the uncertainty can provide some insights. Thermal noise due to the presence of an object and systematic noise due to discrete pixels in images always exist. They lead to the uncertainty of *p*. By defining *p* *p* sin^{2} *θ* and rewriting Eq. 11, we can derive the uncertainty of *p* through error propagation [8]:

$$\frac{\delta p}{p}=\frac{\delta p}{p}=\frac{1}{\mid D\mid}\sqrt{{(\delta ({S}_{2}{S}_{3}))}^{2}{h}_{12}^{2}+{(\delta ({S}_{1}{S}_{2}))}^{2}{h}_{23}^{2}}$$

(12)

where *h _{ij}* is defined as

$${h}_{ij}\equiv {\int}_{p/{R}_{i}^{2}}^{p/{R}_{j}^{2}}\frac{dx}{{x}^{2}}{J}_{0}(x)$$

(13)

with *i*, *j* = 1, 2, 3 and

$$\begin{array}{l}D\equiv ({S}_{1}{S}_{2})\frac{{J}_{0}({\phi}_{3})}{{\phi}_{3}}+({S}_{2}{S}_{3})\frac{{J}_{0}({\phi}_{1})}{{\phi}_{1}}+({S}_{3}{S}_{1})\frac{{J}_{0}({\phi}_{2})}{{\phi}_{2}}\\ =\pi \ell {\rho}_{0}p\left({h}_{12}\frac{{J}_{0}({\phi}_{3})}{{\phi}_{3}}+{h}_{23}\frac{{J}_{0}({\phi}_{1})}{{\phi}_{1}}+{h}_{31}\frac{{J}_{0}({\phi}_{2})}{{\phi}_{2}}\right)\end{array}$$

(14)

Here
${\phi}_{i}\equiv p/{R}_{i}^{2}$ is chosen to be between 0 and *π* for all *i*. In the derivation of Eq. 12, we have already assumed that *R*_{3} < *R*_{2} < *R*_{1} such that the uncertainty from the annular ring region between *R*_{2} and *R*_{3} is uncorrelated with the uncertainty from the area between *R*_{1} and *R*_{2}. The uncertainty from each annular ring consists of the thermal noise and the systematic noise. Both noise sources are also uncorrelated. The former can be approximated by
$\sigma \ell \sqrt{\mathrm{\Delta}x\mathrm{\Delta}y\pi \mid {R}_{i}^{2}{R}_{j}^{2}\mid}$ where *σ* is the standard deviation of the thermal noise in the image and Δ*x*Δ*y* is the image in-plane resolution [4]. The systematic noise is represented by *ε _{ij}* calculated from

Equation 12 can be rewritten as

$$\frac{\delta p}{p}=\frac{\sqrt{{\scriptstyle \frac{\mathrm{\Delta}x\mathrm{\Delta}y}{\pi p{\text{SNR}}^{2}}}\left(\left|{\scriptstyle \frac{1}{{\phi}_{2}}}{\scriptstyle \frac{1}{{\phi}_{3}}}\right|{h}_{12}^{2}+\left|{\scriptstyle \frac{1}{{\phi}_{1}}}{\scriptstyle \frac{1}{{\phi}_{2}}}\right|{h}_{23}^{2}\right)+({\epsilon}_{12}^{2}+{\epsilon}_{23}^{2}){h}_{12}^{2}{h}_{23}^{2}}}{\left|{h}_{12}{\scriptstyle \frac{{J}_{0}({\phi}_{3})}{{\phi}_{3}}}+{h}_{23}{\scriptstyle \frac{{J}_{0}({\phi}_{1})}{{\phi}_{1}}}+{h}_{31}{\scriptstyle \frac{{J}_{0}({\phi}_{2})}{{\phi}_{2}}}\right|}$$

(15)

where SNR *ρ*_{0}/*σ* is the signal to noise ratio of the magnitude image. In order to determine the optimal combination of _{1}, _{2}, and _{3}, for simplicity, we assume *ε*_{12} = *ε*_{23} = 0 in Eq. 15 and numerically search for the minimum of *δp*/*p* with 0 < * _{i}* ≤

The susceptibility can be calculated from the known effective magnetic moment, if the radius of the object can be determined. Two approaches are considered to determine the radius of the object when the object has no spin density. The first approach is to utilize a spin echo sequence. The second approach is to perform a gradient echo sequence at a short echo time.

If a spin echo sequence is used to image the cross section of the object, the MR signal from an arbitrarily uniform area *A* that includes the object is:

$${S}_{SE}={\rho}_{0,SE}\ell (A-\pi {a}^{2})$$

(16)

where *ρ*_{0,}* _{SE}* is the effective spin density (a constant) of the spin echo images. We have neglected the dephasing effect during the sampling time. With two arbitrary areas

$$\delta (\pi {a}^{2})=\frac{\sqrt{\mathrm{\Delta}x\mathrm{\Delta}y}}{\text{SNR}}\sqrt{{A}_{2}+\frac{{({A}_{2}-\pi {a}^{2})}^{2}}{{A}_{1}{A}_{2}}}$$

(17)

Note that SNR in Eq. 17 refers to the signal to noise ratio of the spin echo image. The noise between *A*_{2} and *A*_{1} − *A*_{2} is uncorrelated. In order to minimize the uncertainty, it is obvious from Eq. 17 that the area *A*_{1}*A*_{2} should be chosen as large as possible before heterogeneity is encountered in the images, while the area *A*_{2} should be chosen as small as possible but larger than *πa*^{2}. In fact, *A*_{2} has to be larger than the area of the distorted object in the image so the cross section of the object can be estimated more accurately.

In the gradient echo approach, an additional pseudo circle with radius *R* may be used to solve Eq. 5 as *g* becomes the only unknown after both p and *ρ*_{0} are found. However, Fig. 4a demonstrates that the integral in Eq. 5 oscillates and quickly approaches an asymptotic value when |*g*| or echo time increases (with a fixed *p*). This is due to the *J*_{0} function and 1/*x*^{2} in the integrand. In order to obtain a unique solution of *g* from Eq. 5, the echo time *T _{E}* needs to be chosen short enough such that |

Images from three different gel phantoms are analyzed in this paper. The first set of images come from the phantom images in [4]. A narrow but long air cylinder with a diameter of roughly 1.6 mm appears in the images. Images acquired from a coronal 3D gradient echo sequence (Fig. 3) and a spin echo sequence are available for our re-analyses. The imaging parameters of the gradient echo sequence were: *T _{E}* = 5 ms and 20 ms,

The second set of images is from a phantom that consisted of a hollow straw in the gel. As one end of the straw was sealed, the hollow straw created another air cylinder in the gel phantom. The diameter of the straw was 5.24 ± 0.01 mm. This phantom was prepared with 107 g of the gel powder in 1200 ml distilled water. The straw was placed vertically in the phantom when the gel was in its liquid form. Coronal images of a multiple echo 3D gradient echo sequence and a spin echo sequence were acquired. The imaging parameters of the gradient echo sequence were: *T _{E}* = 5 ms, 10 ms, and 15 ms,

The distortion of objects in images (Fig. 5) due to the susceptibility effect was studied by varying the read bandwidth. These studies were preformed on another straw phantom made on a different day. The read bandwidths of the gradient echo sequence at *T _{E}* = 10 ms were varied from 70 Hz/pixel, 220 Hz/pixel, 350 Hz/pixel, to 700 Hz/pixel. Other imaging parameters remained the same as those listed in the above paragraph.

In order to calculate the correct systematic noise (or uncertainty) due to discrete pixels, both image parameters and measured magnetic moments are used. In addition, the susceptibility difference between the air and gel is assumed to be 9.4 ppm [12]. The simulations and the addition of Gaussian noise were described in [4]. It is worth mentioning again that proper rotations of data after Fourier transformation are required. This is due to different definitions of the Fourier transformation in different computer software, but we follow the definition provided in [1].

In all simulations, we assume that *θ* = *π*/2, Δ*χ*= 9.4 ppm, and *B*_{0} = 1.5 T. The image resolution is set to 1 mm ×1 mm (or 1 pixel ×1 pixel). Other parameters are provided with corresponding simulated result in the next Section.

The one circle approach appears to be the easiest method in the identification of the object center. The radius of the pseudo circle cannot be too small or too large. In the former case, the thermal noise and discrete pixels will contribute significant uncertainty to the process. We would recommend a radius of at least three pixels. In the latter case, the overall MR signal will be dominated from the pixels without magnetic moment information and become insensitive to identifying the center. Equation 36 in Appendix C is consistent with these considerations. Based on the phase images, a circle whose circumference intersects with phase values around ±2 radians along the vertical and horizontal axes is a reasonable choice while the phase of 2.6 radians would be the optimal choice. In the approach of using two concentric circles, we suggest using radii that differ by at least one pixel such that enough pixels are used in the analysis. Based on our experience in simulation and phantom studies, the center of the object identified by any of our three procedures often differs from the actual center by 0.1 to 0.3 pixel. In addition, all the results from the one circle approach agree with the theoretical criterion (Eq. 36 in Appendix C). We have also tried the approach by Sedlacik et al. [5] on some of our simulated images and have identified the same centers determined by our one-circle method. However, as the “subvoxel shift approach by [5] involves a 2D Fourier transformation of every possible center, while our method only requires a sum of complex numbers from the neighboring subpixels, our method takes less than half the computing time of their approach.

Simulations of Eq. 15 reveal that the uncertainty of the magnetic moment will be at a minimum when the three phase values are smallest, roughly 1 radian, and 2.5–3 radians. This result means that one circle should be chosen as large as possible while the other two circles should be chosen when their radii are roughly the lengths from the center of the object to the desired phase values. All circles should be outside the phase aliasing region. Due to the left handed setting by MRI manufacturers, note that the signs in all phase images of Fig. 3 are opposite to the signs in Eq. 1 (right handed system). Practically, as the uncertainty in the phase image is inversely proportional to SNR in the magnitude image [1], the radius of the largest circle can be chosen where the phase value is not smaller than 0.1 radian. Furthermore, as the MR signal is discretized and the center of the object may be different from the center of a pixel, the phase value measured from MR images may not be exactly equal to any of the above desired phase values and may not be symmetric around the object. In addition, phase values close to ±*π* radians may not be available from the phase profiles due to the partial volume effect. Simulated phase profiles shown in Fig. 6 reflect these practical problems.

An infinitely long cylinder with radius 0.8 pixel is simulated without the thermal noise. The center of the object is purposely shifted to the 128.9th pixel in simulations. Our method identified a center that is only 0.1 pixel away from the simulated **...**

In any case, the phase profiles are only used as references for choosing the radii of the three circles for the measurements of the magnetic moment. It is not necessary to identify the exact phase values for our measurements. In fact, any choice of three circles can lead to a solution of the magnetic moment except that the uncertainty may be larger than desired (see below). Nevertheless, in order to better estimate the uncertainty (Eq. 15) due to each radius used in our analyses, we have averaged four radii and their corresponding phase values taken from a vertical and a horizontal phase profile through the center of the object.

Table 1a lists the magnetic moments obtained from our previous data [4]. As the radii of the three circles are chosen as described above, the uncertainties calculated from Eq. 15 are less than 10%. As a comparison, when the radii of the largest circles are reduced, the uncertainties of the measured magnetic moments are increased but are still within 16% (Table 1b).

(a) Magnetic moments of six slices from the existing images at echo times of 5 ms and 20 ms. The first column lists the slice number, with the smaller number representing the top slice. The second and the seventh columns list the values of the three radii **...**

In theory, the effective magnetic moment is proportional to the echo time. The results of the magnetic moment obtained from different echo times clearly show that relationship within uncertainties in Table 1. With the known radius of the object 0.8 mm and an assumed Δ*χ*= 9.4 ppm, the theoretical magnetic moments at *T _{E}* = 5 ms and 20 ms are 6.04 radian·mm

Two slices of the second phantom (straw phantom) at three echo times are analyzed and their results are shown in Table 2. As the diameter of the straw is 5.24 ± 0.01 mm, the theoretical magnetic moments at *T _{E}* = 5 ms, 10 ms and 15 ms are 64.7 radian·mm

Table 3 shows magnetic moments of the straw at four different read bandwidths from the same slice. All the measured magnetic moments agree with each other within uncertainties. Therefore, our measurements of the magnetic moment from a gradient echo sequence do not seem to be affected by the distortion effect even at a relatively low bandwidth. The measured moments in both Tables 2 and and33 are slightly different from the theoretical values as the presence of the gel powder used in the phantoms can shift the susceptibility away from that of pure water. As we choose the radii outside the phase aliasing regions, we also effectively avoid the distortion effect.

When the object has no spin density, the cross section of the object and its uncertainty are calculated based on Eq. 16 and Eq. 17. Three slices of the spin echo images from our previous study are analyzed as they correspond to slice 10, 19, and 37 of the gradient echo images. The measured cross sections are shown in Table 4. The result from slice 19 agrees well with the theoretical value of the object cross section, 2.01 mm^{2}. In general, the results from Table 4 also indicate a collapse at the bottom of the phantom and an enlargement of the object at the top of the phantom.

Three cross sections *A*_{measured} of the narrow air cylinder are analyzed from the spin echo images. The areas of *A*_{1} and *A*_{2} are shown in the number of pixels while the values of *A*_{measured} and *A*_{predicted} are shown in units of mm^{2}. The uncertainties of the **...**

For the straw phantom, one slice of the spin echo images at two different read bandwidths is analyzed. The SNR is 33:1 for images acquired with bandwidth 90 Hz/pixel and is 11:1 for images with bandwidth 590 Hz/pixel. The measured cross sections are 21.4 ± 0.3 mm^{2} and 22.2 ± 0.8 mm^{2} at bandwidth 90 Hz/pixel and 590 Hz/pixel, respectively. In the former case, *A*_{1} = 1148 mm^{2} and *A*_{2} = 120 mm^{2}, while in the latter case *A*_{1} = 821 mm^{2} and *A*_{2} = 69 mm^{2}. The theoretical value of the cross section is 21.6 mm^{2} if the cross section is a perfect disk. As the straw has slightly deformed inside the gel phantom, its cross section may not be a perfect disk but has a smaller cross section than that of the perfect disk, given the same circumference. These results show good agreement between the measurements and the theoretical value within uncertainties predicted by Eq. 17. The results also indicate that the geometric distortion does not significantly affect the quantification.

In addition to our efforts in the previous subsection, we also explore the possibility of quantifying susceptibility from gradient echo images. We investigate this possibility with simulations. If all variables and parameters are known, we can calculate the MR signal in each case. Such an MR signal is a number and is displayed as a horizontal straight line in Fig. 4a or 4b. On the other hand, if the magnetic moment *p* is known but the susceptibility and object radius are not known, then we can simulate the MR signal as a function of *g*, as shown as a curve in Fig. 4a or 4b. The intersections of the curve and the line indicate the possible solutions of the susceptibility. The result shown in Fig. 4a indicates that the minimal value of |Δ*χ*| (with a fixed *p*) can be extracted at a given echo time. In this example, the minimal value of |Δ*χ*| is roughly 1.5 ppm, much less than 9.4 ppm used for the image simulations. If the echo time is shortened, as shown in Fig. 4b, then *g* or susceptibility is likely to be uniquely determined. If the thermal noise is added as one arbitrary unit (compared to *ρ*_{0}, which is 10 units) and an image resolution is assumed as 1 mm^{2}, then the uncertainty of the signal *S* within a radius *R* (3 mm) is roughly 5.3 units (
$\sqrt{\pi {R}^{2}}/\sqrt{\mathrm{\Delta}x\mathrm{\Delta}y}$). This means that the signal *S* is within one standard deviation of the asymptotic value shown in Fig. 4a when *g* is larger than roughly 3.5. On the other hand, with a much shorter echo time (which may be impractical), the susceptibility and radius of the object can be resolved from Eq. 5 and Fig. 4b. However, this result implies that the spin echo approach presented in the above subsection seems better.

The examples shown in Fig. 4a and 4b indicate that one cannot accurately determine the susceptibility and object size through a curve fitting of data from long echo times, where *g* becomes too large. If sufficient noise had been included in the simulations, both Fig. 4c and 4d demonstrate that the volume fraction and susceptibility each can only be roughly determined at orders of magnitude through curve fitting. However, when an object has enough MR signal, the curve fitting method appears applicable to determine the susceptibility and volume of the object [5]. This is because Eq. 4 only contains the susceptibility information (i.e., *g*), rather than the magnetic moment (i.e., *p*).

Three common problems in MRI have been overcome in our method of obtaining the magnetic moment: partial volume effect, dephasing effect (signal loss shown in Fig. 6a and 6c), and the phase aliasing effect (shown in Fig. 6d). In our examples, the partial volume effect is overshadowed by the dephasing effect. These are practical reasons why our method applied on small objects is much better than the conventional least squares fitting method or the method of measuring the relaxation time ${T}_{2}^{\ast}$. Some results of how bad the least squares fitting method has performed are shown in [4].

As shown in Fig. 6, even when the object has no signal, the actual center of the object may not be in the pixel with the lowest magnitude signal. This fact implies that the center of the object has to be determined based on the information outside the object, as suggested by our procedures here or by Sedlacik et al. [5]. However, the approach by [5] will take much longer time than ours.

After the center of the object is identified, two phase profiles through the object center along the vertical and horizontal directions are used for the selections of the radii of three circles. All three circles should be larger than the phase aliasing area around the object but the smallest circle should be as small as possible (i.e., as close to the aliased phase as possible). However, as shown in Fig. 6, the maximum phase value may not be close to ±*π* radians but a choice of a slightly larger radius of the smallest circle does not seem to increase the uncertainty much. Due to the discrete pixels, a function can be used to interpolate the phase values such that more accurate phase values can be used for the selections of three radii at non-integer pixels. For this purpose, we have tried a linear function and a function based on Eq. 1 and we have found no difference of which function to use. This is clearly supported by the fact that Eq. 11 does not rely on the knowledge of phase values, which only serve as references of how to choose the radii. The radii of non-integer pixels are an important feature in our method. For example, if we consider two radii at 2.3 and 2.6 pixels, both radii refer to the same pixel along the vertical and horizontal directions. However, as the signal inside a circle can change even with a slight change of its radius, Eq. 11 offers a much more accurate way to solve the magnetic moment. Practically, a uniform region may appear only around the neighborhood of the object. Even with non-optimal choices of radii, our method appears to be effective as most of the uncertainties listed in Table 1b are still within 10%.

The uncertainty *δp*/*p* derived from Eq. 15 is inversely proportional to SNR and
$\sqrt{p}$ (or
$\sqrt{{T}_{E}}$). This is clearly shown in Table 1. Although Eq. 15 implies that the longer the echo time, the smaller the uncertainty, it is important to remember that the SNR of the material around the object will also reduce when the echo time increases. Thus, if the echo time is longer than the relaxation time *T*_{2} of the surrounding material, the uncertainty will likely increase.

The systematic uncertainty and the uncertainty due to the thermal noise in the measurements of magnetic moments are comparable in most of our phantom studies. Although all sources of uncertainties are considered uncorrelated in Eq. 15, actually, the systematic uncertainties from two annular rings *S*_{1}*S*_{2} and *S*_{2}*S*_{3} may be partially correlated as they share the same boundary (2*πR*_{2}). However, we have neglected that correlation in this research.

Given the uncertainties and magnetic moments shown in Table 1, the internal collapsing of the air cylinder in the first phantom seems the only explanation. The measurements of the cross sections from the spin echo images also support this finding (Table 4). This collapsing is not obvious from the visual inspection of images. Assuming Δ*χ*= 9.4 ppm, the radius of the air cylinder appears to change from 0.8 mm (slice 19) to 0.67 mm (slice 46). With an image resolution of 1 mm ×1 mm, such a subpixel change of the object radius can be clearly distinguished from our magnetic moment measurements.

The presence of the object susceptibility can lead to distortion artifacts and dephasing effects in images. They occur in all images but in the gradient echo images the dephasing effect at the echo dominates the artifacts unless the read bandwidth is very low. A usual way to minimize the distortion artifacts is to increase the read bandwidth in a sequence. However this change will lead to a lower SNR in images and thus a larger uncertainty (Eq. 15 and Table 4). Our phantom studies indicate that the geometric distortion does not seem to affect the measurements of magnetic moment as shown in Table 3. In the spin echo images, as the dephasing effect due to the susceptibility difference can be enhanced during the sampling period, the measured volume of the object should be treated as the maximal volume. This is supported by the results shown in Table 4. The geometric distortion also does not seem to affect the volume measurements of the objects from the spin echo images.

Although the absolute susceptibility of gel is assumed as the same as the absolute susceptibility of water, −9 ppm, in our analyses, the actual susceptibility of the gel may be slightly different than that of water.

Further examination of Eq. 17 indicates some interesting aspects and the limitation of the spin echo approach. If the cross section of the object occupies more than one pixel in the image, then the uncertainty estimated from Eq. 17 does not depend on the image resolution but it depends on the fields of view. This is because SNR is proportional to
$\sqrt{{L}_{x}{L}_{y}\mathrm{\Delta}x\mathrm{\Delta}y}$ where *L _{x}* and

Practically, a spin echo sequence is routinely performed in each clinical diagnosis. As the magnetic moment can be determined from one gradient echo sequence described in Section 2.3, knowing the radius of the object will automatically lead to the value of susceptibility. The only exception is when *θ* = 0. In this special case, only the size of the object can be determined. When a gradient echo image and a spin echo image are analyzed, no registration of the images is required in our approach.

In the gradient echo approach, the susceptibility and volume cannot be resolved individually when the susceptibility is larger than a certain value with a fixed magnetic moment (see Fig. 4a). This is consistent with the theory of electromagnetism, as in the far field the magnetic moment of an object is the leading term [6]. In general, our experience shows that the gradient echo approach for determining object susceptibility is not favorable, as the desired echo time to minimize distortions and signal losses may need to be very short and thus require a very high read bandwidth. For this reason, the spin echo approach is a better method.

Our method can accurately quantify the effective magnetic moment of a narrow cylinder from MR images. It has overcome three general problems in MRI and is applicable to small objects in a locally uniform background. However, the application of our method is not limited to small objects. The effective magnetic moment of an object can be obtained from standard gradient echo images while the volume of the object may be derived from typical spin echo images. We have demonstrated the feasibility of our approaches through phantom studies. In addition, a phantom study has revealed that a subpixel change of the object volume can be distinguished from the change of magnetic moment using our method in MRI.

This work was supported in part by the Department of Radiology at Wayne State University and NIH HL062983-04A2. The authors would like to thank Mr. Zahid Latif, R.T., and Mr. Yang Xuan for technical assistance in phantom imaging. The authors would also like to acknowledge the use of Mathematica. The author Y.-C. N. Cheng would like to thank Dr. William Condit for his stimulating discussions.

We first define *β* *ρ*^{2} = (*ρ*cos*ψ* −*x*_{0})^{2}+(*ρ*sin*ψ*−*y*_{0})^{2} and *α* (*ρ*cos*ψ*−*x*_{0})^{2}− (*ρ*sin*ψ* −*y*_{0})^{2} such that

$$cos2\psi =\frac{{(\rho cos\psi -{x}_{0})}^{2}-{(\rho sin\psi -{y}_{0})}^{2}}{{\rho}^{2}}\equiv \frac{\alpha}{\beta}$$

(18)

We also define

$$f\equiv {e}^{-\mathit{ip}cos2\psi /{\rho}^{2}}={e}^{-ip\alpha /{\beta}^{2}}$$

(19)

such that the MR signal *S* is

$$S=\ell {\rho}_{0}{\int}_{{R}_{1}}^{{R}_{2}}d\rho \phantom{\rule{0.16667em}{0ex}}\rho {\int}_{0}^{2\pi}d\psi \phantom{\rule{0.16667em}{0ex}}f$$

(20)

The first derivatives of *f* with respect to *x*_{0} and *y*_{0} are

$$\frac{\partial f}{\partial {x}_{0}}=-\mathit{ipf}\frac{2}{{\beta}^{3}}({x}_{0}-\rho cos\psi )[3{(\rho sin\psi -{y}_{0})}^{2}-{(\rho cos\psi -{x}_{0})}^{2}]$$

(21)

$$\frac{\partial f}{\partial {y}_{0}}=\mathit{ipf}\frac{2}{{\beta}^{3}}({y}_{0}-\rho sin\psi )[3{(\rho cos\psi {x}_{0})}^{2}-{(\rho sin\psi {y}_{0})}^{2}]$$

(22)

The second derivatives of *f* with respect to *x*_{0} and *y*_{0} at *x*_{0} = *y*_{0} = 0 are

$${\frac{{\partial}^{2}f}{\partial {x}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=-2{p}^{2}f\frac{1}{{\rho}^{6}}(1+cos6\psi )-6\mathit{ipf}\frac{1}{{\rho}^{4}}cos4\psi $$

(23)

$${\frac{{\partial}^{2}f}{\partial {y}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=-2{p}^{2}f\frac{1}{{\rho}^{6}}(1-cos6\psi )+6\mathit{ipf}\frac{1}{{\rho}^{4}}cos4\psi $$

(24)

As

$${\int}_{0}^{2\pi}d\psi f=2\pi {J}_{0}(p/{\rho}^{2})$$

(25)

$${\int}_{0}^{2\pi}d\psi f\xb7cos4\psi =-2\pi {J}_{2}(p/{\rho}^{2})$$

(26)

$${\int}_{0}^{2\pi}d\psi f\xb7cos6\psi =2\pi i{J}_{3}(p/{\rho}^{2})$$

(27)

it becomes clear that

$$Im{\frac{{\partial}^{2}S}{\partial {x}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=-Im{\frac{{\partial}^{2}S}{\partial {y}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}$$

(28)

and

$$Re{\frac{{\partial}^{2}S}{\partial {x}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=Re{\frac{{\partial}^{2}S}{\partial {y}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=2\pi \ell {\rho}_{0}\left(\frac{p}{{R}_{2}^{2}}{J}_{1}\left(\frac{p}{{R}_{2}^{2}}\right)-\frac{p}{{R}_{1}^{2}}{J}_{1}\left(\frac{p}{{R}_{1}^{2}}\right)\right)$$

(29)

As *R*_{2} > *R*_{1} and the maximum of *xJ*_{1}(*x*) occurs at *x* ≈ 2.4, the first root of *J*_{0}(*x*), proper choices of *R*_{1} and *R*_{2} can lead to a negative value of Eq. 29. Therefore, the object center can be identified by maximizing the real part of the signal from an annular ring.

The calculations used in this Section require the results derived in the Appendix A and Eq. 9. At *x*_{0} = *y*_{0} = 0, *r*(*ψ*) = *a* and the relevant MR signal is

$$S=\ell {\rho}_{0}{\int}_{0}^{2\pi}d\psi {\int}_{a}^{R}d\rho \phantom{\rule{0.16667em}{0ex}}\rho \phantom{\rule{0.16667em}{0ex}}f$$

(30)

which is identical to Eq. 20 if the upper and lower limits of the *ρ* integral are properly replaced.

The second derivatives of Eq. 7 are the results derived from Appendix A with proper changes of the integral limits, plus the following term:

$${\frac{{\partial}^{2}S}{\partial {x}_{0}^{2}}|}_{\mathit{extra}}=\frac{\partial}{\partial {x}_{0}}\left({\frac{\partial S}{\partial {x}_{0}}|}_{\mathit{extra}}\right)-\ell {\rho}_{0}{\int}_{0}^{2\pi}{d\psi r(\psi )\frac{\partial r(\psi )}{\partial {x}_{0}}\frac{\partial f}{\partial {x}_{0}}|}_{\rho =r(\psi )}$$

(31)

where

$${\frac{\partial S}{\partial {x}_{0}}|}_{\mathit{extra}}=-\ell {\rho}_{0}{\int}_{0}^{2\pi}{d\psi r(\psi )\frac{\partial r(\psi )}{\partial {x}_{0}}f|}_{\rho =r(\psi )}$$

(32)

The second derivative with respect to *y*_{0} has the identical format with *x*_{0} replaced by *y*_{0}.

After tedious derivations, we obtain

$$Im{\frac{{\partial}^{2}S}{\partial {x}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=-Im{\frac{{\partial}^{2}S}{\partial {y}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}$$

(33)

and

$$Re{\frac{{\partial}^{2}S}{\partial {x}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=Re{\frac{{\partial}^{2}S}{\partial {y}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}=2\pi \ell {\rho}_{0}\left(\frac{p}{{R}^{2}}\right){J}_{1}\left(\frac{p}{{R}^{2}}\right)>0$$

(34)

when 0 < |*p*/*R*^{2}| < 3.8, the second root of *J*_{1}(*x*). That the result of Eq. 34 only depends on *R* rather than both *R* and *a* indicates that this procedure of finding the object center is still valid even if the center of the circle *R* is outside the object. From a different point of view, when an object (such as a nanoparticle) is very small, one can identify the center of the object by replacing the small object by a much larger object but with the identical magnetic moment *p*. This is because when the value of |*g*| is larger than roughly 5, one cannot unambiguously determine the volume of the object based on the gradient echo signal (see Fig. 4a).

Even though the center of the object can be determined by the method described in the Appendix B, the center cannot be accurately located due to the presence of noise in images. Consider the Taylor expansion of the MR signal

$$S({x}_{0},{y}_{0})=S{\mid}_{{x}_{0}={y}_{0}=0}+\frac{1}{2}{x}_{0}^{2}{\frac{{\partial}^{2}S}{\partial {x}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}+\frac{1}{2}{y}_{0}^{2}{\frac{{\partial}^{2}S}{\partial {y}_{0}^{2}}|}_{{x}_{0}={y}_{0}=0}+(\text{higher}\phantom{\rule{0.16667em}{0ex}}\text{order}\phantom{\rule{0.16667em}{0ex}}\text{terms})$$

(35)

where the first order terms vanish, as described in the main text. The center of the object cannot be accurately determined if the sum of the second order terms is less than the noise term and even if the higher order terms are neglected. Using the result of Eq. 34 and considering only the thermal noise, $\sigma \ell \sqrt{\mathrm{\Delta}x\mathrm{\Delta}y\pi {R}^{2}}$, the following criterion can be established

$$\frac{{x}_{0}^{2}+{y}_{0}^{2}}{\mathrm{\Delta}x\mathrm{\Delta}y}\le \frac{\sqrt{{\scriptstyle \frac{\pi {R}^{2}}{\mathrm{\Delta}x\mathrm{\Delta}y}}}}{\pi \text{SNR}\left({\scriptstyle \frac{p}{{R}^{2}}}\right){J}_{1}\left({\scriptstyle \frac{p}{{R}^{2}}}\right)}=\frac{\sqrt{{\scriptstyle \frac{p}{\mathrm{\Delta}x\mathrm{\Delta}y}}}}{\sqrt{\pi}\text{SNR}{\left({\scriptstyle \frac{p}{{R}^{2}}}\right)}^{3/2}{J}_{1}\left({\scriptstyle \frac{p}{{R}^{2}}}\right)}$$

(36)

The maximum of *x*^{1.5}*J*_{1}(*x*) occurs at *x* ≈ 2.6 and the maximum is roughly 2. At *x* = 1.5 and 2, *x*^{1.5} *J*_{1}(*x*) = 1.02 and 1.63, respectively. With a given magnetic moment *p* in an image, Eq. 36 clearly defines how accurate the center can be determined.

1. Haacke EM, Brown RW, Thompson MR, Venkatesan R. Magnetic Resonance Imaging: Physical Principles and Sequence Design. John Wiley & Sons; 1999.

2. Weisskoff RM, Kiihne S. MRI susceptometry: image-based measurement of absolute susceptibility of MR contrast agent and human blood. Magn Reson Med. 1992;24:375–383. [PubMed]

3. Lai S, Hopkins AL, Haacke EM, Li D, Wasserman BA, Buckley P, Friedman L, Meltzer H, Hedera P, Friedland R. Identification of vascular structures as a major source of signal contrast in high resolution 2D and 3D functional activation imaging of the motor cortex at 1.5 T: Preliminary results. Magn Reson Med. 1993;30:387–392. [PubMed]

4. Cheng YCN, Hsieh CY, Neelavalli J, Liu Q, Dawood MS, Haacke EM. A complex sum method of quantifying susceptibilities in cylindrical objects: The first step toward quantitative diagnosis of small objects in MRI. Magn Reson Imaging. 2007;25(8):1171–1180. [PubMed]

5. Sedlacik J, Rauscher A, Reichenbach JR. Obtaining blood oxygenation levels from MR signal behavior in the presence of single venous vessels. Magn Reson Med. 2007;58:1035–1044. [PubMed]

6. Jackson JD. Classical Electrodynamics. New York: John Wiley & Sons; 1999.

7. Robson P, Hall L. Identifying particles in industrial systems using MRI susceptibility artifacts. American Institute of Chemical Engineers. 2005;51(6):1633–1640.

8. Bevington PR, Robinson DK. Data Reduction and Error Analysis for the Physical Sciences. New York: McGraw-Hill; 1992.

9. Press WH, Flannery BP, Teukolsky SA, Vetterling WT. Numerical Recipes in C. New York: Cambridge University Press; 1992.

10. Hsieh C-Y, Cheng Y-CN, Neelavalli J, Haacke EM. An improved approach of quantifying magnetic moments of small in-vivo objects. Proceedings of the International Society for Magnetic Resonance in Medicine; 2007. p. 2596.

11. Haacke EM, Xu Y, Cheng YCN, Reichenbach J. Susceptibility weighted imaging (SWI) Magn Reson Med. 2004;52(3):612–618. [PubMed]

12. Lide DR, editor. Handbook of Chemistry and Physics. 87. New York: CRC Press; 2006–2007.

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |