Fuzzy anatomical connectedness of the brain using single and multiple fibre orientations estimated from diffusion MRI

Computerized Medical Imaging and Graphics 34 (2010) 504–513

Contents lists available at ScienceDirect

Computerized Medical Imaging and Graphics journal homepage: www.elsevier.com/locate/compmedimag

Fuzzy anatomical connectedness of the brain using single and multiple fibre orientations estimated from diffusion MRI Stamatios N. Sotiropoulos a , Li Bai b , Christopher R. Tench a,∗ a b

Division of Clinical Neurology, University Hospital, University of Nottingham, UK School of Computer Science, University of Nottingham, UK

a r t i c l e

i n f o

Article history: Received 7 April 2009 Received in revised form 22 July 2009 Accepted 20 August 2009 Keywords: Tractography White matter DTI Q-ball Distributed connectivity Fibre crossing ODF

a b s t r a c t A new fuzzy algorithm for assessing white matter connectivity in the brain using diffusion-weighted magnetic resonance images is presented. The proposed method considers anatomical paths as chains of linked neighbouring voxels. Links between neighbours are assigned weights using the respective fibre orientation estimates. By checking all possible paths between any two voxels, a connectedness value is assigned, representative of the weakest link of the strongest path connecting the voxel pair. Multiple orientations within a voxel can be incorporated, thus allowing the utilization of fibre crossing information, while fibre branching is inherently considered. Under the assumption that paths connected strongly to a seed will exhibit adequate orientational coherence, fuzzy connectedness values offer a relative measure of path feasibility. The algorithm is validated using simulations and results are shown on diffusion tensor and Q-ball images. Crown Copyright © 2009 Published by Elsevier Ltd. All rights reserved.

1. Introduction In brain white matter (WM), water molecules diffuse faster along than across the underlying neuronal tracts [1]. These preferred diffusion orientations can be calculated using a set of diffusion-weighted (DW) magnetic resonance images [2] and can be used as estimates of WM fibre orientations in each image voxel [1]. By utilizing the orientation estimates through the whole brain, tractography methods are able to reconstruct WM fibre bundles and offer anatomical connectivity information non-invasively and in vivo [3–6]. Early work on tract reconstruction propagated curves, known as streamlines, that follow the fibre orientation estimates [3–6]. Despite its success delineating many WM tracts in the human brain [4,7], streamline tractography has two major limitations: (a) it does not account for fibre branching and (b) there is no intrinsic way to characterize the likelihood of a putative tract. Hard binary relationships are utilized in streamline tractography to assess the connections, with two voxels being considered either connected or not connected. However, medical images are by nature fuzzy

∗ Corresponding author at: Division of Clinical Neurology, B Floor, Medical School, University Hospital, University of Nottingham, Nottingham NG7 2UH, UK. Tel.: +44 115 823 1192; fax: +44 115 970 9738. E-mail addresses: [email protected] (S.N. Sotiropoulos), [email protected] (L. Bai), [email protected] (C.R. Tench).

[8] and DW images are no exception. Experimental noise, hardware limitations, limited spatial resolution, and partial volume artefacts are some of the factors that contribute to the fuzziness of the images. Therefore, a fuzzy framework seems more appropriate. The benefit of using fuzzy logic in tract segmentation has been presented in [9]. In this study, we present a fuzzy algorithm for assessing anatomical connectivity. Different methods have been proposed to overcome the limitations of streamline tractography. Probabilistic techniques [10–16] treat the fibre orientation as a random variable and define a distribution for it. In [10] the distribution is predicted using simulations, in [11,14] it is estimated using Monte-Carlo Markov Chain (MCMC) methods and in [16] it is sampled using particle filtering. A random walk model [12], non-parametric bootstrap techniques [13], and direct Bayesian inference [15] have also been utilized to characterize the fibre orientation uncertainty. Streamlines are then generated in a Monte-Carlo fashion and in each propagation step a sample is drawn randomly from the local orientation distribution. Many streamlines are generated from a seed point and the fraction of streamlines that pass through a voxel is defined as its index of connectivity to the seed. A drawback of these techniques is the inevitable reduction of the connectivity index with distance from the seed [10,11]. Moreover, the repetitive streamline generation increases execution time, while the connectivity values can depend on the total number of streamlines launched. All the previous methods, deterministic and probabilistic, determine white matter bundles using locally greedy criteria, i.e.

0895-6111/$ – see front matter. Crown Copyright © 2009 Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.compmedimag.2009.08.006

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

tracking sequentially through orientation estimates in adjacent voxels. Another way of tract reconstruction is to identify the optimal path between two voxels of interest, according to some global criterion, rather than identifying paths arising from a single voxel. These global approaches are less sensitive to local artefacts introduced by various noise sources [17]. Front propagation techniques [18–22] are an example. They identify a best path from a seed to all other brain voxels by evolving a surface from this seed. The surface front evolves faster along the fibre orientation estimates of the underlying tissue. Surface evolution is performed using the fast marching algorithm in [18–20], an empirical front propagation method in [21] and an anisotropic wavefront evolution algorithm in [22]. A distributed index of connectivity can then be defined using either functions of the front arrival times in each voxel or the agreement between the surface trajectory and the underlying fibre orientation estimates. Distributed indices of connectivity are also calculated by graph-based tractography techniques [23,24]. These methods utilize graph theory to exhaustively search the image for the strongest path connecting any two voxels of interest. In this study, we introduce a new global tractography method using a fuzzy framework for assessing brain anatomical connectivity. We adapt fuzzy connectedness [8,25], an algorithm initially proposed for image segmentation, to the tractography problem. Anatomical paths are considered as chains of linked neighbouring voxels. A connectedness value is assigned between any voxel and a seed of interest by checking all the possible paths that connect them on the discrete image grid and finding the strongest. The connectedness value reflects the weakest link along this strongest path. It quantifies the orientational coherence along a path and provides a relative measure for path feasibility, under the condition that fibre orientations change smoothly along an anatomical path. Therefore, we expect that paths exhibiting high orientational coherence will acquire high connectedness values relative to the background. Compared to other tractography methods, our algorithm provides converged connectivity values that do not drop systematically with the distance from the seed, in a single iteration and for all image voxels. Furthermore, path propagation is performed with relatively high angular resolution and combined in a single step with connectivity assignment. Algorithm implementation with dynamic programming results in fast execution. The new method inherently accounts for fibre branching, and fibre crossing information can be incorporated, extending the preliminary work presented in [26,27]. To reduce the chance of false positives when propagating through crossing regions, multiple links are considered – when appropriate – between neighbouring voxels to account for the existence of multiple tracts. Our method can be combined with any DW magnetic resonance imaging (MRI) reconstruction method that provides N ≥ 1 fibre orientation estimates. We present results in simulated data, as well as, in human data obtained using diffusion tensor imaging (DTI) [28] and Q-ball imaging [29]. The former allows estimation of a single fibre orientation in each voxel, while the latter estimates multiple orientations when partial volume is present.

