PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
IEEE Trans Med Imaging. Author manuscript; available in PMC 2010 July 28.
Published in final edited form as:
PMCID: PMC2896997
NIHMSID: NIHMS187341

FRATS: Functional Regression Analysis of DTI Tract Statistics

Hongtu Zhu*
Department of Biostatistics, and Biomedical Research Imaging Center, University of North Carolina at Chapel Hill, Chapel Hill, NC 27599, USA
Martin Styner
Department of Psychiatry, University of North Carolina at Chapel Hill, Chapel Hill, NC 27599, USA
Niansheng Tang
Department of Statistics, Yunnan University, Kunming, Yunnan, P.R.China
Zhexing Liu
Department of Psychiatry, University of North Carolina at Chapel Hill, Chapel Hill, NC 27599, USA
Weili Lin
Department of Radiology, University of North Carolina at Chapel Hill, Chapel Hill, NC 27599, USA

Abstract

Diffusion tensor imaging (DTI) provides important information on the structure of white matter fiber bundles as well as detailed tissue properties along these fiber bundles in vivo. This paper presents a functional regression framework, called FRATS, for the analysis of multiple diffusion properties along fiber bundle as functions in an infinite dimensional space and their association with a set of covariates of interest, such as age, diagnostic status and gender, in real applications. The functional regression framework consists of four integrated components: the local polynomial kernel method for smoothing multiple diffusion properties along individual fiber bundles, a functional linear model for characterizing the association between fiber bundle diffusion properties and a set of covariates, a global test statistic for testing hypotheses of interest, and a resampling method for approximating the p-value of the global test statistic. The proposed methodology is applied to characterizing the development of five diffusion properties including fractional anisotropy, mean diffusivity, and the three eigenvalues of diffusion tensor along the splenium of the corpus callosum tract and the right internal capsule tract in a clinical study of neurodevelopment. Significant age and gestational age effects on the five diffusion properties were found in both tracts. The resulting analysis pipeline can be used for understanding normal brain development, the neural bases of neuropsychiatric disorders, and the joint effects of environmental and genetic factors on white matter fiber bundles.

Keywords: Diffusion tensor imaging, fiber bundle, functional regression, global test statistic, registration

I. Introduction

Water molecules in the human brain diffuse preferentially along the major axis of white matter fiber bundles. Thus, Diffusion Tensor Imaging (DTI), which can track the effective diffusion of water in the human brain in vivo, is used to accurately map in vivo the structure and orientation of fiber tracts in the white matter of the brain [1], [2]. The directional dependence of water diffusion in each voxel can be characterized by a matrix, called a diffusion tensor (DT), and the degree of diffusivity can be quantified by using the three eigenvalue-eigenvector pairs of DT and its related parameters, such as fractional anistropy (FA) [3], [4], [5], [6]. There has been a wealth of neuroimaging studies using a set of water diffusion related parameters including FA, mean diffusivity (MD), and the three eigenvalues of diffusion tensor as a marker for white matter tract maturation and integrity found in [7], [8], [9], [10], and many others. Analyzing these water diffusion related parameters across subjects requires specific statistical methods for the group analysis of diffusion imaging data.

In the current literature, three major approaches to the group analysis of diffusion imaging data are region-of-interest (ROI) analysis, voxel based analysis, and fiber tract based analysis [11], [12], [13]. The region-of-interest (ROI) method used in some neuroimaging studies [14], [15] primarily averages diffusion properties in some manually drawn ROIs for each subject and then creates a single statistic per ROI [13]. Some major drawbacks of ROI analysis are the difficulty in identifying meaningful ROIs, particularly the long curved structures common in fiber tracts, the instability of statistical results obtained from ROI analysis, and the partial volume effect in relative large ROIs. The ROI method is also based on a stringent assumption that diffusion properties in all voxels of the same ROI are essentially homogeneous, which is largely false for DTI.

Voxel based analysis is used more commonly than ROI analysis in neuroimaging studies [16], [17], [18], [19]. The first step involves fitting a statistical model to the smoothed and registered diffusion properties imaging data from multiple subjects at each voxel to generate a parametric map of test statistics (or p-values). The second step is to correct for the multiple comparisons across the many voxels of the imaging volume [20], [21], [22]. It suffers from issues of alignment quality and the arbitrary choice of smoothing extent [23], [20], [11], [24]. As shown in Jones et al. (2005), the final results of voxel-based analysis can strongly depend on the amount of smoothing in the smoothed diffusion imaging data.

With the drawbacks mentioned of ROI and voxel based analysis, there is a growing interest in the DTI literature in developing fiber tract based analysis of diffusion properties, such as eigenvalues and fractional anisotropy (FA) values [11], [12], [25], [26]. For instance, Smith and his coauthors develop a tract-based spatial statistics framework for constructing local diffusion properties along the white matter skeleton followed by performing pointwise hypothesis tests on the skeleton [11]. Yushkevich and coauthors propose a model-based framework for the statistical analysis of diffusion properties on the medial manifolds of fiber tracts followed by testing pointwise hypotheses on the medial manifolds [25]. Similar to the work proposed in this paper, Goodlett and coauthors propose to use a functional data analysis method to compare a univariate diffusion property, such as fractional anisotropy, across two (or more) populations for a single hypothesis test per tract [26]. Their method is limited to a univariate diffusion property and cannot control for other covariates of interest, such as age, gender and behavioral variables. Moreover, the permutation test and the Hotelling T2 statistic used in [26] may be statistically over-sensitive, because only the data in the reduced principal component analysis space are permuted. Statistically, such a procedure, which ignores substantial noise in the original data, can lead to misleading results.

What these three methods do not account for is the comparison of fiber bundle diffusion properties across groups and the development of fiber bundle diffusion properties along time, while controlling for other covariates of interest, such as gender [16], [14], [11], [17], [18], [19]. Making these comparisons requires a regression modeling framework for the analysis of fiber bundle diffusion properties and a set of covariates of interest, such as age, diagnostic status and gender. This paper presents a functional regression analysis of DTI tract statistics, called FRATS, for modeling the relationship between fiber bundle diffusion properties and covariates of interest. Compared with [26] and other existing literature, there are four methodological contributions in this paper. First, the local polynomial kernel method is used to regularize multiple diffusion properties along individual fiber bundles. Second, a functional linear model is developed to characterize the association between fiber bundle diffusion properties and any covariate of interest. Third, a global test statistic is proposed for testing hypothesis of interest. Fourth, a resampling method is developed for estimating the p-value of the global test statistic. The proposed methodology is used to characterize the development of five diffusion properties along the splenium and right internal capsule tracts in a clinical study of neurodevelopment in the University of North Carolina at Chapel Hill.

