An in vivo MRI Template Set for Morphometry, Tissue Segmentation, and fMRI Localization in Rats.

Over the last decade, several papers have focused on the construction of highly detailed mouse high field magnetic resonance image (MRI) templates via non-linear registration to unbiased reference spaces, allowing for a variety of neuroimaging applications such as robust morphometric analyses. However, work in rats has only provided medium field MRI averages based on linear registration to biased spaces with the sole purpose of approximate functional MRI (fMRI) localization. This precludes any morphometric analysis in spite of the need of exploring in detail the neuroanatomical substrates of diseases in a recent advent of rat models. In this paper we present a new in vivo rat T2 MRI template set, comprising average images of both intensity and shape, obtained via non-linear registration. Also, unlike previous rat template sets, we include white and gray matter probabilistic segmentations, expanding its use to those applications demanding prior-based tissue segmentation, e.g., statistical parametric mapping (SPM) voxel-based morphometry. We also provide a preliminary digitalization of latest Paxinos and Watson atlas for anatomical and functional interpretations within the cerebral cortex. We confirmed that, like with previous templates, forepaw and hindpaw fMRI activations can be correctly localized in the expected atlas structure. To exemplify the use of our new MRI template set, were reported the volumes of brain tissues and cortical structures and probed their relationships with ontogenetic development. Other in vivo applications in the near future can be tensor-, deformation-, or voxel-based morphometry, morphological connectivity, and diffusion tensor-based anatomical connectivity. Our template set, freely available through the SPM extension website, could be an important tool for future longitudinal and/or functional extensive preclinical studies.

Introduction

A magnetic resonance image (MRI) “template set” comprises a reference or template image, which can be an image of average intensity and/or shape, and a set of either discrete or probabilistic images representing tissue types, all in the reference space. An undeniable landmark in this area is the statistical parametric mapping (SPM) software package (Friston et al., 2010), which provides a human average template and tissue priors of white matter, gray matter, and cerebrospinal fluid, as well as methods for spatial normalization and segmentation. The construction of this human template set can be found elsewhere (Mazziotta et al., 1995a,b, 2001) and its main applications are functional MRI (fMRI) localization, voxel-based morphometry (Ashburner and Friston, 2000) SPM (Friston et al., 2010), among others. To enrich the deliverables of the template set, various discrete digitalizations with finer parcellation of gray (Tzourio-Mazoyer et al., 2002; Rasser et al., 2008) and white matter (Hua et al., 2008; Mori et al., 2008), as well as their probabilistic versions (Shattuck et al., 2008), have been proposed. Combined with spatial normalization and tissue segmentation, these digital parcellations afford further morphometric applications, such as the statistical comparison of volumes and areas of regions-of-interest (ROI), deformation- (Chung et al., 2003) or tensor- (Kipps et al., 2005) based morphometry, ROI-based morphological connectivity (Sanabria-Diaz et al., 2010); anatomical connectivity between gray matter ROI-based diffusion tensor imaging (DTI) tractography (Iturria-Medina et al., 2007; Iturria-Medina et al., 2008); and localization of fMRI activation within ROIs (Tzourio-Mazoyer et al., 2002; Maldjian et al., 2003), which have been also useful for pre-surgical planning (Bizzi, 2009) or even providing seeds for DTI-tractography (Saur et al., 2010).

The aforementioned human neuroimaging applications are also required for rodents, since they constitute the animal models of preclinical research for human diseases. The importance of building mouse templates and structure digitalizations for morphometry is demonstrated in several papers. In contrast to manual 2D morphometry, SPM5 voxel-based morphometry in the R6/2 transgenic mouse detected differences between genotypes (Sawiak et al., 2009). Other phenotyping examples are the reduced cerebellar volume of the inferior colliculus and the olfactory bulbs reported for cdf mutant mice (Bock et al., 2006); the morphologic abnormalities found in various areas of the Reeler mouse brain (Badea et al., 2007b) and in the striatum of the YAC128 Huntington disease mouse model (Lerch et al., 2008); and the morphologic differences between the 129S1/SvImJ, C57/Bl6, and CD1 mouse strains (Chen et al., 2006). An excellent review on mice “templating” can be found in (Dorr et al., 2008), with more recent efforts in (Ma et al., 2008; Aggarwal et al., 2009; Chuang et al., 2011; Hawrylycz et al., 2011).

However, in spite of the recent uprising of rat models of epilepsy (Tenney et al., 2003; Nersesyan et al., 2004; David et al., 2008; Englot et al., 2008), Parkinson (Van Camp et al., 2003; Pelled et al., 2005; Hou et al., 2010), autism (Pletnikov et al., 2001; Osumi, 2009), Alzheimer (Ganten, 1998; Yang et al., 2011), stroke (Li et al., 2000; Aoki et al., 2003), spinal cord injury (Ramu et al., 2006; Endo et al., 2007; Takagi et al., 2009; Ghosh et al., 2010), depression (Huang et al., 2011), among others, there is less focus on the construction of rat templates and either discrete or probabilistic digitalizations of tissues or structures. In fact, although there have been different approaches (Toga et al., 1995; Leergaard et al., 2003; Keilholz et al., 2006; Hjornevik et al., 2007), only two papers have provided rat MRI template sets with average reference images, tissue segmentations, or discrete parcellations (Schweinhardt et al., 2003; Schwarz et al., 2006). Both sets however present various drawbacks:

They are based on medium field MRI (4.7 T), though recent animal scanners are 7.0 T or higher. In fact, the most updated mouse templates come even from 9.4 T systems.

Their average templates are built after linear registration of individual MRIs to a common space. As recently defended by (Lu et al., 2010), that policy responded to the sole goal of qualitative evaluations of roughly localized fMRI activations in stereotaxic spaces: BOLD elicited by forepaw stimulation in (Schweinhardt et al., 2003) and functional perfusion MRI after pharmacological administration in (Schwarz et al., 2006). In contrast, mouse templates have been built by means of hierarchical, multi-resolution non-linear registration schemes in order to achieve the highest possible spatial correspondence between all images of the sample. This means that the maximum possible information about the variability of the sample is stored in the deformation fields (Ashburner and Friston, 2000; Bookstein, 2001; Mechelli et al., 2005). Provided that the conditions for the warps are accomplished – i.e., the transformation is a diffeomorphism or a regularized enough to guarantee that the Jacobian is positive at every point of the region of interest – the higher the spatial match of the registered images the finer the details of the average image. Thus, the spatial normalization of any individual MRI can be achieved with the highest possible parameterization providing an optimized statistical performance for the discrimination of pathological individuals or samples.