2. Methods 2.1. Fuzzy connectedness framework Given a three-dimensional digital image, we define a local fuzzy relation between any two neighbouring voxels i and j of the image called affinity (i, j) (: Z3 × Z3 → [0,1]) [25]. The affinity is determined by the similarity of some image intensity-based features at voxels i and j and it is zero between non-neighbouring voxels. The fuzzy connectedness algorithm (FCA) utilizes the local affinity values to assign strengths FC(a,b) of a global fuzzy relation

505

called fuzzy connectedness, between any pair of voxels a and b. This assignment is performed by checking all the possible paths that connect the two voxels on the discrete image grid and finding the strongest. There are many possible paths that connect a and b. Each path can be considered a chain of voxels, with successive elements being neighbouring. In this framework, the smallest affinity  along a path, i.e. the weakest link of the chain, determines the path strength. Considering all possible paths connecting a and b, the strength of the strongest will be the connectedness FC(a,b). Therefore, after selecting a seed voxel or region of interest (ROI) and given an affinity function, a connectedness value FC(a,s) between the seed s and every other image voxel a can be computed. In the following sections, we will refer to this connectedness value as FC(a) for simplicity, since there is always a seed s. The fuzzy connectedness algorithm can be implemented using dynamic programming, which results in execution at almost interactive speeds [30]. 2.2. Fuzzy connectedness tractography To apply FCA to WM tractography, two main aspects of the algorithm are considered and modified; the affinity function and the path generation procedure. The affinity is determined by the similarity of the fibre orientation estimates between neighbouring voxels, as derived from the intensities of the DW images. More details on the affinity function are presented in the next section. Regarding path generation, FCA takes into account all possible links that are legal voxel connections on the discrete image grid. However, this does not guarantee that the paths are anatomically realistic. Being an image segmentation algorithm, FCA tends to grow both laterally and longitudinally; ideally, we would like to restrict propagation along the WM fibre orientations. Furthermore, 180◦ turns are not inherently avoided; such turns are possible given that diffusion MRI provides orientations and not directions. A modification is needed to force forward path propagation. We, therefore, introduce a memory property for each voxel. Once a voxel a is identified as part of a path, its previous voxel c in the path chain, through which a is connected to the seed in the strongest way, is stored. Knowing the vector rca connecting voxel c to a, we can constrain propagation only to directions n, so that rca ·n > 0, prohibiting lateral or backward path propagation. Seed voxels, i.e. voxels without memory, are treated differently. All neighbours, both forward and backward, are considered for the seeds to initiate propagation in both directions. The fuzzy connectedness tractography (FCT) algorithm is described in detail in Appendix A. In the following sections, different aspects of the algorithm are presented. 2.3. Affinity and connectedness based on fibre orientation estimates Assuming that anatomical paths are smooth curves, fibre orientation estimates along a path are expected to vary smoothly. Therefore, we assess the affinity between neighbouring voxels i and j, using a coherence condition on the orientation estimates obtained from diffusion MRI. Our method can utilize N fibre orientations per voxel, thus allowing fibre crossing information to be incorporated. The affinity function we use is a modified version of the symmetric voxel linking function suggested in [19,31]. It is defined between the lth orientation of voxel i and the kth orientation of j, as: lk (i, j) = kl (j, i) =

1 C · (1 − min(|el (i) · n(i, j)|, |ek (j) · n(i, j)|, |el (i) · ek (j)|))

(1)

where el (i) is the lth fibre orientation vector at voxel i, n(i,j) is the unit vector connecting the centres of voxels i and j and C a

506

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

normalization constant to keep affinities in the [0,1] interval. C is chosen as the maximum of the non-normalized affinity values across the whole 3D image. The above function is peaked when all three vectors are collinear, i.e. when the lth orientation vector of voxel i points to the centre of voxel j and the kth orientation vector of voxel j points to the centre of voxel i. In this case, the affinity as presented in Eq. (1) becomes infinite, which is handled in practice by multiplying the min function by a number slightly less than 1. If L (1 ≤ L ≤ N) orientations exist in voxel i and K (1 ≤ K ≤ N) orientations in voxel j, the indices l and k range from 1 to L and 1 to K respectively, giving L × K affinity values between the two voxels. To account for the multiple affinity values between two voxels, we treat each of the orientations per voxel as a separate element in the fuzzy connectedness framework. Then, the structural element of the algorithm becomes (i, l), the lth orientation of voxel i, rather than voxel i itself. Therefore, FCT does not compute a scalar connectedness value FC(i) for a voxel i, but a vector FC(i) with length L; lth vector entry FCl (i) corresponds to the FC value of the lth orientation. The final FC value for the specific voxel i will be given by the maximum of these L values. 2.4. Increasing angular resolution Affinities are calculated between each image voxel i and its immediate neighbours. The neighbourhood size effectively determines the number of candidate propagation directions n. A straightforward implementation uses the 3 × 3 × 3 neighbourhood, i.e. 26 neighbours in 3D. We can increase the angular resolution of path propagation by increasing the neighbourhood size to 5 × 5 × 5, including next nearest neighbours, as shown in Fig. 1. This gives 124 neighbours in 3D. In this case, though, neighbours that are not directly connected to the focal voxel i exist, such as voxel j in Fig. 1b. To avoid discontinuities in paths, we require that at least one nearest neighbour of the focal voxel is included in the path. More specifically, for each next nearest neighbour j (light gray in Fig. 1b) we consider the two nearest neighbours c and d (darker gray in Fig. 1b) that are traversed by the vector connecting the focal voxel i to j. We then choose between c and d, the one with the highest connectedness FCh value to the seed. If this connectedness is smaller than the connectedness value of j and propagation is decided from i to j, we assign FCh = FC(j). Implementation details are included in Appendix A.

