Longitudinal Analysis of Image Time Series with Diffeomorphic Deformations: A Computational Framework Based on Stationary Velocity Fields

We propose and detail a deformation-based morphometry computational framework, called Longitudinal Log-Demons Framework (LLDF), to estimate the longitudinal brain deformations from image data series, transport them in a common space and perform statistical group-wise analyses. It is based on freely available software and tools, and consists of three main steps: (i) Pre-processing, (ii) Position correction, and (iii) Non-linear deformation analysis. It is based on the LCC log-Demons non-linear symmetric diffeomorphic registration algorithm with an additional modulation of the similarity term using a confidence mask to increase the robustness with respect to brain boundary intensity artifacts. The pipeline is exemplified on the longitudinal Open Access Series of Imaging Studies (OASIS) database and all the parameters values are given so that the study can be reproduced. We investigate the group-wise differences between the patients with Alzheimer's disease and the healthy control group, and show that the proposed pipeline increases the sensitivity with no decrease in the specificity of the statistical study done on the longitudinal deformations.

1. Introduction

An important topic in neuroimaging is to analyse the progression of morphological changes in the brain observed in time series of images, in order to model and quantify normal or pathological biological evolutions (Scahill et al., 2002). Deformation-Based Morphometry (DBM) (Ashburner et al., 1998) characterizes the morphological changes of the brain in terms of spatial transformations (here called deformations), estimated by means of non-linear registration. A sub-field of DBM, called Tensor-Based Morphometry (TBM) focuses on the first derivatives of the deformation. Depending on the cross-sectional or longitudinal nature of the dataset used, we can define on one hand cross-sectional DBM and on the other hand longitudinal DBM (Chung et al., 2001) that we will focus on in this article. Longitudinal DBM main steps can be summarized as (i) quantifying the evolution of the morphology of each subject by estimating the individual's longitudinal deformation from the time series of images, and (ii) characterizing how this evolution varies among a sample using a suitable normalization for the individual biological variability. A variety of DBM approaches can be found in the literature (e.g., Davatzikos et al., 2001; Cardenas et al., 2007; Lorenzi et al., 2011; Südmeyer et al., 2012), each of them associated to specific non-linear registration methods, and processing pipelines. The comparison between the different DBM methods is not straightforward: the efficiency of each DBM pipeline is generally demonstrated on different data sets (or different subsets of the same data set) and the tools the processing pipeline is composed of are generally not all available. In the existing DBM pipelines—e.g., SPM (Friston, 2007), FreeSurfer (Reuter et al., 2012), PipeDream¹, Anima²—the multivariate information coming from the three-dimensional deformation is generally not used for the statistical analysis. To do so, one would need to express the three-dimensional deformation of every subject in a common space to compare them. There exists few algorithms that compute this 3D transport (e.g., Lorenzi and Pennec, 2013) and in the absence of this tool, the DBM analysis often becomes a TBM analysis only. Studies are thus generally performed on the Jacobian determinant of the deformation or on the segmented regions of interest—since it is easier to compute these scalar maps in a common space. Moreover, in the developing context of reproducible research that has gained interest over the last years (Nature, 2013; McCormick et al., 2014), a good practice should be for researchers to publish the full details of their methodology: source code, data and parameters.

This is the objective of this article: to gather all the details in the same paper and propose a pipeline for the community, following the examples of Avants et al. (2011) and Ashburner and Ridgway (2013). Our computational framework is a complement to the existing processing pipelines. It enables researchers to replicate and verify their findings with a third party reproducible pipeline, thus enhancing the convincing power of their results. Our pipeline is based on Lorenzi et al. (2011), who proposed a hierarchical framework for the group-wise analysis of time series of images using diffeomorphic deformations parameterized by Stationary Velocity Fields (SVF). We bring a complement to the already existing literature by explicitly detailing all the processing steps required for the longitudinal analysis of neuroimages by relying on freely available tools. In addition to this contribution, we integrate a modification to the non-linear registration algorithm by adding a masking to the similarity term as proposed by Brett et al. (2001) while keeping the symmetry of the formulation. This change increases the robustness of the results with respect to intensity artifacts located in the brain boundaries. The proposed processing pipeline is based on freely available software and tools [the complete list can be found in Appendix (Supplementary Materials)].

The paper is structured as follows: in Section 2, we develop a comprehensive processing pipeline called Longitudinal Log-Demons Framework (LLDF); we present each elementary modules it is based on, and after introducing the mathematical formalism related to DBM, we modify the LCC log-Demons to incorporate a confidence mask. Experimental results show that this contribution leads to increased sensitivity of the statistical study on the longitudinal deformations. In Section 3, we show an illustration of the pipeline on the statistical analysis of longitudinal brain changes in Alzheimer's disease. Because it is freely and easily available for benchmarking, we use the data from the longitudinal Open Access Series of Imaging Studies (OASIS) database (Marcus et al., 2010). We finally conclude and present the perspectives of this work in Section 4.

2. Processing Pipeline for the Analysis of Longitudinal Images

We consider longitudinal observations of MRI scans for a given subject S_i, at N_i time points t₀, t₁,…, t_Ni−1 (all the subjects do not necessarily have the same number N_i of time points). The corresponding images are denoted as $I_0^i$, $I_1^i$,…, $I_{N_i-1}^i$ respectively. The aim of the processing pipeline is to estimate each subject's longitudinal deformation from the image time series, and then transport the deformations in a common space to perform statistical group-wise analyses. The construction of the pipeline is based on elementary modules described in the following paragraphs and it can therefore be divided into three main parts (cf. Figure 1): (1) Pre-processing, (2) Position correction, and (3) Non-linear deformation analysis. The pipeline proposed in this work relies on a number of neuroimaging tools previously proposed and validated by different groups. Our choice was motivated by our personal experience and by the optimal performances obtained in the presented application. We however acknowledge that other tools could have been employed. For this reason, the modular nature of the pipeline allows the replacement of the proposed tools with specific ones, such as in the case of longitudinal analysis in postnatal brain development (cf. Section 2.1.3).

FIGURE 1

Figure 1. Proposed processing pipeline for longitudinal analysis. The pipeline is composed of three major steps. Starting with raw images, we first pre-process them, then correct the spatial position differences to end up with the longitudinal deformations for each subject in the template space. Dotted lines correspond to evaluated transformations whereas plain lines correspond to applied transformations.

2.1. Pre-Processing

In this initial part of the pipeline all the individuals' images are processed independently of the time points. The pre-processing consists of the following chain of elementary steps: (1) Standard reorientation, (2) Field of view reduction and, (3) Intensity non-uniformity correction. Different criteria have been taken into account for choosing the tools and software used to perform these elementary steps. Firstly, we only selected freely available tools part of well-established software—so that the pipeline can be reproduced by anyone—relying on already validated tools. Secondly, to make the pipeline user-friendly, we chose tools that necessitate minimal fine tuning in terms of parameters.

2.1.1. Standard Reorientation

Images from the MRI scanner are not necessarily oriented following the standard orientation defined by the MNI152 (Fonov et al., 2009) template (Figure 2). This misorientation would prevent us from properly processing the images.

FIGURE 2

Figure 2. Pre-processing steps. (A) Reorientation of a subject coronal view. Left: what is displayed initially as the coronal view is the sagittal one. Right: after reorientation it is truly the coronal view that is displayed. (B) Field of View Reduction. Left: the original Field of View (FOV) including the head and neck (red rectangle). Right: after reduction, the cropped FOV does not contain the neck, but only the head. (C) Intensity Inhomogeneity Correction. Left: the image has an intensity non-uniformity. The same tissue class has a lower intensity in the bottom left (red ellipse), and a higher intensity in the bottom right part of the image (green ellipse). Middle: after correction, the intensity appearance of the image is more homogeneous (cf. red and green ellipses). Right: estimated multiplicative field. (D) Skull-stripping. Left: the head with its skull. Middle: the brain after the whole process of skull-stripping and image masking. We see that the resulting image has the same intensity as the original one; this is not the case of the image output by Robex (right image).