Their template spaces are biased, i.e. do not minimize the non-linear normalization deformations of all brains in the studied population. This introduces a bias in the registration toward those brains closest in shape to the template (Guimond, 2000; Kochunov et al., 2001a,b).

Although probabilistic images of brain parenchyma and cerebrospinal fluid were provided in (Schwarz et al., 2006), there is a lack of probabilistic white matter and gray segmentations, used for the popular and worldwide used prior-based methods such as the SPM8 unified segmentation (Ashburner and Friston, 2005) or FSL–FAST. To the knowledge of the authors, the only published priors for rodents are presented in (Sawiak et al., 2009) for mice.

For the aforementioned reasons, both template sets are not the most suitable for studies depending on (1) unbiased normalization methods (such as deformation- or tensor-based analysis, or detailed characterization of anatomical or functional features); (2) normalization combined with tissue segmentation [such as voxel-based morphometry, as done in (Sawiak et al., 2009) for mice]; or (3) both combined with digitalization-based classification (also called “automatic labeling”) of structures of the individual MRI (useful for ROI deformation-based morphometry as done in (Lerch et al., 2008) for mice, structure volumetry, as done in (Badea et al., 2007a) for mice, and possibly ROI-based DTI and morphological connectivity).

Therefore we present a new in vivo template set for rats with the first clear advantage of being based on a 7.0-T scanner. The set comprises a minimal deformation template (MDT), as in (Kovacevic et al., 2005), based on T2 MRIs of 30 male Wistar rats. The non-linear transformations were achieved by means of a multi-resolution, multi-grid, high dimensional registration method provided by a very general and adaptable software: “elastix” (Klein et al., 2010). Our set also contains probabilistic gray matter and white matter segmentations, which are reported for the first time for rats. These were assessed by averaging warped individual segmentations that were obtained from fitting various parametric and non-parametric curves to the histograms of the semi-manually non-brain stripped T2 MRIs. Additionally, we include a digitalization of 96 cortical structures of latest edition of the Paxinos and Watson (P & W) atlas (Paxinos and Watson, 2007), traditionally necessary for the anatomical interpretation of fMRI results since (Schwarz et al., 2006), and now useful for the aforementioned ROI-based analyses, within the cerebral cortex. This digital parcellation was warped to the average template by means of the approximate thin plate splines (TPS) non-linear registration method (Rohr et al., 2001) based on matching 36 pairs of landmarks in the template and atlas spaces. Our template set as well as the documentation about the chosen landmarks and digitalized regions are freely available either through the SPM extensions webpage or at http://www.idac.tohoku.ac.jp/bir/en/db/rb/101028.html.

Among the possible applications for which the template is intended we firstly prove its usefulness for three basic ones: (1) spatial normalization and (2) tissue segmentation of individual MRIs, and (3) automatic labeling of cortical structures in the individual MRI space. Note that in the case of studies not intended for precise anatomical characterizations relying on the structural T2 MRI, like fMRI or DTI, acquisition time must be reduced at expenses of the T2 image resolution, which becomes a secondary aspect, though necessary for these basic applications. Therefore, we also test them with “coarser” MRIs, i.e., those having a lower resolution and much shorter acquisition time than those used for the construction of the template set. Besides, we validated the performance of spatial normalization in localizing fMRI activations elicited by forepaw and hindpaw stimulation within the expected cortical regions, in both hemispheres. The aforementioned three basic applications also served to calculate and report the average and SD of the volumes of the cortical structures defined by the atlas. To illustrate morphometric applications in vivo, we probe volume–body weight correlations, evidencing some features of Wistar brain ontogenesis.

Materials and Methods

Notation

A lower case bold symbol, e.g., f, denotes a column vector. An upper case bold symbol, e.g., A, denotes a matrix whose size is specified by the corresponding subscripts, if they appear. Non-bold symbols denote scalar magnitudes. 1 and 0 are n × m matrices of ones and zeros respectively. ||f||2 is the Frobenius norm of f. The symbol ∧ denotes the logical operator “AND.” The symbol⊗denotes the Kronecker product. The superscript T denotes transposition.

Animal preparation

For the template construction and volumetry

Thirty (30) Wistar rats (Charles River Japan, Yokohama, Japan) were used for the template construction. They were distributed into five groups of six rats for each week of age (6, 7, 8, 9, 10 weeks old, body weight ranges of 214 ± 19, 275 ± 21, 312 ± 9, 331 ± 17, 364 ± 18 g, respectively). All procedures and protocols were performed in agreement with the policies established by the “Animal Care Committee” at Tohoku University, Sendai, Japan. Each rat was initially anesthetized with 5% isoflurane and then secured on a custom-built holder using adhesive tape and a bite bar. A breathing sensor (SA Instruments Inc., NY, USA) was placed under the ventral surface of the rat body. Anesthesia was further kept with isoflurane (at 1 L/min oxygenation), administered via a face mask. A constant breathing rate was maintained at around 50 breaths/min for the entire span of MRI acquisition by manually keeping the concentration of isoflurane in the range of 1.5–2.5%. This empirical rate is suitable for an MRI acquisition with satisfactorily low motion artifacts and little compromise of the scanning time. Core body temperature was maintained at 37.0 ± 1°C by means of a hot water-circulating pad.

For fMRI experiments

An additional 6 week old rat underwent fMRI experiment. The details of the experimental procedures are described in (Sumiyoshi et al., 2011). Briefly, the rat underwent mechanical ventilation (SAR-830AP Ventilator, CWE, Inc., Ardmore, PA, USA) with a combination of α-chloralose anesthesia (40 mg/kg/h, i.v.) and muscle relaxation administration (pancuronium bromide, 2 mg/kg/h i.v.). The physiological variables were maintained at normal limits by small adjustments in ventilation (end-tidal CO2 = 27.5 ± 4 mmHg). A pair of small needle electrodes (NE-224S, Nihon Koden, Tokyo, Japan) was inserted for electrical stimulation.

Magnetic resonance imaging