Fig. 1. Different neighbourhood systems of focal voxel i, 3 × 3 × 3 in (a) and 5 × 5 × 5 in (b). Only neighbours within the same slice as i are presented. Black arrows represent candidate directions for path propagation. Dark gray: nearest neighbours of voxel i, light gray: next nearest neighbours of voxel i.

2.5. Simulations A phantom containing a circular tract (Fig. 2a) and another containing two curved crossing tracts (Fig. 2b) were numerically generated. The tracts were embedded within a diffusively isotropic surround, to simulate gray matter (GM), where fibre orientations are random. The noise-free DW signals in each voxel were simulated as described in [32], using a Gaussian mixture model: Sk = So ·

N 

fm exp(−bg Tk Dm g k )

(2)

m=1

where Sk is the image intensity after the application of the kth diffusion-sensitizing magnetic field gradient, gk is the direction vector of the kth diffusion-sensitizing gradient and b depends on the magnitude and duration of the gradient. So is the image intensity at the same voxel without the application of any diffusionsensitizing gradient, N is the number of tracts simulated and Dm is a 3 × 3 diffusion tensor simulating the diffusion profile along the mth tract, having a volume fraction fm [28]. Sixty-one (k = 61) evenly spaced directions [33], b = 1000 s/mm2 for the former and b = 3000 s/mm2 for the latter phantom were used. The number of tracts was set to N = 1 for all voxels, apart from the voxels contained in the crossing region of the second phantom, where N = 2. The tensors in simulated WM tracts had a fractional anisotropy FA = 0.8

Fig. 2. Fuzzy connectedness (FC) maps obtained from simulations (SNR = 20). Connectedness values are shown colour coded. Seed voxels are indicated with white arrows. (a) From left to right: fibre orientation estimates for the circular phantom, FC map obtained using a 3 × 3 × 3 neighbourhood, FC map obtained using a 5 × 5 × 5 neighbourhood and the same seed, Streamline generated from the same seed and corresponding binary connectivity index (white = connected, black = not connected). The black solid line in the orientations plot defines a cross-section used for quantitative simulations (see Fig. 3). (b) From left to right: fibre orientation estimates for the crossing phantom (only crossing region is shown magnified), FC maps obtained using a 5 × 5 × 5 neighbourhood and three different seed points. In both (a) and (b), orientation estimates are superimposed on diffusion anisotropy values, which are high for WM voxels, low for GM voxels and intermediate for voxels in the crossing region.

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

[1] and an orientation tangent to the desired simulated shape at each point. The tensors in GM voxels had an FA = 0.25 and random orientation. All tensors had a trace of 2.1 × 10−3 mm2 /s. Zero-mean Gaussian noise was added in quadrature, as described in [32] to simulate the Rician nature of MRI noise. Given the simulated DW data, we estimated fibre orientation estimates. For the first phantom we used the diffusion tensor imaging (DTI) model [28] that can estimate N = 1 orientation per voxel. DTI fits an ellipsoid to the data and estimates an orientation as the major axis of this ellipsoid. For the second phantom Q-ball imaging (QBI) [29] was utilized. QBI can estimate more than one orientation per voxel, when appropriate, e.g. when fibre crossing occurs. To achieve that a probability density function for the fibre orientations, called the orientation distribution function (ODF), is computed directly from the data in each voxel. The peaks of the ODF give the orientation estimates. We computed the ODFs using even spherical harmonics up to 4th order as basis [34]. The CAMINO diffusion toolkit was then used [35] to find the ODF peaks. A maximum of N = 2 orientations were calculated in each voxel. FCT was applied to both phantoms to test its performance in single (N = 1) and multiple orientation (N = 2) fields. 2.6. MR acquisition and processing A DTI [28] and a Q-ball [29] whole-brain scan of a healthy male subject were performed; local research ethics committee approval and informed consent were obtained. A single-shot, spin-echo, echo-planar, diffusion-weighted sequence was used (acquisition matrix 112 × 112 with in-plane resolution 2 × 2 mm, interpolated during reconstruction to 224 × 224) in a Philips 3T Achieva clinical imaging system (Philips Medical Systems, Best, The Netherlands). A parallel imaging factor of 2 was used. Six b = 0 s/mm2 images were acquired and averaged. Diffusion weighting was applied in 61 evenly spaced directions [33] with b = 1000 s/mm2 for the DTI scan (TE = 57.6 ms, TR = 10,990 ms) and b = 3000 s/mm2 for the Qball scan (TE = 72 ms, TR = 15,292 ms). 52 slices were acquired with a thickness of 2 mm. Total imaging time was roughly 35 min. Images were corrected for eddy current distortion using FSL’s diffusion toolbox [36]. Brain was extracted using BET [36]. Corrected images were then tri-linearly interpolated along the out-of-slice axis to get isotropic voxels with 1 mm × 1 mm × 1 mm dimensions. Fibre orientations were estimated in each case as described in the simulations section. FCT was performed on both the DTI and Q-ball datasets. Voxels with a small degree of diffusion anisotropy (FA < 0.2) were excluded from the analysis, assuming that these represent GM or CSF. For comparison purposes, we performed streamline tractography using a modified FACT algorithm [5] (FA threshold = 0.2, curvature threshold = 45◦ ), implemented in CAMINO [35]. Q-ball orientation estimates were used, so that streamlines could propagate through crossings. We also implemented the distributed graph-based tractography, as described in [23], using the whole QBI ODFs and performed multi-tensor (N = 2) probabilistic tractography [37] using CAMINO [35]. 3. Results 3.1. Simulations Fig. 2a presents fuzzy connectedness maps obtained using FCT on the circular phantom. The gray scale corresponds to FC values between each voxel and to the seed, which is indicated by the white arrow. For comparison, the binary connectivity index obtained from streamline tractography for the same seed is shown. The distributed nature of FCT compared to streamline tractography is evident. FCT identified strong paths arising from the seed,

507

