Application of MRI-Based Partial-Volume Correction to the Analysis of PET Images of μ-Opioid Receptors Using Statistical Parametric Mapping

Step 1.

The subject-specific MR image was segmented into 3 compartments, GM, white matter (WM), and cerebrospinal fluid (CSF), using probability maps in SPM99. A little radiofrequency (RF) nonuniformity correction was chosen for MR segmentation to compensate for the differences in position of intracranial structures within the MRI head coil (12).

Step 2.

Voxels were assigned into GM, WM, and CSF tissue classes according to the maximum probability encountered for each voxel across the 3 datasets. These 3 tissue classes were subsequently binarized (given a value of 0 for absence of tissue and 1 for presence of tissue).

Step 3.

A mock-PET dataset was simulated from the binarized MR image by weighting the corresponding tissue classes with realistic, empirically determined relative values (GM × 1.0 + WM × 0.3 + CSF × 0.1) and smoothing using a 3-dimensional gaussian kernel approximation of the PSF of the PET scanner. This allows a mock-PET–PET intramodality coregistration, which we have found to be faster and more efficient than MRI–PET coregistration.

Step 4.

The MR image was coregistered to the PET image using the simulated PET data and applying the resultant transformation matrix to the output of step 2. At this point, the binarized GM, WM, and CSF were resliced using nearest-neighbor interpolation to remain in the space (i.e., dimensions and voxel size) of the MR image.

Step 5.

The GM, WM, and CSF binarized images were convolved using a 3-dimensional gaussian kernel corresponding to the PSF of the PET scanner.

Step 6.

The smoothed data were resliced and resampled to the space of the PET image using trilinear interpolation.

Step 7.

The WM and CSF averages were determined by applying a much higher-probability threshold (P > 0.90) to the output from step 1 to classify and binarize only voxels that had a high probability of being the respective tissue. After coregistration and reslicing, these 2 files were smoothed to the resolution of the PET scanner, and another threshold (P > 0.95) was applied to get the subset of WM and CSF voxels that were relatively unaffected by coregistration and smoothing. These WM and CSF voxels were used as a mask to sample the observed PET image and compute a reliable estimate of the average WM and CSF activity.

Step 8.

The scaled WM and CSF compartments obtained after multiplying by the respective average WM and CSF activities were subtracted from the observed PET scan to correct for spill-in of counts from the surrounding medium to GM voxels. Finally, the partial-volume–corrected PET image was obtained by dividing this quantity by the GM recovery coefficients image (GM from the output of step 5) to correct for spill-out of counts from the GM voxels to the surrounding medium. Because the recovery coefficients image also contains low recovery coefficients, a correction mask including only voxels greater than 0.21 (empirically determined value) in the GM recovery coefficients image (where the value varied between 0 and 1) was used to avoid large correction amplification.

Segmentation and Inhomogeneity Correction.

SPM99 uses a clustering algorithm that classifies MR voxels into GM, WM, and CSF on the basis of their intensities and a priori information about their spatial localization, derived from probability maps obtained from many healthy subjects. To minimize partial-volume effects between tissue types, high-resolution and contrast T1-weighted MR images were obtained, and a little RF nonuniformity correction (12) was applied to correct for smooth intensity variations in the MR image. This correction was based on the approximate amount of nonuniformity in our MRI data but does not compensate for systematic biases in the segmentation method that result from a decrease in GM/WM contrast with aging.

