|Home | About | Journals | Submit | Contact Us | Français|
Diffusion-weighted magnetic resonance imaging holds substantial promise as a technique for non-invasive imaging of white matter (WM) axonal projections. For diffusion imaging to be capable of providing new insight into the connectional neuroanatomy of the human brain, it will be necessary to histologically validate the technique against established tracer methods such as horseradish peroxidase and biocytin histochemistry. The macaque monkey provides an ideal model for histological validation of the diffusion imaging method due to the phylogenetic proximity between humans and macaques, the gyrencephalic structure of the macaque cortex, the large body of knowledge on the neuroanatomic connectivity of the macaque brain and the ability to use comparable magnetic resonance acquisition protocols in both species. Recently, it has been shown that high angular resolution diffusion imaging (HARDI) can resolve multiple axon orientations within an individual imaging voxel in human WM. This capability promises to boost the accuracy of tract reconstructions from diffusion imaging. If the macaque is to serve as a model for histological validation of the diffusion tractography method, it will be necessary to show that HARDI can also resolve intravoxel architecture in macaque WM. The present study therefore sought to test whether the technique can resolve intravoxel structure in macaque WM. Using a HARDI method called q-ball imaging (QBI) it was possible to resolve composite intravoxel architecture in a number of anatomic regions. QBI resolved intravoxel structure in, for example, the dorsolateral convexity, the pontine decussation, the pulvinar and temporal subcortical WM. The paper concludes by reviewing remaining challenges for the diffusion tractography project.
Despite continuing improvements in structural magnetic resonance imaging (MRI) technology, individual axons remain beyond the resolving capability of conventional structural MRI methods. Recently, it has been proposed that white matter (WM) axon pathways can be reconstructed non-invasively using an MRI technique called diffusion-weighted MRI (Conturo et al. 1999; Mori et al. 1999; Basser et al. 2000; Koch et al. 2002; Lori et al. 2002; Mori & Van Zijl 2002; Parker et al. 2002; Behrens et al. 2003). Diffusion-weighted MRI, or simply diffusion MRI, measures the molecular mobility of the endogenous water in brain tissue (Basser et al. 1994; Pierpaoli et al. 1996). The reconstruction of WM pathways from diffusion MRI is based on the observation that in cerebral WM, water diffusion is greater in the direction parallel to the fibre bundle than in the transverse direction (Chenevert et al. 1990; Chien et al. 1990; Moseley et al. 1991; Beaulieu 2002). This phenomenon is referred to as diffusion anisotropy. The biophysical basis for diffusion anisotropy in WM is not fully known but may be due to the diffusion barrier presented by the myelin sheath and the cell membrane (Beaulieu 2002).
Investigators have sought to reconstruct WM pathways from diffusion images by computationally tracing the direction of maximal diffusion. Borrowing from a hydrodynamic analogy, the diffusion images provide the flow vector field, and the WM pathways represent flow streamlines. The reconstruction of WM pathways from diffusion imaging is a technique referred to as diffusion tractography.
Diffusion tractography has been applied widely in both research and clinical applications, but the technique has yet to be validated histologically in any systematic fashion. The tractography technique will need to be validated against established tracer methods before the tractography reconstructions can be interpreted with any reasonable confidence. Table 1 lists the advantages and disadvantages of the diffusion MR method compared with conventional histological techniques. In general, there is a trade-off between the invasiveness and model-dependence of a method. For example, diffusion MR tractography is non-invasive but requires a model of the relationship between the diffusion signal and the underlying myeloarchitecture. Invasive tracers such as horseradish peroxidase and biocytin can follow the continuity of an axon without need for a mathematical model of the transport process.
The accuracy of diffusion MR tractography is currently limited by a number of factors including the low spatial resolution of the MR image relative to the tract curvature (Lori et al. 2002), the ambiguity between crossing and bending fibres (Basser et al. 2000; Tuch 2004) and distortions in the diffusion images (Jezzard & Clare 1999). The latter problem hampers accurate registration between diffusion images and conventional T1-weighted structural images from which anatomical labels are conventionally obtained. More fundamentally, the tractography method is based on a presumed correspondence between the measured diffusion function and the underlying myeloarchitecture, although this relationship has never been tested.
The macaque monkey provides an ideal model for histological validation of the tractography method due to the phylogenetic proximity between humans and macaques, and the gyrencephalic structure of macaque cortex. Additionally, there is a large body of existing knowledge on the anatomic connectivity of the macaque brain. In particular, validation of diffusion tractography against the CoCoMac database (Stephan et al. 2001; Kotter 2004; http://www.cocomac.org) is a promising avenue of research. The CoCoMac database contains a catalogue of reported tracing results in the macaque in a format which can be directly compared with the MR.
The macaque also provides an excellent validation model due to the ability to use comparable MR acquisition protocols in both humans and macaques. Parker and colleagues (2002) have previously demonstrated the technical feasibility of acquiring high resolution diffusion tensor imaging (DTI) in the macaque. It has recently been shown that a high angular resolution diffusion imaging (HARDI) method called q-ball imaging (QBI) can resolve multiple intravoxel axon orientations in human WM (Tuch et al. 2003; Tuch 2004). This capability promises to substantially boost the accuracy of the diffusion MR tract reconstructions. Other HARDI methods include persistent angular structure MRI (PASMRI; Jansons & Alexander 2003), mixture model decomposition (Tuch et al. 2002; Parker & Alexander 2003), generalized DTI (Ozarslan & Mareci 2003; Liu et al. 2004), spherical harmonic transformation (Frank 2001) and spherical harmonic deconvolution (Tournier et al. 2004).
The conventional diffusion imaging method for mapping fibre orientation is DTI (Basser et al. 1994; Pierpaoli et al. 1996). DTI can only resolve a single diffusion orientation within each imaging voxel however (Wiegell et al. 2000; Pierpaoli et al. 2001). This limitation is because of the assumption of Gaussian diffusion implicit to the tensor model.
In contrast, QBI poses no assumptions on the underlying diffusion process and can resolve multiple intravoxel diffusion orientations. Rather than measuring the diffusion tensor, QBI measures the diffusion orientation distribution function (ODF). The diffusion ODF ψ(u) describes the probability distribution for a water molecule to displace in a direction u, where u is a unit vector with no polarity, also called a director. In anatomic regions containing multiple axon orientations, the diffusion ODF exhibits multiple diffusion peaks, or multimodal diffusion (MMD; Tuch et al. 2003; Tuch 2004). The QBI method is thus capable of detecting the presence of multiple diffusion directions within an individual voxel.
The resolution of a QBI experiment is defined in diffusion displacement space and is determined by the experimental diffusion wavevector q (Tuch 2004). The diffusion wavevector q is proportional to the duration δ and magnitude g of the diffusion-sensitizing magnetic field gradient; specifically, q=(2π)−1γgδ where γ is the gyromagnetic ratio (Callaghan 1993). In QBI, the diffusion wavevector defines a projection function ρ(r), which describes the radial fall-off of the projection beam in diffusion space. For a standard single wavevector QBI acquisition, the projection function is a Bessel beam, ρ(r)=J0(2πqr) (Tuch 2004). The projection function has a width in the order of microns for a conventional QBI acquisition. The microscopic resolving power of QBI enables the technique to resolve multiple diffusion orientations within an individual imaging voxel. Previously, such intravoxel complexity has confounded diffusion imaging methods, such as DTI, which are only capable of resolving a single diffusion orientation within each imaging voxel (Wiegell et al. 2000; Pierpaoli et al. 2001).
It is important to note that the resolution of a QBI acquisition has dimensions length and not angle. Thus, the diffusion wavevector does not provide a direct indication of the angular resolution. The angular resolution also depends on the q-space sampling density on the sphere. The final resolution will be a function of the sampling wavevector, the sampling density and any smoothing operations applied in post-processing. Developing an imaging framework for QBI which will allow for the description of the angular resolution for a given sampling and reconstruction scheme is an area of active research.
QBI has demonstrated the capability of resolving complex intravoxel WM architecture in the human brain (Tuch et al. 2003; Tuch 2004) but it is not clear if QBI can also resolve intravoxel WM structure in the macaque brain. If the macaque is to serve as a model for histological validation of diffusion MR tractography, it will be necessary to show that QBI can also resolve complex WM architecture in the macaque. The objectives of the present study were therefore twofold: to demonstrate that QBI can resolve subvoxel architecture in macaque WM, and to develop a platform for future validation of diffusion tractography against histological tracer methods and the CoCoMac database.
In macaque WM, QBI resolved intravoxel architecture in a number of anatomic regions. QBI resolved intravoxel WM structure in, for example, the dorsolateral convexity, the pontine decussation, the pulvinar and temporal subcortical WM. The intravoxel fibre crossings were consistent with the known association, projection and callosal fibre architecture of the macaque brain. The paper concludes by discussing remaining challenges for the diffusion tractography project.
MRI data were acquired from one juvenile male rhesus monkey (Macaca mulatta, 3.5kg, ID #M304). The data were acquired on a Siemens Trio 3T MRI scanner located at the Athinoula A. Martinos Center for Biomedical Imaging, Massachusetts General Hospital (Charlestown, Massachusetts). All procedures conformed to Massachusetts General Hospital and National Institutes of Health guidelines for the care and use of laboratory animals (Subcommittee on Research Animal Care protocol #2003N000338).
The animal was anaesthetized with a ketamine and xylazine dosage (10 and 0.5mgkg−1, i.m.). After induction of the anaesthesia, an intravenous line (saphenous vein) was used to give additional ketamine injections to maintain anaesthesia during the scan. The monkey was placed into a magnet compatible stereotactic apparatus (Kopf Instruments, Tujunga, California). Local anaesthetic (lidocaine cream) was applied to the ends of the ear bars and ophthalmic ointment was applied to the eyelids to minimize discomfort induced by the stereotactic apparatus. A heating pad was placed beneath the subject to keep the animal warm during the scan session. A custom-built 13cm diameter flat surface (receive-only) coil was positioned over the skull and was fixed at a dorsal–ventral level of the outer ear canals. Care was taken that the brain was positioned in the middle of the coil.
Whole-head, high resolution structural MRI data were acquired using an MPRAGE sequence (Mugler & Brookeman 1990) with TR/TI/TE=2500/1100/4ms, α=8°, 350×350×350μm resolution, total acquisition time: 16min 2s. Diffusion MRI QBI data were acquired using a twice-refocused balanced echo sequence 90–g1–180–g2–g3–180–g4–acq, where the 180 pulse-pair was positioned to minimized eddy current-induced image distortions (Reese et al. 2003). Two QBI datasets were acquired: a ‘low wavevector’ acquisition with b=4000smm−2, q=537cm−1 and a ‘high wavevector’ acquisition with b=8000smm−2, q=682cm−1. These wavevector settings correspond to diffusion spatial resolutions of, respectively, Δr=7.13 and 5.61μm (Rayleigh criterion). Both acquisitions used the same slice prescription: 52 contiguous coronal slices (1.6mm thick, 0mm gap), 102.4×102.4mm FOV, 64×64 acquisition matrix, 1.6×1.6mm in-plane resolution. The sequence parameters for the low wavevector acquisition were TR/TE=7000/128ms, b=4000smm−2, g(max)=26mTm−1, δ(eff)=48ms, q=537cm−1, and for the high wavevector acquisition were TR/TE=10000/152ms, b=8000smm−2, g(max)=26mTm−1, δ(eff)=61ms, q=682cm−1.
For each acquisition, 10 T2 and 252 diffusion-weighted images were acquired. The 252 diffusion gradient directions were obtained from the vertices of a fivefold regularly tessellated icosahedron projected onto the sphere. Figure 1a,b shows the q-space sampling scheme. The total acquisition time was 30min 34s for the low wavevector acquisition and 43min 40s for the high wavevector acquisition.
To correct for drift between the two diffusion acquisitions in the phase-encode direction, the high wavevector volume was registered to the low wavevector volume using the FLIRT (Jenkinson & Smith 2001; Jenkinson et al. 2002) tool from FSL (http://www.fmrib.ox.ac.uk/fsl). The registration process entailed solving the transformation from the high wavevector average T2 volume to the low wavevector average T2 volume, allowing only for displacements in the phase-encode direction. The normalized mutual information cost function was used. The transform was then applied to the high wavevector volumes, and the volumes were resampled using nearest neighbour interpolation. Nearest neighbour interpolation was used, instead of sinc interpolation, to prevent MMD from arising artificially from interpolation between unimodal diffusion in neighbouring voxels. A brain mask was defined by thresholding the average of the low wavevector diffusion volumes. The diffusion signal from the high wavevector and low wavevector datasets were then masked, normalized by their respective T2 volumes, and averaged.
The diffusion ODF maps were reconstructed using the algorithm described in Tuch (2004). The diffusion ODF ψ(u) for each voxel was reconstructed by applying the matrix Funk–Radon transform (FRT) to the averaged diffusion data. The FRT was implemented using spherical radial basis function interpolation with a spherical Gaussian kernel (σ=8°). The ODF was reconstruction onto n=752 reconstruction points. The reconstruction points were obtained from the vertices of a fivefold regularly tessellated dodecahedron projected onto the sphere. The same set of points was used as the basis function centres. Each ODF was then smoothed using a spherical Gaussian kernel with σ=3°. The generalized fractional anisotropy (GFA) metric (Tuch 2004) and the peak diffusion direction u* were computed for each voxel. The GFA metric is defined as GFA=s.d.(ψ)/r.m.s.(ψ) where s.d. is the standard deviation and r.m.s. is the root mean square. The peak diffusion direction is defined as u*=argmaxψ(u).
Figure 1b shows the theoretical projection functions ρ(r) for the low, high and combined wavevector acquisitions. The projection function provides a measure of the spatial resolution of the QBI acquisition. The projection function gives the radial fall-off of the projection beam in diffusion displacement space (Tuch 2004). The combined low/high wavevector acquisition gives a diffusion displacement resolution of Δr=6.33μm (Rayleigh criterion). Figure 1c shows the raw data and the combined data for a subset of the diffusion gradient sampling directions.
The diffusion ODF maps were rendered as fields of colour-coded spherical polar plots. For visualization purposes, the ODFs were minimum–maximum normalized and then scaled by GFA. This rescaling emphasizes the angular contrast within individual ODFs, while reducing the visual contribution from ODFs in regions of low anisotropy (Tuch 2004). The equation for the rendering mapping is given in the caption for figure 3.
Anatomical labels were assigned based on the Martin & Bowden (2000) atlas and using the NeuroNames (Bowden & Martin 1995; Martin & Bowden 2000) vocabulary where applicable. The identified structures are summarized in table 2.
Colour-coded peak diffusion orientation maps of the macaque are shown in figure 2. The peak diffusion orientation maps show the direction of maximal diffusion at each voxel location. The maps delineate major WM structures including the internal capsule (ic), the splenium of the corpus callosum (scc), the body of the corpus callosum (bcc), the fornix (fx), the posterior part of the anterior commissure (acp), the middle cerebellar peduncle (mcp), the decussation of the superior cerebellar peduncle (xscp), the optic radiation (or) and the optic tract (opt). Detailed diffusion ODF maps from select regions are shown in figures 3–7.
Figure 3 (central panel) shows a diffusion ODF map of precentral gyrus (PrG) and the dorsolateral convexity which includes, respectively, projections to the motor, parietal and temporal cortices. MMD can be seen at the intersection (figure 3, arrow 1) between the projections of the scc (red) and the superior–inferior directed fibres (blue–green) connecting to the motor cortex, and where the projections segregate into the superior parietal lobule (SPL) and postcentral gyrus (PoG; figure 3, arrow 2). MMD is also seen in the WM beneath the supramarginal gyrus (SMG; figure 3, arrow 3), including a lateral (red) component along the gyral axis and a transverse anterior–posterior component (bright green). The parietal lobule has dense reciprocal connections with visual association, cingulate, premotor and prefrontal cortical areas (Petrides & Pandya 1984). The intersection of intrinsic and extrinsic parietal lobule fibres is demonstrated in the ODF map (figure 3, arrow 2). Extrinsic fibres from this area of the parietal lobule, specifically the SPL and dorsal bank of the intraparietal sulcus (itps), contribute to the superior longitudinal fasciculus (slf), fronto-occipital fasciculus and bcc.
Figure 4 shows the thalamocortical projections from the pulvinar (Pul). MMD is seen where the projections from the Pul (red) turn from a medial–lateral orientation (red) to a superior–inferior orientation (arrow) to project to the parietal cortex. The Pul sends projection fibres to the somatosensory and visual association areas in the parietal lobe (Schmahmann & Pandya 1990). The fibres can be seen passing the or.
MMD can also be seen in subcortical WM. Figure 5 shows MMD in the subcortical WM underlying superior temporal gyrus (STG). Similar to SMG (figure 3, arrow 3), a medial–lateral component (red) directed along the gyral axis is apparent, as well as a transverse component in the anterior–posterior (green) direction. This anterior–posterior component is most likely due to temporal–prefrontal projections passing through STG (Seltzer & Pandya 1989).
Figure 6 shows a diffusion ODF map of the midbrain (MB) and the pontine decussation. While the T1-weighted image shows no contrast between the different MB structures, the diffusion ODF map clearly reveals the superior colliculus (SC), the inferior colliculus (IC) and the central grey substance of MB (CGMB; blue). The mcp (green) projecting to the cerebellum (Cb) is also apparent. The oblique-oriented fibres of the trapezoid body (tz; blue–purple) are seen at the level of the pons (Pons), inferior to the CGMB. Intravoxel fibre crossing is seen where the fibres from the commissure for the superior colliculus (csc; red) and the commissure for the inferior colliculus (cic; red) intersect the fibres emanating from the CGMB (figure 6, arrow 1). Additionally, more ventro-caudal horizontally crossing fibres are seen at the level of the Pons (figure 6, arrow 2). The lateral component here may reflect connections between the CGMBs. At the level of the csc and cic, the lateral architecture is dominant, and at the level of CGMB, the superior–inferior and lateral fibres contribute equally.
Figure 7 shows fibres crossing in the centrum semiovale (semioval centre, sovc). The diffusion ODF map shows crossing between the corona radiata (cr; blue), slf (green) and the transcallosal projections from bcc (red). In Figure 7, arrow 1 shows the crossing between the cr and the transcallosal fibres projecting to PrG. Arrow 2 shows the crossing between transcallosal projections and the slf.
QBI resolved composite subvoxel architecture in macaque WM. Intravoxel structure was resolved in a number of anatomic regions including the centrum semiovale, the pontine decussation, the pulvinar, and temporal subcortical WM. The intravoxel fibre crossings were consistent with the known association, projection and callosal fibre architecture of the macaque brain.
Notwithstanding QBI's ability to resolve intravoxel architecture, there still remain a number of obstacles to accurate reconstruction of WM pathways with QBI. One considerable obstacle is that a number of different fibre arrangements are consistent with an observed diffusion function. For example, for a diffusion ODF with two peaks, the two peaks may be due to two fibres crossing, a single fibre bundle segregating into two clusters, two fibres ‘kissing’ (Basser et al. 2000) or a single fibre bending. This ambiguity is due to the reflective symmetry of the ODF, that is, the diffusion probability in any given direction is the same in the opposite direction.
This topological ambiguity is apparent in the diffusion ODF maps. For example, in figure 4 it is unclear if the observed MMD in the Pul is due to crossing lateral and ascending fibres or because a single fibre population has rapidly turned from a medial–lateral orientation to a superior–inferior orientation. It is worth noting that bending fibre bundles can only produce MMD if the fibre curvature is variable, since a fibre bundle with constant curvature would simply exhibit a single broad peak. Alternatively, MMD could also arise from multiple bending fibres with different orientations.
Figure 3 (arrow 2) also shows the degeneracy between bending and crossing fibres. It is not clear if the MMD arises from the projections segregating into SPL and PoG, or from intrinsic connections within and between the gyri (Pandya & Seltzer 1982; Seltzer & Pandya 1986). Similarly, for the MMD seen in figure 3 (arrow 3) and figure 5, it is not clear if the transverse diffusion component seen in SMG and STG, respectively, is due to WM fibres bending rapidly to insert into the gyral wall or due to WM association fibres passing through the gyrus (Seltzer & Pandya 1989). The complexity of the WM architecture at the gyral wall poses a substantial challenge for diffusion tractography. Two open questions include whether this topological ambiguity can be resolved with diffusion imaging, and what the implications of this ambiguity are for tractography reconstructions.
The interpretation of diffusion tract solutions as anatomical tracts is based on the assumption that the diffusion ODF is closely related to the underlying fibre ODF. However, the relationship between the diffusion and fibre ODFs is not entirely clear. The relationship is obscured by the lack of knowledge about the physiological mechanism for diffusion anisotropy in WM (Beaulieu 2002). Clarifying the relationship between the diffusion and fibre ODFs, and incorporating this relationship into tractography algorithms, would help support the interpretation of the diffusion tract reconstructions as fibre pathways instead of diffusion pathways.
Given the above concerns, validation of the diffusion connectivity method against established tracer techniques will be essential. Histological validation will be particularly necessary for the three-dimensional tract reconstructions derived from diffusion imaging. The present study has demonstrated the feasibility of acquiring QBI in the macaque, and that the technique can resolve subvoxel WM architecture. Future work will focus on validating the diffusion connectivity maps against histological tracer methods.
The authors are grateful to Thomas Benner and Timothy Reese for developing the diffusion MRI protocol, and to Larry Wald for constructing the MRI coil. The authors acknowledge support from NINDS NS46532, NCRR RR14075, GlaxoSmithKline, the Athinoula A. Martinos Foundation, the Mental Illness and Neuroscience Discovery (MIND) Institute, the National Alliance for Medical Imaging Computing (NIBIB EB005149), which is funded through the NIH Roadmap for Medical Research, the Human Frontiers Science Program (RGY 14/2002) and Deutsche Forschungsgemeinschaft (GRK320 and KO 1560/6-2). The authors would also like to thank the organizers of the 2004 Brain Connectivity Workshop for providing a stimulating collaborative environment.
One contribution of 21 to a Theme Issue ‘Multimodal neuroimaging of brain connectivity’.