All MRI data were acquired using a 7.0-T Bruker PharmaScan system (Bruker Biospin, Ettlingen, Germany) with a 38-mm-diameter birdcage coil. Prior to all MRI experiments, we firstly performed a global magnetic field shimming inside the coil, which was completed in a ROI around gray matter (∼300 μL) by using a point resolved spectroscopic protocol (Sanganahalli et al., 2009). The line width (full width at half maximum) at the end of the shimming procedure ranged from 15 to 20 Hz in the ROI.

High resolution/long duration T2 MRI for the template construction and volumetry

High resolution T2-weighted images were obtained using a respiratory-gated 2D Turbo RARE sequence with fat suppression under the following parameters: TR = 10971 ms, TEeff = 30 ms, RARE factor = 4, effective spectral bandwidth = 100 kHz, flip angle = 90°, field of view = 32 mm × 32 mm, matrix size = 256 × 256, in-plane resolution = 125 μm × 125 μm, number of slices = 128, slice thickness = 0.3 mm, slice gap = 0 mm, and number of averages = 16. The total scanning time for T2-weighted imaging was around 3 h, depending on the respiration rate of the rat.

Low resolution/short duration T2 MRI for the fMRI experiments

In order to perform the necessary image processing (spatial registration and tissue segmentation), the fMRI experiment was accompanied by a low resolution T2 MRI with exactly the same parameters as the high quality one, excepting for the number of averages (10), the slice thickness (0.5 mm) and the number of slices (54). These parameters were chosen to reduce the scanning time of this image to around 50 min, making possible its inclusion in fMRI experiments.

Functional MRI experiments

We separately performed an fMRI experiment for each of the four types of electrical stimulations, i.e., forepaw and hindpaw, both left and right. A block-design paradigm (10 blocks) was employed for each fMRI experiment. Each block comprised 30 s of electrical stimulation followed by 40 s of the resting condition. The electrical pulses (3 Hz, 3.0 mA, and 0.3 ms width) were produced by a generator (SEN-3401, Nihon Koden) and an isolator (SS-203J, Nihon Koden).

The fMRI images were obtained by using a gradient-echo echo-planar imaging with fat suppression under the following parameters: TR = 2000 ms, TE = 15 ms, effective spectral bandwidth = 250 kHz, flip angle = 30°, field of view = 25 mm × 14 mm, matrix size = 125 × 70, in-plane resolution = 200 μm × 200 μm, number of slices = 7, slice thickness = 1.5 mm, slice gap = 0 mm, number of volumes = 370, and dummy scans = 4. The total scanning time for the fMRI experiment was 12 min and 20 s.

Construction of the minimal deformation template

The intensity of each of the 30 high resolution MRIs was corrected for intensity inhomogeneities due to a non-ideal receiver coil profile (Ashburner, 2002). Then, similar to a procedure described in (Kovacevic et al., 2005), an intensity “equalization” of all corrected MRI was performed. This operation (described in Intensity “Equalization” of the MRIs in Appendix) simultaneously shifts the position and scales the width of the peak of the histogram of each MRI to be around those corresponding to a randomly chosen MRI.

The corrected high resolution MRIs were segmented into “brain only” MRIs by removing all voxels not belonging to the brain parenchyma, cerebellum, and ventricles. The binary brain mask was extracted using a modified version of the method described by (Uberti et al., 2009) for the mouse head. In this semi-automatic method, based on the level set theory, a 2D closed contour in each slice expands until reaching the edge of the brain, where an automatic stopping criterion is fulfilled. We found that this method shows imperfect results, e.g. the occasional inclusion of part of the trigeminal nerve in the mask. Therefore we decided to stop every contour manually, after visually inspecting it at each expanding step. We demanded a brain mask as more precise as possible for the non-linear registration procedures, but avoiding the excessive burden of a fully manual delineation.

The MDT was created with a multi-resolution, multi-grid registration scheme similar to (Kovacevic et al., 2005). After choosing the brain only MRI of the rat having the median body weight, all other were registered to its space trough affine registration. The affine transformations, say M i = 1…30, mapping the reference to the individual spaces were averaged and inverted. The resulting transformation, M0, was post multiplied by each individual transformation to obtain the transformations, M0M i = 1…30, to the best representative affine space. These were applied to obtain the first linear average template. Using elastix (Klein et al., 2010), a non-linear transformation mapping each point in template space to the i-th individual space, say f(x), was estimated by optimizing a cost function based on the normalized mutual information between the template and the i-th individual brain only MRI. The non-linearly transformed images were averaged to obtain the second template. This procedure was iterated five times, using a multi-resolution scheme starting from a downsampling/blurring factor of 8 to 1 and randomly sampled grid spacing from 8 to 4. The final concatenated transformations were inverted and averaged, namely Finally, all transformations were concatenated and applied to the original brain only MRIs and whole head corrected MRIs to obtain the final corresponding templates:

where T(x) represents either the MDT or the whole head MRI template and I(x) is the i-th individual MRI.

The MDT image extends from −8.29 to 8.58 mm left to right, from −15.66 to 6.71 mm caudal to ventral and from −4.56 to 8.06 mm ventral to dorsal; while for the whole head template these values are −15.29 to 16.58 mm, −21.41 to 16.84 mm, and −18.07 to 13.81 mm, respectively. The origin (anterior commissure) and orientation of the coordinate system is very similar to the convention established in (Hawrylycz et al., 2011), i.e., the “Waxholm space.” Voxel size (isotropic 1.25 mm3) is the actual value multiplied by 10 to make images usable in standard human software like SPM, i.e., 1 mm means 100 μm. The orientation of the images is left–right (neurological).

Assessment of the probabilistic tissue segmentations

We provide probabilistic images of the gray matter, white matter, and cerebrospinal fluid of the Wistar rat in the P & W space. These types of images are considered as spatial priors in SPM segmentation (Ashburner and Friston, 2005; Friston et al., 2010). To obtain these priors, the native “brain only” MRIs were firstly segmented into tissue types by using a simple histogram analysis. Second, these native segmentations were warped to the MDT space using the transformations obtained in the subsection above to finally average them through individuals.

Figure 1B depicts the histogram of the individual “brain only” MRI shown in Figure 1A. Despite the evident contrast between white and gray matter in this image, suggesting the existence of at least two populations of voxels with intensities around two different intensity peaks, only a single global peak is clearly seen in Figure 1B. This is possibly due to noise and partial volume effects. Therefore, before any further analysis of the histogram, we applied an anisotropic filter (Perona and Malik, 1990) to the image. The effect of this filter is well known: its denoising effect narrows the distribution of the existing intensity populations while penalizing the smoothing operation in those regions where sharp edges are present (high gradient) avoiding additional partial volumes caused by the smoothing.