In both young and elderly subjects, the MR segmentation algorithm performed well as assessed visually by the high spatial concordance of gray nuclei and cortical structures (i.e., small cortical sulci) between the actual MR image and the GM tissue class as well as other tissue classes. The mismatch in resolution between our MRI data (0.93 × 0.93 × 1.5 mm) and the SPM template (2 × 2 × 2 mm) used for MR segmentation did not introduce any noticeable segmentation inaccuracies as a result of partial-volume averaging. More formal assessment of the quality of segmentation using simulated data (12) has shown that this segmentation method is robust and accurate for high-quality T1-weighted MR images, such as our MRI data. However, some voxels were misclassified as GM at the interface of WM and ventricles. A small amount of nonbrain tissue was also misclassified as GM but was removed almost entirely using morphologic operations in SPM99, by which rendered surface brain is extracted. To improve tissue classification, segmented data for which the probability varied between 0 and 1 were binarized by finding for each voxel its maximum probability among the 3 classes of tissue and assigning this voxel to the appropriate class. In contrast to another study, in which only voxels with P > 0.33 (13) in each of the 3 tissue classes were kept, no additional threshold needed to be applied to our data because no differences were noted between data for which a threshold of P > 0.33 was applied and data for which it was not.

MR–PET Coregistration.

To register MRI segments with PET data, a simulated PET image was first obtained after convolving a weighted sum of the 3 MRI brain compartments to the PSF of the actual PET data. Parameters describing the rigid-body registration between simulated and actual PET images were estimated, and the transformation was applied to MRI brain compartments. The minimization of the squared differences between the simulated PET image and the actual PET image used the Gauss–Newton optimization and in all cases converged rapidly, giving high registration accuracy that was comparable to the level of accuracy achieved by automated MRI–PET (14) or MRI–SPECT with SPM registration (15) based on a similar methodology. To preserve voxel intensities and the ratio of brain to nonbrain present in the original high-resolution MRI data, nearest-neighbor interpolation was used to reslice the registered MR image to MRI space.

Smoothing to Obtain Recovery Coefficients (Spill-out Correction).

Data were convolved with a slightly overestimated PET PSF and resliced to the PET space using trilinear interpolation to compensate for MRI data downsampling. The measured PET PSF was shown to present a 5%–22% variation within 10 cm of the center of the field of view, where the subject’s head is usually confined, which may translate into less than a 7% error in GM recovery coefficients based on PVC simulations (16). These simulations have also shown that errors due to resolution mismatch between convolved MRI data and PET are less severe when scanner resolution is overestimated, as in our study, than when scanner resolution is underestimated, indicating that conservative estimates of PET resolution are accompanied by a more accurate PVC.

Spill-in Correction.

Prior PET experiments after naloxone blockade (11) have shown low and uniform activity in the WM, consistent with the low levels of μ-OR reported by autoradiographic studies (17) and the assumption of WM homogeneity of this PVC algorithm. Simulations found that a 10% variation in WM uniformity will result in less than a 1% error in PVC if WM homogeneity is assumed.

The rationale for this step was to obtain an average WM value from WM voxels that were spatially remote from GM to avoid spill-in of GM activity. For this purpose, a high-probability threshold (P > 0.90) was applied to the WM segments to classify and binarize only voxels that have a very high probability of being WM voxels. After coregistration and reslicing, these subsets of voxels were smoothed to the resolution of the PET scanner and another threshold (P > 0.95) was applied to get the subset of WM voxels relatively unaffected by coregistration and smoothing. These voxels were located at greater than 2 times the PET resolution from any GM voxel in the actual PET or MR image and were consequently relatively unaffected by GM spill-in. These WM voxels are used as a mask to sample the observed PET image and compute a reliable estimate of the average WM activity. Because, in a normal-sized ventricular system, spill-in from GM is likely to contaminate all CSF voxels, CSF radioactivity can be measured invasively by lumbar puncture or noninvasively with PET only in subjects with an enlarged ventricular system. In this study, the assumption of near-zero CSF activity at equilibrium after intravenous administration of ¹¹C-carfentanil is consistent with the low CSF radioactivity that we measured recently using a high-resolution PET scanner in subjects with an enlarged ventricular system (unpublished data, 2003).

Correction Mask.