II. Methods

This paper introduces a set of statistical tools for the group analysis of fiber bundle diffusion properties while controlling for a set of covariates of interest. Comparing diffusion properties in populations of DTIs requires imaging methods for establishing correspondence among regions of anatomy in DTIs. For this purpose, we use the DTI atlas building methods proposed in [26], which are included here for completion. The focus of this paper is to develop the functional regression model pipeline, called FRATS, for assessing the association between fiber bundle diffusion properties and covariates of interest.

A. DTI Atlas Building

Before applying FRATS, a quality control pipeline is applied to each of the subjects' diffusion weighted images (DWIs) to check if the DWIs contain large slice brightness artifacts, intra-gradient Venetian blind artifacts, and motion artifacts using DTIPrep (www.ia.unc.edu/dev). Skull stripping is then performed automatically via a brain mask computed from both the non-diffusion weighted baseline image and isotropic diffusion weighted image using a free software itkEMS (www.ia.unc.edu/dev/download/itkems). Masking out the non-brain tissues stabilizes the subsequent DTI estimation and registration step. Diffusion tensors are estimated for each subject from a series of cleaned diffusion weighted images using a weighted least square tensor estimation [27], [2].

After DTI estimation, registration occurs and its accuracy is crucial for identifying meaningful group differences. We use a nonlinear fluid deformation based high-dimensional, unbiased atlas computation method to carry out a large deformation non-linear registration [28]. The atlas building procedure is initialized by an affine registration and followed by a nonlinear registration of a feature image which is sensitive to the geometry of white matter. Further reference of the DTI atlas building procedure has been described in [26]. Then, we warp each of the tensor images into the unbiased space and reorient them using the finite strain approximation proposed by [29] to get the aligned DTIs. After averaging all the registered DTIs, a study specific unbiased DTI atlas is created. The tensor warping and averaging are performed using a novel Log-Euclidean geometry [30].

Major fiber bundles are tracked in the atlas space within 3D Slicer (www.slicer.org). First, ROIs are drawn with the Editor module and then fibers are tracked using a stream line method with the Tractography module. Usually, the fibers directly generated in this way contain some fibers not belonging to the bundle that we are interested in due to either the noise in DTI data or not proper ROI delineation. We use FiberViewer to visualize and clean the fibers. Tractography of atlas fiber bundles are more reliable than individual ones because of the improved signal-to-noise ratio of the atlas.

With the fiber bundles in atlas space, each subject's DTI data is transformed into the atlas space. When transforming the diffusion properties from native space to the atlas space, tri-linear interpolation is used to interpolate diffusion properties from the native space. Finally, we get a set of individual tracts with corresponding geometry but varying diffusion properties. Within FiberViewer, a cutting plane is first placed perpendicular to the fibers at a location where the fibers are organized neatly to define the origin of arc length. The diffusion properties along fiber tracts are then measured and plotted as functions of arc length. Using FiberViewer, these plots from each subject are generated. This allows us to conduct a group comparison analysis of the along-bundle tensor properties. Except for 3D Slicer, all DTI tools and itkEMS are part of the open source UNC NeuroLib software repository (https://www.ia.unc.edu/dev/).

B. FRATS: Nonparametric Model

After spatial normalization of tensor images, we propose to use a functional regression model to analyze diffusion properties along the same fiber bundle from multiple subjects. FRATS consists of four major components: a regularization of individual tracts, a functional linear model, a global test statistic for hypothesis testing, and a resampling method for estimating the p-value of the global test statistic. A schematic overview of FRATS is in Figure 1. We describe each of these components in detail below.

Fig. 1
A schematic overview of FRATS: a nonparametric model for regularizing individual tracts, a functional linear model, a global test statistic for hypothesis testing, and a resampling method for estimating the p-value of the global test statistic.

To regularize individual tracts, we develop a nonparametric model for multiple diffusion properties (e.g., FA) along the same fiber bundle as a smooth function of arc length from each subject as follows. For the i-th subject, we consider a m × 1 vector of diffusion properties, denoted by yi,j = (yij,1, (...), yij,m)T, and its associated arc length sj for the j-th location grid point on the fiber bundle for j = 1, (...),L0 and i = 1, (...), n, where L0 and n denote the numbers of grid points and subjects, respectively. The nonparametric model is given by

yi,j=fi(sj)+ϵi,j,
(1)

where fi(s) = (fi,1(s), (...), fi,m(s))T is an m × 1 vector of continuous functions with second-order continuous derivative, E[yi,j|fi(sj)] = fi(sj), and Cov[yi,j|fi(sj)] = Σ(sj). Our aim is to reconstruct the continuous function vector fi(s) based on observed fiber bundle diffusion properties.

To simultaneously construct all individual function vectors fi(s), we develop an adaptive local polynomial kernel smoothing technique [31], [32], [33], [34], [35], [36]. Specifically, using Taylor's expansion, we can expand fi(sj) at s to obtain

fi(sj)=fi(s)+f.i(s)(sjs)=Aizj,
(2)

where zj = (1, sjs)T and Ai=[fi(s)f.i(s)] is an m × 2 matrix with f.i(s)=dfi(s)ds. We develop an algorithm to estimate Ai as follows.

Step (1.1). Step (1.1) is to construct an initial estimate of fi(s) for each i. Let ai;k be the k-th row of Ai and K(·) be a kernel function, such as the Gaussian, Epanechnikov, or uniform kernels [31]. Without loss of generality, we assume that the kernel function K(·) is always a symmetric probability density function. Throughout the paper, we take K(t) = 0.75(1 − t2)1(|t| ≤ 1), where 1(·) is an indicator function. For each k and a fixed bandwidth hk, we estimate ai;k by minimizing an objective function given by

j=1L0(yij,kai;kTzj)2Khk(sjs),
(3)

where Khk(·) = K(·/hk)hk is a rescaled kernel function. With some calculation, it can be shown that

a^i;k=[j=1L0zjKhk(sjs)zjT]1j=1L0Khk(sjs)zjyij,k.
(4)

Let e1,2 = (1, 0)T. Then, since ai;k=(fi,k(s),f.i,k(s))T, fi,k(s) can be estimated by

f^i,k(s)=e1,2Ta^i;k=j=1L0K~hk0(sjs,s)yij,k,
(5)