along with other weaker connections. WM voxels with coherent fibre orientations were distinguished from the non-coherent GM background. The effect of increasing angular resolution using more neighbours is also shown in Fig. 2a. As expected, a larger neighbourhood assisted the algorithm in capturing more curvature and increased the connectedness values. The ability of the method to go through fibre crossings and utilize multiple fibre orientations per voxel is shown in Fig. 2b. When seeding within the individual tracts, paths in the appropriate tract were identified. When seeding in the crossing region, paths belonging to both tracts were found. To test the behaviour of FCT against noise, 100 Monte-Carlo simulations were performed using the circular phantom at a given SNR level. A linear ROI, indicated in Fig. 2a (left subfigure) with a black solid line, defined a cross-section of the phantom. The central part of this cross-section (7 voxels) was within the WM circular ring, while its edges (5 voxels each) were within the GM region. The fuzzy connectedness values FC of the voxels along this ROI were studied across the 100 simulations. The seed voxel was the same across simulations, at the location indicated by the white arrow in Fig. 2a. Fig. 3 shows the mean and the standard deviation of the FC values at each of these discrete voxel locations, using the 3 × 3 × 3 and the 5 × 5 × 5 neighbourhoods. The noise-free connectedness values along the ROI are shown by a dashed line, which is continuous just for visualization purposes. The SNR in these simulations was 15, representative of the noise conditions in a DW-MRI acquisition. The plots show the robustness of the calculated values against noise and their good agreement with the noise-free ones. Given that the simulated tract is very curved and that the affinity we use (Eq. (1)) penalizes curvature, it is expected that the connectedness of the stronger paths will be smaller than one. Deviation from this ideal value is reduced in the case of the 5 × 5 × 5 neighbourhood (Fig. 3b). Compared to the values obtained with a 3 × 3 × 3 neighbourhood (Fig. 3a), the connectedness of the voxels belonging to the circular ring is almost doubled. Furthermore, the contrast between the more central and stronger paths and the surrounding weaker paths is enhanced with the 5 × 5 × 5 neighbourhood. Using both neighbourhood sizes, however, FCT differentiated regions that are coherently connected to the seed from the non-coherent ones. 3.2. In vivo DW images The performance of FCT in different brain regions is shown in Figs. 4–7. A 5 × 5 × 5 neighbourhood was used in all cases. Typical execution times for a seed ROI were in the order of 1 min for DTI FCT and 2 min for Q-ball FCT on a 3.2 GHz PC. A horizontal ROI within the body of the corpus callosum seeded the algorithm in Fig. 4. In Fig. 4a, coronal maximum intensity projections (MIP) of the raw connectedness values FC are shown when the DTI and the Q-ball orientation estimates are used, respectively (the intensity of a voxel on a MIP is the maximum intensity at this voxel location across all the slices parallel to the projection plane). For comparison, streamlines generated within the multi-fibre field of the Q-ball reconstruction are presented. The colour-coded connectedness maps had higher intensities in regions that were relatively more likely to belong to the corpus callosum. The distributed nature of the algorithm is evident, since all image voxels are assigned a connectedness value. In Fig. 4b and c, thresholds were set to keep the top connected voxels and remove the relatively homogeneous FC values of the background. The thresholds were determined from the histograms of the FC values, where sharp discontinuities indicated the transition from the background to the top connected voxels (Fig. 5). These thresholds corresponded roughly to the 97.5th percentile of the distribution of the FC values for DTI FCT (Fig. 4b) and to the 95th percentile for Q-ball FCT (Fig. 4c), reflecting the higher coverage of the callosal tracts by the latter method. A MIP

508

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

Fig. 3. Mean (circle) and standard deviation (horizontal bars) of fuzzy connectedness values across 100 simulations (SNR = 15), at voxel locations defined by the linear ROI of Fig. 2a (black solid line). Voxel counting starts from the top left corner of the phantom and increases towards the phantom centre. FC values were obtained using (a) a 3 × 3 × 3 neighbourhood, (b) a 5 × 5 × 5 neighbourhood. Noise-free FC values are indicated by the dashed line.

of thresholded connectedness values, along with the FC values on different coronal slices is presented in each case. DTI FCT (Fig. 4b) captured the medial portion of the corpus callosum, as expected when using the DTI orientations [7]. However, more branching was captured compared to the medial portion resolved from streamline tractography (Fig. 4a). Q-ball FCT (Fig. 4c) went through the crossing at the centrum semiovale and the results captured the lateral projections to the cortex, apart from the medial callosal paths. We should point out that the resolved crossings in Fig. 4c were not perfectly symmetric in the two hemispheres and greater fanning was evident in the left hemisphere. This was caused by asymmetric Q-ball reconstruction in the contralateral regions of the centrum semiovale, possibly due to noise. A similar organization of tractography maps is presented in Fig. 6, with the seed being a coronal ROI in the cingulum. In this case, streamlines (Fig. 6a) were stopped at the level of the splenium of the corpus callosum and did not continue towards the parahippocampal gyrus, probably due to noise artefacts. The distributed nature of FCT, however, allowed the algorithm to continue, with Q-ball FCT (Fig. 6c) having a slightly higher connectedness through that point than DTI FCT (Fig. 6b). Given that there are not many fibre crossings along this tract, only small differences can be observed when single

and multiple orientations per voxel are used (Fig. 6b and c). Q-ball FCT (Fig. 6c) assigned higher values in the subgenual area, where partial volume artefacts distort the DTI orientations. Furthermore, the body of the cingulum was slightly thicker (Fig. 6c) compared to the DTI FCT result (Fig. 6b), since the apparent (due to limited resolution) crossing between the cingulum and the corpus callosum is not resolved by DTI. For the same reason, paths with relatively low connectedness arising from the main body of the cingulum towards the cortex that were evident with DTI FCT were not observed with Q-ball FCT. MIPs of thresholded FC values in other brain regions are shown in Fig. 7. In all cases DTI orientations were utilized, apart from Fig. 7a and c, where Q-ball estimates were used. In Fig. 7a, a saggittal ROI within the medial portion of the superior longitudinal fasciculus was the seed and only Q-ball orientations supported the existence of paths arising from the parietal cortex and ending up in motor cortical regions, in agreement with [14]. In Fig. 7b a single seed voxel was selected in the inferior part of the fornix, close to the amygdaloid complex and in Fig. 7c, a horizontal seed ROI at the upper pons was used. In Fig. 7d, a single voxel in the splenium of the corpus callosum seeded the algorithm. Finally, a seed placed at the apex of Meyer’s loop was found to be strongly connected to