Finally, a correction mask including only voxels greater than 0.21 in the GM recovery coefficients image (in which values vary between 0 and 1) was used to avoid large correction amplification due to low recovery coefficients at the image periphery. This threshold, empirically determined, resulted in the best tradeoff in sensitivity and specificity in identifying GM voxels because it preserved GM voxels of interest but removed misclassified voxels at interfaces between GM or WM and CSF where the recovery coefficient was small.

Spatial Normalization and Smoothing.

First, 0- to 90-min average images were created to generate an image with high anatomic detail. This image was then normalized to an ¹¹C-carfentanil template (average, 0–90 min) in a standard space (Montreal Neurological Institute [MNI]) having voxel size equal to 2 × 2 × 2 mm, and the transformation matrix was applied to the 34- to 82-min average image of the non-PVC and PVC data used for subsequent analysis (18). Stereotactic normalization used a 12-parameter affine transformation followed by nonlinear deformations using a linear combination of 3-dimensional discrete cosine transform basis functions (4 × 5 × 4) and 12 nonlinear iterations with medium nonlinear regularization (19). Parametric maps of μ-OR binding were calculated using a tissue-ratio measure ([regional counts/mean occipital counts] − 1) of the spatially normalized average 34- to 82-min data, where the mean occipital counts represented the nonspecific binding. The tissue-ratio measure is proportional to the binding potential B_max/K_D, where B_max represents the total receptor density and K_D is the affinity constant. The tissue-ratio measure is independent of changes in tracer transport resulting from alterations in regional CBF (20). Data were smoothed using a 12-mm FWHM gaussian smoothing kernel to a final isotropic resolution of approximately 14 mm.

Statistical Analysis and Results.

The relationship between μ-OR binding and age was determined at each voxel using the general linear model (21). Analysis of covariance (ANCOVA) was used to remove the confounding effect (nuisance covariate) of intersubject differences in global μ-OR binding (21) to obtain adjusted μ-OR binding values at each voxel. Regional effects of age were assessed using linear regression of μ-OR binding with age at each voxel, which resulted in a t statistic at each voxel and an associated P value. Correction of the P value for the multiple comparisons of voxels in the total volume S used random-field-theory formulas and took into account the high correlation between contiguous voxels due to the image reconstruction process and smoothing (22). Connected sets or clusters of voxels were considered significant if they survived a height threshold of P ≤ 0.001 (uncorrected) and their spatial extent k (size) at this threshold had a volume greater than 49 voxels, the spatial resolution of the PET scanner (i.e., PET resolution/voxel volume = 7 × 7 × 8 mm/2 × 2 × 2 mm).

VOI Analysis of the Age Regression on μ-OR Binding.

The mean adjusted value computed over all voxels in significant clusters shown by PVC was obtained for each subject and compared with mean adjusted values obtained by sampling exactly the same spatial location on non-PVC data (after lowering the height threshold). A test of homogeneity of slopes of the linear regression of age on adjusted μ-OR binding for PVC and non-PVC data was obtained. One-factor ANCOVA with non-PVC and PVC data as grouping nominal variables, age as a covariate, and μ-OR as a dependent variable was used to assess statistical differences in regression slopes. In addition, the correlation coefficients between age and non-PVC data and between age and PVC data for each cluster were assessed for significant statistical differences using the method of Cohen and Cohen (23).

Rate of Change with Age in the Absolute Amount of GM.

To assess how regional differences in GM volume in native space (before spatial normalization to MNI space) may influence PVC, spatially normalized GM segments were multiplied (modulated) by the Jacobean determinants calculated from the spatial normalization process. This produced a measure of absolute amount (density) of GM in native space (12,24), in contrast to nonmodulated data, for which a measure of relative amount (concentration) of GM is obtained. These images were then smoothed with a 12-mm gaussian kernel. Statistical analysis used ANCOVA to remove global differences in GM density, and contrasts [1 0 0] and [−1 0 0] were used to assess, respectively, positive and negative correlation of GM density with age. VOIs were placed on significant clusters, and the percentage change in GM density per year was computed.