Figure 1

(A) Individual low resolution MRI of a rat and (B) its intensity histogram. Note that there is a global peak only, though it seems there is an unclear second peak.

The histogram of each filtered “brain only” MRI was fitted with four Gaussian curves, representing non-brain tissue, white matter, gray matter, and cerebrospinal fluid. Additionally three non-parametric curves, representing partial volume populations, fitted the residual part of the histogram not fitted by the Gaussians. The latter curves were obtained with a roughly empirical method: we subtracted the Gaussian curves from the histogram and zeroed the negative values, dividing the result into the three partial volume regions bounded by the centers of the Gaussians. Each voxel was labeled with seven probabilities, Pnb, Pnb/WM, PWM, PWM/GM, PGM, PGM/CSF, and PGM/CSF, corresponding to the values, at the intensity of the voxel, of the seven fitted functions normalized to sum 1. Based on purely empirical criteria, these probabilistic classifications were reduced to three tissue types, namely:

where μ = arg max P(i), m = arg mediani P(i), H(i) is the Heaviside step function and f (i; μ) is a smooth function of the intensity i, decaying from its maximum at μ until zero at m and m; y and z being the tissue types at both sides of tissue x in the histogram. This procedure is inspired by (Kovacevic et al., 2002) for segmenting human brains. Notwithstanding, they are merely empirical segmentations, they served as spatial priors in SPM8 (without the normalization step) to obtain the conclusive and robust probabilistic tissue segmentations of white matter, gray matter, and cerebrospinal fluid of the individual MRI. These final segmentations were binarized by means of the maximum probabilistic (MaxProb) criterion, registered to the MDT space and averaged. Finally, these averages were manually corrected by an expert to constitute the probabilistic priors of our new template set.

The MaxProb criterion is:

Note that the probabilistic segmentations might be the logical choice for the averaging. We found however they yield an underestimation of the probabilistic gray matter prior, even after the proper normalization of the priors to sum 1. This consequently produced an incorrect too low probability gray matter map from the SPM automatic segmentation of the individual MRIs. Applying the MaxProb criterion overcame this situation.

Digitalization of the Paxinos and Watson atlas

The sixth edition of the P & W stereotaxic atlas was originally provided in PDF format (Paxinos and Watson, 2007). We selected a total of 48 cortical structures on each hemisphere (total: 96) that were manually delineated with Adobe Illustrator in the coronal slices of this PDF book. The acronyms used by Paxinos and Watson appear in Table 1. We obtained a BITMAP image of each coronal slice containing a set of closed and filled contours of different colors representing the cortical structures present in the slice. Ad hoc software was developed to align and interpolate these bitmaps to produce a digital volume in the P & W space (isotropic voxel size of 0.12 mm) representing a digital parcellation image with the 96 selected cortical structures.

Table 1

Mean and SD of the volume of the anatomical structures of our sample.

Structure name	Sample mean ± SD (mm3)		Structure name	Sample mean ± SD (mm3)
	Right	Left		Right	Left
AID	8.10 ± 1.52	8.36 ± 1.08	PtPR	0.56 ± 0.09	0.59 ± 0.09
AIP	3.72 ± 0.62	4.52 ± 0.57	RSD	10.51 ± 0.83	9.01 ± 0.69
AIV	1.06 ± 0.17	1.01 ± 0.11	RSGb	1.31 ± 0.19	1.22 ± 0.23
APir	5.76 ± 1.45	3.19 ± 1.99	RSGc	1.90 ± 0.29	2.10 ± 0.37
Au1	10.19 ± 1.53	12.06 ± 1.30	S1	0.76 ± 0.13	0.84 ± 0.14
AUD	3.75 ± 0.51	3.84 ± 0.43	S1BF	14.27 ± 1.65	16.80 ± 2.50
AuV	2.97 ± 0.83	3.56 ± 0.82	S1DZ	5.06 ± 0.73	6.11 ± 0.85
Cg1	13.39 ± 1.45	12.58 ± 1.52	S1DZ0	2.17 ± 0.32	2.07 ± 0.23
Cg2	2.91 ± 0.48	3.65 ± 0.55	S1FL	8.55 ± 0.99	8.66 ± 1.01
DI	5.10 ± 0.60	5.06 ± 0.43	S1HL	4.78 ± 0.54	4.44 ± 1.14
DIEnt	2.69 ± 0.30	3.76 ± 0.42	S1J	4.25 ± 0.66	4.65 ± 0.52
DLEnt	8.21 ± 0.68	8.90 ± 1.01	S1Sh	0.22 ± 0.05	0.25 ± 0.06
DLO	3.75 ± 1.15	3.32 ± 0.95	S1Tr	1.50 ± 0.23	1.38 ± 0.32
Ect	7.01 ± 1.27	7.03 ± 1.12	S1ULp	9.46 ± 1.03	11.60 ± 1.26
Fr3	4.77 ± 0.82	5.63 ± 0.70	S2	4.45 ± 0.87	5.35 ± 1.00
GI	4.69 ± 0.54	4.92 ± 0.63	TeA	2.79 ± 0.41	3.28 ± 0.52
GIDI	0.20 ± 0.06	0.19 ± 0.06	V1	4.41 ± 0.55	4.68 ± 0.72
LPtA	3.77 ± 0.43	3.57 ± 0.61	V1B	9.46 ± 1.15	10.80 ± 1.08
M1	16.82 ± 2.36	16.67 ± 2.76	V1M	5.28 ± 0.84	3.85 ± 0.70
M2	8.27 ± 1.26	5.59 ± 1.19	V2L	3.59 ± 0.46	3.71 ± 0.61
MEnt	1.92 ± 0.41	1.74 ± 0.34	V2ML	2.34 ± 0.43	1.70 ± 0.53
MPtA	1.25 ± 0.29	1.07 ± 0.31	V2MM	0.84 ± 0.17	0.63 ± 0.19
PRh	0.74 ± 0.14	0.72 ± 0.12	VIEnt	0.62 ± 0.12	0.65 ± 0.22
PtPC	0.09 ± 0.03	0.08 ± 0.03	Cortex	222.28 ± 9.96	227.69 ± 10.66
PtPD	2.10 ± 0.25	2.26 ± 0.32	Hemisphere	862.80 ± 40.81	902.12 ± 45.33
Whole brain	1764.92 ± 85.57		White matter	1110.75 ± 44.27
Parenchyma	1729.92 ± 83.71		Cerebrospinal fluid	99.11 ± 16.95
Gray matter	619.17 ± 46.34