where K~hk0(,) are the empirical equivalent kernels [31]. Thus, f^i,k=(f^i,k(s1),,f^i,k(sL0))T=Si,kyi:,k, where yi:,k = (yi1,k, (...), yiL0,k)T and Si,k is the smoother matrix for the k-th measurement of the i-th subject. For each k, we pool the data from all n subjects and select the optimal bandwidth hk, denoted by h^k,opt(1), by minimizing the generalized cross-validation score given by

GCVk(hk)=(n)1i=1nj=1L0[yij,kf^i,k(sj)]21L01tr(Si,k).
(6)

Based on the optimal h^k,opt(1), we can estimate f^i,k(s) for all i.

Step (1.2). Step (1.2) is to construct an estimator of the covariance matrix Σ(sj) at sj. Specifically, we consider the unbiased sample covariance matrix at sj given by

Σ^(sj)=(nm)1i=1n[yi,jf^i(sj)]2,
(7)

where f^i(s)=(f^i,1(S),,f^i,m(s))T and a[multiply sign in circle]2 = aaT for any vector a. It can be shown that Σ^(sj) converges to the true Σ(sj) in probability as both n and L0 go to infinity.

Step (1.3). Step (1.3) is to compute an adaptive estimator of fi(s) for all i using the initial results from Steps (1.1) and (1.2). For all k and a fixed bandwidth h, we estimate Ai by minimizing an objective function given by

j=1L0(yijAizj)TΣ^(sj)1(yijAizj)Kh(sjs).
(8)

Let Zj=block diagonal(zjT,zjT,,zjT) be an m × 2m matrix and Bi=(ai;1T,,ai;mT)T. It can be shown that

B^i=[j=1L0Kh(sjs)ZjTΣ^(sj)1Zj]1j=1L0Kh(sjs)ZjTΣ^(sj)1yij,
(9)

which leads to a new estimator of fi,k(s), denoted by f^i,k(s)sec for each i and k. Let S~i,k be the smoother matrix for the k-th measurement of the i-th subject such that f^i,ksec=(f^i,k(s1)sec,,f^i,k(sL0)sec)T=S~i,kyi:,k. We pool the data from all n subjects and m measurements and select the optimal bandwidth h, denoted by h^opt by minimizing the generalized cross-validation score, denoted by GCV(h), given by

i=1nj=1L0[yijf^i(sj)sec]TΣ^(sj)1[yijf^i(sj)sec]nn(mL0)1k=1mtr(S~i,k),
(10)

where f^i(s)sec=(f^i,1(s)sec,,f^i,m(s)sec)T. Based on the optimal h^opt, we can estimate f^i,k(s)sec for all i and k. Similar to the arguments in [35], it can be shown that when Σ(s) varies across s, f^i,k(s)sec based on the optimal h^opt is more accurate than f^i,k(s) obtained from Step (1.1).

Step (1.4). Step (1.4) is to estimate the mean function f(s) and the covariance function Γ(s, t) of fi(s). Specifically, following [36], [34], we can estimate f(s) and Γ(s, t) by using their empirical counterparts of the estimated f^i(s)sec as

f^(s)=n1i=1nf^i(s)secand
(11)

Γ^(s,t)=n1i=1n[f^i(s)secf^(s)][f^i(t)secf^(t)]T
(12)

The diagonal of Γ(s, s) reflects the variance of fi(s) at the location s.

C. FRATS: Functional Linear Model

We develop a functional linear model to characterize the relationship between all diffusion properties along fiber tracts and a set of covariates of interest, such as age, group, and gender. We assume that

fi(s)=B(s)xi+ηi(s),i=1,,n,
(13)

where B(s) is a m × p matrix of functions of s, xi is a p × 1 vector of covariates of interest, and ηi(s) satisfies Ei(s)|xi] = 0 and Cov[ηi(s), ηi(t)|xi] = Γη(s, t). B(s) characterizes the association between fiber bundle diffusion properties and the covariates of interest xi.

As an illustration, in our clinical study on early brain development, we are interested in studying the evolution of the three eigenvalues λi of diffusion tensor (λ1 ≥ λ2 ≥ λ3) along two selected fiber tracts in 128 healthy pediatric subjects (Figs. 2 (a) and 5 (a)). See clinical data and results section for details. We consider a functional linear model of three eigenvalues along a specific tract as follows:

fi,1(s)=(β11(s),β12(s),β13(s),β14(s))xi+ηi1(s),fi,2(s)=(β21(s),β22(s),β23(s),β24(s))xi+ηi2(s),fi,3(s)=(β31(s),β32(s),β33(s),β34(s))xi+ηi3(s),
(14)

where fi,k(·) equals the smoothed λk curve from the i-th subject for k = 1, 2, 3, and xi = (1, gi, Gagei, agei)T, in which gi and Gagei denote gender and gestational age, respectively. Moreover, gestational age is the age of a newborn infant, both gestational age and age are continuous covariates, and gender is a discrete covariate. In this case, m = 3, B(s) = (βjk(s)) is a 3 × 4 matrix, and xi = (1, gi, Gagei, agei)T.

Fig. 2
Splenium tract and diffusion properties along the splenium tract: (a) the splenium tract extracted from the tensor atlas with color representing mean FA value; (b) FA; (c) MD; (d) λ1; (e) λ2; (f) λ3. The diffusion properties in ...
Fig. 5
Results from the analysis of FA and MD on the splenium tract: the −log10(p) values of test statistics Sn(sj) for testing gender effect in panel (a), gestational age effect in panel (b), and age effect in panel (c) on FA; the −log10(p) ...

We develop an estimation algorithm to estimate B(s) and Γη(s, t) as follows.

Step (2.1). Step (2.1) is to estimate B(s). Let Bk(s) be the kth row of B(s). Then, based on the estimated function vectors f^i(s), we calculate the least-squares estimator of B(s), denoted by B^(s), by minimizing an objective function given by

i=1n[f^i(s)secB(s)xi]T[f^i(s)secB(s)xi].
(15)

Specifically, the least-squares estimator of Bk(s), denoted by B^k(s), is given by

B^k(s)T=(i=1nxi2)1i=1nxif^i,k(s)sec
(16)

for k = 1, (...), m.

Step (2.2). Step (2.2) is to estimate Γη(s, t). Let η^i(s)=f^i(s)secB^(s)xi. Then, the covariance matrix Γη(s, t) can be estimated by

Γ^η(s,t)=(nm)1i=1nη^i(s)η^i(t)T.
(17)

The covariance matrix Γη(s, t) characterizes the variation of ηi(s), which is different from Γ(s, t).

To carry out statistical inference on B(s), we need some asymptotic results for B^(s). Let vec(B(s)) = (B1(s)T, (...),Bm(s)T)T. Theoretically, following the arguments in [36], it can be shown that under some regularity conditions,