Fig. 4. Fuzzy connectedness maps obtained from a horizontal ROI within the corpus callosum. (a) From left to right: coronal maximum intensity projection (MIP) of the raw FC values using the DTI orientation estimates, coronal maximum intensity projection of the raw FC values using the Q-ball orientation estimates, streamlines generated using the Q-ball orientation estimates, (b) MIP of thresholded FC values using the DTI estimates, along with coronal slices, (c) MIP of thresholded FC values using the Q-ball estimates, along with coronal slices. In all cases, the FC maps are superimposed on diffusion fractional anisotropy (FA) images. Seed locations are indicated by white arrows.

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

Fig. 5. Histogram-based detection of background connectedness threshold. The histogram of the FC values obtained with Q-ball FCT when seeding in the corpus callosum (Fig. 4) is shown. A very sharp transition in the histogram, pointed by the black arrow, was used to identify the background connectedness value.

509

the optic radiation fanning to the occipital cortex and the lateral geniculate nucleus of the thalamus (Fig. 7e). Apart from connectedness maps, paths generated by FCT can be reconstructed using the memory information. The strongest path connecting any voxel to a seed can be found in a backward fashion. Starting from an end point, we could recursively find the previous voxels and fibre orientations within those voxels that led to the seed. By setting a threshold to the FC values, we could plot the paths that corresponded to the top connected voxels. Fig. 8 presents paths belonging to different brainstem tracts reconstructed this way. During the backward propagation, a moving average filter was utilized to smooth the trajectories. The tracts of interest were the middle, superior and inferior cerebellar peduncle, the medial lemniscus and the corticospinal tract. FCT was applied individually for each of the five tracts, using one seed voxel per tract, and only the paths arising from the top 1% connected voxels were plotted (for the medial lemniscus and superior cerebellar peduncle a smaller threshold was used due to the smaller size of these tracts). These paths agree well with a-priori anatomical knowledge of the brainstem [7] and verify the anatomically realistic propagation of our algorithm. FCT results were compared against the results of two other quantitative tractography algorithms, the latest graph-based algorithm [23] and the multi-tensor (N = 2) probabilistic tractography [37]. The graph-based method is a global tractography approach

Fig. 6. Fuzzy connectedness maps obtained from a coronal ROI within the cingulum. (a) From left to right: sagittal maximum intensity projection (MIP) of the raw FC values using the DTI orientation estimates, sagittal maximum intensity projection of the raw FC values using the Q-ball orientation estimates, streamlines generated using the Q-ball estimates, (b) MIP of thresholded FC values using the DTI estimates, along with horizontal slices, (c) MIP of thresholded FC values using the Q-ball estimates, along with horizontal slices. In all cases the FC maps are superimposed on diffusion fractional anisotropy (FA) images. Seed locations are indicated by white arrows.

Fig. 7. Fuzzy connectedness maps obtained using various seeds in the brain. Maximum intensity projections of thresholded FC values are shown in each case. The seeds are placed (a) at the superior longitudinal fasciculus, (b) at the inferior part of the fornix, (c) at the upper pons within the pyramidal tract, (d) at the splenium of the corpus callosum, and (e) at the apex of Meyer’s loop. DTI orientations are used in all cases, apart from (a) and (c) where Q-ball estimates are utilized. Seed locations are indicated by white arrows.

510

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

hemisphere. On the other hand, FCT propagated through the crossing bilaterally and assigned relatively high values only to corpus callosum paths. This is indicative of the algorithm’s ability to propagate correctly through crossing configurations. We should point out that apart from using the same data for all three algorithms, FCT and GT also utilized exactly the same diffusion ODFs. However, FCT is benefited from the utilization of multiple affinities (i.e. links) between neighbours in crossing regions that reflect the existence of multiple fibre populations. 4. Discussion

Fig. 8. Reconstruction of paths in the brainstem using FCT. Different colours are used for paths belonging to different WM tracts. Five different seed voxels were used at the middle cerebellar peduncle (red), at the right inferior cerebellar peduncle (orange), at the right superior cerebellar peduncle (green), at the right medial lemniscus (yellow) and at the right corticospinal tract (blue). The seed voxels are indicated with a sphere. Only the top paths connected to each of the seeds are shown. DTI orientation estimates were used. The right thalamus rendered and in white colour and a horizontal slice of the FA map are shown for anatomical reference. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)

similar to FCT, so a comparison between them is more direct. Multitensor probabilistic tractography is a local tracking approach, but it was employed due to the popularity of probabilistic streamline techniques [17]. All methods utilized the Q-ball data and a horizontal ROI within the body of the corpus callosum as a seed. Maximum intensity projections of the FC values, of the path strengths obtained from graph-based tractography (GT) and of the index of connectivity from probabilistic tractography are presented in Fig. 9. We can observe that GT (Fig. 9b) resolves lateral portions of the corpus callosum (indicated by black arrows), but the strengths assigned to them are relatively low. Furthermore, paths belonging to the crossing corona radiata are assigned a high strength, meaning that propagation through the crossing at the centrum semiovale has given rise to false positive connections. Multi-tensor probabilistic tractography (Fig. 9c) improves on GT results, since the lateral callosal tracts on the left hemisphere can be differentiated from the background. However, a higher index of connectivity is assigned to the crossing corona radiata tracts on the right hemisphere compared to the right callosal tracts. Furthermore, the inherent reduction of connectivity with distance from the seed due to the repetitive Monte-Carlo sampling is evident, especially in the right