The structure names follow the same acronym used in the sixth edition of the P and W atlas.

This digital parcellation was registered and resampled with nearest neighbor interpolation to the template by using the approximate TPS registration method (Rohr et al., 2001). This is a modified version of the widely used standard (interpolating) TPS method (Bookstein, 1989), which estimates a smooth spatial transformation that maps the coordinates of a set of selected landmarks in a target space to their homologous in a source space. Unlike the standard TPS method, the approximate TPS method takes into account the uncertainty of the coordinates of the landmarks so the map between both sets is not exact. A total of 36 landmarks were defined in the P & W space (target) and their corresponding pairs were manually located in the template. Therefore, landmark coordinates are affected by human error making the use of the approximate TPS registration method more adequate than the standard one. The way this error is assigned to each landmark and implemented, taking into account the dimensions and anisotropy of the voxel, is described in Section “Details of the Implementation of the Approximate Thin Plate Spline Registration Method” in Appendix. The landmarks were located using T2-CARB software, a version of IMAGIC for rats.

The voxels of the warped atlas belonging to white matter, cerebrospinal fluid, or another non-brain part of the head, where removed. Further minor manual corrections were made with MRIcron.

SPM normalization, segmentation, and automatic labeling of structures, using SPM8

Using the tissue priors, the SPM8 unified segmentation method was tested in both the high and low resolution MRIs. This method provides individual tissue segmentation into gray matter, white matter and cerebrospinal fluid, and the forward and inverse transformations to the template space. The automatic labeling of the native MRI into different gray matter structures was also tested by applying the inverse transformation to digitalized atlas. This classification was further refined so that only those voxels belonging to an individual binary gray matter mask remained. By means of the hidden random Markov field theory we excluded isolated points and filled small holes in the segmentation (Cuadra et al., 2005). This has been usually done with IBASPM (an extension toolbox of SPM5; Alemán-Gómez et al., 2006) in several applications (Alemán-Gómez et al., 2006; Iturria-Medina et al., 2007; Sotero et al., 2007; Sanabria-Diaz et al., 2010).

Functional MRI test

After slice timing and inter-scan motion correction, spatial smoothing was performed to the EPI volumes using a Gaussian kernel of 0.8 mm full width at half maximum. Single-subject analysis was performed with a critical T-value for each voxel calculated for a significance level of p < 0.05 (FWE corrected). Then the average EPI was registered to the low resolution T2 MRI with the standard SPM8 12 parameters registration method, based on the normalized mutual information criterion. The resulting transformation was applied to the fMRI results to correct for geometrical distortions provoked by B0 inhomogeneities. We verified that this approach performed similarly to the correction method based on gradient field maps described in (Anderson et al., 2001; Hutton et al., 2002). After spatially normalizing the anatomical low resolution T2 MRI to the MDT space, the transformation was applied to the corrected fMRI result. The anatomical evaluation of fMRI localization could be done by overlaying the significant activations on the MDT and digital parcellation image.

Volumetry

Firstly, each individual MRI was warped using “elastix” to the MDT space by means of a highly parameterized, but still diffeomorphic (Jacobians > 0), transformation. The latter was then applied to the individual SPM8 gray matter, white matter, and cerebrospinal fluid segmentations. The volume of a desire cortical region is the sum of the Jacobians of the transformation through the voxels belonging simultaneously to the region (according to the digital parcellation image) and to the MaxProb normalized gray matter segmentation, times template voxel size. Likewise, the volumes of whole brain, gray matter, white matter, brain parenchyma, and cerebrospinal fluid were also calculated.

We then performed a simple survey about Wistar rat ontogenetic development by correlating the volumes of the structures with a developmental variable. We found that the performance of a volume–age regression model is limited due to the very low time accuracy in which age is usually recorded in many laboratories, including ours (1 week). Even in this lapse, rats undoubtedly exhibit dramatic increases in size and body weight, and consequently the brain could change as well. Based on empirical rationale, body weight, which is measured using a much higher resolution than age (1 g vs. 1 week), might be a proxy of age in the regression models in the studied postnatal period (P42–P70). In fact, we found a very strong and significant linear relationship between both variables for the age range of our sample [weight = (35.67 g/week) × age + 13.32, pslope ≈ 0, but pintercept > 0.1]. Therefore, we sought for volume–body weight relations by means of a robust linear regression method (DuMouchel and O’Brien, 1989), consisting of the estimation of the slope β and intercept η in:

This equation was fitted for the volume of each cortical structure, cortical hemisphere, whole brain, brain parenchyma, and cerebrospinal fluid.

It should be noted that we used the same individual MRIs for both the template construction and volume calculations. This could favorably bias our reports with respect to new MRIs from now on since the template is the minimum deformation space of the sample. However we only intend to probe the possible applications of the template set, encouraging rigorous experimental designs for future work.

Results

Minimum deformation template and digitalized parcellation image

Figures 2A,B show the histograms of the MRIs before and after the equalization respectively. Although not demonstrated in this paper, the high variability of the intensity peaks in Figure 2A suggests an undesirable effect when averaging without equalization. The first set of peaks in Figure 2B corresponds to non-brain tissue (excepting fluid), while the last two correspond to brain and cerebrospinal fluid. Note that the optimization procedure described in Section “Intensity ‘Equalization’ of the MRIs” in Appendix successfully “equalized” the last two peaks, as desired.

Figure 2

Histograms of all intensity corrected MRIs (A) before and (B) after the equalization, showing that the widths and positions of the peaks are transformed to very similar values in each individual.

Figures 3A,B depict the initial T2 MRI average template, based on linear registration only (namely modifying Eq. 1 into ), and the MDT, respectively. As previously demonstrated by Kovacevic et al. (2005) for mice, the MDT possesses a higher (though slighter) degree of detail. Moreover, our template appears to be less noisy and with better resolution than those presented in (Schweinhardt et al., 2003; Schwarz et al., 2006), since it was based a MRI system with a higher main magnetic field (7.0 vs. 4.7 T) and improved T2 protocol. The digitalized parcellation image is overlaid on several coronal slices of the MDT in Figure 4 to provide an idea about the spatial match between both images.