n[vec(B^(s))vec(B(s))+O(h^opt2)]LGP(0,ΓB),
(18)

where →L denotes the convergence in distribution and GP(0, ΓB) is a Gaussian process with mean zero and covariance function ΓB(s,t)=Γη(s,t)ΩX1, in which [multiply sign in circle] denotes the Kronecker product of two matrices and ΩX=limnn1i=1nxi2. The asymptotic result in (18) ensures the validity of constructing the global test statistic detailed below.

D. FRATS: Global Test Statistic

We develop a global test statistic to test linear hypotheses of B(s) in order to answer various scientific questions involving a comparison of fiber bundle diffusion properties along fiber bundles across two (or more) diagnostic groups and the development of fiber bundle diffusion properties along time. We can formulate these questions as linear hypotheses of B(s) as follows:

H0:Cvec(B(s))=b0(s)for allsvs.
(19)

H1:Cvec(B(s))b0(s),
(20)

where C is a r × mp matrix of full row rank and b0(s) is a given r × 1 vector of functions.

As an illustration, we are interested in testing the age effect on the evolution of the three eigenvalues of the diffusion tenors along the two selected fiber tracts in 128 healthy pediatric subjects in our clinical study on early brain development. Statistically, for model (1), testing the age effect can be formularized as follows:

H0:β14(s)=β24(s)=β34(s)=0for allsvs.

H1:β14(s)0orβ24(s)0orβ34(s)0.

In this case, we have

C=(000100000000000000010000000000000001)

and b0(s) [equivalent] (0, 0, 0)T for all s. The use of multiple diffusion properties in model (13) allows us to compare different functions in B(s) associated with different diffusion properties. For instance, suppose that we want to test

H0:β14(s)=β24(s)for allsvs.H1β14(s)β24(s).

In this case, we have

C=(000100010000)

and b0(s) [equivalent] (0) for all s.

We test the null hypothesis H0 : Cvec(B(s)) = b0(s) using a global test statistic Sn defined by

Sn=n0F0d(s)T[C(Γ^η(s,s)Ω^X1)CT]1d(s)ds,
(21)

where Ω^X=n1i=1nxi2,d(s)=Cvec(B^(s))b0(s), and F0 is the whole arc length of a specific fiber bundle. In order to use Sn as a test statistic, we need an asymptotic result. Specifically, similar to the arguments in [36], we can show that under some conditions and H0, nd(s) and Sn converge weakly to N(0,C(Γη(s,s)ΩX1)CT) and a weighted χ2 distribution, respectively, as n → ∞.

In addition, at a given grid point sj on a specific tract, we can also test the local null hypothesis H0(sj) : Cvec(B(sj)) = b0(sj) using a local test statistic Sn(sj) defined by

Sn(sj)=nd(sj)T[C(Γ^η(sj,sj)Ω^X1)CT]1d(sj).
(22)

Under some conditions and H0(sj), nd(sj) and Sn(sj) converge weakly to N(0,C(Γη(sj,sj)ΩX1)CT) and a weighted χ2 distribution with r degrees of freedom, respectively, as n → ∞.

E. FRATS: Resampling Method

We develop a resampling method (or wild bootstrap method) to approximate the p-value of Sn [37], [38]. The resampling method has four key steps as follows.

Step (3.1): Fit the functional linear model f^i(s)=B(s)xi+ηi(s), i = 1, (...), n, under the null hypothesis H0, which yields B^(s) and η^i(s)=f^i(s)B^(s)xi.

Step (3.2): Generate a random sample τi(g) from a N(0, 1) generator for i = 1, (...), n and then construct f^i(s)(g)=B^(s)xi+τi(g)η^i(s). Then, based on f^i(s)(g), we calculate

B^k(s)(g)T=(i=1nxi2)1i=1nxif^i,k(s)(g)
(23)

for k = 1, (...), m, where B^k(s)(g)T and f^i,k(s)(g) are, respectively, the kth row of B^(s)(g) and f^i(s)(g). Finally, let d(s)(g)=Cvec(B^(s)(g))b0(s), we compute Sn(g)=n0F0d(s)(g)T[C(Γ^η(s,s)Ω^X1)CT]1d(s)(g)ds and Sn(sj)(g)=nd(sj)(g)T[C(Γ^η(sj,sj)Ω^X1)CT]1d(sj)(g) for j = 1, (...), L0.

Step (3.3): Repeat Step (3.2) G times to obtain {Sn,max(g)=max1jLoSn(sj)(g):g=1,,G} and calculate

p(sj)=G1g=1G1(Sn,max(g)Sn(sj))

for each sj. The p(sj) is the corrected p-value at the location sj.

Step (3.4): Repeat Step (3.2) G times to obtain {Sn(g):g=1,,G} and calculate

p=G1g=1G1(Sn(g)Sn).

If p is smaller than a pre-specified value α, say 0.05, then we reject the null hypothesis H0.

We note several advantages of using the resampling method in the above test procedure. Computationally, the above procedure only requires the repeated calculation of B^k(g) and Sn(g). Thus, because computing B^k(g) is very fast, the proposed test procedure is computationally efficient. The proposed resampling method also performs better than parametric bootstrap. Specifically, the parametric bootstrap requires parametric assumptions for ηi(s). Moreover, for functional linear model (13), permutation method is not directly applicable without very strong assumption such as complete exchangeability.

F. Monte Carlo Simulations

We conducted a set of Monte Carlo simulations to evaluate the Type I and II error rates of the global test statistic Sn based on the resampling method. We simulated FA and MD measures along the right internal capsule tract (Fig. 7(a)) obtained from our clinical data according to

(FAi(sj),MDi(sj))T=(fi,1(sj),fi,2(sj))T+ϵi,j,fi,1(s)=(β11(s),β12(s),β13(s),β14(s))xi+ηi1(s),fi,2(s)=(β21(s),β22(s),β23(s),β24(s))xi+ηi2(s),

where xi = (1, gi, Gagei, agei)T, ηi(s) = (ηi1(s), ηi2(s))T is a 2 × 1 vector of Gaussian process with zero mean and covariance matrix Γη(s, t) and ϵi,j is a 2×1 vector of Gaussian random variables with zero mean and covariance matrix Σ(sj). To mimic imaging data in practice, we used FA and MD measures along the right internal capsule tract from the 128 subjects in our clinical data to estimate Σ^(sj) of Σ(sj) via equation (7), B^(s) of B(s) via equation (16) and Γ^η(s,t) via equation (17). Then, except for (β12(s), β22(s)) for all s, we fixed all other parameters at their values obtained from our clinical data, whereas we assumed (β12(s),β22(s))=c(β^12(s),β^22(s)), where c is a scalar specified below and (β^12(s),β^22(s)) were estimators obtained from our clinical data.