We have presented a new fuzzy framework to assess white matter anatomical connectivity in the brain. Fuzzy connectedness tractography is a global tracking approach that utilizes any number of fibre orientations per voxel to produce distributed connectedness maps, so that every image voxel is eventually connected to the seed through some path. A connectedness value for a voxel characterizes how feasible the path connecting it to the seed is, given an affinity function. This value is computed by trying all possible valid connections between a voxel and the seed and finding the strongest. The smallest affinity, i.e. the weakest link, along this strongest path will be the connectedness of the voxel to the seed. The affinity employed in this study reflects the anatomical assumption of slowly varying fibre orientations along WM paths. Thus, the computed connectedness values provide a relative measure for path feasibility that depends on the orientational coherence along a path. The algorithm discriminates between voxels that exhibit orientational coherence to the seed and the background. Apart from the fuzzy nature of the algorithm, our algorithm expands other distributed tractography methods [18–24] in two ways: (a) FCT propagates paths using a relatively high angular resolution and (b) it can utilize any number of fibre orientations per voxel by treating each orientation as a different voxel instance. It can thus be combined with any DW reconstruction technique that provides orientation estimates, such as multi-tensor models [38], spherical deconvolution [39], persistent angular structure [40] and diffusion spectrum imaging [41]. In this study, we used Q-ball imaging as a working example to estimate multiple orientations per voxel. Q-ball imaging has been combined with other types of distributed tractography, such as fast marching [18] and graph-based tractography [23]. However, both previous studies define a single affinity value between any two voxels, utilizing the whole Q-ball ODF, regardless of the number of fibre populations being present in each voxel. In the case of a fibre crossing, this may allow propagation towards all crossing orientations and artificially elevate the connectivity index of crossing tracts over paths belonging to the tract of interest, as shown in Fig. 9b and explained in detail in [42]. Our approach, defines – if needed – many affinities between two voxels and propagates more strongly along the appropriate orientation in the case of crossings.

Fig. 9. (a) Coronal maximum intensity projection of the raw FC values using the Q-ball FCT and a horizontal seed ROI in the body of the corpus callosum. (b) Coronal maximum intensity projection of the path strengths calculated using graph-based tractography and the Q-ball ODFs to define graph weights. (c) Coronal maximum intensity projection of the index of connectivity calculated using two-tensor probabilistic tractography. Colour coding is qualitatively the same in all cases with dark red being low and white being high. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

FCT produced reasonable connectedness maps in many different WM regions (Figs. 4–7). Given that validation of a tractography algorithm is in general very difficult and impossible on a single subject basis, we performed advanced streamline tractography using Q-ball orientation estimates and the same seed ROIs; this at least shows that the FCT results are commensurate with the, now familiar from the literature, streamline results for our data [4,7]. Connectedness maps expanded the connectivity picture from streamline tractography (Figs. 4–6), exhibiting greater branching. Even when only the DTI orientations were used, the results seem more distributed than the Q-ball streamlines, as shown in the cingulum (Fig. 6). We should point out, that our algorithm operates on the discrete image grid and connects voxels centre-to-centre, in order to identify the strongest path between two voxels, according to a global criterion. Streamline tractography can move anywhere within a voxel using locally greedy criteria. This can also contribute to the differences observed between the two techniques. Furthermore, when multiple orientations are available in a single voxel, streamline tractography selects one of the available orientations. This might be a good option for clear cases of fibre crossings, but it may cause an underestimation of the underlying structure when fibre branching occurs. FCT will consider all available orientations and assign a connectedness value to each of them, depending on how coherently connected they are to the seed. Therefore, it will capture more branching. Compared to probabilistic streamline tractography approaches [10–15], FCT assigns a fuzzy rather than a probabilistic index of connectivity. The latter gives the number of times a voxel has been transversed by a tract over repetitive streamline propagations in a Monte-Carlo fashion and is inherently sensitive to the distance from the seed (Fig. 9c). The former quantifies the fuzziness of paths and provides a relative measure of connectivity, given an affinity function. The affinity reflects the anatomical assumption of smoothness in WM paths and quantifies the coherence of adjacent fibre orientations. It can be influenced by factors that affect the estimation of orientations, such as limited spatial resolution, noise, diffusion profile reconstruction limitations and partial volume artefacts, and these factors can contribute to the observed spread of the FC values. Clearly, fuzzy and probabilistic indices are different and their interpretation in terms of absolute anatomical connectivity remains an open question for further investigation. Nevertheless, the fuzzy framework presented here extends the binary framework of streamline tractography, where only two outcomes are allowed (“not connected” vs. “connected”), and relates to it more naturally rather than to probabilistic approaches. Similar to all tractography algorithms, FCT can be influenced by unresolved partial volume artefacts. For example, false positives are evident in Fig. 7d, where paths of the inferior longitudinal fasciculus are shown connected to Meyer’s loop. However, in most cases these artefacts were characterized by relatively lower connectedness values. We should point out the difference between FCT and a fuzzy segmentation algorithm, such as the one presented in [9]. Our target (actually the target of any tractography algorithm) was to reconstruct anatomically realistic paths and scale them in terms of their connectedness to a seed rather than delineate and segment WM tracts. There can be voxels belonging to the same tract that are not physically connected, for example voxels along two parallel fibre bundles. Using the results of tractography and carefully chosen seeds, tract segmentation can be performed as a post-processing step. Our approach is governed by certain assumptions, which should be pointed out. For obtaining orientation estimates from Q-ball images, we have used the local ODF maxima, an approach followed in many recent studies (e.g. [43]). This approach has been shown to

511

introduce bias in the estimated orientations, due to the finite width of each peak [44]. However, obtaining orientation estimates from the peaks of the ODF has been utilized as an example. As pointed out before, FCT will work with any method that provides fibre orientations, such as the tensor decomposition approach described in [44] that does not bias the estimates. Furthermore, orientational uncertainty, which is the basis for probabilistic tractography [10–15] is not explicitly taken into account in this study. Some of the uncertainty is indirectly considered as fuzziness, but voxel similarity is assessed using strictly the fibre orientation estimates. The algorithm is independent of the affinity function, therefore more advanced affinities that incorporate uncertainty can be used. We are moving towards this direction, as shown recently in [42]. The discreteness of the algorithm both in the spatial and the orientational domain may potentially cause problems, especially when very sharp turns are present within the range of a very few voxels. To overcome these situations, we used an interpolated image grid and we also increased the angular resolution of path propagation (Fig. 1). As shown, the results were satisfactory in many WM regions. A future extension of the current framework will consider more than one location within a voxel, so that not only centre–centre connections exist. A further remedy will be the incorporation of curvature in the affinity function, which currently penalizes deviations from straight-line connections. Identifying the background connectedness intensity and keeping only the relevant paths is a challenging issue in distributed tractography techniques. In this study, we used a histogram-driven threshold to choose the top connected voxels to the seed. More advanced methods that involve hypothesis testing [23,45] may be used. The fact that the connectedness values are based on the weakest link along a path can make them sensitive to noise. However, unless a large neighbourhood is contaminated, the algorithm can go around noise artefacts due to its distributed nature and global tracking criteria it uses. Furthermore, variations of the fuzzy connectedness algorithm exist with increased robustness against noise. Scale-based fuzzy connectedness [46] is an example and scale-based FCT can be a direct extension of our algorithm. In summary, we have presented a fuzzy distributed algorithm to assess white matter connectivity using orientation estimates obtained from diffusion MRI. Fuzzy connectedness tractography can use any number of fibre orientations per voxel and can therefore be combined with any model-based or model-free reconstruction technique that provides such estimates. In this work, we illustrated its performance using two popular methods, diffusion tensor and Q-ball imaging. The algorithm checks all valid paths on the discrete image grid that connect two voxels of interest and assigns a connectedness value between them, representative of the strength of the strongest path. Under the assumption that WM tracts exhibit orientational coherence, it can differentiate between regions that are highly connected to a seed and paths that are not. Its flexibility in the way neighbouring voxel similarity is assessed, allows the incorporation of different affinity functions and image modalities that may even be application-specific.

