MAGNETIC RESONANCE WATER SELF-DIFFUSION TENSOR ENCODING OPTIMIZATION METHODS FOR FULL BRAIN ACQUISITION

Water diffusion tensor magnetic resonance imaging (DT-MRI) is a non-invasive and sensitive modality that is becoming increasingly popular in diagnostic radiology. DT-MRI provides in vivo directional information about the organization and microdynamics of deep brain tissue that is not available by other MRI relaxationbased methods. The DT-MRI experiment involves a host of imaging and diffusion parameters that influence the efficiency (signal-to-noise ratio per unit time), accuracy, and specificity of the information sought. These parameters may include typical imaging parameters such as TE, TR, slice thickness, sampling rate, etc. The DTI relevant parameter space includes pulse duration, separation, direction, number of directions (Ne), order, sign and strength of the diffusion encoding gradient pulses. The goal of this work is to present and compare different tensor encoding strategies used to obtain the DT-MRI information for the whole brain. In this paper an evaluation of tensor encoding advantage is presented using a multi-dimensional non-parametric Bootstrap resampling method. This work also explores the relationship between different tensor encoding schemes using the analytical encoding approach. This work shows that the minimum energy optimization approach can produce uniformly distributed tensor encoding that are comparable to the icosahedral sets. The minimum condition encoding sets are not uniformly distributed and are shown to be suboptimal and related to a commonly used heuristic tensor encoding set. This work shows that the icosahedral set is the only uniformly distributed set with Ne = 6. At equal imaging time, the Bootstrap experiments show that optimal tensor encoding sets can have 6 < Ne < 24.


INTRODUCTION
Recent advances in magnetic resonance based diffusion tensor imaging (DT-MRI) provide important in vivo information about water diffusion properties (direction, magnitude and anisotropy) of biological tissues (Basser et al., 1994;Le Bihan, 1995).A typical DT-MRI experiment is commonly acquired using a spin-echo prepared single-shot echo-planar sequence with bipolar diffusion encoding gradient pulses applied along the three physical gradient channels (Taylor and Bushell, 1985).The measured diffusion sensitized signal is a superposition of spatially dependent incoherent phases, acquired due to random translational Brownian motion during the spatially dependent diffusion encoding magnetic field (Stejskal and Tanner, 1965;Hoehn-Berlag, 1995).The diffusion magnitude and anisotropy may reflect the state of the cellular membrane permeability, myelination and axonal integrity, compartmentalization, and hindranceboth intrinsic and geometric-to mobility (Norris 2001;Cercignani and Horsfield, 2001;Gass et al., 2001;Hoehn-Berlag, 1995).The DT-MRI technique is quantitative and has been compared with anatomical and histochemical findings (Makris et al., 1997;Hsu et al., 1998;Axer et al., 2000;Wall et al., 2000).The DT-MRI information is important in mapping deep white matter architecture (Mori et al., 2001;Melhem et al., 2002) and can be useful in linking the geometric contributions to tissue ionic conductance (Haueisen et al., 2002).Although there are different acquisition strategies and interpretation flavors such as q-space and multi-tensors, the single tensor model approach is adopted here.It is assumed that the imaging parameters are selected such that an acceptable image spatial resolution and SNR has been achieved.A recent analysis by Basser (2002) has shown that the single tensor model is a practical option in a clinical setup if the acquisition parameters are carefully selected.This paper reviews the fundamental aspects of the single tensor encoding and estimation problem.The tensor encoding theory is used to define figures-of-merit that can be used in the design and evaluation of different encoding sets.A Bootstrap statistical procedure is used to compare the performance of different diffusion tensor encoding strategies.This paper also presents some results of full human brain DT-MRI maps acquired at 3T, using optimally designed tensor encoding sets.