Fig. 7
Right internal capsule tract and diffusion properties along the tract: (a) the right internal capsule tract extracted from the tensor atlas with color representing mean FA value; (b) FA; (c) MD; (d) λ1; (e) λ2; (f) λ3. The diffusion ...

We tested the hypotheses H0 : β12(s) = β22(s) = 0 for all s along the right internal capsule tract against H1 : β12(s) ≠ 0 or β22(s) ≠ 0 for at least a s on the tract. We assumed c = 0 to assess the Type I error rates for the global test statistic, and then we assumed c = 0.1, 0.2, 0.3, and 0.4 to examine the Type II error rates for Sn. In both cases,

C=(0100000000000100)andb0(s)(00)

for all s. We set n = 128 and 64. For n = 128, we used the same values of age, gestational age, and gender from all 128 subjects in our clinical data. For n = 64, we randomly chose 32 males and 32 females from our clinical data and used their values of age, gestational age, and gender to simulate the values of FA and MD along the right internal capsule tract.

For each simulation, the significance levels were set at α = 5% and 1%, and 1,000 replications were used to estimate the rejection rates. For a fixed α, if the Type I rejection rate is smaller than α, then the test is conservative, whereas if the Type I rejection rate is greater than α, then the test is anticonservative, or liberal.

G. Clinical Data

This study was approved by the Institute Review Board of the University of North Carolina at Chapel Hill. A total of 128 healthy full-term subjects (75M and 53F) during the first year of age were recruited and written informed consent was obtained from their parents before imaging acquisition. The subjects were taken from a larger study designated to investigating early brain development at our institution. A total of 128 neonates with a mean postnatal age of 25±17.9 days (range: 9 to 144 days) were included in this study. Efforts were made to ensure the subjects slept comfortably inside the MR scanner, and thus none of the subjects was sedated during imaging session. All subjects were fed and calmed to sleep on a warm blanket with proper ear protection.

All images were acquired on a 3T Allegra head only MR system (Siemens Medical Inc., Erlangen, Germany) with a maximal gradient strength of 40 mT/m and a maximal slew rate of 400 mT/(m·msec). A single shot EPI DTI sequence (TR/TE=5400/73 msec) with eddy current compensation was used to obtain DTI images. Diffusion gradients with a b-value of 1000 s/mm2 were applied in six non-collinear directions, (1,0,1), (−1,0,1), (0,1,1), (0,1,−1), (1,1,0), and (−1,1,0). The reference scan (b=0) was also acquired for the construction of diffusion tensor matrices. Contiguous slices with slice thickness of 2mm were scanned to cover the whole brain. The voxel resolution was isotropic 2mm, and the in-plane field of view was 256mm in both directions. A total of five scans were acquired and averaged to improve the signal-to-noise ratio of the images.

A weighted least square estimation method was used to construct the diffusion tensors [27], [2]. Then the image processing steps in the DTI atlas building were used to process all 128 DTI data sets and compute diffusion properties along all fiber tracts of interest. For the sake of space, we chose two tracts of interest including the splenium of the corpus callosum tract and the right internal capsule tract (Figs. 2(a) and 7(a)) and then computed fractional anistropy (FA), mean diffusivity (MD), and the three eigenvalues of the diffusion tensors, denoted by λ1 ≥ λ2 ≥ λ3, at each grid point on both tracts for each of the 128 subjects. FA denotes the inhomogeneous extent of local barriers to water diffusion, while MD measures the averaged magnitude of local water diffusion. The three eigenvalues of diffusion tensor may, respectively, reflect the magnitude of water diffusivity along and perpendicular to the long axis of white matter fibers [39].

We applied FRATS to the joint analysis of FA and MD values along the splenium tract as follows. We fitted the functional linear model (13) to these smoothed FA and MD functions from all 128 subjects, in which xi = (1, gi, Gagei, agei)T and m = 2, and then we used equation (16) to estimate the function of regression coefficient vector β^(s). Secondly, we constructed the global test statistic Sn to test the effects of all the three covariates for FA alone, MD alone, and joint FA and MD, respectively, and performed hypothesis testing on the whole splenium tract. The p-value of Sn was approximated using the resampling method with G = 10,000. Based on the significant results for MD, we further analyzed the three eigenvalues of diffusion tensors along the splenium tract using the same steps.

Similar to the analysis of the splenium tract, we also applied FRAST to the analysis of FA and MD values along the right internal capsule tract. Then, we fitted the functional linear model (13) to these smoothed FA and MD functions from all 128 subjects, in which xi = (1, gi, Gagei, agei)T and m = 2. We performed formal hypothesis testings in order to examine the effects of gender, gestational age, and age on FA and MD values along the right internal capsule tract. We further analyzed the three eigenvalues of diffusion tensors along the right internal capsule tract.

III. Results

A. Monte Carlo Simulations

Overall, the rejection rates for Sn based on the resampling method were accurate for all sample sizes (n = 64, or 128) at both significance levels (α = 0.01 or 0.05) (Figs. 3(a) and (b)). Consistent with our expectations, the statistical power for rejecting the null hypothesis increased with the sample size.

Fig. 3
Simulation study: Type I and Type II error rates. Rejection rates of Sn based on the resampling method are calculated at five different values of c for sample sizes of 64, 128 subjects at the 5% (green) and 1% (red) significance levels: (a) n = 128; (b) ...

B. Clinical Data: Splenium Tract

We used the adaptive local polynomial kernel smoothing technique to simultaneously smooth FA and MD function vectors for each of 128 subjects, which leads to the smoothed curves (Figs. 4(a) and (b)). We computed the sample covariance matrix and obtained correlation coefficient between FA and MD at each grid point. We observed negative correlations at most grid points (Fig. 4 (c)). Based on the smoothed FA and MD functions, we estimated their mean functions and their covariance functions using equations (11) and (12) (Figs. 4 (d) and (e)). Then, we fitted the functional linear model (13) to the smoothed FA and MD functions from all 128 subjects, in which xi = (1, gi, Gagei, agei)T and m = 2. The elements of β^(s) corresponding to FA are shown in Fig. 4(f).

Fig. 4
Results from the analysis of FA and MD on the splenium tract: reconstructed curves f^i(s)sec for FA in panel (a) and MD in panel (b); (c) estimated correlation between FA and MD along the tract; estimated covariance matrices Γ^(s,t) for FA in ...