Acknowledgements This study was funded by the European Commission Fp6 Marie Curie Action Programme (MEST-CT-2005-021170) under the CMIAG project. We would like to thank Dr Paul S. Morgan from the Division of Academic Radiology at the University of Nottingham, UK for generously providing us the diffusion images. Also Prof. Jayaram K. Udupa from the Department of Radiology at the University of Pennsylvania, USA for stimulating discussions.

512

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513

Appendix A. The FCT algorithm is based on the label-setting fuzzy connectedness algorithm presented in [30]. A priority queue is used as an auxiliary data structure. Each queue element (i, l) comprises of the coordinates of voxel i and an index l identifying which of the L orientations of voxel i the element refers to. For each of the queue elements, a fuzzy connectedness value FCl (i) is stored and elements are sorted so that the one with the maximum connectedness value is always the first element of the queue. Neighbouring elements are sorted, so that nearest neighbours are considered first and next nearest neighbours follow when affinities are calculated. The algorithm is presented as comprised of two functions:

References [1] Basser PJ. Inferring microstructural features and the physiological state of tissues from diffusion-weighted images. NMR Biomed 1995;8:333–44. [2] Turner R, Le Bihan D, Maier J, Vavrek R, Hedges LK, Pekar J. Echo-planar imaging of intravoxel incoherent motion. Radiology 1990;177:407–14. [3] Basser PJ, Pajevic S, Pierpaoli C, Duda J, Aldroubi A. In vivo fiber tractography using DT-MRI data. Magn Reson Med 2000;44:625–32. [4] Catani M, Howard RJ, Pajevic S, Jones DK. Virtual in vivo interactive dissection of white matter fasciculi in the human brain. Neuroimage 2002;17:77–94. [5] Mori S, Crain BJ, Chacko VP, van Zijl PC. Three-dimensional tracking of axonal projections in the brain by magnetic resonance imaging. Ann Neurol 1999;45:265–9. [6] Conturo TE, Lori NF, Cull TS, Akbudak E, Snyder AZ, Shimony JS, et al. Tracking neuronal fiber pathways in the living human brain. Proc Natl Acad Sci USA 1999;96:10422–7. [7] Mori S, Wakana S, Van Zijl PC, Nagae-Poetscher LM. MRI atlas of human white matter. Elsevier; 2005.

[8] Udupa JK, Saha PK. Fuzzy connectedness and image segmentation. Proc IEEE 2003;91:1649–69. [9] Awate SP, Zhang H, Gee JC. A fuzzy, nonparametric segmentation framework for DTI and MRI analysis: with applications to DTI-tract extraction. IEEE Trans Med Imaging 2007;26:1525–36. [10] Parker GJ, Haroon HA, Wheeler-Kingshott CA. A framework for a streamlinebased probabilistic index of connectivity (PICo) using a structural interpretation of MRI diffusion measurements. J Magn Reson Imaging 2003;18:242–54. [11] Behrens TE, Woolrich MW, Jenkinson M, Johansen-Berg H, Nunes RG, Clare S, et al. Characterization and propagation of uncertainty in diffusion-weighted MR imaging. Magn Reson Med 2003;50:1077–88. [12] Hagmann P, Thiran JP, Jonasson L, Vandergheynst P, Clarke S, Maeder P, et al. DTI mapping of human brain connectivity: statistical fibre tracking and virtual dissection. Neuroimage 2003;19:545–54. [13] Jones DK, Pierpaoli C. Confidence mapping in diffusion tensor magnetic resonance imaging tractography using a bootstrap approach. Magn Reson Med 2005;53:1143–9. [14] Behrens TE, Berg HJ, Jbabdi S, Rushworth MF, Woolrich MW. Probabilistic diffusion tractography with multiple fibre orientations: What can we gain? Neuroimage 2007;34:144–55. [15] Friman O, Farneback G, Westin CF. A Bayesian approach for stochastic white matter tractography. IEEE Trans Med Imaging 2006;25:965–78. [16] Bjornemo M, Brun A, Kikinis R, Westin C-F. Regularized stochastic white matter tractography using diffusion tensor MRI. In: Proceedings of the international conference on medical image computing and computer-assisted intervention (MICCAI). 2002. p. 435–42. [17] Behrens TEJ, Jbabdi S. MR diffusion tractography. In: Johansen-Berg H, Behrens TE, editors. Diffusion MRI: from quantitative measurement to in vivo neuroanatomy. Elsevier; 2009. p. 333–51. [18] Campbell JS, Siddiqi K, Rymar VV, Sadikot AF, Pike GB. Flow-based fiber tracking with diffusion tensor and q-ball data: validation and comparison to principal diffusion direction techniques. Neuroimage 2005;27:725–36. [19] Parker GJ, Wheeler-Kingshott CA, Barker GJ. Estimating distributed anatomical connectivity using fast marching methods and diffusion tensor imaging. IEEE Trans Med Imaging 2002;21:505–12. [20] Staempfli P, Jaermann T, Crelier GR, Kollias S, Valavanis A, Boesiger P. Resolving fiber crossing using advanced fast marching tractography based on diffusion tensor imaging. Neuroimage 2006;30:110–20. [21] Tournier JD, Calamante F, Gadian DG, Connelly A. Diffusion-weighted magnetic resonance imaging fibre tracking using a front evolution algorithm. Neuroimage 2003;20:276–88. [22] Jackowski M, Kao CY, Qiu M, Constable RT, Staib LH. White matter tractography by anisotropic wavefront evolution and diffusion tensor imaging. Med Image Anal 2005;9:427–40. [23] Iturria-Medina Y, Canales-Rodriguez EJ, Melie-Garcia L, Valdes-Hernandez PA, Martinez-Montes E, Aleman-Gomez Y, et al. Characterizing brain anatomical connections using diffusion weighted MRI and graph theory. Neuroimage 2007;36:645–60. [24] Zalesky A. DT-MRI fiber tracking: a shortest paths approach. IEEE Trans Med Imaging 2008;27:1458–71. [25] Udupa JK, Samarasekera S. Fuzzy connectedness and object definition: theory, algorithms, and applications in image segmentation. Graph Models Image Process 1996;58:246–61. [26] Sotiropoulos SN, Tench CR, Bai L. In-vivo brain anatomical connectivity using diffusion magnetic resonance imaging and fuzzy connectedness. In: IEEE International Conference on BioInformatics and BioEngineering (BIBE). Athens, Greece, 2008. p. 1–8. [27] Sotiropoulos SN, Tench CR, Bai L. Fuzzy Anatomical Connectedness using Diffusion MRI: An Approach to Tractography of the Brain. In: ISMRM Meeting, Toronto, Canada, 2008. p. 1846. [28] Pierpaoli C, Jezzard P, Basser PJ, Barnett A, Di Chiro G. Diffusion tensor MR imaging of the human brain. Radiology 1996;201:637–48. [29] Tuch DS. Q-ball imaging. Magn Reson Med 2004;52:1358–72. [30] Nyul LG, Falcao AX, Udupa JK. Fuzzy-connected 3D image segmentation at interactive speeds. Graph Models 2003;64:259–81. [31] Poupon C, Clark CA, Frouin V, Regis J, Bloch I, Le Bihan D, et al. Regularization of diffusion-based direction maps for the tracking of brain white matter fascicles. Neuroimage 2000;12:184–95. [32] Kingsley PB. Introduction to diffusion tensor imaging mathematics. Part III. Tensor calculation, noise, simulations, and optimization. Concepts Magnet Reson A 2006;28A:155–79. [33] Cook PA, Symms M, Boulby PA, Alexander DC. Optimal acquisition orders of diffusion-weighted MRI measurements. J Magn Reson Imaging 2007;25:1051–8. [34] Hess CP, Mukherjee P, Han ET, Xu D, Vigneron DB. Q-ball reconstruction of multimodal fiber orientations using the spherical harmonic basis. Magn Reson Med 2006;56:104–17. [35] Cook PA, Bai Y, Nedjati-Gilani S, Seunarine KK, Hall MG, Parker GJ, et al. Camino: open-source diffusion-MRI reconstruction and processing. In: ISMRM Meeting. 2006. p. 2759. [36] Smith SM, Jenkinson M, Woolrich MW, Beckmann CF, Behrens TE, JohansenBerg H, et al. Advances in functional and structural MR image analysis and implementation as FSL. Neuroimage 2004;23(Suppl. 1):S208–19. [37] Parker GJ, Alexander DC. Probabilistic Monte Carlo based mapping of cerebral connections utilising whole-brain crossing fibre information. Inf Process Med Imaging 2003;18:684–95.