DIFFUSION TENSOR PARAMETER SPACE
The symmetric positive definite diffusion tensor (DT) has six independent elements, which can be diagonalized into three eigenvalues and three orthonormal eigenvectors (Fig. 1); therefore, a minimum of six independent encoded measurements, in addition to a reference image, is required to decode the diffusion tensor (Basser et al., 1995).In general, the total number of DT-MRI images, N T , is composed , where N e (s) is the number of encoding directions each averaged N d (s) times at some diffusion weighting step (s).The diffusion weighting b-factor contribution for a bipolar gradient pulse is given by b = (δGγ) 2 (∆-δ/3) (Stejskal and Tanner, 1965), where γ is the nuclear spin gyromagnetic ratio and G is the amplitude of the diffusion gradient pulses.The two-step DT-MRI experiment is acquired by varying the direction, sign, order, strength (G), number of reference (N ref ; b~0) and diffusion weighted images (N d ), at fixed imaging parameters: repetition time (T R ), effective echo time (T E ), diffusion pulse duration (δ) and separation time (∆), slice thickness, and sampling parameters, etc. (Fig. 2).The total imaging time, T, for a two b-factor DT-MRI experiment (N s = 2) can be expressed as T = (N ref + N d N e )T R .The error associated with the set of DT-MRI encoded measurements is a function of the tensor gradient encoding scheme (Conturo et al., 1995) and the partitioning of the total diffusion imaging time, N T (Xing et al., 1997;Cercignani and Horsfield, 1999;Jones et al., 1999).