From our model, we performed formal hypothesis testings in order to examine the effects of gender, gestational age, and age on FA and MD values along the splenium tract. We first considered tests based on FA and MD, respectively, and performed hypothesis testing at each grid point along the splenium tract (Figs. 5). For FA alone, no significant effect of any covariate was found, even though the −log10(p) value of Sn(s) for age at a single grid point was slightly greater than 2 (Fig. 5 (c)). For MD alone, no significant gender effect was found (Fig. 5 (d)), whereas the effects of age and gestational age were significant at the head and tail of the splenium tract (Figs. 5 (e) and (f)).

For FA alone, no significant effects of all covariates were found in Table 1, which agrees with our previous findings based on the local test statistics Sn(s). This indicates that the degree of anisotropy on the splenium tract does not differ significantly between male and female groups and does not change significantly with gestational age and age. For MD alone, we observed a significant association between gestational age and MD and a marginal significance for age. This indicates a significant change of the degrees of diffusivity, not the degree of anisotropy, along the splenium tract. For joint FA and MD on the splenium tract, no significant gender effect was found, whereas the effects of both age and gestational age were significant.

TABLE I
The p-values of Sn for testing the effects of gender, gestational age (Gage), and age on the splenium and right internal capsule tracts

Along the splenium tract, we observed positive correlations among the three eigenvalues (Fig. 6(d)). Particularly, large correlations between λ2 and λ3 along the splenium tract indicate the small differences between λ2 and λ3 on the splenium tract (Figs. 2 (e) and (f)), whereas relatively weak correlations between λ1 and λ3 indicate large differences between λ1 and λ3 along the splenium tract (Fig. 2(d)). Then, we fitted the functional linear model (13) to the smoothed three eigenvalues functions from all 128 subjects, in which xi = (1, gi, Gagei, agei)T and m = 3.

Fig. 6
Results from the analysis of the three eigenvalues of diffusion tensor on the splenium tract: the −log10(p) values of test statistics Sn(sj) for testing gender effect in panel (a), gestational age effect in panel (b), and age effect in panel (c) ...

We examined association between all covariates of interest and all eigenvalues along the splenium tract by performing formal hypothesis testing at each grid point. No significant gender effect was found for all the three eigenvalues (Fig. 6 (a)), whereas significant gestational age effect was observed at the head and tail of the splenium tract (Fig. 6 (b)). We observed a similar association pattern between gestational age and the three eigenvalues of DT (Fig. 6 (b)). We picked a grid point with a significant p value of Sn(s) and observed that the three eigenvalues decrease with gestational age. For age, a similar association pattern was observed for λ2 and λ3, whereas λ1 was significantly associated with age at the head of the splenium tract (Fig. 6 (c)).

C. Clinical Data: Right Internal Capsule Tract

We simultaneously smoothed the FA and MD function vectors of all 128 subjects and then estimated their mean functions (Figs. 8 (a) and (d)) and their covariance functions (Figs. 8 (b) and (e)). Then, we fitted the functional linear model (13) to the smoothed FA and MD functions from all 128 subjects, in which xi = (1, gi, Gagei, agei)T and m = 2. We considered tests based on FA alone, MD alone, and (FA, MD), respectively, and performed hypothesis testing at each grid point along the right internal capsule tract (Figs. 8(c) and (f)). No significant gender effect was found for all tests (Table 1). For FA alone, significant gestational age and age effects were observed for the right internal capsule tract (Figs. 8 (c) and (f)). We observed negative correlations between FA and MD at most grid points (Fig. 8 (g)). We picked a grid point with the significant p value of Sn(s) at the significance level 5% and observed the increasing pattern of FA with the standardized gestational age (Fig. 8 (i)). For MD alone, significant gestational age and age effects were found in the middle of the right internal capsule tract (Fig. 8 (f)). We picked a grid point with the significant p value of Sn(s) and observed the decreasing pattern of MD with the standardized gestational age (Fig. 8 (h)).

Fig. 8
Results obtained from the analysis of FA and MD on the right internal capsule tract: the estimated mean function and covariance function for FA in panels (a) and (b), respectively; the estimated mean function and covariance function for MD in panels (d) ...

We also constructed the global test statistic Sn to assess the covariates of interest based on FA alone, MD alone, and (FA, MD), respectively, and performed hypothesis testing on the right internal capsule tract. The p-value of Sn was approximated using the resampling method with G = 10,000. The resulting p-values are reported in Table 1. No significant gender effect was found for all scenarios, which indicates no significant difference between male and female groups. We observed significant gestational age and age effects for both FA and MD. Our results agree with previous DTI findings. Current DTI studies including neonates have revealed that FA and MD respectively increases and decreases with age. In white matter, neonates have significantly lower anisotropy values and significantly higher MD values compared to adults [40], [41].

We also observed positive correlations among the three eigenvalues (Fig. 9(b)). We observed large correlation between λ2 and λ3 and relatively weak correlation between λ1 and λ3. This agrees with the small differences between λ2 and λ3 along the right internal capsule tract (Figs. 7 (d)–(f)) and the relatively large differences between λ1 and λ3 in the middle of the right internal capsule tract (Figs. 7(d)–(f)). Then, we fitted the functional linear model (13) to the smoothed three eigenvalues functions from all 128 subjects, in which xi = (1, gi, Gagei, agei)T and m = 3. We performed formal hypothesis testings by examining the association between all covariates of interest and the three eigenvalues along the right internal capsule tract at each grid point. No significant gender effect was found for all the three eigenvalues (Fig. 9 (d)). Gestational age was strongly associated with all the three eigenvalues (Fig. 9 (e)). The scatter plot at a selected grid point shows the decreasing pattern for all eigenvalues with the gestational age. We observed a moderate age effect on all the three eigenvalues (Fig. 9 (f)).

Fig. 9
Results from the analysis of the three eigenvalues of diffusion tensor on the right internal capsule tract: (a) the estimated mean functions for λ1 (red), λ2 (blue), and λ3 (blue); (b) correlations among λ1, λ2 ...

IV. Conclusion and Discussion

The contributions of our work are two folds. Technically, we have presented a novel functional regression framework, called FRATS, for analyzing multiple fiber bundle diffusion properties with a set of covariates of interest, such as age, diagnostic status and gender, in real applications. FRATS can be used to understand normal brain development, the neural bases of neuropsychiatric disorders, and the joint effects of environmental and genetic factors on white matter fiber bundles. From the application end, FRATS is applied to characterizing the development of diffusion properties on fiber bundles in a clinical study of neurodevelopment. Our approach is able to reveal the complex inhomogeneous spatiotemporal maturation patterns as the apparent changes in FA, MD, and the three eigenvalues of DT. Specifically, our results suggest that white matter maturation patterns are different in different white matter regions.