We thus use FSL - fslreorient2std (Jenkinson et al., 2012), to reorient each image to match the standard orientation. Starting with I, the image acquired by the scanner, this tool applies rotations of 0, 90, 180, or 270 degrees around the image axes to get I^std, the reoriented output image. Notice that this reorientation only changes the header and does not perform any interpolation.

2.1.2. Reduction of the Field of View

Brain scans can sometimes include the neck or the shoulders (cf. Figure 2), and analysing the whole image would increase the image processing time and lead to increased errors due to intensity artifacts. Therefore, it is preferable to reduce the Field of View (FOV) of the image to include the head only.

For this purpose, we use FSL–robustfov (Jenkinson et al., 2012): given an image I, comprising the head and the neck, it automatically crops the neck and other regions outside the head by re-sizing the height of the image, starting at the top of the skull, to a default size of 170 mm so that we finally obtain I^head, the image containing the head only.

In some rare cases (in another study not reported here, one case out of 120), this automatic tool might provide a wrong result, leaving an important part of the neck in the image or cropping the head. In that case, one can still manually set the correct height of the head.

2.1.3. Intensity Inhomogeneity Correction

One of the most common artifact in MRI scans is the shading one: an intensity non-uniformity for voxels of the same tissue class (cf. Figure 2). Therefore, each MR image I undergoes an intensity non-uniformity correction using ANTs–N4BiasFieldCorrection (Tustison et al., 2010; Avants et al., 2011) to obtain the corrected image I^Hom. This algorithm improves the N3 Intensity Inhomogeneity correction (Sled et al., 1998) and is based on the assumption that there exists a smooth, slowly varying multiplicative field F corrupting the image intensities: I = I^Hom×F. In the specific case of early brain development where heterogeneous myelination occurs, the default correction algorithm might be insufficient and a dedicated correction method could be used following Prastawa et al. (2004) example. The choice of the most appropriate algorithm is let to the user. In any case, the Local Correlation Criteria (similar to ANTS Cross-correlation Avants et al., 2011) we use for the non-linear registration in Section 2.3.2 is robust to local intensity bias and is potentially able to cope with an incomplete inhomogeneity correction.

2.1.4. Skull-Stripping

It is often necessary (e.g., in Section 2.2.1) to process the brain without its surrounding skull. For this reason, the pipeline includes a skull-stripping step (also called non-brain removal tool). We selected Robex (Iglesias et al., 2011) for the robustness of its results with no parameter fine tuning: Iglesias et al. (2011) showed it generally performs better than six other popular algorithms (BET, Smith, 2002; BSE, Shattuck et al., 2001; FreeSurfer³; AFNI⁴; BridgeBurner, Mikheev et al., 2008; and GCUT, Mahapatra, 2012). Our experiments were in agreement with this affirmation: when using Robex on our datasets, we no longer had large parts of the skull remaining which was sometimes the case when using FSL–BET with the default parameters. Inputing I, the image with the brain and its surrounding skull, Robex outputs I^robex and I^mask, the skull stripped brain and the corresponding region mask respectively. In fact, Robex applies an additional intensity inhomogeneity correction and thus modifies the intensity of the output image I^robex. Therefore, one has to use the output mask I^mask and mask the original image I to obtain I^brain, the image with the brain only (cf. Figure 2).

2.2. Position Correction

Contrary to the previous section, the images are now treated depending on the subject (and time point). This module consists of two combined steps: (1) Longitudinal rigid registration, and (2) Affine spatial normalization. We first present these modules before explaining how we combine them.

2.2.1. Longitudinal Rigid Registration

For a single subject, the acquisition at different time points is usually not performed with the same position of the head in the scanner. This creates a global rigid (six degrees of freedom) misalignment of each subject data series. Since the aim of this work is to model the subtle local longitudinal brain changes, we need to account for this source of variability that generally exceeds the longitudinal variability. Taking the baseline I₀ as the reference position, we rigidly align the follow-up images I₁,…, I_N−1 to the baseline I₀, using the rigid transformations ${\varphi }_R^1$,…,${\varphi }_R^{N-1}$, to obtain the rigidly aligned image $I_1^{al}$,…, $I_{N-1}^{al}$ (cf. Figure 3).

FIGURE 3

Figure 3. Position correction steps. (A) Rigid registration of subject images. The image on the left is the follow-up image of a subject, the baseline (used as the reference) being the image in the middle. The image on the right is the subject image after rigid alignment. (B) Affine normalization of a subject image. Left: subject image. Middle: the MNI152 template. The subject image and the template differ in size and orientation. Right: result of the affine normalization.

We choose to use FSL–FLIRT (Jenkinson and Smith, 2001; Jenkinson et al., 2002) for the linear registration as it is the benchmark linear registration framework used in the influential work of Klein et al. (2009) for the comparison of several state-of-the-art non-linear registration algorithms. The different steps of the rigid registration step are described in Algorithm 1. We note that despite the optimization in two steps, only one single rigid transformation is applied. Composing the transformations from the—whole head and skull-stripped head—intra-subject rigid registrations minimizes the potential resampling artifacts introduced by the repeated resampling of the data (during the different rigid registration steps). Lastly, we use B-splines as the interpolation method (more accurate than the standard tri-linear interpolation, Parker et al., 1983) and the normalized correlation as the cost function.

ALGORITHM 1

Algorithm 1. Longitudinal Rigid Registration between 2 Images

2.2.2. Affine Spatial Normalization

Each brain differs in size and shape. In preparation for the group analysis and in order to align each subject anatomy in a common reference space, we normalize each subject head (shape and pose) to the MNI152 reference space using an affine (twelve degrees of freedom) transformation. Practically, the brain normalization consists in resampling each subject baseline image I₀ in a common standard space S_MNI (MNI152 space) using an affine transformation ϕ_A computed with FSL–FLIRT to obtain the normalized image $I_0^{MNI}$ (see Figure 3). We use B-splines as the interpolation method and the normalized correlation as the cost function.

2.2.3. Combined Longitudinal Rigid Registration and Spatial Normalization

In the spirit of Section 2.2.1, we avoid as much as possible the potential resampling artifacts by composing the two spatial transformations ϕ_R and ϕ_A from the previous steps. The baseline I₀, is spatially normalized to the MNI152 space using ϕ_A (cf. Section 2.2.2). Concerning the follow-up images I_j, we apply the composition of ϕ_A and ${\varphi }_R^j$ to I_j. Since $I_j^{al}$ and I₀ are already rigidly aligned the transformations that map both of them to the template S_MNI are the same.

2.3. Non-Linear Deformation Analysis

After the correction of the images in position and intensity, we can estimate the residual longitudinal morphological differences using non-linear registration. For this non-linear registration step, all the subjects are processed independently in order to compute each individual longitudinal deformation (expressed in every subject anatomy but with the same coordinate space). The final step is done in three stages: (1) Estimation of the subject-specific longitudinal deformation trajectory using the previously computed longitudinal deformations, (2) Study-specific template creation, and (3) Transport of the subject-specific longitudinal deformation trajectory in the template (cf. Figure 6). Before going further, we introduce the mathematical formalism related to Deformation Based Morphometry.

2.3.1. Mathematical Formalism for Deformation-Based Morphometry

The longitudinal evolution of a point x of the brain between the initial biological time point t₀ = 0 and the biological time t₁ is defined by the deformation ϕ that maps the initial position x(t₀) to the position x(t₁):

$\varphi :{ℝ}^n\times ℝ\to {ℝ}^n$

$(x,t)\mapsto x(t)=\varphi (x,t)$

In neuroimaging, the preservation of the brain topology is important; it can be obtained under the large deformation diffeomorphic setting (Joshi and Miller, 2000; Beg et al., 2005). In this framework, we define the transformations φ that belong to the group ${G}$ of diffeomorphisms: differentiable bijections with differentiable inverse. The transformations are parameterized by the flow of time-dependent velocity vector fields v(x, s) (with the parametrization time s ∈ [0, 1]) specified by the following ordinary differential equation:

$\frac{\partial \phi (x,s)}{\partial s}=v(\phi (x,s),s),$

with φ(x, 0) = Id(x) (identity transformation). The resulting deformation ϕ, mapping x(t₀) to x(t₁) is given by the flow at s = 1: ϕ(x, t₁) = φ(x, 1). In the spirit of the log-Euclidean framework, Arsigny et al. (2006) proposed to restrict to the one-parameter subgroup of diffeomorphisms where the velocity vectors are stationary (i.e., constant over the parametrization time s): v(x, s) = v(x). In this case, the transformation ϕ(x, t₁) is encoded by the stationary velocity field (SVF) v(x) via the Lie group exponential map: ϕ(x, t₁) = exp(v(x)); the exponential map is defined as the flow of the stationary ordinary differential equation:

$\frac{\partial \phi (x,s)}{\partial s}=v(\phi (x,s)),$

with φ(x, 0) = Id(x) and s∈[0, 1].

2.3.2. Non-Linear Symmetric Diffeomorphic Registration with Confidence Mask

We estimate the subtle longitudinal changes using symmetric non-linear diffeomorphic registration. The diffeomorphic deformations are parameterized using Stationary Velocity Fields (SVF), providing us with a rich mathematical and computational setting (see Arsigny et al., 2006; Vercauteren et al., 2008; Lorenzi et al., 2011).

To non-linearly register I_i to I_j, we estimate the Stationary Velocity Field v_i-j (cf. Figure 4) via an alternate minimization of the following log-Demons energy with respect to v_i-j and the auxiliary SVF v_c (Cachier et al., 2003). Instead of minimizing a global energy, a correspondence field v_c is introduced, so that two simple, fast, and more efficient minimization steps are performed, respectively for E_Sim and E_Reg. In the first step, E_Sim is minimized using a gradient descent method, whereas in the second step E_Reg can be solved explicitly as the Gaussian convolution of v_c when the regularization term is chosen adequately:

In this formula, σ_i is the parameter linked to the noise in the image, σ_x is linked to the uncertainty of the matching in the correspondence term, σ_T is the regularization weight, Sim is the similarity criterion, Reg the regularization term, and Corr is the correspondence term that links v_i-j to v_c. The LCC log-Demons (Lorenzi et al., 2013) uses ρ the Local Correlation Coefficient (LCC) similarity metric (Cachier et al., 2003) since it is robust to local intensity artifacts:

$\rho (I_i,I_j)={\displaystyle {\int }_{\Omega }\frac{\overline{I_i{I}_j}}{\sqrt{\overline{I_i^2}\overline{I_j^2}}}}\text{ with }\overline{I}=G_{\sigma }*I,$

where G_σ is the Gaussian smoothing operator with a kernel size of σ and Ω is the image domain.

FIGURE 4

Figure 4. Comparison of three non-linear diffeomorphic registration methods. First and second column: we see the intensity bias affecting the source and target images. (A) Registration of the head with no confidence mask: strong deformation fields are estimated in the skull and meninges that diffuse to the outer cortex region and bias the results (cf. red ellipse where a non-realistic expansion of 38% is found). (B) Registration of the skull-stripped images (no confidence mask): the use of the skull-stripped images biases the result at the level of the outer cortex (cf. red ellipses) where non-existing high value deformations are found due to the high intensity gradient. In fact, skull-stripping imposes the outside brain intensity to be zero creating a high intensity gradient that biases the registration results (the update δv_i-j is directly proportional to the image gradient). (C) Registration of the head with confidence mask: the registration using the confidence mask enables us to estimate realistic transformations in the outer cortex.

Therefore, by considering the symmetric resampling $I_i^{\prime }=I_i\circ exp(\frac{v_c}{2})$ and $I_j^{\prime }=I_j\circ exp(-\frac{v_c}{2})$, the first term of Equation (1) can be written as:

$\begin{array}{l}\mathrm{\text{Sim}}(v_c,I_i,I_j)={\rho }^2(I_i^{\prime },I_j^{\prime })={\rho \prime }^2(v_c,I_i,I_j)\\ =\frac{{\left[\overline{I_i\circ exp(\frac{v_c}{2}).I_j\circ exp(-\frac{v_c}{2})}\right]}^2}{\overline{{\left[I_i\circ exp(\frac{v_c}{2})\right]}^2}.\overline{{\left[I_j\circ exp(-\frac{v_c}{2})\right]}^2}}\mathrm{\text{.}}\end{array}$

If we define the update field δv_i-j through the zeroth order term of the Baker-Campbell-Hausdorff (BCH) formula (Bossa et al., 2007):

$\delta v_{i\mathrm{\text{-}}j}=log(exp(-v_{i\mathrm{\text{-}}j})\circ exp(v_c))\approx -v_{i\mathrm{\text{-}}j}+v_c,$

then in the first part of the alternate optimization of Equation (1), E_Sim has to be minimized with respect to δv_i-j :

$E_{Sim}(\delta v_{i\mathrm{\text{-}}j},I_i,I_j)=-\frac{1}{{\sigma }_i^2}{\rho {}^{\prime }}^2(\delta v_{i\mathrm{\text{-}}j},I_i,I_j)+\frac{1}{{\sigma }_x^2}||\delta v_{i\mathrm{\text{-}}j}|{|}^2,$

with $\mathrm{\text{Corr}}(v_{i\mathrm{\text{-}}j},v_c)=||log(exp(-v_{i\mathrm{\text{-}}j})\circ exp(v_c))|{|}^2=||\delta v_{i\mathrm{\text{-}}j}|{|}^2$.

In the second part of the optimization, E_Reg should be minimized with respect to v_i-j :

$E_{Reg}(v_{i\mathrm{\text{-}}j},v_c)=\frac{1}{{\sigma }_x^2}||log(exp(-v_{i\mathrm{\text{-}}j})\circ exp(v_c))|{|}^2+\frac{1}{{\sigma }_T^2}Reg(v_{i\mathrm{\text{-}}j})\mathrm{\text{.}}$

The registered images generally comprise the brain and its surrounding skull which can lead to corrupted results. In fact, the resulting deformation field generally exhibits high values in the region of the meninges and the skull that diffuse through regularization in the outer cortex (see Figure 4), potentially yielding to misleading discoveries.

One solution is to only register the brain tissues and the cerebrospinal fluid (CSF) obtained through skull-stripping. However, this solution may be prone to errors (small parts of the outer cortex could be cropped) and puts the outside brain intensity to zero creating a high intensity gradient that biases the registration results (as shown on Figure 4), since the update δv_i-j is directly proportional to the image gradient.

Therefore, we modified the LCC log-Demons algorithm to incorporate the use of a confidence mask as proposed by Brett et al. (2001), and first introduced in the Demons algorithm by Stefanescu et al. (2004). We consider that we do not want to align the structures outside the brain (skull, meninges,…). Therefore, the voxels outside the brain should have no influence in the similarity minimization step. We define a probabilistic mask ω(x) such that its value is ω(x) = 1 for a voxel inside the brain, ω(x) = 0 outside, and in-between depending on the confidence we have for the voxel. The new log-Demons energy to minimize is:

$\begin{array}{l}E(v_{i\mathrm{\text{-}}j},v_c,I_i,I_j)=\omega \frac{1}{{\sigma }_i^2}\mathrm{\text{Sim}}(v_c,I_i,I_j)+\frac{1}{{\sigma }_x^2}\mathrm{\text{Corr}}(v_{i\mathrm{\text{-}}j},v_c)\\ +\frac{1}{{\sigma }_T^2}\mathrm{\text{Reg}}(v_{i\mathrm{\text{-}}j})\mathrm{\text{.}}\end{array}$

Thus, only the first part of the minimization (E_Sim) is modified and we still get a closed-form solution leading to an effective computational scheme for the optimization of E_Sim [cf. demonstration in Appendix (Supplementary Materials)]:

$\delta v_{i\text{-}j}=\{\begin{array}{l}-\frac{2\Lambda }{||\Lambda |{|}^2+\frac{1}{\omega }\frac{4}{{\rho }^2}\frac{{\sigma }_i^2}{{\sigma }_x^2}}\text{, if }\omega >0\hfill \\ 0,\text{ if}\omega =0\hfill \end{array}$

with

$\begin{array}{r}\Lambda =\frac{G_{\sigma }*(I_i\nabla I_j^T)}{G_{\sigma }*(I_i{I}_j)}-\frac{G_{\sigma }*(I_j\nabla I_i^T)}{G_{\sigma }*(I_i{I}_j)}+\frac{G_{\sigma }*(I_i\nabla I_i^T)}{G_{\sigma }*(I_i^2)}\\ -\frac{G_{\sigma }*(I_j\nabla I_j^T)}{G_{\sigma }*(I_j^2)}\mathrm{\text{.}}\hfill & (2)\end{array}$

In order to keep a symmetric formulation of the registration, the probabilistic mask ω is defined using two masks. The first one is the brain mask M of the moving image and the second one is the brain mask F of the fixed image. The mask ω is then defined as the average of the symmetric resampling of the two brain masks in the halfway space:

$\omega =\frac{1}{2}\left[M\circ exp(\frac{v_c}{2})+F\circ exp(-\frac{v_c}{2})\right]$

Hence, the registration problem is still defined on the whole image domain but the update is weighted differently depending on the confidence on the brain areas. In our experiments, we defined the initial brain masks (for both fixed and moving images) as binary masks.

2.3.3. Estimation of the Subject-Specific Longitudinal Trajectory via Fully Symmetric SVF Regression

Given the previously estimated series of longitudinal deformations ϕ_i-j = exp(v_i-j) with 0 ≤ i < j ≤ N − 1 for a subject, we then model the subject-specific longitudinal deformation trajectory $\widehat{\varphi }$ as :

$\widehat{\varphi }(x,t)=exp(t·\widehat{v}(x))\mathrm{\text{with}}t\in ℝ,$

where $\widehat{v}$ is the best fit of a fully symmetric linear model in time—through the origin—of the series of SVFs v_i-j :

$\begin{array}{l}\widehat{v}=\underset{v}{\text{arg}\text{min}}{\displaystyle \sum _{0\le i<j\le N-1}\Vert (}{t}_j-t_i)v-v_{i\text{-}j}{\Vert }^2\\ =\frac{{\displaystyle {\sum }_{0\le i<j\le N-1}(t_j-t_i)v_{i-j}}}{{\displaystyle {\sum }_{0\le i<j\le N-1}{(t_j-t_i)}^2}}\text{.}\end{array}$

This model uses all the possible combinations of SVFs v_i-j between the different time points while using the symmetry of the pairwise registration (v_i-j = −v_j-i) to simplify the problem. $\widehat{v}$ and $\widehat{\varphi }(t=1)=exp(\widehat{v})$ represent the subject-specific evolution trajectory over a year. One should note that a linear model of the longitudinal SVFs does not lead to a linear model of the deformations. For up to three time points, our experience showed that a linear model in time is sufficient to explain the data. A higher-order model could be used for a higher number of time points at the cost of increasing the statistical complexity.

2.3.4. Unbiased Study-Specific Template Construction

In order to compare all the subject-specific longitudinal deformation trajectories, we need to have these deformations normalized in the same common reference anatomy called study-specific template A. Although each subject brain is normalized to the standard space (cf. Section 2.2.2), the affine alignment is not sufficient to compensate for the local anatomical differences (there is no voxel-to-voxel correspondence yet between the different anatomies). Among the available methods for the template construction, we chose to use the method from Guimond et al. (2000) consisting in the iterative averaging of intensities and deformations. This iterative process is described in Algorithm 2 and illustrated on 136 subjects (Figure 5). In the following experiments, the iterative algorithm was stopped at the seventh iteration. At a given iteration there are two successive image resamplings due to the application of two deformations; this can bias the centering of the template. To ensure it is centered, we minimize the number of image resamplings at a given iteration by using a zeroth order term of the BCH: $log(exp(v_k^i)\circ exp(-{\stackrel{-}{v}}_k))\approx v_k^i-{\stackrel{-}{v}}_k$. Moreover, a good practice for the selection of the initialization image for A₀ is to manually choose a subject image that is roughly centered with respect to the considered sample in order to avoid being blocked in a local minimum. In practice, we checked that changing the reference image for A₀ changed the final template A by only a negligible amount as shown on Figure 5.

ALGORITHM 2

Algorithm 2. Creation of an Unbiased Template A

FIGURE 5

Figure 5. Top: Iterative template Construction. Example of the construction of the template (green frame) of a study of 136 subjects, 9 subjects are displayed. Red frame: the reference subject (OAS2_0017) used for the initialization. Bottom: Influence of the reference subject used to initialize the study-specific template. We built a second study-specific template by initializing it with a different subject (OAS2_0077). We computed $\frac{1}{V}\sum _i|Templat{e}_{OAS2\text{\_}0077}(x_i)-Templat{e}_{OAS2\text{\_}0017}(x_i)|$ over the brain mask at each iteration. Although the initial reference images are dissimilar, we obtain two very similar templates at the end.

Here again the non-linear registrations are performed using our modified LCC log-Demons algorithm with confidence mask (we used the subjects images masks), in order to estimate the study-specific template while being robust to the artifacts on the brain boundaries.

Another point concerns the choice of the time point at which the template is created. There is no golden rule and the choice of the time point is usually let to the user. As for us, we use the images I₀ at the first time point t₀ to create the template.

2.3.5. Parallel Transport of the Subject-Specific Longitudinal Stationary Velocity Field

Now that a common brain anatomical image is defined, we need to express each subject-specific longitudinal deformation trajectory $\widehat{\varphi }$ in the template anatomy to be able to compare them. To do so, we use the parallel transport computed with the Pole ladder (Lorenzi and Pennec, 2013) of the subject-specific longitudinal SVF trajectory $\widehat{v}$ along the inter-subject SVF w₀ parameterizing the cross-sectional transformation ψ₀ = exp(w₀) that maps I₀ to A (cf. Figure 6).

FIGURE 6

Figure 6. Illustration of the parallel transport of the study transformations to the study-specific template. After each subject longitudinal SVF transport, the mean transformation ${\stackrel{-}{\varphi }}_T$ is computed by taking the exponential of the average of all the transported subject-specific longitudinal SVFs ${\Pi }_{w_0^i}({\widehat{v}}_T^i)$.

The result is ${\widehat{v}}_T={\Pi }_{w_0}(\widehat{v})$, the subject-specific longitudinal SVF trajectory normalized in the template space. We can then compute the subject-specific longitudinal deformation trajectory in the template space ${\widehat{\varphi }}_T=exp({\widehat{v}}_T)$. The different steps necessary for the parallel transport are described in Algorithm 3. It is then possible to perform a statistical analysis on these transported subject-specific longitudinal stationary velocity fields ${\widehat{v}}_T$ as shown in Section 3.

ALGORITHM 3

Algorithm 3. Pole Ladder for the Parallel Transport of the Longitudinal Stationary Velocity Field

3. Application to the Analysis of the Longitudinal Changes in Alzheimer's Disease

The aim of this section is to show an application of the proposed processing pipeline. We focused our illustration on Alzheimer's disease, a neuro-degenerative disease that causes dramatic changes in the brain anatomy over time. We use the OASIS database (Marcus et al., 2010).

3.1. OASIS Database

The clinical cohort considered in this study is composed of 64 patients diagnosed with very mild to moderate Alzheimer's disease, and 72 healthy individuals. For these subjects, 2 to 5 longitudinal brain acquisitions (T1 Magnetic Resonance Imaging) were available, corresponding to a follow-up time t_0-j = t_j−t₀ of 0.5 to 6.9 years. Further information can be found in Appendix (Supplementary Materials).

3.2. Methods and Results

After applying the processing pipeline to the database (the parameters used for the different steps are summarized Table 1), we obtain the transported subject-specific longitudinal deformation trajectories ${\widehat{\varphi }}_T^i(t)=exp(t·{\widehat{v}}_T^i)$ for each subject i in the study-specific template: we thus get 72 subject-specific longitudinal SVFs ${\widehat{v}}_T^i$ for the healthy controls and 64 for the patients with Alzheimer's disease.

TABLE 1

Table 1. Parameters used for each module of the longitudinal study.

Concerning the non-linear registration parameters for the LCC log-Demons, the optimal parameters we propose here would of course be different for another study, and we recommend to fine-tune in priority the amount of regularization (-b) and the number of iterations (-a). As for the SVF exponentiation (and log-Jacobian), all the computations were performed using an Euler forward integration scheme (option -z 1 in SVFLogJacobian tool).

Before discussing the results of the group-wise comparisons of the longitudinal evolutions, let us focus on an illustrative result concerning a single subject (OAS2_0002). We computed the log-Jacobian map—which quantifies the relative volume changes associated to the longitudinal deformation—for the SVF v_0-2 of the longitudinal evolution between t₀ and t₂; the result can be seen on Figure 7. We can observe the expansion of the ventricles and more particularly in the temporal horn of the lateral ventricles, as well as the contraction in the hippocampi. Moreover, there exists an artifact outside the brain (left hand edge of the follow-up image on Figure 7). The use of the non-linear registration with confidence mask enables us to avoid any artifactual volume change in our log-Jacobian map and therefore provides more stable results. This is illustrated in Figure 4, where we compare the deformation found with and without the use of the confidence mask; we see on the left hand of the image (red circle on image A.) that this kind of artifact can locally bias the estimation of longitudinal deformations when the mask is not explicitly accounted for (Ashburner and Ridgway, 2013).

FIGURE 7

Figure 7. Log-Jacobian map for the subject OAS2_0002. We computed the log-Jacobian map—which represents the relative change of volume—for the SVF of the longitudinal evolution between t₀ and t₂. We can observe an expansion in the ventricles and more particularly in the temporal horn of the lateral ventricles and a contraction in the hippocampi. Moreover, although there is an artifact outside the brain (left hand edge of the follow-up image at t₂), the use of the non-linear registration with confidence mask enables us to avoid any artifactual volume change in our log-Jacobian map.

Concerning the groups study, we consider the subject-specific deformations over a year (t = 1) so that we study the SVFs ${\widehat{v}}_T^i$. It is then possible to visualize the mean volume changes during 1 year for each group of patients with Alzheimer's disease and healthy controls. After computing the average SVF for the non-demented group and the Alzheimer's one, we compute the associated log-Jacobian maps⁵ (cf. Figure 8), and compare the modeled group-wise evolutions. We can see that the main expansion region is located in the lateral ventricles with higher values for the Alzheimer's patients group than for the healthy control one. Moreover, for the patients with Alzheimer's disease we can see an expansion in the temporal horn of the lateral ventricles that does not exist in the control group. Finally, the atrophy is higher for the Alzheimer's patients and mainly located in several parts of the white matter, in the thalamus and in the hippocampi whereas there is no visible contraction in the hippocampi or in the thalamus for the control group. These results are coherent with the findings reported in the literature (Braak and Braak, 1991; Fox et al., 1996; Jack et al., 2004; Schott et al., 2005; de Jong et al., 2008).

FIGURE 8

Figure 8. Template for the 136 OASIS subjects at t₀ and log-Jacobian maps (1 year evolution) of the patients with Alzheimer's disease and the healthy control group. The main expansion region concerns the lateral ventricles where the Alzheimer's patients exhibit higher values when compared to the healthy subjects. Moreover, for the patients with Alzheimer's disease we can see an expansion in the temporal horn of the lateral ventricles that does not exist in the non-demented control group. Finally, the atrophy is higher for the Alzheimer's patients and mainly located in several parts of the white matter, in the thalamus and in the hippocampi whereas there is no visible contraction in the thalamus or the hippocampi for the healthy group.

3.2.1. Two-Sample t-Test: Alzheimer's Patients vs. Healthy Controls

We now statistically investigate the group-wise differences between the modeled longitudinal evolutions of the Alzheimer's patients group and the healthy control group by using a voxel-wise two-sample t-test on the log-Jacobian maps. For illustrative purposes, we show here a standard univariate analysis on a scalar map, but the use of the parallel transport in our pipeline enables us to do statistics directly on the subject-specific SVFs as shown in Section 3.2.3. The null hypothesis is that there exists no difference between the mean of the two groups. We used SPM8 (see Friston, 2007) for this test and corrected for multiple testing using the Family-Wise Error rate (FWE) with a corrected p-value of 0.05 in order to control for the same level of specificity. The t-test was limited to the brain mask. The result map with the thresholded t-values can be seen on Figure 9. The statistically different volume changes occur in the lateral ventricles, more particularly in the temporal horn, and also in the thalamus.

FIGURE 9

Figure 9. Corrected t-statistic map for the volume changes differences between the patients with Alzheimer's disease and the healthy control group (for the 3 registration methods) on one slice. The three results present similar patterns with statistical differences in the ventricular region, more particularly in the temporal horn of the lateral ventricles, and also in the thalamus. The volume of the regions of statistical significant differences are 10.4, 16.5, and 17.5 cm³ for respectively “Pipeline Skull-stripped” (C), “Pipeline Head” (B), and “LLDF” (A). Moreover, the t-values are higher with the “LLDF” than with the two other methods. (Correction for multiple testing using the Family-Wise Error rate with a corrected p-value of 0.05).

3.2.2. Reliability of the LCC Log-Demons with a Confidence Mask

We tested the reliability of the implemented LCC log-Demons registration with a confidence mask. We compared it with the original LCC log-Demons applied to full head images or skull-stripped images. We therefore ran three similar processing pipelines where the only difference was the non-linear registration method used; the processing pipeline using the LCC log-Demons with a confidence mask is denoted as LLDF, the one using the registration of the whole head is called Pipeline Head, and the pipeline registering skull-stripped images is denoted as Pipeline Skull-stripped. Similarly to Section 3.2.1, we investigated the differences between the Alzheimer's patients group and the healthy control group in each case and compared the obtained results to see which method has the highest statistical sensitivity to find volume changes between the two groups.

The three corrected t-maps are presented Figure 9⁶. The three results present similar patterns with most of the statistical differences in the ventricular region and more particularly in the temporal horn of the lateral ventricles. Other statistical differences can be found in the thalamus. The volume of the regions of statistical significant differences are 10.4, 16.5, and 17.5 cm³ for respectively “Pipeline Skull-stripped”, “Pipeline Head”, and “LLDF”. Moreover, the t-values are higher with the “LLDF” than with the two other methods. In average on the same statistical region (the smallest region, obtained by computing the intersection of the three statistically significant regions), we obtain an absolute t-value of 6.13 with “LLDF” against 5.98 with “Pipeline Head”, and 5.69 with “Pipeline Skull-stripped”. This increase of the t-values can be explained by the increased group difference for “LLDF” compared to the group differences of the other two methods and not by a reduction of the variance. On the same statistical region, we observe a relative increase of 23.4% with respect to “Pipeline Head” and of 23.7% with respect to “Pipeline Skull-stripped”. Therefore, the LLDF pipeline enables us to have an increased statistical sensitivity with no decrease of the specificity.

3.2.3. Illustration of a DBM Analysis: Hoteling's Two-Sample T²-Test

Finally, we illustrate the main advantage of the LLDF: by using the parallel transport in our pipeline it is then possible to perform statistics directly on the subject-specific longitudinal trajectories. We therefore perform a multivariate Hoteling's two-sample T²-test to show the group-wise differences between the modeled subject-specific longitudinal trajectories of the Alzheimer's patients group and the healthy control group—obtained using the confidence mask. The null hypothesis is that there exists no difference between the mean of the two groups. We corrected for multiple testing using 5000 permutations and we limited the test to the brain mask. The resulting T²-map thresholded for a corrected p-value of 0.05 can be seen on Figure 10⁷.

FIGURE 10

Figure 10. Top: Group longitudinal trajectories for the patients with Alzheimer's disease and the healthy control group (obtained with the LLDF method). We can see that the mean trajectory for the demented group has a higher magnitude than the control one. Bottom: Corrected T²-map for the longitudinal trajectories differences between the patients with Alzheimer's disease and the healthy control group (for the LLDF method) on one slice. The statistical differences between the demented and the control groups are located in the lateral ventricles, in the temporal horn of the ventricles, in the hippocampi, and in the caudate nuclei. The volume of the regions of statistical significant differences is 41.0 cm³. (The Hoteling's T²-test was corrected for multiple testing using 5000 permutations and the map is thresholded for a corrected p-value of 0.05).

We can see that the statistical differences between the demented and the control groups are located in the lateral ventricles, in the temporal horn of the ventricles, in the hippocampi and, in the caudate nuclei. The volume of the regions of statistical significant differences is larger than the one found using the univariate test: 41.0 cm³. The observed differences in the statistically significant regions between the univariate t-test (cf. Section 3.2.1) and the multivariate Hoteling's T²-test can be explained by the fact that in the first case the study is restricted to the volumetry only whereas in the second case it focuses on the displacement field—which in addition to the volumetry also includes translations and rotations. With this difference in mind, we can say that the patterns found in the two tests are coherent. For example concerning the caudate nuclei, although there is no statistically significant difference in the volume changes between the patients with Alzheimer's and the healthy subjects, there exist statistically significant differences in the displacements of the caudate nuclei between the two groups.

4. Conclusion and Discussion

We proposed and detailed a new processing pipeline⁸ for the longitudinal analysis of image data series. It is based on freely available software and tools so that anyone can reproduce our study, use this pipeline to replicate and verify findings conducted with other pipelines or use it to perform new studies. Moreover, we also implemented a masking of the similarity term in the non-linear registration (with a formulation that ensures symmetry). It enhances the robustness of the registration results with respect to intensity artifacts in the boundary of the brain, thereby increasing the sensitivity of the statistical studies done on the longitudinal deformations. We finally showed on an open-access database that the results obtained with this pipeline are consistent with the findings from the literature.

The use of the parallel transport in our pipeline enables us to perform both standard univariate analysis on a scalar map and also statistics directly on the SVFs as illustrated by the multivariate Hoteling's T²-test. Therefore, changes other than the ones linked to volumetry (like rotations or translations of the brain structures) could be studied. Concerning the confidence mask, initializing it with probabilistic masks of the fixed and moving image (instead of binary ones) could be used to take into account the uncertainty linked to the skull-stripping at the brain boundaries. However, in our experiment the use of binary masks was sufficient to increase the sensitivity of the statistical group-wise analysis while not decreasing the specificity. Intensities artifacts inside the brain such as prominent blood vessels could also be incorporated in the confidence mask if a blood vessels segmentation was available.

The most important issue for the longitudinal processing pipeline is related to the asymmetry biases (Ridgway et al., 2015) that need to be avoided in the processing. Two types of asymmetries can be distinguished. The first one, described in Reuter and Fischl (2011) and Yushkevich et al. (2010), is introduced by the resampling of all the follow-up images except the baseline. In our case, all the images (including the baseline image I₀ at t = 0) are resampled once and only once in the common reference space. In the case of the follow-up images, the transformation used to resample the image is the combination of a rigid and an affine transformation (cf. Section 2.2.3), whereas in the case of the baseline image, we use the subject to reference space affine transformation only. This aspect of the pipeline has some similarity to that of Rohrer et al. (2013)—where again, some repeated interpolations are avoided, while other interpolations are symmetric by virtue of being in MNI space rather than in the native baseline space. It could be possible to go one step further and to avoid any explicit interpolation by initializing the non-linear registration (in the LCC log-Demons software) with the combined affine/rigid transformation using the software parameter:—initial-linear-transform. However, this would still imply an implicit internal resampling and in this case we would no longer follow the assumption made in LDDMM and the SVF framework that all the field tends toward zero when we get away from the center of the image (i.e., beyond the borders of the image). In practice, we observe edge-effects and a proper way to deal with the problem should be to revise the LCC log-Demons algorithm in order to explicitly handle the two transformations separately and make sure that the criterion (and the discretization) would be affine invariant.

The second type of bias is related to the non-centrality of the time point where the subject longitudinal deformations are computed (also referenced as favoring a particular time point). Several non-stationary velocity fields-based methods (LDDMM) have taken great care of that (Avants et al., 2011; Niethammer et al., 2011; Ashburner and Ridgway, 2013). In these methods, the initial velocity (or equivalently the momentum map) is different at different time points along a geodesic. In that case, for more than two time points, it is necessary to choose a time point for the subject-specific template, and this time point is generally the average (or median) of the observed time points. The momentum maps (from the template to all the time points) can then be compared in the template reference space only. In the stationary velocity field framework, the velocity field is—by definition—stationary. Thus, the SVF resulting from the registration is the same all along the trajectory: it is not expressed in material coordinates at a specific time point but in Eulerian coordinates which are not attached to a given time point. Therefore, in the symmetric LCC log-Demons any subject time point can be chosen to perform the pairwise registrations without needing a subject-specific template. Moreover, the annualized log-Jacobian map is valid for all time points even if its value for a material point changes with time along its trajectory. Finally, even if each registration is fundamentally pairwise, the effect of the multiple time points is taken care of using the fully symmetric linear model in time described in Section 2.3.3. This model uses all the possible combinations of SVFs in order to avoid favoring any specific time point. Notice that this approach is sub-optimal with unbalanced data where large variations exist in the number of time points N_i between the subjects. This can be corrected using methods like the one described in Guillaume et al. (2014). However, in the study presented here, only 13 subjects out of 136 had more than three time points. The majority had two or three time points which did not unbalance the data too much.

Apart from the bias, one can wonder what would be the best method between LDDMM and the SVF framework. At first sight, LDDMM might appear as a better theoretical model for an elastic mechanical deformation since it is based on the conservation of the Hamiltonian. However, it is not completely clear that the longitudinal evolution of a brain (intra-subject) is an elastic deformation that conserves the energy. Moreover, in practice Lorenzi and Pennec (2013) showed that for the longitudinal registration the differences between the two methods are very subtle and the stationary velocity field framework can be used.

Author Contributions

MH has worked on the design, implementation, and experimental study of the pipeline. ML, NA, and XP have co-supervised him in the design and revision of the work, and they gave their final approval of the version to be published.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This work was partially funded by the European Research Council through the ERC Advanced Grant MedYMA 2011-291080 (on Biophysical Modeling and Analysis of Dynamic Medical Images).

OASIS grant numbers: P50 AG05681, P01 AG03991, R01 AG021910, P20 MH071616, U24 RR021382.

Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fnins.2016.00236

Footnotes

1. ^http://sourceforge.net/projects/neuropipedream/

2. ^https://github.com/Inria-Visages/Anima-Public/wiki

3. ^http://surfer.nmr.mgh.harvard.edu

4. ^http://afni.nimh.nih.gov/

5. ^The log-Jacobian maps (for OAS2_0002 and the different groups) are available on NeuroVault (Gorgolewski et al., 2015) at http://neurovault.org/collections/YBADDEIH/.

6. ^The t-maps as well as the group difference and estimated variance maps for the three methods are available at http://neurovault.org/collections/YBADDEIH/.

7. ^The T²-map is available at http://neurovault.org/collections/YBADDEIH/.

8. ^The whole pipeline will be released as a complement of the already available LCC log-Demons software.

References

Arsigny, V., Commowick, O., Pennec, X., and Ayache, N. (2006). “A log-euclidean framework for statistics on diffeomorphisms,” in Medical Image Computing and Computer-Assisted Intervention – MICCAI 2006, Vol. 4190 of Lecture Notes in Computer Science, eds R. Larsen, M. Nielsen, and J. Sporring (Berlin; Heidelberg: Springer), 924–931.

Ashburner, J., Hutton, C., Frackowiak, R., Johnsrude, I., Price, C., and Friston, K. (1998). Identifying global anatomical differences: deformation-based morphometry. Hum. Brain Mapp. 6, 348–357.

PubMed Abstract | Google Scholar

Ashburner, J., and Ridgway, G. R. (2013). Symmetric diffeomorphic modelling of longitudinal structural MRI. Front. Neurosci. 6:197. doi: 10.3389/fnins.2012.00197

PubMed Abstract | CrossRef Full Text

Avants, B. B., Tustison, N. J., Song, G., Cook, P. A., Klein, A., and Gee, J. C. (2011). A reproducible evaluation of ANTs similarity metric performance in brain image registration. Neuroimage 54, 2033–2044. doi: 10.1016/j.neuroimage.2010.09.025

PubMed Abstract | CrossRef Full Text | Google Scholar

Beg, M., Miller, M., Trouve, A., and Younes, L. (2005). Computing large deformation metric mappings via geodesic flows of diffeomorphisms. Int. J. Comput. Vis. 61, 139–157. doi: 10.1023/B:VISI.0000043755.93987.aa

CrossRef Full Text | Google Scholar

Bossa, M., Hernandez, M., and Olmos, S. (2007). Contributions to 3D diffeomorphic atlas estimation: Application to brain images. Med. Image Comput. Comput. Assist. Interv. 10, 667–674. doi: 10.1007/978-3-540-75757-3_81

PubMed Abstract | CrossRef Full Text | Google Scholar

Braak, H., and Braak, E. (1991). Alzheimer's disease affects limbic nuclei of the thalamus. Acta Neuropathol. 81, 261–268.

PubMed Abstract | Google Scholar

Brett, M., Leff, A. P., Rorden, C., and Ashburner, J. (2001). Spatial normalization of brain images with focal lesions using cost function masking. Neuroimage 14, 486–500. doi: 10.1006/nimg.2001.0845

PubMed Abstract | CrossRef Full Text | Google Scholar

Cachier, P., Bardinet, E., Dormont, D., Pennec, X., and Ayache, N. (2003). Iconic feature based nonrigid registration: the PASHA algorithm. Comput. Vis. Image Understand. 89, 272–298. doi: 10.1016/S1077-3142(03)00002-X

CrossRef Full Text | Google Scholar

Cardenas, V. A., Studholme, C., Gazdzinski, S., Durazzo, T. C., and Meyerhoff, D. J. (2007). Deformation-based morphometry of brain changes in alcohol dependence and abstinence. Neuroimage 34, 879–887. doi: 10.1016/j.neuroimage.2006.10.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Chung, M., Worsley, K., Paus, T., Cherif, C., Collins, D., Giedd, J., et al. (2001). A unified statistical approach to deformation-based morphometry. Neuroimage 14, 595–606. doi: 10.1006/nimg.2001.0862

PubMed Abstract | CrossRef Full Text | Google Scholar

Davatzikos, C., Genc, A., Xu, D., and Resnick, S. M. (2001). Voxel-based morphometry using the RAVENS maps: methods and validation using simulated longitudinal atrophy. Neuroimage 14, 1361–1369. doi: 10.1006/nimg.2001.0937

PubMed Abstract | CrossRef Full Text | Google Scholar

de Jong, L. W., van der Hiele, K., Veer, I. M., Houwing, J. J., Westendorp, R. G., Bollen, E. L., et al. (2008). Strongly reduced volumes of putamen and thalamus in Alzheimer's disease: an MRI study. Brain 131, 3277–3285. doi: 10.1093/brain/awn278

PubMed Abstract | CrossRef Full Text | Google Scholar

Fonov, V., Evans, A., McKinstry, R., Almli, C., and Collins, D. (2009). Unbiased nonlinear average age-appropriate brain templates from birth to adulthood. Neuroimage 47(Suppl.1), S102. doi: 10.1016/S1053-8119(09)70884-5

CrossRef Full Text

Fox, N. C., Warrington, E. K., Freeborough, P. A., Hartikainen, P., Kennedy, A. M., Stevens, J. M., et al. (1996). Presymptomatic hippocampal atrophy in Alzheimer's disease. Brain 119, 2001–2007. doi: 10.1093/brain/119.6.2001

PubMed Abstract | CrossRef Full Text | Google Scholar

Friston, K. (2007). “Statistical parametric mapping,” in Statistical Parametric Mapping, eds K. Friston, J. Ashburner, S. Kiebel, T. Nichols, and W. Penny (London: Academic Press), 10–31. doi: 10.1016/B978-012372560-8/50002-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Gorgolewski, K. J., Varoquaux, G., Rivera, G., Schwartz, Y., Ghosh, S. S., Maumet, C., et al. (2015). NeuroVault.org: a web-based repository for collecting and sharing unthresholded statistical maps of the human brain. Front. Neuroinform. 9:8. doi: 10.3389/fninf.2015.00008

PubMed Abstract | CrossRef Full Text | Google Scholar

Guillaume, B., Hua, X., Thompson, P. M., Waldorp, L., and Nichols, T. E. (2014). Fast and accurate modelling of longitudinal and repeated measures neuroimaging data. Neuroimage 94, 287–302. doi: 10.1016/j.neuroimage.2014.03.029

PubMed Abstract | CrossRef Full Text | Google Scholar

Guimond, A., Meunier, J., and Thirion, J.-P. (2000). Average brain models: a convergence study. Comput. Vis. Image Understand. 77, 192–210. doi: 10.1006/cviu.1999.0815

CrossRef Full Text | Google Scholar

Iglesias, J., Liu, C.-Y., Thompson, P., and Tu, Z. (2011). Robust brain extraction across datasets and comparison with publicly available methods. Med. Imaging IEEE Trans. 30, 1617–1634. doi: 10.1109/TMI.2011.2138152

PubMed Abstract | CrossRef Full Text | Google Scholar

Jack, C., Shiung, M., Gunter, J., O'Brien, P., Weigand, S., Knopman, D., et al. (2004). Comparison of different MRI brain atrophy rate measures with clinical disease progression in AD. Neurology 62, 591–600. doi: 10.1212/01.wnl.0000110315.26026.ef

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, M., Bannister, P., Brady, M., and Smith, S. (2002). Improved optimization for the robust and accurate linear registration and motion correction of brain images. Neuroimage 17, 825–841. doi: 10.1006/nimg.2002.1132

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, M., Beckmann, C. F., Behrens, T. E., Woolrich, M. W., and Smith, S. M. (2012). FSL. Neuroimage 62, 782–790. doi: 10.1016/j.neuroimage.2011.09.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, M., and Smith, S. (2001). A global optimisation method for robust affine registration of brain images. Med. Image Anal. 5, 143–156. doi: 10.1016/S1361-8415(01)00036-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Joshi, S., and Miller, M. (2000). Landmark matching via large deformation diffeomorphisms. Image Proces. IEEE Trans. 9, 1357–1370. doi: 10.1109/83.855431

PubMed Abstract | CrossRef Full Text | Google Scholar

Klein, A., Andersson, J., Ardekani, B. A., Ashburner, J., Avants, B., Chiang, M.-C., et al. (2009). Evaluation of 14 nonlinear deformation algorithms applied to human brain MRI registration. Neuroimage 46, 786–802. doi: 10.1016/j.neuroimage.2008.12.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Lorenzi, M., Ayache, N., Frisoni, G., and Pennec, X. (2011). “Mapping the effects of Aβ₁₋₄₂ levels on the longitudinal changes in healthy aging: hierarchical modeling based on stationary velocity fields,” in Medical Image Computing and Computer-Assisted Intervention – MICCAI 2011, Vol. 6892 of Lecture Notes in Computer Science, eds G. Fichtinger, A. Martel, and T. Peters (Berlin; Heidelberg: Springer), 663–670.

Lorenzi, M., Ayache, N., Frisoni, G., and Pennec, X. (2013). LCC-demons: A robust and accurate symmetric diffeomorphic registration algorithm. Neuroimage 81, 470–483. doi: 10.1016/j.neuroimage.2013.04.114

PubMed Abstract | CrossRef Full Text | Google Scholar

Lorenzi, M., and Pennec, X. (2013). Geodesics, parallel transport & one-parameter subgroups for diffeomorphic image registration. Int. J. Comput. Vis. 105, 111–127. doi: 10.1007/s11263-012-0598-4

CrossRef Full Text | Google Scholar

Mahapatra, D. (2012). Skull stripping of neonatal brain MRI: using prior shape information with graph cuts. J. Digit. Imaging 25, 802–814. doi: 10.1007/s10278-012-9460-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Marcus, D. S., Fotenos, A. F., Csernansky, J. G., Morris, J. C., and Buckner, R. L. (2010). Open access series of imaging studies: longitudinal MRI data in nondemented and demented older adults. J. Cogn. Neurosci. 22, 2677–2684. doi: 10.1162/jocn.2009.21407

PubMed Abstract | CrossRef Full Text | Google Scholar

McCormick, M. M., Liu, X., Ibanez, L., Jomier, J., and Marion, C. (2014). ITK: enabling reproducible research and open science. Front. Neuroinform. 8:13. doi: 10.3389/fninf.2014.00013

PubMed Abstract | CrossRef Full Text | Google Scholar

Mikheev, A., Nevsky, G., Govindan, S., Grossman, R., and Rusinek, H. (2008). Fully automatic segmentation of the brain from T1-weighted MRI using bridge burner algorithm. J. Magn. Reson. Imaging 27, 1235–1241. doi: 10.1002/jmri.21372

PubMed Abstract | CrossRef Full Text | Google Scholar

Nature (2013). Announcement: reducing our irreproducibility. Nature 496:398. doi: 10.1038/496398a

CrossRef Full Text

Niethammer, M., Huang, Y., and Vialard, F.-X. (2011). “Geodesic Regression for Image Time-Series,” in Medical Image Computing and Computer-Assisted Intervention – MICCAI 2011: 14th International Conference, Toronto, Canada, September 18-22, 2011, Proceedings, Part II (Berlin; Heidelberg: Springer).

Google Scholar

Parker, J., Kenyon, R. V., and Troxel, D. (1983). Comparison of interpolating methods for image resampling. Med. Imaging IEEE Trans. 2, 31–39.

PubMed Abstract | Google Scholar

Prastawa, M., Gilmore, J., Lin, W., and Gerig, G. (2004). “Automatic segmentation of neonatal brain MRI,” in Medical Image Computing and Computer-Assisted Intervention – MICCAI 2004, Vol. 3216 of Lecture Notes in Computer Science, eds C. Barillot, D. Haynor, and P. Hellier (Berlin; Heidelberg: Springer), 10–17.

Reuter, M., and Fischl, B. (2011). Avoiding asymmetry-induced bias in longitudinal image processing. Neuroimage 57, 19–21. doi: 10.1016/j.neuroimage.2011.02.076

PubMed Abstract | CrossRef Full Text | Google Scholar

Reuter, M., Schmansky, N. J., Rosas, H. D., and Fischl, B. (2012). Within-subject template estimation for unbiased longitudinal image analysis. Neuroimage 61, 1402–1418. doi: 10.1016/j.neuroimage.2012.02.084

PubMed Abstract | CrossRef Full Text | Google Scholar

Ridgway, G., Leung, K., and Ashburner, J. (2015). Computing brain change over time. Brain Mapp. 1, 417–428. doi: 10.1016/B978-0-12-397025-1.00313-4

CrossRef Full Text | Google Scholar

Rohrer, J. D., Caso, F., Mahoney, C., Henry, M., Rosen, H. J., Rabinovici, G., et al. (2013). Patterns of longitudinal brain atrophy in the logopenic variant of primary progressive aphasia. Brain Lang. 127, 121–126. doi: 10.1016/j.bandl.2012.12.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Scahill, R., Schott, J., Stevens, J., Rossor, M., and Fox, N. (2002). Mapping the evolution of regional atrophy in alzheimer's disease: unbiased analysis of fluid-registered serial MRI. Proc. Natl. Acad. Sci. U.S.A. 99, 4703–4707. doi: 10.1073/pnas.052587399

PubMed Abstract | CrossRef Full Text

Schott, J. M., Price, S. L., Frost, C., Whitwell, J. L., Rossor, M. N., and Fox, N. C. (2005). Measuring atrophy in Alzheimer disease: a serial MRI study over 6 and 12 months. Neurology 65, 119–124. doi: 10.1212/01.wnl.0000167542.89697.0f

PubMed Abstract | CrossRef Full Text | Google Scholar

Shattuck, D. W., Sandor-Leahy, S. R., Schaper, K. A., Rottenberg, D. A., and Leahy, R. M. (2001). Magnetic resonance image tissue classification using a partial volume model. Neuroimage 13, 856–876. doi: 10.1006/nimg.2000.0730

PubMed Abstract | CrossRef Full Text | Google Scholar

Sled, J., Zijdenbos, A., and Evans, A. (1998). A nonparametric method for automatic correction of intensity nonuniformity in MRI data. Med. Imaging IEEE Trans. 17, 87–97. doi: 10.1109/42.668698

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, S. M. (2002). Fast robust automated brain extraction. Hum. Brain Mapp. 17, 143–155. doi: 10.1002/hbm.10062

PubMed Abstract | CrossRef Full Text | Google Scholar

Stefanescu, R., Commowick, O., Malandain, G., Bondiau, P.-Y., Ayache, N., and Pennec, X. (2004). “Non-rigid atlas to subject registration with pathologies for conformal brain radiotherapy,” in Medical Image Computing and Computer-Assisted Intervention – MICCAI 2004, Vol. 3216 of Lecture Notes in Computer Science, eds C. Barillot, D. Haynor, and P. Hellier (Berlin; Heidelberg: Springer), 704–711.

Südmeyer, M., Pieperhoff, P., Ferrea, S., Krause, H., Groiss, S. J., Elben, S., et al. (2012). Longitudinal deformation-based morphometry reveals spatio-temporal dynamics of brain volume changes in patients with corticobasal syndrome. PLoS ONE 7:e41873. doi: 10.1371/journal.pone.0041873

PubMed Abstract | CrossRef Full Text | Google Scholar

Tustison, N. J., Avants, B. B., Cook, P. A., Zheng, Y., Egan, A., Yushkevich, P. A., et al. (2010). N4ITK: Improved N3 bias correction. IEEE Trans. Med. Imaging 29, 1310–1320. doi: 10.1109/TMI.2010.2046908

PubMed Abstract | CrossRef Full Text | Google Scholar

Vercauteren, T., Pennec, X., Perchant, A., and Ayache, N. (2008). “Symmetric log-domain diffeomorphic registration: a Demons-based approach,” in Medical Image Computing and Computer-Assisted Intervention – MICCAI 2008, Vol. 5241 of Lecture Notes in Computer Science, eds D. Metaxas, L. Axel, G. Fichtinger, and G. Sz¨kely (Berlin; Heidelberg: Springer), 754–761.

Yushkevich, P. A., Avants, B. B., Das, S. R., Pluta, J., Altinay, M., and Craige, C. (2010). Bias in estimation of hippocampal atrophy using deformation-based morphometry arises from asymmetric global normalization: an illustration in ADNI 3 T MRI data. Neuroimage 50, 434–445. doi: 10.1016/j.neuroimage.2009.12.007

PubMed Abstract | CrossRef Full Text | Google Scholar