SELF-DIFFUSION TENSOR ENCODING THEORY
In this section we describe the single diffusion tensor estimation problem.The encoding theory is essential in the optimization of the encoding parameters.Additionally, it is convenient to recognize and use the mathematical analogy of the DT-MRI theory with the theory developed to characterize the chemical shift tensor (Alderman et al., 1990) and the thermal expansivity tensor (Nye, 1957).Since the tensor is diagonally symmetric, a minimum of six noncollinear encoding directions is needed to obtain the six independent elements of the diffusion tensor, D (Basser et al., 1995 Alternatively, the unique tensor elements can be represented be represented by the column vector (Nye, 1957) [ ] The single tensor model that relates the unknown tensor D with the encoding diffusion unit vector, £(k) = [g x (k) g y (k) g z (k)] t , the diffusion weighting bfactor, b, the reference intensity image, I(0), and the measured diffusion weighted images, I(k), can be expressed as To enhance the SNR and reduce phase fluctuations, the reference and diffusion weighted images are repeated and magnitude averaged N ref and N d times respectively.The scalar quantity can be viewed as the projection of the diffusion tensor onto the measurement encoding unit direction g.D m (k) can also be expressed as a matrix product where Define the column vector Ó, which is formed from the measured data, with components where k = 1...N e .The relationship between Ó and d can be expressed in matrix notation as where H is the design encoding matrix, which is formed from the row vectors ©(k) t .The least-squares optimization and the singular value decomposition methods (Hasan et al., 2001a) can now be used to estimate the diffusivity vector, d , where the superscript -1 refers to matrix pseudo-inversion:

ESTIMATION OF THE TENSOR MAPS
The estimated diffusion tensor elements, The mean apparent diffusivity coefficient is defined as ADC = µ(D) = I 1 /3.In addition, define the fractional tensor anisotropy index, FA (Pierpaoli and Basser, 1999).
It has also been demonstrated by Hasan et al. (2001b) that the three principal invariants can be used in the computation of the eigenvalues and the orthonormal eigenvectors of the tensor and in relating different tensor anisotropy measures.

ERROR ESTIMATION OF THE DIFFUSIVITIES AND ENCODING QUALITY FIGURES-OF-MERIT
The covariance matrix of the estimated diffusivities (Peebles, 1993;Rust, 2001) is given by: Eq. 7 predicts that the cov( Y ) depends on the b- factor, the I(0) and I(k) SNR levels.The log-linear functional relationship also couples the cov( d ) to the estimated unknown values d (Anderson, 2001).In the case where the noise is independent, we can write The structure of the encoding matrix is important in performing the inversion procedure as well as defining figures-of-merit, such as condition number, energy, variance, etc., to optimize the spatial distribution of the encoding vectors.We define the total tensor variance as and the encoding matrix condition number as Borrowing from the analogy with the electrostatic repulsion problem of N e identical unit charged pairs we compute the interaction energy, E, of N e identical pair charges positioned at £ on the surface of a unit sphere (Conturo et al., 1995;Jones et al., 1999;Hasan et al., 2001).
Note that since the unit charges are constrained to the surface of a unit sphere, ||£(k)||=1, we obtain the value for the relative distance matrix element R kr = || £(k) -£(r) || = 2sin(γ kr /2), where γ kr is the angle between the two unit vectors £(k) and £(r).To attain an equilibrium state or configuration with minimum energy, the angle between any two encoding vertices of some encoding set has to be as large as possible.

REVIEW OF OPTIMIZATION AND SPATIAL DISTRIBUTION METHODS
In general, the distribution of tensor shapes and orientations are not known a priori for a given tissue of interest.Consequently, it is normally assumed that the optimum encoding set will have uniformly distributed (UD) directions over a sphere.Several optimization methods have been reported to distribute the encoding N e vertices on the diffusion gradient encoding sphere (Hasan et al., 2001a).These methods can be classified into three main approaches: heuristic, natural (polyhedral), and numerical.The heuristic sets use combinations of the cubic vertices; these sets use the sides (orthogonal: x or y or z), corners (tetrahedral: x and y and z) (Conturo et al., 1995;1996) and edge bisectors of the cube (x and y, x and z, y and z) (Basser and Pierpaoli, 1998).Heuristic sets can provide N e = 6, 7, 10, or 13 directions; these sets are not generally uniformly distributed in space.The exception is a combination of orthogonal and tetrahedral encoding with N e = 7 directions (Shimony et al., 1999;Hasan et al., 2001a).The analytical solution for uniformly distributing unit vectors on the sphere is an unsolved problem in mathematics for an arbitrary number N e required for tensor encoding (see references in Hasan et al., 2001a).The solution exists only for the icosahedron (20-faces, 12-vertices, 30 edges) and its dual the dodecahedron (12-faces, 20-vertices).The icosahedral vertex, face and edge centroids can thus provide 6,10 and 15 directions.Face triangularization can also provide more discrete sets with N e = 5n 2 + 1, n = 1,2,3,.…(Hasan et al., 2001a).
For an arbitrary number of encoding directions N e > 6, the spatial distribution problem can be solved by numerically optimizing and ordering tiles of equal solid angles on a unit sphere (Hasan et al., 2001a).The optimization cost function can be minimized using the constrained least squares method subject to the scalar figure of merits discussed earlier.A comparison of numerical optimization criteria, including minimum condition number (MC) (Skare et al., 2000), minimum tensor variance (MV) (Papadikas et al., 1999), minimum force (MF) (Jones et al., 1999) and minimum energy (ME), has been investigated previously (Hasan et al., 2001a).The minimum force and energy approaches establish an analogy with the electrostatic interaction problem of N e identically charged pair dipoles constrained on the surface of a sphere.The minimum energy approach proved to be the fastest to converge, and it was able to give the natural icosahedral set for N e = 6 (Hasan et al., 2001a).The numerical sets generated have to be further balanced to minimize the linear cross term couplings.The balancing adjusts the sign of each encoding diffusion gradient vector and seeks to minimize the magnitude of the net sum

DIFFUSION TENSOR IMAGE ACQUISITION EXPERIMENTS: PARAMETERS AND SUBJECTS
Full human brain DT-MRI measurements were acquired with General Electric SIGNA (Waukesha, WI) 1.5 Tesla and 3 Tesla MRI scanners.The data were acquired in accordance with IRB approved protocols, and signed consent forms were obtained.
The imaging and diffusion tensor parameters were T R /T E /δ/∆ = 9000/65/16/21 ms with G max = 40 mT/m at 3Τ.The sequence was a spin-echo prepared EPI (see Fig. 2) that utilized a ramp and partial k y sampling.

IMAGE AND DTI PROCESSING
The DW magnitude images (after complex kspace 2d-fft construction) were averaged by the scanner.The data were then corrected for motion and distortion using the Automated Image Registration Package (AIR) (Wood et al., 1998).The registration template, for each diffusion weighted image, was selected based on a diffusion weighted image that shows the least distortion sensitivity to the diffusion along the phase encoding direction.The diffusion tensor decoding, anisotropy computation and visualization of the raw and diffusion tensor model were constructed according to the formalism in the theory using a custom developed IDL-MATLAB DTI toolbox.

BOOTSTRAP EVALUATION EXPERIMENT
A single section DTI data set was acquired for Bootstrap statistical analysis.Data were acquired by repeating the measurements for different encoding schemes according to N T = N ref + N d N e .The Bootstrap is a computer-intensive, non-parametric statistical method that is done by randomly resampling the original data pool (with replacement) for the Bootstrap number of trials, N B > 1000 (Efron and Tibshirani, 1993).The advantage of the Bootstrap analysis over Monte Carlo simulations is that the latter assumes that the measurement noise distribution is normally distributed; whereas, the former works on experimentally measured data with no a priori assumptions on the noise distribution (Pajevic and Basser, 1999).Here the comparison was done at equal imaging time (Hasan et al., 2000) for different tensor encoding schemes and different signal-to-noise levels.The Bootstrap method is used to provide the distribution of any DTI metric on a single voxel basis as well as the variance of the DTI scalar and vectorial maps.The Bootstrap iteration was performed efficiently, using the analytical diagonalization approach (Hasan et al., 2001b).To reduce memory usage, the Bootstrap variance in some scalar function f is estimated using the real-time expression: , where E(f) refers to the expected value <f>.

THE RELATIONSHIP BETWEEN TENSOR ENCODING SCHEMES
To study the relationship between different encoding schemes, the analytical encoding method suggested by Hasan et al. (2001a) was used.The method constructs the simplest encoding vector matrix for N e = 6 from cyclic permutations of the normalized vector [ ] can be constructed from the rows of the generating encoding set where x is a variable that allows analytical computation.The H matrix is given by (Eq. 6) This analytical and general encoding matrix gives three familiar cases.For x=1 we obtain the cube edge bisectors described by Basser and Pierpaoli (1998); this encoding scheme is referred to here as the double oblique gradient scheme (ODG).The principal icosahedral encoding scheme is obtained for , with condition number = 10 2 .
The minimum condition set (Skare et al., 2000) is obtained with x = mc 5 ± 21 x = 2 , with condition number = 7 2 .The icosahedral family is the only uniformly distributed set with N e = 6.The minimum energy spatial distribution algorithm will also provide uniformly distributed sets for arbitrary N e including N e = 6.In Fig. 3, we plot four scalar functions versus x for the minimum energy, variance, condition number and minimum angle between any two pairs constrained to the unit encoding sphere.Notice that 1) the minimum energy and minimum variance scheme is provided by the icosahedral set, 2) the encoding sets at x = 1 and x = x mc have the same tensor variance and are biased as illustrated by the minimum angle plot, 3) the minimum condition criterion does not yield uniformly distributed encoding directions, nor does the commonly used (ODG) encoding scheme for x = 1, and 4) the icosahedral solution at x = τ is expected to be more stable to small perturbations or errors in the design diffusion gradient directions.