There are still limitations that need to be addressed. First, we have focused on the analysis of a set of water diffusion related parameters based on diffusion tensor image rather than those based on high angular resolution diffusion image (HARDI), because diffusion tensor image is commonly used for characterizing major white matter fiber bundles. In the future, we will apply FRATS to water diffusion related parameters obtained from HARDI [42], [43], [44]. Second, we have focused on tensor derived univariate measures. It is also interesting to extend FRATS to principal directions and full diffusion tensors on fiber bundles [45], [46], [47], [48], [49]. Third, we have focused on analyzing fiber bundle diffusion properties with a set of covariates. The proposed methodology can be readily extended to more complex fiber structures, such as the medial manifolds of fiber tracts [25]. Fourth, we have focused on cross-sectional studies. Extending FRATS to longitudinal studies and family studies requires further research [50].

We have developed FRATS for statistically analyzing multiple diffusion properties along fiber bundle and assessing their association with a set of covariates of interest in real applications. FRATS integrates four statistical tools for formally testing hypothesis of interest and carrying out statistical inference in real applications. The proposed methodology is demonstrated in a clinical study of neurodevelopment. In this study, significant age and gestational age effects on multiple diffusion properties were examined and localized in two representative tracts. We expect that this novel statistical tool will lead to new findings in our clinical applications.

Acknowledgment

The authors would like to thank the editor, an associate editor and three referees for valuable suggestions.

This work was supported in part by NSF grant BCS-08-26844 and NIH grants UL1-RR025747-01, MH086633, and AG033387 to Dr. Zhu, NIH grants MH064065, HD053000, and MH070890 to Dr. Gilmore, NIH grant R01NS055754 to Dr. Lin, Lilly Research Laboratories, the UNC NDRC HD 03110, and NIH Roadmap Grant U54 EB005149-01, NAMIC to Dr. Styner.

References