S.N. Sotiropoulos et al. / Computerized Medical Imaging and Graphics 34 (2010) 504–513 [38] Sotiropoulos SN, Bai L, Morgan PS, Auer DP, Constantinescu CS, Tench CR. A regularized two-tensor model fit to low angular resolution diffusion images using basis directions. J Magn Reson Imaging 2008;28:199–209. [39] Tournier JD, Calamante F, Gadian DG, Connelly A. Direct estimation of the fiber orientation density function from diffusion-weighted MRI data using spherical deconvolution. Neuroimage 2004;23:1176–85. [40] Jansons KM, Alexander DC. Persistent angular structure: new insights from diffusion magnetic resonance imaging data. Inverse Probl 2003;19:1031–46. [41] Wedeen VJ, Hagmann P, Tseng WY, Reese TG, Weisskoff RM. Mapping complex tissue architecture with diffusion spectrum magnetic resonance imaging. Magn Reson Med 2005;54:1377–86. [42] Sotiropoulos SN, Tench CR, Morgan PS, Bai L. Robust graph-based tracking through crossing fibre configurations. In: IEEE International Symposium on Biomedical Imaging (ISBI). 2009. p. 1394–7. [43] Haroon HA, Morris DM, Embleton KV, Alexander DC, Parker GJ. Using the modelbased residual bootstrap to quantify uncertainty in fiber orientations from Qball analysis. IEEE Trans Med Imaging 2009;28:535–50. [44] Schultz T, Seidel HP. Estimating crossing fibers: a tensor decomposition approach. IEEE Trans Vis Comput Graph 2008;14:1635–42. [45] Morris DM, Embleton KV, Parker GJ. Probabilistic fibre tracking: differentiation of connections from chance events. Neuroimage 2008;42:1329–39. [46] Saha PK, Udupa JK, Odhner D. Scale-based fuzzy connected image segmentation: theory, algorithms, and validation. Comput Vis Image Underst 2000;77:145–74.

513

Stamatios N. Sotiropoulos received his M.Sc. in Biomedical Engineering from the University of Minnesota-Twin Cities, USA in 2005 and his B.Sc. (Hons) in Electronics and Computer Engineering from the Technical University of Crete, Greece in 2003. Since 2006 he has been a Marie Curie fellow at the University of Nottingham, UK, working towards his Ph.D. in diffusion-weighted MRI of the brain and white matter tractography. Li Bai has a B.Sc. and a M.Sc. in Mathematics, and a Ph.D. in Computer Science from the University of Nottingham, UK. She is currently associate Professor at the School of Computer Science, University of Nottingham. Her main research interests are in the area of pattern recognition and computer vision, including face recognition, advanced interfaces for games, video based personal navigation, and medical image analysis. She is a referee for numerous international journals and funding councils. Christopher R. Tench received his B.Sc. (Hons) in physics with theoretical physics from the University of Manchester, UK in 1996. He received his Ph.D. in neuroimage analysis in 2003 from the University of Nottingham, UK while working as a research associate. Since then he has continued as a research fellow, in the division of clinical neurology, at the University of Nottingham. His main research interests are in applications of MRI to multiple sclerosis, including diffusion and spinal cord MRI.