ENCODING SCHEME MINIMUM ANGLE TEST
Fig. 4 compares the minimum angle computed from the vector dot product of all pairs in a family of isotropic icosahedral encoding sets constrained to the unit encoding sphere, as well as other heuristic approaches.The icosahedral family is used to provide the theoretical values.Note that the observations made in the previous section are still valid.In addition, the minimum energy encoding sets (ME) are indistinguishable from the icosahedral for N e = 6, 21, and 46.The Conventional-7 (x, y, z, xy, xz, yz, xyz), Jones30, and Force10 encoding schemes are described by Jones (1999).The decahedral, DSM-6 and DSM-10 are minimum condition sets, as generated by Skare (2000).The ORTH scheme uses x, y, z, xy, xz, and yz cubic vertices.The scaled tetraorth scheme is used in Shimony (2000) and Nye (1957).The minimum angle is a sensitive test to characterize the bias and sensitivity of tensor encoding sets.The angle has to be as large as possible when an optimal UD configuration is attained.Fig. 5 shows a triangularization of the minimum energy and icosahedral sets.Notice that the ME procedure can generate isotropically distributed encoding sets with arbitrary N e ; whereas, the icosahedral based sets are limited to the discrete values of N e = 5n 2 + 1, n = 1,2,3,4….(Hasan et al., 2001a).