[1] Basser PJ, Mattiello J, LeBihan D. Mr diffusion tensor spectroscopy and imaging. Biophysical Journal. 1994;66:259–267. [PubMed]
[2] Basser PJ, Mattiello J, LeBihan D. Estimation of the effective self- diffusion tensor from the nmr spin echo. Journal of Magnetic Resonance Ser. B. 1994;103:247–254. [PubMed]
[3] Pierpaoli C, Basser PJ. Toward a quantitative assessment of diffusion anisotropy. Magn Reson Med. 1996;36:893–906. [PubMed]
[4] Hasan KM, Basser PJ, Parker DL, Alexander AL. Analytical computation of the eigenvalues and eigenvectors in dt-mri. J Magn Reson. 2001;152:41–47. [PubMed]
[5] Hasan KM, Narayana PA. Computation of the fractional anisotropy and mean diffusivity maps without tensor decoding and diagonalization: theoretical analysis and validation. Magn Reson Med. 2003;50:589–598. [PubMed]
[6] Zhu HT, Xu D, Amir R, Hao X, Zhang H, Alayar K, Ravi B, Peterson B. A statistical framework for the classification of tensor morphologies in diffusion tensor images. Magnetic Resonance Imaging. 2006;24:569–582. [PMC free article] [PubMed]
[7] Moseley M. Diffusion tensor imaging and aging-a review. NMR Biomed. 2002;15:553–560. [PubMed]
[8] Mukherjee P, McKinstry RC. Diffusion tensor imaging and tractography of human brain development. Neuroimaging Clin N Am. 2006;16:19–43. [PubMed]
[9] Cascio CJ, Gerig G, Piven J. Diffusion tensor imaging: application to the study of the developing brain. J Am Acad Child Adolesc Psychiatry. 2007;46:213–223. [PubMed]
[10] Rollins NK. Clinical applications of diffusion tensor imaging and tractography in children. Pediatr Radiol. 2007;37:769–780. [PubMed]
[11] Smith SM, Jenkinson M, Johansen-Berg H, Rueckert D, Nichols TE, Mackay CE, Watkins KE, Ciccarelli O, Cader MZ, Matthews PM, Behrens TE. Tractbased spatial statistics: voxelwise analysis of multi-subject diffusion data. NeuroImage. 2006;31:1487–1505. [PubMed]
[12] O'Donell LJ, Westin CF, Golby AJ. Tract-based morphometry for white matter group analysis. NeuroImage. 2009;45:832–844. [PMC free article] [PubMed]
[13] Snook L, Plewes C, Beaulieu C. Voxel based versus region of interest analysis in diffusion tensor imaging of neurodevelopment. NeuroImage. 2007;34:243–252. [PubMed]
[14] Bonekamp D, Nagae LM, Degaonkar M, Matson M, Abdalla WMA, Barker PB, Mori S, Horsk A. Diffusion tensor imaging in children and adolescents: reproducibility, hemispheric, and age-related differences. NeuroImage. 2008;34:733–742. [PMC free article] [PubMed]
[15] Gilmore JH, Smith LC, Wolfe HM, Hertzberg BS, Smith JK, Chescheir NC, Evans DD, Kang C, Hamer RM, Lin W, Gerig G. Prenatal mild ventriculomegaly predicts abnormal development of the neonatal brain. Biol Psychiatry. 2008;64:1069–1076. [PMC free article] [PubMed]
[16] Chen YS, An HY, Zhu HT, Stone T, Smith JK, Hall C, Bullitt E, Shen DG, Lin WL. White matter abnormalities revealed by diffusion tensor imaging in non-demented and demented hiv+ patients. NeuroImage. 2009;47:1154–1162. [PubMed]
[17] Focke NK, Yogarajah M, Bonelli SB, Bartlett PA, Symms MR, Duncan JS. Voxel-based diffusion tensor imaging in patients with mesial temporal lobe epilepsy and hippocampal sclerosis. NeuroImage. 2008;40:728–737. [PubMed]
[18] Cmara E, Bodammer N, Rodrguez-Fornells A, Tempelmann C. Age-related water diffusion changes in human brain: A voxel-based approach. NeuroImage. 2007;34:1588–1599. [PubMed]
[19] Snook L, Paulson LA, Roy D, Phillips L, Beaulieu C. Diffusion tensor imaging of neurodevelopment in children and young adults. NeuroImage. 2005;26:1164–1173. [PubMed]
[20] Ashburner J, Friston KJ. Voxel-based morphometry: the methods. NeuroImage. 2000;11:805–821. [PubMed]
[21] Wager TD, Keller MC, Lacey SC, Jonides J. Increased sensitivity in neuroimaging analyses using robust regression. Neuroimage. 2005;26:99–113. [PubMed]
[22] Worsley KJ, Taylor JE, Tomaiuolo F, Lerch J. Unified univariate and multivariate random field theory. NeuroImage. 2004;23:189–195. [PubMed]
[23] Hecke WV, Sijbers J, Backer SD, Poot D, Parizel PM, Leemans A. On the construction of a ground truth framework for evaluating voxel-based diffusion tensor mri analysis methods. NeuroImage. 2009;46:692–707. [PubMed]
[24] Jones DK, Symms MR, Cercignani M, Howard RJ. The effect of filter size on vbm analyses of dt-mri data. NeuroImage. 2005;26:546–554. [PubMed]
[25] Yushkevich PA, Zhang H, Simon TJ, Gee JC. Structure-specific statistical mapping of white matter tracts. Neuroimage. 2008;41:448–461. [PMC free article] [PubMed]
[26] Goodlett CB, Fletcher PT, Gilmore JH, Gerig G. Group analysis of dti fiber tract statistics with application to neurodevelopment. NeuroImage. 2009;45:S133–S142. [PMC free article] [PubMed]
[27] Zhu HT, Zhang HP, Ibrahim JG, Peterson BG. Statistical analysis of diffusion tensors in diffusion-weighted magnetic resonance image data (with discussion) Journal of the American Statistical Association. 2007;102:1085–1102.
[28] Joshi S, Davis B, Jomier M, Gerig G. Unbiased diffeomorphic atlas construction for computational anatomy. Neuroimage. 2004;23:S151–S160. [PubMed]
[29] Alexander DC, Pierpaoli C, Basser PJ, Gee JC. Spatial transformations of diffusion tensor magnetic resonance images. IEEE Transactions on Medical Imaging. 2001;20:1131–1139. [PubMed]
[30] Arsigny V, Fillard P, Pennec X, Ayache N. Log-euclidean metrics for fast and simple calculus on diffusion tensors. Magnetic Resonance in Medicine. 2006;56:411–421. [PubMed]
[31] Fan J, Gijbels I. Local Polynomial Modelling and Its Applications. Chapman and Hall; London: 1996.
[32] Wand MP, Jones MC. Kernel Smoothing. Chapman and Hall; London: 1995.
[33] Wu HL, Zhang JT. Nonparametric Regression Methods for Longitudinal Data Analysis. John Wiley & Sons, Inc.; Hoboken, New Jersey: 2006.
[34] Ramsay JO, Silverman BW. Functional Data Analysis. Springer-Verlag; New York: 2005.
[35] Welsh AH, Yee TW. Local regression for vector responses. Journal of Statistical Planning and Inference. 2006;136:3007–3031.
[36] Zhang J, Chen J. Statistical inference for functional data. The Annals of Statistics. 2007;35:1052–1079.
[37] Zhu HT, Ibrahim JG, Tang N, Rowe DB, Hao X, Bansal R, Peterson BS. A statistical analysis of brain morphology using wild bootstrapping. IEEE Trans Med Imaging. 2007;26:954–966. [PMC free article] [PubMed]
[38] Lin DY. An efficient monte carlo approach to assessing statistical significance in genomic studies. Bioinformatics. 2005;6:781–787. [PubMed]
[39] Song SK, Sun SW, Ju WK, Lin SJ, Cross AH, Neufeld AH. Diffusion tensor imaging detects and differentiates axon and myelin degeneration in mouse optic nerve after retinal ischemia. Neuroimage. 2003;20:1714–1722. [PubMed]
[40] Neil JJ, Shiran SI, McKinstry RC, Schefft GL, Snyder AZ, Almli CR, Akbudak E, Aronovitz JA, Miller JP, Lee BC, Conturo TE. Normal brain in human newborns: apparent diffusion coefficient and diffusion anisotropy measured by using diffusion tensor mr imaging. Radiology. 1998;209:57–66. [PubMed]
[41] Zhai G, Lin W, Wilber KP, Gerig G, Gilmore JH. Comparisons of regional white matter diffusion in healthy neonates and adults performed with a 3.0-t head-only mr imaging unit. Radiology. 2003;229:673–681. [PubMed]
[42] Lenglet C, Campbell JSW, Descoteaux M, Haro G, Savadjiev P, Wassermann D, Anwander A, Deriche R, Pike GB, Sapiro G, Siddiqi K, Thompson PM. Mathematical methods for diffusion mri processing. NeuroImage. 2009;45:S111–S122. [PMC free article] [PubMed]
[43] Leow AD, Zhu S, Zhan L, McMahon K, de Zubicaray GI, Meredith M, Wright MJ, Toga AW, Thompson PM. The tensor distribution function. Magnetic Resonance in Medicine. 2009;61:205–214. [PMC free article] [PubMed]
[44] Tuch DS, Reese TG, Wiegell MR, Makris N, Belliveau JW, Wedeen VJ. High angular resolution diffusion imaging reveals intravoxel white matter fiber heterogeneity. Magnetic Resonance in Medicine. 2002;48:577–582. [PubMed]
[45] Schwartzman A. PhD thesis. Stanford University; California: Jul, 2006. Random Ellipsoids and False Discovery Rates: Statistics for Diffusion Tensor Imaging Data.
[46] Lepore N, Brun CA, Chou Y, Chiang M, Dutton RA, Hayashi KM, Luders E, Lopez OL, Aizenstein HJ, Toga AW, Becker JT, Thompson PM. Generalized tensor-based morphometry of hiv/aids using multivariate statistics on deformation tensors. IEEE Transactions in Medical Imaging. 2008;27:129–141. [PMC free article] [PubMed]
[47] Schwartzman A, Dougherty RF, Taylor JE. Cross-subject comparison of principal diffusion direction maps. Magn Reson Med. 2005;53:1423–1431. [PubMed]
[48] Zhu HT, Cheng YS, Ibrahim JG, Li YM, Hall C, Lin WL. Intrinsic regression models for positive definitive matrices with applications in diffusion tensor images. Journal of the American Statistical Association. 2009;104:1203–1212. [PMC free article] [PubMed]
[49] Whitcher B, Wisco JJ, Hadjikhani N, Tuch DS. Statistical group comparison of diffusion tensors via multivariate hypothesis testing. Magnetic Resonance in Medicine. 2007;57:1065–1074. [PMC free article] [PubMed]
[50] Fang Y, Wang Y. Testing for familial aggregation of functional traits. Statistics in Medicine. 2010 [PMC free article] [PubMed]