Figure 3

An axial slice of the (A) average MRI after linear registration to the “average” affine space and the (B) minimal deformation template.

Figure 4

Digitalized Paxinos and Watson atlas, overlaid on the minimal deformation template after being warped to the template space using the approximate thin plate splines (TPS) method, plus further manual corrections.

Probabilistic tissue segmentations in the MDT space

The effect of the anisotropic filter on the individual high resolution “brain only” MRI can be seen in Figure 5A. This image and the resulting histogram in Figure 5B can be directly compared with Figures 1A,B, respectively. Instead of the single peak shown in Figure 1B, two very distinguishable peaks are now present in the histogram of Figure 5B, corresponding to white matter and gray matter populations respectively. The Gaussians and non-parametric curves representing different tissue types and partial volumes are shown in this figure. The final individual probabilistic segmentations of gray matter, white matter, and cerebrospinal fluid are overlaid on the individual “brain only” MRI in Figure 5C.

Figure 5

(A) Filtered low resolution MRI shown in Figure 1 (right) and illustrative picture of the effect of the anisotropic filter on the image in the region enclosed by the yellow box (left): the image is more denoised along the low changes of intensities; the smoothing kernel is a tensor calculated with the gradient of the image at each point. (B) Histogram showing two identifiable peaks corresponding to the white matter and gray matter intensity populations (lilac curve). Four Gaussian distributions were fitted (cyan: non-brain, green: white matter, red: gray matter and blue: cerebrospinal fluid) and three non-parametric curves representing partial volumes. (C) Final tissue segmentations (same color code as the fitted Gaussian curves), after SPM8 post processing of the preliminary segmentations obtained from the histogram analysis.

The probabilistic segmentations (the so-called spatial priors in SPM) of tissue types are shown with colors in Figure 6. As a quantitative measure of the overlapping of the different tissue types between the template and the digitalized atlas, we found that, among the voxels classified as part of the cortex by the digitalized atlas, 6.0% and 5.7% are classified by the corresponding MaxProb segmentations as white matter and cerebrospinal fluid, respectively.

Figure 6

Axial slices of the gray matter (red), white matter (green), and cerebrospinal fluid (blue) probabilistic segmentations.

Our gray and white matter tissue priors resemble those obtained by (Sawiak et al., 2009) for mice, the sole paper in the literature assessing these segmentations for rodents. In that paper, the three tissue classes were obtained directly from the average template image and smooth variability was provided by applying a 250-μm Gaussian filter. In contrast, we exploited the actual variability of our sample by averaging all individual MaxProb binary segmentations, achieving tissue probabilities.

SPM8 unified segmentation results

For the sake of brevity, we shall only present the results for the low resolution MRI. The spatially normalized low resolution MRI is shown together with the average template in Figure 7 for three orthogonal views. Figures 8A,B show the individual tissue segmentations and automatic labeling, respectively, overlaid on a set of coronal slices of the individual low resolution MRI in its native space.

Figure 7

Registered orthogonal views of the average MRI template (right) and spatially normalized low resolution MRI (left). In order to provide a comprehensive view, we constructed this picture with the “Check Reg” button of SPM8, which is the usual tool of SPM users for the visualization of registration results.

Figure 8

(A) Tissue segmentations (same color code as 6) and (B) atlasing overlaid on the individual coronal slices of the low resolution MRI in the native space.

We have described a histogram-based analysis for obtaining the preliminary individual segmentations. We clarify that this is not a practical technique for automatically segmenting an individual MRI since histogram construction requires extraction of the brain mask. In contrast, the SPM segmentation methods (or another Bayesian technique), based on spatial priors, can directly deal with the individual raw whole head MRI.

The mean and SD of the brain structure volumes are shown in Table 1. A way to at least evaluate the performance of brain volume estimations is by means of the similarity index, originally defined as the Dice’s coefficient (Dice, 1945) and used in (Zijdenbos et al., 1994) for Neuroimaging. This index was calculated between the brain mask estimated via MaxProb SPM segmentations and the one extracted with the semi-automatic method described in the Section “Materials and Methods,” which we consider as the ground truth. For our 30 rats, the average SI is very high (around 0.97 ± 0.00), if we take into account that SI takes values from 0 to 1 corresponding to 0 and 100% of overlap, respectively. where the subscripts M and A denote MaxProb and Atlas brain masks, respectively.

The reported volume means and deviations might constitute tentative references for future work in detecting altered volumes of structures in rat models of diseases or atrophy provoked by injuries.

Table 2 contains the regions presenting a significant slope (p < 0.01) in the regression model (4) as well as the putative expected percent of volume change within the range of weights of our sample. This percent is calculated as: % = (β × Δw/〈V〉) × 100, where Δw is the range of weights and 〈V〉 is the mean volume of our sample. The average of the slope of the regression is 1.042 ± 0.228 mm3/g for the entire brain. A simple and intuitive linear dimensional analysis gives us an idea about the plausibility of those values. Table 1 of the P & W book (Paxinos and Watson, 2007) presented, for different body weights of similar rats, the distances between the Bregma to the Interaural Plane in both the anterior–posterior (AP) and dorso-ventral (DV) directions. The slope of the regression between (AP*DV)3/2 (assuming an isometric scaling between brain linear measures and brain volume) and the corresponding rat body weights is very significant, and equal to 1.234 ± 0.117 mm3/g, which is within the range of slope values (1.042 ± 0.228 mm3/g) between brain volume and body weight in our sample.

Table 2

Structures with a significant slope in the regression of volume vs. body weight.

Structure	% of volume change in 6–10 weeks	p-Value (* < 0.001)	Plot X axis: weight (g); Y axis: volume (mm3)
Gray matter	7.53 ± 2.03	0.0009*
White matter	16.50 ± 4.31	0.0007*
Cerebrospinal fluid	34.64 ± 10.27	0.0022
Parenchyma	10.52 ± 2.76	0.0007*
Right hemisphere	11.28 ± 2.53	0.0001*
Left hemisphere	12.05 ± 2.66	0.0001*
Whole brain	11.74 ± 2.57	0.0001*