APPLICATION OF UD SETS AND VALIDATION OF ISOTROPIC PRINCIPAL VECTOR IN THE HUMAN BRAIN
In this section we seek to study the distribution of the principal eigenvector in a normal human brain DT-MRI data set.A full brain DT-MRI set was obtained from a healthy male volunteer.A balanced ME23 encoding scheme, acquired at 1.5T, was used.The relative angle between the principal eigenvector and each of the 23 directions was computed.In addition, the same procedure was applied on the vertex directions obtained from the triangularization of the principal icosahedron with a total number of vertices almost matching the number of voxels in the normal human brain, 10n 2 + 2 for n = 135.Fig. 6 shows that the normal subject's brain principal vector with FA > 0.2 is approximately distributed in all directions.Notice the resemblance with the predicted theoretical relationship for a uniformly distributed vector field: dN(α)/dα ~ sin(α) (Mardia and Jupp, 2000).Fig. 6.A test of the isotropic distribution of brain fibers hypothesis on a normal subject (A,B).The expected angular distribution for the same tensor encoding scheme (ME23 @1.5T) and a uniformly distributed vector distribution on a unit sphere (C,D).

BOOTSTRAP VALIDATION OF ISOTROPIC TENSOR ENCODING ADVANTAGE N E _ 6
Fig. 8 illustrates the application of the multidimensional Bootstrap statistical method in the study of tensor metric (e.g., ADC, FA) dependence on the DTI parameter space.The figure shows the expectation values of ADC, FA and their corresponding errors for equivalent imaging times.The standard error and total variance maps of the ME6 (ICOSA6, N e = 6), ME12 and ME24 are noticeably smaller than the heuristic ORTH, ODG and ME48 schemes.Note that an optimal isotropic tensor encoding scheme may have 6 < N e < 24 for the quoted image parameters.
In Fig. 9, the single voxel bootstrap distributions of ADC and FA are shown for two low and high anisotropy representative voxels at two different SNR levels for the icosahedral encoding scheme.The lowest SNR level corresponds to an SNR value of about 40, N T = N ref + N d N e = 55.Note that 1) the distribution is approximately normally distributed, 2) the anisotropy is generally overestimated at low SNR and low anisotropy and 3) since the three eigenvalue images are combined to obtain the ADC we expect the ADC to be less biased to noise than the tensor fractional anisotropy index, FA.In Fig. 9, the single voxel angular dispersion of the principal eigenvector, e 1 , is shown for different encoding schemes, low and high anisotropy voxels and two SNR levels.The reference vector, e ref , is computed from the highest SNR possible for each encoding set.The analysis is done as a function of the tensor encoding scheme, SNR and anisotropy.Constant total acquisition times, T E , T R , b, I(0), etc. were maintained.The encoding sets compared have N e = 6, 12, 24 and 48.