The second column shows the expected percent of volume growth of each structure within the range of weights in our sample, which might be considered a proxy of age within the period from 6 to 10 weeks. The fourth column depicts body weight (g) – volume (mm.

Although the volume of the entire hemispheres, brain, parenchyma, white matter, and cerebrospinal fluid are significantly correlated with body weight, the cerebral cortex is not. This result resembles those found in (Maheswaran et al., 2009a,b) for wild type mice, where brain, cerebellum, and hippocampus significantly and positively correlated with age, while neocortex did not. Moreover (Chuang et al., 2011) has recently fitted a sigmoid function to C57BL/6 mouse volume development from embryonic stages (E0–E18) to postnatal ages (P0–P80 and older). Inspecting the figures of that paper, we believe that indeed a linear model is suitable for the postnatal age period (P42–P70) used in our paper. We also observe that while the brain, hippocampus, and cerebellum again present a clear tendency to grow with age, cortex volume is likely to remain constant during this period. As postulated in (Chuang et al., 2011), the constancy of cortical volume appear to be at expenses of the growth of other brain structures, say white matter.

On the other hand, no substructures constituting the cortex resulted to be significantly correlated (according to p < 0.01) with body weight. Although this appears to be reasonable upon the fact of cortical volume constancy, it is in disagreement, for example, with (Riddle et al., 1992), where a growth of the somatosensory cortex was demonstrated using succinic dehydrogenase histochemistry and digital image analysis. This might be due to the age period used in this paper (P42–P70) in contrast to the wider range used in (Riddle et al., 1992; two groups: P7 and P70–P84).

ROI-based anatomical interpretation of fMRI

Figure 9 shows, on the average template, the fMRI activations caused by forepaw and hindpaw stimulations. As expected (Spenger et al., 2000), they are located within the contralateral forelimb (S1FL) and hindlimb (S1HL) regions of the primary somatosensory cortex, respectively. To illustrate the accuracy of this localization we calculated the percent of voxels within the corresponding anatomical region which have a T-value above 6, 7, 8, and 9. Respectively, these values are 83.1, 89.5, 94.2, and 98.9% for the left forepaw; 91.8, 98.4, 99.6, and 100.0% for the right forepaw; 85.3, 95.3, 100.0, 100.0% for the left hindpaw; and 81.8, 97.3, 100.0, and 100.0% for the right hindpaw. This means that the center and most parts of the significantly activated regions are within the corresponding expected anatomical areas.

Figure 9

The fMRI activations induced by forepaw and hindpaw electrical stimulations. The left (green) and right (violet) forepaw activations are located within the contralateral forelimb regions of the primary somatosensory cortex, denoted as S1FL in the P and W atlas (orange and red regions respectively). Likewise, the left (blue) and right (yellow) activations are within the contralateral hindlimb regions, denoted as S1HL (yellow and cyan regions respectively).

Our fMRI forepaw results resemble those obtained in (Schweinhardt et al., 2003) and (Lu et al., 2010). Like (Schweinhardt et al., 2003), we obtained fewer activated regions than (Lu et al., 2010). Note that in the latter, other somatosensory related areas (both ventrolateral and ventral posterolateral thalamic nuclei, S1DZ and S1ULp) were significantly activated, probably in relation to the use of an iron-oxide contrast agent, which enhances the signal to noise ratio of the image. Moreover and controversially, barrel field “S1BF” was also activated in (Lu et al., 2010).

Discussion

In this paper, we constructed a new MRI template set of Wistar rat comprising in the same space: (1) a MDT, intensity, and shape average of 30 T2-weighted MRIs, (2) probabilistic images of gray matter, white matter, cerebrospinal fluid, and brain, and (3) a digitalization of 96 cortical structures of the last version of the P & W atlas. We aim to provide a template set for the most basic (and interdependent) neuroimaging applications, namely unbiased spatial normalization, tissue (in special gray/white matter) segmentation and automatic labeling of cortical structures of individual rat structural T2 MRIs, since they are mandatory for the majority of neuroimaging studies. Thus, we validated them not only with both individual high resolution T2-weighted MRIs (acquired in 3 h) but with coarser ones, acquired in the relatively short time of 50 min. The latter case is of particular importance for certain long lasting experiments, such as fMRI or any type of diffusion weighted imaging (DWI), where the time must be gained at expenses of the structural MRI.

Comparison with previous rat MRI templates

Table 3 summarizes the main characteristics of our template set, the procedures for its construction and the tested applications, facing them against the two previous attempts. In this table, bold letters highlight the more evidently advantageous features. The major disadvantages of our paper are logistic: a medium sample size (though we used Wistar, a strain of rat never used before for a template construction) and the lack of the highest MR main magnetic fields available, e.g., 9.4 T.

Table 3

Comparison between three approaches to obtain a template set for rats.

Paper	N/strain	MRI intensity corrections	Individual MRI registration previous average template	Template space	Probabilistic tissue segmentations	P & W atlas edition	Atlas-template registration	Neuroimaging applications
	MR system main field					N
Schweinhardt et al. (2003)	5/in vivo SD	Histogram equalization (not described)	Affine, based on minimizing the sum of the weighted distances between 49 landmarks pairs in the individual and atlas space	Paxinos and Watson	–	Second (Paxinos and Watson, 1989)	–	SPM spatial normalization. Right and left forepaw BOLD fMRI localization on the primary forelimb area of the somatosensory cortex (S1FL), visually evaluated (tested)
	4.7 T
Schwarz et al. (2006)	97/in vivo SD	Receiver coil correction and square root to reduce inter-variability and enhance tissue separation	PCA-based rigid to P & W atlas, matching brains and atlas outer surfaces	Arbitrary	Brain, brain parenchyma, and CSF	Fourth (Paxinos and Watson, 1998)	–	SPM spatial normalization. Pharmacological-induced perfusion MRI (phMRI) localization with volume of interest (VOI) evaluation of both localization and spatial extent (tested)
	4.7 T					468
Current	30/in vivo Wistar	Receiver coil correction, histogram equalization (1D linear registration) and, for segmentation, anisotropic filtering to enhance tissue separation	Highly non-linear, based on a multi-resolution, multi-grid scheme	Minimal deformation	Brain, gray matter, white matter, and CSF	Sixth and latest (Paxinos and Watson, 2007)	ATPS based on weighted landmarks	SPM spatial normalization, SPM tissue segmentation and atlasing in low resolution MRI. Right and left of both forepaw and hindpaw fMRI localization with quantitative ROI evaluation and volumetry (all tested). ROI–ROI DTI-based connectivity (future work) morphological connectivity (future work)
	7.0 T					96

Sprague-Dawley; P and W, Paxinos and Watson; ATPS, approximate thin plate splines.

Intended to the sole purpose of a rough localization functional features, both previous template sets were built by means of linear registration methods to biased spaces. Indeed, several studies have been using these templates for localization of BOLD fMRI (Ramu et al., 2006; Endo et al., 2007; Reyt et al., 2010; Sumiyoshi et al., 2011), position emission tomography (PET; Fujita et al., 2005; Casteels et al., 2006; Frumberg et al., 2007; Hjornevik et al., 2008; Sung et al., 2009), pharmacological MRI (phMRI; Gozzi et al., 2007, 2008; Schwarz et al., 2007; Dodd et al., 2010), or manganese-enhanced MRI (MEMRI; Eschenko et al., 2010). In contrast, our registration method sought for the MDT of our sample, avoiding a bias in registration methods toward those heads most similar to the template (Guimond, 2000; Kochunov et al., 2001a,b). Besides, the use of highly parameterized non-linear deformations allowed for the achievement of a high detailed MDT, storing the most part of the sample variability in the deformation fields. Thus, our template set is not only useful for functional localization but also for unbiased, robust and detailed anatomical studies based on non-linear registrations. In addition we provide, for the first time in rats, gray and white matter probabilistic segmentations. This allows for the application of Bayesian segmentation methods, such as FAST–FSL or SPM8 unified segmentation (Ashburner and Friston, 2005).

Only (Schwarz et al., 2006) and our paper have provided digitalization of the P & W structures, being the former digitalization overwhelmingly superior to ours (468 vs. 96). Increasing the number of digitalized structures must therefore be a mandatory part of our main future work. However, we have used the latest version of the P & W atlas (sixth) while they used the fourth one and we have established a non-linear spatial correspondence between template and atlas.

Exemplified applications

In our paper, we proved that the fMRI activations elicited by forepaw and hindpaw stimulations in both hemispheres can be located within the expected anatomical regions of the P & W digitalization. Note that hindpaw stimulation was not tested previously, probably due to the difficulty in obtaining a robust and distinguishable activation.

As an example, we also tested the feasibility of using the template set for a volumetric analysis of the rat brain. In fact, to the knowledge of the authors, there is a scarcity of volumetric studies on rodent brains and, particularly, no such studies on rats. We proposed to calculate the volumes of the subdivided cortical structures of the atlasing. Thus, we were able to explore the possible relationship between these volumes and body weight, as a proxy of age, for our sample of male Wistar rats, opening a window for the study of the ontogenesis of the rat brain. Whole brain, whole gray and white matter, parenchyma, and cerebrospinal fluid, but not whole cortex or its subdivisions, are correlated with body weight and thus very likely with age within the period of our sample (P42–P70). Therefore, we suggest that either body weight or age should be taken as a covariate in any regression analysis between volumes and other variables of interest. Also, the normal developing curves of brain structures and substructures can be used as phenotypic signatures in rat disease models.

Comparison with mouse templates

As mentioned before, like all mice templates since the seminal paper of (Kovacevic et al., 2005), our registration method was non-linear to the minimum deformation space and we provide gray and white matter segmentations like (Sawiak et al., 2009). Mouse templates are also accompanied by 3D discrete images of labeled structures, such as (Dorr et al., 2008; 62 structures) and (Chen et al., 2006; 42 structures), or probabilistic approaches (Ma et al., 2005; Badea et al., 2007a). Like in our template set, the annotation of structures in those templates is usually based on standard references, e.g. (Paxinos and Franklin, 2000). Besides, we have defined a coordinate system which is inspired in a very recent effort for standardization of mice “atlasing,” i.e., the so-called Waxholm space (Hawrylycz et al., 2011).

Our template set is in vivo. This is the most contrastive feature when compared to all in vitro mouse templates in the literature, excepting the single case of (Ma et al., 2008). The disadvantage of in vitro templates is the loss of representativeness of in vivo samples, due to ventricle collapse, myelin water fraction changes, and manipulation related distortions due to the skull stripping procedure despite various papers proposed to overcome the last problem by leaving the brain in situ (inside the cranium). While in vitro template have a very high resolution, in vivo templates are more suitable for in vivo applications, such as longitudinal studies of the same animals and fMRI.

Applications of our template set

Apart from precise localization of functional features like (Schwarz et al., 2006), identifying phenotypes like (Bock et al., 2006; Badea et al., 2007b), or the quantification of volume and surface of structures like (Badea et al., 2007a), the automatic labeling might be used for ROI to ROI DTI-based tractography. Seed and target ROIs can be obtained from the automatic labeling, whereas tracts can be restricted to go through the parenchyma or even white matter. An immediate future work must be the use of this tractography for assessing robust and systematical measures of ROI–ROI anatomical connectivity, as has been proposed elsewhere (Iturria-Medina et al., 2007, Iturria-Medina et al., 2008). We perceive a future interest in using this technique to probe the anatomical substrates of various fMRI-proven connections between either ipsilateral or contralateral areas (Lu et al., 2007; Pawela et al., 2008; Bifone et al., 2010; Magnuson et al., 2010). The feasibility of combining different neuroimaging techniques to assess connectivity, such as fMRI, DTI, and the recent so-called “morphological connectivity,” which is based in either volume or area ROI/ROI correlations (Sanabria-Diaz et al., 2010), opens a wide window for understanding rat brain anatomical and functional networking.

As was already demonstrated in primates and humans, neuroimaging characterization of rats will aid the understanding of the anatomical and physiological substrates of different brain pathologies and contribute to unraveling the ontogenetic and phylogenetic development of mammals.

Practical advices for neuroimagers

Our template set can be used directly with any neuroimaging software for spatial registration, segmentation, and automatic labeling of structures. We only alert the reader about helping the registration with a previous alignment of the MRI according to the template orientation and setting the origin in the anterior commissure. Our template set and some helpful codes are publicly available through the SPM extension webpage http://www.fil.ion.ucl.ac.uk/spm/ext/, or at http://www.idac.tohoku.ac.jp/bir/en/db/rb/101028.html.

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.