DISCUSSION AND CONCLUSIONS
Diffusion tensor imaging is a sensitive modality that is very promising in radiological and psychiatric applications (Le Bihan et al., 2001;Melhem et al., 2002).A full brain DT-MRI exam can be obtained in less than 10 min using the single shot echo-planar method (see Fig. 10).An efficient, accurate and specific DTI experiment seeks to partition the total imaging time into the different parameter space variables to obtain high quality images with minimal artifacts.The data sufficiency condition for tensor decoding requires at least seven independent measurements.Without a priori knowledge about the full brain fiber distribution, the encoding directions are required to be uniformly distributed over the encoding unit sphere in order to reduce the encoding directional bias.This work has presented the fundamental aspects of the MRI based single diffusion tensor encoding and estimation theory.This work discussed some criteria (figures-of-merit) for the evaluation and optimization of single tensor encoding directions assuming a two-diffusion weighting scheme.The icosahedral and minimum energy sets have been found to be optimal for single tensor encoding when the tensor distribution is uniformly distributed.This work also related the minimum condition set with N e = 6 to a commonly used tensor encoding scheme (ODG).The minimum condition criterion does not lead to uniformly distributed tensor encoding sets.
The Bootstrap non-parametric statistical method was used at equal imaging time to study the tensor encoding advantage with more than the least required optimal six icosahedral vertices.The Bootstrap method is less restrictive than the Monte Carlo approach and is more useful in obtaining the distribution of any DT-MRI related map.This work shows that minimum energy sets with 6 ^ N e ^ 24 are able to reduce the total tensor variance and the variance in the estimated ADC and FA maps.The observations made on the behavior of the heuristic ORTH and ODG and minimum variance encoding schemes were confirmed with the Bootstrap method and are consistent with an earlier study using Monte Carlo simulations (Hasan et al., 2001a).The Bootstrap evaluation results in this work, although limited to one section and one subject, are useful in providing feedback on the protocol design.The Bootstrap method presents a framework to test the effect of other important DTI parameters, that have not been addressed directly in this paper, such as N ref and the b-factor on the accuracy of the estimated tensor (Xing et al., 1997;Cercignani and Horsfield, 1999;Jones et al., 1999).The SNR-anisotropyencoding-angular dispersion results are useful also in optimizing and regularizing fiber tracking algorithms.Tensor encoding with uniformly distributed sets seems to be more justified for whole brain studies, as it has been demonstrated that the human brain fiber orientations can sweep all possible orientations.In addition, tensor encoding with the sufficient icosahedral set with N e = 6 and optimal sets with 6 < N e < 24 seems to reduce the dependency on systematic errors per encoding direction.For full brain DT-MRI data acquisition at 1.5 T, this work recommends using optimally distributed tensor encoding sets with 6 ^ N e < 24.Further Bootstrap studies with more subjects are needed to evaluate the theoretical predictions and the experimental observations presented in this work at 3 T.

Fig. 1 .
Fig. 1.The diffusion tensor ellipsoid, the local tissue and physical gradient coordinate systems.Note that the positive definite and symmetric diffusion tensor can be diagonalized (De i =λ i e i ) into three eigenvalues and three eigenvectors.

Fig. 2 .
Fig. 2. Schematic of the spin-echo bipolar gradient diffusion tensor imaging pulse sequence and the DT-MRI parameter space.
= D yx , D xz = D zx and D yz = D zy .
d , can be used to estimate the magnitude, anisotropy and orientation of the local diffusion tensor.Define the three rotationally independent principal tensor invariants, I 1 = Trace(D), I 3 = Det(D), I 2 = I 3 Trace(D -1 ) and the Frobenius norm ||D|| which is the root of The balanced numerical sets can then be called UD tensor encoding isotropic sets.

Fig. 3 .
Fig. 3. Scalar functions obtained from the analytical encoding matrix.The energy, variance, condition number and minimum angle as function of the variable, x..

Fig. 4 .
Fig. 4. A summary of the minimum angle subtended between any two pairs of different heuristic, minimum energy and icosahedral tensor encoding schemes.The minimum angle should be as large as possible to warrant uniform distribution.