Multimodal partial volume correction: Application to [ 11 C]PIB PET/MRI myelin imaging in multiple sclerosis

Abstract

The [¹¹C]PIB PET tracer, originally developed for amyloid imaging, has been recently repurposed to quantify demyelination and remyelination in multiple sclerosis (MS). Myelin PET imaging, however, is limited by its low resolution that deteriorates the quantification accuracy of white matter (WM) lesions. Here, we introduce a novel partial volume correction (PVC) method called Multiresolution–Multimodal Resolution-Recovery (MM-RR), which uses the wavelet transform and a synergistic statistical model to exploit MRI structural images to improve the resolution of [¹¹C]PIB PET myelin imaging. MM-RR performance was tested on a phantom acquisition and in a dataset comprising [¹¹C]PIB PET and MR T1- and T2-weighted images of 8 healthy controls and 20 MS patients. For the control group, the MM-RR PET images showed an average increase of 5.7% in WM uptake while the grey-matter (GM) uptake remained constant, resulting in +31% WM/GM contrast. Furthermore, MM-RR PET binding maps correlated significantly with the mRNA expressions of the most represented proteins in the myelin sheath (R²= 0.57 ± 0.09). In the patient group, MM-RR PET images showed sharper lesion contours and significant improvement in normal-appearing tissue/WM-lesion contrast compared to standard PET (contrast improvement > +40%). These results were consistent with MM-RR performances in phantom experiments.

Keywords

Partial volume correction PET/MR wavelet myelin multiple sclerosis

Introduction

Multiple sclerosis (MS) is a complex neurological disorder that represents the first non-traumatic cause of neurological disability among young adults.¹ MS pathophysiology is characterized by an auto-immune aggression of myelin sheaths resulting in demyelinated lesions and axonal degeneration. Myelin is, however, a dynamic tissue and there is strong evidence in animal models that myelin repair is an efficient process which may follow myelin insult.² However, little is known about the dynamics of re-myelination in patients with MS. Sensitive imaging tools are now required to measure myelin dynamics in-vivo to investigate spontaneous remyelination in MS patients and to act as biomarkers for new pro-myelinating therapies.

Positron emission tomography (PET) with [¹¹C]PIB was originally developed to image amyloid deposition in neurodegenerative disorders and dementia³ and has been recently repurposed for myelin imaging in vivo in humans^4–6 (Figure 1). The essential steps to measure de-myelination and re-myelination consist of lesions detection, quantification of myelin content and comparison between baseline and follow-up scans. In this perspective, the influence of partial volume effect (PVE) is relevant when detecting and quantifying myelin changes in lesions, considering the typical poor image resolution of PET and the quantification bias resulting from activity spill-in.

Figure 1.

Representation of [¹¹C]PIB PET myelin signal. (a) Axon structure with visible alternation of myelin lipid bilayer and myelin water content; (b) detailed representation of myelin structure where the lipid bilayers and the myelin basic proteins (MBP) in between are highlighted; (c) schematic representation of [¹¹C]PIB PET myelin signal: [¹¹C]PIB, in analogy to its interaction to amyloid fibrils, might have multiple interactions with the myelin structure: on one side it might get trapped into β-sheet structures of myelin proteins such as MBP,^41–43 and on the other side it might be highly soluble in the myelin associated lipid bilayer.⁴⁴

The literature provides a wide range of partial volume correction (PVC) techniques to address the PVE in PET imaging.⁷ One group is represented by region-of-interest (ROI) based methods, which are, however, limited by the assumption of radiotracer homogeneous distribution within each ROI.^7,8 This condition is clearly inapplicable to MS lesions due to their heterogeneity.^9,10 Alternatively, voxel-based techniques like partition-based¹¹ and multiresolution methods^12,13 are also available. A distinct final approach to PVC consists of incorporating a model of the system point spread function (PSF) into the reconstruction algorithm.^14,15

In this work, we introduce the Multiresolution–Multimodal Resolution-Recovery (MM-RR) method to achieve improved myelin quantification in [¹¹C]PIB PET brain images by performing effective PVC using T1-weighted MR sequences. MM-RR stems from previous works on brain PET/MR data^13,16–18 and whole-body PET/CT data.¹⁹ Similarly to these approaches, the relation between functional and anatomical images is exploited in a synergistic fashion for a realistic noise-controlled recovering of PET image resolution. However, the new algorithm takes a step further by modelling the differential sensitivity and specificity to myelin of the two modalities without the requirement of any anatomical segmentation which would not be applicable in the instance of the MS population. The ultimate clinical aim is to accurately quantify lesions uptake according to their degree of myelin loss and repair.

Hence, we focused the validation of the approach on clinically relevant metrics¹⁹ and specifically we wanted to demonstrate that MM-RR has the potential to: (1) significantly improve white/grey matter contrast in normal [¹¹C]PIB PET data without quantitative distortion of grey matter values; (2) enhance the contrast between normal appearing white matter (NAWM) and T2 positive perilesional and lesion signals, black holes and Gadolinium positive lesions; (3) ameliorate the partial volume effect in small lesions and (4) reduce the partial volume effect around the ventricular spaces. In addition to the clinical investigation, the method was further validated with two additional analyses. Firstly, the performance of the MM-RR method was compared against standard PVC methods in phantom experiments. Then, similarly to Veronese et al.,²⁰ MM-RR myelin maps from healthy controls were tested against a brain-wide mRNA atlas of myelin-associated proteins that acted as an independent predictor of myelin density.

Materials and methods

The MM-RR method

The MM-RR algorithm stems from previously developed methodology on partial volume correction^13,19 and denoising²¹ where the structural information (CT or MRI) was exploited to improve the image quality of the functional counterpart (PET). As in these techniques, MM-RR is based on the wavelet decomposition of both functional and structural images.

In the existing techniques, the high-frequency wavelet coefficients are transferred from the high-resolution anatomical image to the PET after appropriate statistical weighting. However, prior to the wavelet decomposition, the anatomical image is segmented by means of tissue classification based on image intensity levels or by using anatomical atlases, under the assumption of radiotracer homogeneous distribution within each ROI. This is inapplicable in case of demyelinating diseases like MS, due to the heterogeneity of myelin content between and within lesions.^9,10

In the formulation presented here, the image resolution recovery depends on the ad hoc modelling of the relationship between low-frequency information of PET and T1-weighted MR that is then used to insert high-frequency MR information into the PET image. No a priori brain tissue segmentation is performed and no modelling steps are required to remove MRI components that are not related to myelin. The MM-RR algorithm can be summarized as follows:

Given a T1-weighted MR and [¹¹C]PIB PET binding map (co-registered to the same space), the two individual images undergo a 3D wavelet decomposition. A model is fitted to the scatter plot defining the relationship between MR and PET signal intensities at the lowest resolution scales (low-resolution modelling) (Figure 2(a)). This relationship relies on the assumption of association between myelin-bound water (as measured by T1-weighted MR) and basic myelin protein (as measured by [¹¹C]PIB PET) and can be modelled by a linear regression through the origin with positive slope. Note that this relationship is proposed only for T1 relaxation and should be re-evaluated for use with other MRI modalities.

The linear model is calculated in the same fashion between the MR and PET coefficients for all the remaining resolution levels of the 3D wavelet transforms (high-resolution modelling). In this instance, given the loss of resolution of the PET data, the linear slope differs from the one at step A. (Figure 2(b));

For each resolution level of the 3D PET wavelet transform, new wavelet coefficients are calculated, rescaled from the difference between the models fitted in A) and B) and corrected to account for MR-PET differences in term of information content and image noise (Figure 2(c));

The new PET wavelet coefficients are mapped back to the image space through inverse wavelet transform, generating the high-resolution [¹¹C]PIB PET image (Figure 2(d)).

Figure 2.

Graphical representation of the MM-RR algorithm. (a) wavelet decomposition and model definition of the relationship between MR and PET signals at low resolution scales (Low-resolution modelling); (b) model definition of the relationship between MR and PET high-resolution coefficients (High-resolution modelling); (c) definition of new wavelet coefficients from the models defined in (a) and (b) after appropriate scaling and de-noising; (d) inverse wavelet transform of the final wavelet coefficients originates the improved PET image. To note that the wavelet coefficients in XY plots are unitless.

The detailed mathematical derivation of the method is reported in the following paragraph.

The MM-RR pipeline

Wavelet decomposition MM-RR algorithm

The wavelet transform (WT) decomposes a 1D signal $f (x)$ into its high- and low- frequency components though a filter bank as follows

f (x) = \sum_{k} \sum_{j - 1}^{J} d_{j} (k) \cdot ψ_{j, k} (x) + \sum_{k} C_{j} (k) \cdot φ_{J, k} (x)

(1)

where j is usually referred to as the decomposition/resolution level, ψ is a high band-pass function and φ is a low-pass scaling function, while

d_{j} (k)

and

C_{j} (k)

are their resulting coefficients (respectively wavelet and low resolution coefficients).

In tomographic imaging applications, the input signal is no longer mono dimensional, and therefore its components have to be represented by three-dimensional (3D) wavelet and scaling functions. In this case, the transform generates eight components from a 3D input signal while applying the one-dimensional decomposition (high- and low- pass filters) successively along the three $x, y, z$ directions. The low-resolution coefficients resulting from the application of the low-pass filter to all $x, y, z$ directions represent the input signal for the second level wavelet decomposition $j = 2$ . A graphical representation of the 3D wavelet decomposition is reported in Supplementary Material – Figure 1.

In this work, as in previously proposed PVC techniques,^13,19 we preferred the dual-tree complex wavelet transform (DT-CWT)²² to the more widespread used dyadic wavelet transform (DWT)²³ in order to achieve better directionality selectivity, shift invariance and perfect reconstruction. The wavelet decomposition has been carried out using the 3D Complex Dual-Tree Discrete Wavelet Transform function from the Matlab Wavelet package developed by the Polytechnic University, Brooklyn, NY (http://taco.poly.edu/WaveletSoftware/). The setting parameters comprise a maximum decomposition scale level of 2 and the selection of Farras and Kingsbury Q class filters.

Scaling factors

Initially, both functional (PET) and anatomical (MRI) images are decomposed by means of the DT-CWT into several resolution elements up to resolution level $j = 2$ and the wavelet ( $W^{PET}$ , $W^{MRI}$ ) and low resolution ( $L^{PET}$ , $L^{MRI}$ ) coefficients collected into separate matrixes.

The algorithm requires two scaling factors to account for the resolution difference between the two imaging modalities (inter-modality resolution and global coefficient scaling factors) and a weighting factor accounting for the difference between subsequent wavelet resolution levels (intra-modality resolution scaling factor).

The inter-modality resolution scaling factor R compensates for the difference in resolution between the two imaging modalities. It is important to account for this difference since the wavelet transform decomposes the image into several resolution levels. In case the two images being transformed do not have the same initial resolution, then the correspondent decomposition levels will be out of phase compromising the rest of the algorithm. The magnitude of wavelet coefficients depends on the spatial resolution of the image, therefore the scaling coefficient is obtained as the ratio of the coefficients of the original anatomical image with a smoother version which is degraded to PET scanner resolution through a 3D Gaussian filtering in image space.

R^{L} = \frac{\sum L^{MRI}}{\sum L^{MRIs}} R^{W} = \frac{\sum W^{MRI}}{\sum W^{MRIs}}

(2)

The global coefficient scaling factor G compensates for the difference in intensity between PET and MRI coefficients. Note that in equation (3), the PET coefficients are multiplied with the inter-modality resolution scaling factor R for consistency purpose.

G^{L} = \frac{R^{L} \cdot \sum L^{PET}}{\sum L^{MRI}} G^{W} = \frac{R^{W} \cdot \sum W^{PET}}{\sum W^{MRI}}

(3)

The PET and MRI wavelets and low-resolution coefficients are then respectively multiplied within each iterative step by the inter-modality resolution and the global coefficient scaling factors that remain unvaried during the whole iterative procedure outlined in the next section.

MM-RR algorithm

The first step of the algorithm consists of creating a correlation model between the functional and structural information taking into account only the low-resolution coefficients (Figure 2(a)). This step is also called low-resolution modelling and corresponds to

L^{model} = α \cdot G^{L} \cdot L^{MRI} + β

(4)

{residuals}_{L} = G^{L} \cdot L^{MRI} - L^{model}

(5)

where α and β are the slope and intercept of MR and PET regression at low-resolution scales (

L^{model}

), stored together with the

{residuals}_{L}

before moving on with the analysis of the wavelet coefficients.

In the following iterative steps, each set of wavelet coefficients related to a specific direction and resolution level is analysed separately. For each set, a correlation analysis is performed between the functional and structural high-resolution wavelet coefficients $(W^{MRI}, W^{PET})$ . The coefficients that fall in the second and fourth quadrant of the Cartesian graph (Figure 2(b)) are discarded before calculating the correlation model. These coefficients with negative correlation refers to image component that are present in only one image modality thus resulting from image noise.

W^{model} = α \cdot G^{W} \cdot W^{MRI} + β

(6)

Once the linear model is computed, the residuals are calculated and stored (Figure 2(b)).

{residuals}_{W} = G^{W} \cdot W^{MRI} - W^{model}

(7)

The residuals are then added to the correlation model $(L^{model})$ obtained from the low-resolution coefficients at the first step to obtain the final coefficients of the improved PET image (Figure 2(c)). The initial model containing information on low-resolution uniform activity distribution is integrated with its missing high-resolution information for the specific set of coefficients described by the residuals

W^{new_PET} = L^{model} + {residuals}_{W}

(8)

To note that in equation (8), the correlation model obtained with low-resolution coefficients of equation (4) is used as foundation in the computation of the coefficients for the new high resolution PET. Since the magnitude of the coefficients increases as the decomposition proceeds to finer resolution levels, an intra-modality resolution scaling factor $Λ$ accounting for this difference is needed.

We used a robust measure of the noise variance computed for each resolution level ( $σ_{L}$ and $σ_{W}$ ) to account for this scaling factor. Precisely, we computed the median absolute deviation MAD which is then divided by 0.6745 for calibration purposes.^24,25

σ_{L} = MAD {{residual}_{L}} / 0.6745

(9)

σ_{W} = MAD {{residual}_{W}} / 0.6745

(10)

Once the measure of the noise variance are known, the intra-modality resolution scaling factor $Λ$ is computed as their ratio between different resolution levels as in equation (11) and integrated in the final model of equation (8) to obtain equation (12)

Λ = σ_{W} / σ_{L}

(11)

W^{new_PET} = L^{model} + Λ \cdot {residuals}_{W}

(12)

Once the new wavelet coefficients $(W^{new_PET})$ are calculated iteratively for all the wavelet quadrants, an inverse wavelet transform is performed to obtain the new high resolution PET image (Figure 2(d)).

Phantom experiments

Given the impossibility to generate simulated PET/MR myelin imaging data accounting for the complex biology of brain tissues, the evaluation of the MM-RR method was done using the PET/CT phantom data generated in Grecchi et al.¹⁹ where there is a similar positive association between the PET fluorine signal and the CT measured density.

Briefly, the NEMA IEC body phantom with six spheres of different volumes (range from 0.5 cm³ to 26.5 cm³) was considered. Compartments were filled with both iodinated contrast media (CM) Omnipaque300TM (300 mg/ml organic Iodine) and radioactive tracer 18F-Fluoride. Two different experiments changing the layout of CT and PET contrasts were performed, accounting for possible mismatches between functional and anatomical images (Figure 3(a)). Full experimental details are reported in the original reference.¹⁹

Figure 3.

Phantom experiment. (a) [¹⁸F]Fluoride PET/CT transaxial images of three different phantom experiment acquisitions (one for each line). Alongside the CT images (1st column) are three different type of functional images: scanner PET images (2nd column), images resulting from the inclusion of a PSF model into the reconstruction (3rd column), images resulting after the application of the SFS resolution recovery algorithm (4th column), and images resulting after the application of the MMRR resolution recovery algorithm (5th column). (b–d) Root mean square error and noise analysis. For each sphere (S1–S6), four values corresponding to images obtained with different modalities are reported: Scanner PET (black circle), PET with PSF reconstruction (grey triangle), PET corrected with SFS-RR algorithm (white diamond), and PET corrected with MM-RR algorithm (black cross). (a) Root mean square error for the six spheres obtained as an average between the two phantom experiments. (b–c) Contrast to noise ratio computed for each sphere against a uniform region in the phantom background. Representative phantom figures are reported next to the graphs for reference purposes.

The resulting images were then processed with the MM-RR algorithm, with the CT image used instead of the MRI as the high-resolution low-specificity reference image for the algorithm. In addition, the same data were reconstructed with the inclusion of a Point Spread Function (PSF) model into the standard OSEM algorithm (GE Q.Core VuePoint FX-S, henceforth called PET-PSF). For comparative purposes, PVC with the SFS-RR algorithm was also included.¹⁹

Consistently with previous work,¹⁹ the performance of PVC methods was expressed in term of root mean square error (RMSE) and contrast to noise ratio (CNR).

Dataset

The MM-RR technique was tested with two measured datasets from the same study comprising [¹¹C]PIB PET and structural T1-weighted MR images. For the first part of the study, eight healthy control subjects (three male, five female, age: 31.6 ± 6.3 years) were enrolled.²⁰ Inclusion criteria consisted in an age between 18 and 55 years and the absence of any known neurological or psychiatric condition.

The second set comprising 20 relapsing-remitting MS patients age- and gender-matched to healthy controls (7 male, 13 female, age: 32.3 ± 5.6 years) with at least one gadolinium enhancing (Gd+) lesion over 6 mm in diameter on magnetic resonance imaging (MRI) at study entry.⁶

For both studies, ethical approval was granted by the ethics committee of the Pite-Salpetriere Hospital (Approval No. P080503) and informed consent was obtained from all participants. This study was conducted according to the Declaration of Helsinki.

Positron emission tomography and magnetic resonance imaging

All PET scans were performed on the brain dedicated PET research tomograph ECAT-HRRT (CPS Innovations, Knoxville, TN, USA). This high resolution scanner achieves an intra-slice spatial resolution of ∼2.5 mm full width at half maximum, with 25 cm and 31.2 cm of axial and transaxial field of view.²⁶

The 90-min emission scan was initiated coincident with a 1-min intravenous bolus injection of [¹¹C]PIB (mean 358 ± 34 MBq). There were no differences between healthy controls and patients in term of injected radioactivity or specific activity.⁶ Images were reconstructed using the 3D ordinary Poisson ordered subset expectation maximization algorithm (POSEM)²⁷ implementing point spread function (PSF) modelling up to 10 iterations. An additional smoothing filter implementing the PSF was applied to the reconstructed images to remove image inhomogeneity and spiky artefacts.²⁸

The resulting dynamic PET images consisted of 25 frames of data (6 × 1, 6 × 2, 4 × 3, 6 × 5, 3 × 10 minutes) with a voxel size of 1.22 mm × 1.22 mm × 1.22 mm. Inter-frame subject motion correction was applied by realigning each PET frame to a common reference space through a procedure similar to those reported by Montgomery et al.²⁹ Data were corrected for carbon-11 decay. A reference region was determined using previously published methodology developed specifically for [¹¹C]PIB that uses a supervised clustering approach to determine a set of grey matter voxels with kinetics that are the closest to a predefined set of normative dynamics.³⁰ The average time-activity curve of these voxels was used as input for the Logan plot to derive parametric maps of volume of distribution (DVR).³¹ We refer to Veronese et al.²⁰ for a detailed description of the PET [¹¹C]PIB dynamic quantification.

In addition to the PET acquisition, MR images of all subjects were collected using a 3T Siemens system (Siemens, Erlangen, Germany; TRIO 32-channel TIM system). Specifically, 3D T1-weighted MPRAGE, T2-weighted (T2-w) and 3D FLAIR sequences were performed for all subjects, while patients had additional pre- and post- gadolinium T1 spin-echo sequences. Only the T1-weighted images were used for PVC using MM-RR. Lesions were segmented on T2-w images with reference to co-registered FLAIR scans.

For each subject, T2- and T1- weighted images were co-registered to the PET space using a rigid transformation using Flirt (http://fsl.fmrib.ox.ac.uk/fsl/) maintaining PET data in its original space and avoiding further loss of resolution due to interpolation. Using VBM8 toolbox (hhttp://dbm.neuro.uni-jena.de/vbm/download/) and FIRST (http://fsl.fmrib.ox.ac.uk/fsl/), image voxels were automatically classified as grey matter (GM), white matter (WM) and cerebrospinal fluid (CSF) using a tissue probability threshold of 90%. The quality of all the image segmentations and coregistrations was assessed via visual inspection by an experienced neuroimaging scientist.

Additional regions of interest (ROIs) were manually defined on patient images by an experienced neurologist using the co-registered information from the T2-weighted and T1-weighted images. The ROIs were classified as follows: (1) GM; (2) normal-appearing WM; (3) perilesional WM; (4) lesion WM; (5) black holes and (6) Gadolinium enhancing (Gd+) lesions. Only lesions above 2.5 mm (corresponding to the resolution of the PET-camera) were retained as effective ROI.

Image analysis and evaluation

Qualitative and quantitative evaluations were carried out comparing the original data, here labelled as “Scanner PET”, images with the output of the resolution recovery procedure, here named “MM-RR PET”, in order to evaluate the performance of the MM-RR algorithm.

In the first instance, we analysed the images of the healthy subjects to evaluate resolution recovery and to control for any bias that may have been introduced. For each subject, we calculated the average [¹¹C]PIB DVR in three ROIs (GM, WM and whole brain respectively) and the mean and standard deviation of all subjects was compared between Scanner and MM-RR PET. As an additional check, we calculated the [¹¹C]PIB DVR as a function of the distance from the cerebral spinal fluid for both WM and GM ROIs as a loss of signal near the ventricles was evident in the original data as a result of PVE.²⁰ This is very relevant for the quantification of all lesions close to these boundaries. Finally, for a quantitative evaluation of resolution recovery, a contrast analysis was performed taking into consideration the GM/WM contrast in small spherical regions of 2.5, 5, 7.5 and 10 mm diameter. For each subject, a pair of random regions was drawn on the GM and WM segmented on the T1-weighted images. Both ROIs were then moved to the PET images (Scanner and MM-RR) and the mean of voxel estimates within each region was computed ( ${WM}_{ROI}$ , ${GM}_{ROI}$ ) and the mean local contrast calculated as

The procedure was repeated up to 500 times for each sphere diameter size and mean and variability of mcl % were compared between different ROI sizes as well as between Scanner and MM-RR PET.

For the patient group, we compared the average [¹¹C]PIB DVR with the same statistical methodology used for the healthy subjects group focusing, however, on clinically relevant ROIs, specifically the normal-appearing WM, perilesional area, lesions, black holes and Gadolinium enhancing (Gd+) lesions. These ROIs were chosen consistently with previous studies.⁶ Additionally, we performed a correlation analysis between the average [¹¹C]PIB DVR and the lesional volumes to check whether the effect of the PVC depended on the size of the lesion, as one would have expected. To perform this type of analysis, we grouped all the lesions of all the patients according to their volumes. Then, for each lesions size, we calculated the average [¹¹C]PIB DVR of all patients’ lesions of that specific volume. The analysis was limited to WM lesions, to control for the effect of tissue surrounding contrast on PVC.

MM-RR PET myelin mapping: correlation with myelin associated protein mRNA expressions

In order to evaluate the validity of the myelin maps produced by MM-RR method, [¹¹C]PIB PET maps derived from the healthy control group (both MM-RR and Scanner PET) were tested against the brain mRNA expression atlas of a set of myelin-associated proteins contained in the Allen Human Brain Atlas.³² We have previously shown that these mRNA maps are highly predictive of protein levels in-vivo measured with PET for all those transcripts that do not undergo significant post-transcriptional modifications,³³ a condition that is met for the myelin system.²⁰ The comparison between imaging and gene expressions was performed by correlation analysis at the voxel level in stereotaxic MNI (Montreal Neurological Institute) coordinate accordingly to MENGA pipeline³⁴ (https://www.nitrc.org/projects/menga/). The proteins considered for the analysis were selected based on their prevalence in the myelin sheath³⁵: myelin-associated glycoprotein (MAG), myelin basic protein (MBP), proteolipid protein-1 (PLP1), 2′,3′-Cyclic-nucleotide 3′-phosphodiesterase (CNP), myelin oligodendrocyte glycoprotein (MOG), and oligodendrocyte basic protein (MOBP). To exclude for spurious correlations, [¹¹C]PIB imaging data were also compared with a series of non-myelin proteins including dopamine D2 receptor (DRD2), 5HT_1A serotonin receptor (HTR1A), brain-derived neurotrophic factor (BDNF), and Aquaporin 4 (AQP4).

Results

Phantom data

The performance of the PVC methods is summarized in Figure 3.

As expected, the smaller the sphere, the bigger the bias in the activity estimation, regardless of the PVC method used.

Given the lack of segmented anatomical information, as expected the performance of MM-RR is inferior to the one of SFS-RR (RMSE relative difference +10% ± 11%; CNR relative difference −9% ± 13%) but comparable with the PSF reconstruction method (RMSE relative difference 0% ± 4%; CNR relative difference +6% ± 7%) across spheres. Compared to scanner PET, the images refer to a representative healthy control. The MM-RR PET images show sharper contours compared to the Scanner PET and better WM definition.

The qualitative improvement of the MM-RR [¹¹C]PIB DVR maps is confirmed by the quantitative analysis. Figure 4(b) reports the comparison of the three tissue ROIs (whole brain, GM and WM) in the healthy control dataset. Each box collects the mean value for the specific ROI of all the healthy subjects. The average DVR in the WM is 5.7% higher for the MM-RR PET compared to the Scanner PET (paired t-test p < 0.001). The average DVR in the GM is comparable between Scanner PET and MM-RR PET (paired t-test p = 0.72) confirming the lack of quantitative distortion in the latter. Moreover, the resolution recovery does not come at the cost of increased noise as the standard deviations of the two imaging modalities are comparable for all the three ROIs (F-test p > 0.35 for all ROIs).

Figure 4.

Scanner and MM-RR PET in a representative healthy control. (a) Images: Left column MRI images; central column Scanner PET; right column MM-RR PET. PET colorbars indicate DVR estimates (unitless). (b) [¹¹C]PIB DVR values: The box-plot diagrams show the median (middle line) and range of DVR values for the healthy control group in three different ROIs (from left to right: whole brain, grey matter, white matter). The table reports the numerical value of mean and standard deviation for each box. **indicates p-value<0.001

A quantitative evaluation of the improved image quality achieved with the MM-RR algorithm is shown by the WM/GM contrast analysis (Supplementary Material – Figure 2). The MM-RR PET images show an average increase of 30% in WM/GM contrast compared to Scanner PET, corresponding to an average [¹¹C]PIB DVR WM/GM ratio greater than 15%. These results are independent from ROI size; please note that the highest signal improvement is reported for the smallest ROI (diameter 2.5 mm) with a 46% increase in WM/GM contrast compared to Scanner PET.

A quantitative measure of the PVC that can be achieved with MM-RR is shown in Figure 5, where the DVR for WM and GM is displayed as a function of the distance from the CSF. In GM, [¹¹C]PIB DVR estimates are comparable within the first 3 mm (mean relative difference of 3% ± 2%), after which MM-RR PET returns lower values compared to Scanner PET (mean relative difference −9% ± 1% at 5 mm). To note that the further from the CSF, the lower the number of voxels, especially in the GM, which introduces an additional element of variability. For the WM, instead, the highest discrepancy between Scanner and MM-RR [¹¹C]PIB DVR estimates is found within 5 mm of distance from CSF (up to 31% in mean relative difference), while the two images become comparable as the distance increases (mean relative difference +2% ± 1%). This suggests that the recovered resolution concurs in correcting for spill-out effect in the WM border area. This corresponds to both WM/GM and WM/CSF boundaries. This is of relevance given that neurodegeneration around the ventricles is an early process in the pathophysiology of MS.

Figure 5.

The graph reports the [¹¹C]PIB DVR values measured for white matter (white markers) and grey matter (grey markers) as function of the distance from the cerebral spinal fluid (CSF). Lines with circle markers refer to the Scanner PET and the lines with triangle markers refer to the MM-RR PET. Light grey markers on grey matter lines refer to significant lower fraction of grey matter tissue compared to the sample closer to the CSF.

Correlation with myelin associated protein mRNA expressions

MM-RR did not affect the level of cross-correlation between imaging and genomic data (Figure 6(a)). High-resolution maps returned the same levels of correlation with the mRNA expressions of myelin-associated proteins of standard resolution maps (paired t-test p < 0.05; mean relative difference: −1% ± 4%) (Figure 6(b)). Genomic vs. imaging cross-correlation was significant for all tested cases (p < 0.05), ranging from R²= 0.48 to R²= 0.67. On the opposite, for both resolution levels, imaging vs. genomic cross-correlation was not significant when non-myelin-associated proteins were considered (R²< 0.26 for all subjects and all genes) (Figure 6(b)).

Figure 6.

Comparison between [¹¹C]PET myelin imaging and genomics. Panel A compares imaging vs. genomic cross-correlation between Scanner and MM-RR PET images. All subjects and proteins are reported. Blue circles refer to myelin-associated proteins (i.e. MAG, MBP, PLP1, CNP, MOG and MOBP) while red circles to non-myelin-associated proteins (i.e. DRD2, HTR1A, BDNF and AQP4). Panel B compares mean±SD imaging vs. genomic cross-correlation between subjects for MM-RR and Scanner PET resolution, respectively. Blue bars refer to myelin-associated proteins. Red bars refer to non-myelin-associated proteins. Bright bars refer MM_RR resolution PET. Light bars refer to Scanner resolution PET. All the correlations above significance threshold line result statistically significant (p value<0.05).

Application to [¹¹C]PIB PET myelin imaging in patients with MS

Figure 7 reports the impact of the MM-RR application to the MS patient dataset. Figure 7(a) provides an example of the improved lesion detectability achieved with the MM-RR algorithm compared to the Scanner PET. The images refer to a representative MS patient. As for the healthy control in Figure 4(a), the MM-RR PET images show sharper contours compared to the Scanner PET. Additionally, a noticeable improvement in lesion detectability can be appreciated in the two representative lesions (zoomed in the red circle). Additional specific example for Gadolinium enhanced lesion is reported in Figure 6 – Supplementary material.

Figure 7.

Axial and sagittal views from a representative MS patient. (a) Images: From left to right: MRI images; Scanner PET; MM-RR PET; zooms of the area delimited in the MRI image with a red circle. PET colorbars indicate DVR estimates (unitless). (b) [¹¹C]PIB DVR values: The box-plot diagrams show the median (middle line) and range of DVR for the healthy control group in three different ROIs (from left to right: grey matter GM, normal-appearing white matter NAWM, Perilesions, T2-w lesions, Black Holes, Gad enhancing lesions). The table reports the numerical value of mean and standard deviation for each box. **indicates p-value<0.001. *indicates p-value<0.05

The quantitative impact of the MM-RR application on patient images is reported in Figure 7(b). As previously done for healthy controls (Figure 4(b)), each box collects the mean value for a specific clinically relevant ROI for all the patients. The average [¹¹C]PIB DVR in the lesional tissue (T2-w lesions, black holes and Gadolinium enhancing lesions) is significantly lower for the MM-RR PET (mean difference 4–12%, p < 0.05). On the other hand, there is a consistent increase of the DVR in the perilesional area, normal-appearing WM and GM (paired t-test p < 0.001) up to 5%. Compared to Scanner PET, the average normal-appearing WM/lesion contrast improves in high-resolution PET: +34% for T2-w lesions, +39% for black holes and +101% for Gadolinium enhancing lesions. As for the analysis of the healthy control dataset, the standard deviation of the two imaging modalities is comparable for all the ROIs (F-test p > 0.87 for all ROIs).

These results are then consistent with the resolution analysis carried out on the collected lesions of varying size (Supplementary Material – Figure 3). Lesions on the MM-RR PET tend to have lower [¹¹C]PIB DVR than the Scanner PET ones but this difference is more evident on small lesions, while for big lesions, the two modalities produce very similar DVRs.

Application to [¹⁸F]FDOPA PET imaging

The robustness to anatomical and functional mismatch of PVC wavelet-based approached has been previously demonstrated.¹⁹ However, an additional test was performed to ensure that no bias or artefacts are introduced when the structural and functional signals are not related. We applied the MM-RR algorithm to a sample [¹⁸F]FDOPA PET/MR scan of a healthy subject (www.nitrc.org/projects/spmtemplates/) where the anatomical information had no relation with the PET signal located in the basal ganglia. As predicted, no artefacts are introduced after the application of MM-RR (Supplementary Material – Figure 4).

Discussion

MM-RR has been developed to give a mathematical account of the better specificity to myelin of PET at low resolution as well as the good precision of MRI at finer resolution scales. The application to clinical data showed that MM-RR led to improved WM/GM contrast compared to Scanner PET imaging. In MS patients, MM-RR allowed a better signal characterization of [¹¹C]PIB PET in lesions. One main advantage of MM-RR methods over more conventional approaches is the straightforward applicability without prior knowledge or pre-processing. Standard PVC techniques often require ad hoc measurements of the scanner PSF or tissue classification and/or segmentation, which cannot be easily performed on such heterogeneous lesion distribution.

Our results on the healthy control dataset showed a qualitative significant improvement especially in the delineation of WM structures and an average 30% increase in contrast between white and grey matter at the expenses of no quantification bias or noise increase. The signal recovery and good noise control assessed in the healthy controls were then reflected into patient images that demonstrated lesions with sharper and better-defined contours and the resulting accurate quantification of [¹¹C]PIB DVR in different lesion types and perilesional areas.

Additional validations of the method included application to phantom data and comparison of MM-RR myelin maps with genomic datasets. In the first case, MM-RR returned comparable PVC results to PSF reconstruction method but worse performances in term of bias and image contrast than the atlas-based approach implemented in SFS which is expected given the stronger and more informed prior implemented by SFS compared to MM-RR. Please note that in this work the resolution recovery algorithm was applied to clinical data that were reconstructed with PSF inclusion and POESM settings that were optimized according to our previous work.⁵ The application of resolution recovery methods to parametric maps of lower quality may degrade the final quality of the images.

The reliability of the MM-RR method was confirmed by the genomic analysis as MM-RR returned high-resolution [¹¹C]PIB maps without altering the degree of correlation with the myelin-associated protein gene expressions.

The effectiveness of the MM-RR relies on the validity of some important assumptions. First of all the presence of a shared myelin content between [¹¹C]PIB PET and T1w MR imaging. This relationship allows the identification of the low-resolution model in wavelet space, which is at the basis of the fusion between the two modalities. This assumption limits the MM-RR extendibility to those cases when the information content of both anatomical and functional images is mismatched. At the same time, MM-RR requires the accurate spatial coregistration of both modalities. In case of misalignment, MM-RR method can lead to image artefacts (Supplementary Material – Figure 5). This result was quite expected as all PVC techniques which rely on anatomical information are naturally sensitive to mis-registration errors. Accurate image co-registration should be achieved prior PVC application. The choice of PVC technique should depend on application and image co-registration reliability.³⁶

In order to assess the robustness of the algorithm towards possible artefacts due to mismatch between anatomical and functional images, we additionally tested the algorithm using [¹⁸F]FDOPA PET and T1 MR images from healthy subjects. Since T1 sequences are not sensitive to dopamine, one would not expect significant changes in the final PET images as the analysis demonstrated. It follows that MM-RR method can only be applied when functional and structural data share an association to the biological target. In this investigation, the structural information was provided by T1-w images; however, this MRI sequence is not the most specific to myelin content. Future work will consider the use of more myelin-sensitive MR imaging sequences such as the magnetic transfer imaging (MTR)^37,38 or the multi-component driven equilibrium single pulse observation of T1 (mc-DESPOT).^39,40

Conclusion

We have developed a robust multimodal methodology for the quantitative resolution recovery of brain [¹¹C]PIB PET data, specifically designed for myelin imaging, and tested on a set of clinical healthy controls and MS patients data. The technique produces images with significantly improved quantitative properties and visual quality and it is of straightforward and rapid application. The work presented here is highly relevant for the newly introduced PET-MRI scanners as it presents one of the first synergistic approaches combining PET and MRI data for the same target.

Footnotes

Funding

The author(s) declared the following potential conflicts of interest with respect to the research, authorship, and/or publication of this article: EG was supported by an EPSRC scholarship (EP/K502868/1); FET and MV were supported by the MRC PET Methodology Program grant (G1100809/1) and by the Wellcome Trust “Strategic Award: Inflammation in AD and MDD”. BB was supported by the ARSEP postdoctoral fellowship.

Acknowledgments

The author would like to thank Emilie Poirion, PhD student in Catherine Lubetzki and Bruno Stankoff’s team for her support.

Declaration of conflicting interests

The author(s) declared the following potential conflicts of interest with respect to the research, authorship, and/or publication of this article: Bruno Stankoff reports receiving consulting and lecture fees from Biogen Idec, Novartis, Merck Serono, Sanofi-Aventis, Teva-Pharma and research support from Genzyme. Other authors have no conflict of interest to declare. Federico Turkheimer reports receiving lecture fees from Eli-Lilly.

Authors’ contributions

EG, MV, BB, BS and FET made substantial contribution to the study design. EG, BB and BS made substantial contribution to the data collection. EG, MV, DGL, MB and FET made substantial contribution to the data analysis. EG, MV, DGL, MB and FET made substantial contribution to the data interpretation. EG, MV and FET drafted the manuscript. EG, MV, BB, DGL, MB, BS and FET revised critically the article. All authors approved the last version.

Supplementary material

Supplementary material for this paper can be found at the journal website:

References

Ramagopalan

Dobson

Meier

et al.

Multiple sclerosis: risk factors, prodromes, and potential causal pathways. Lancet Neurol 2010; 9: 727–739.

Franklin

RJM

Gallo

. The translational biology of remyelination: Past, present, and future. Glia 2014; 62: 1905–1915.

Klunk

Engler

Nordberg

et al.

Imaging brain amyloid in Alzheimer's disease with Pittsburgh Compound-B. Ann Neurol 2004; 55: 306–319.

Stankoff

Freeman

Aigrot

et al.

Imaging central nervous system myelin by positron emission tomography in multiple sclerosis using [methyl-11C]-2-(4′-methylaminophenyl)-6-hydroxybenzothiazole. Ann Neurol 2011; 69: 673–680.

Veronese

Bodini

Garcia-Lorenzo

et al.

Quantification of [(11)C]PIB PET for imaging myelin in the human brain: a test-retest reproducibility study in high-resolution research tomography. J Cereb Blood Flow Metab 2015; 35: 1771–1782.

Bodini

Veronese

García-Lorenzo

et al.

Dynamic imaging of individual remyelination profiles in multiple sclerosis. Ann Neurol 2016; 79: 726–738.

Rousset

Rahmim

Alavi

et al.

Partial volume correction strategies in PET. PET Clin 2007; 2: 235–249.

Aston

JAD

Cunningham

Asselin

M-C

et al.

Positron emission tomography partial volume correction&colon; estimation and algorithms. J Cereb Blood Flow Metab 2002; 22: 1019–1034.

Metz

Weigand

Popescu

et al.

Pathologic heterogeneity persists in early active multiple sclerosis lesions. Ann Neurol 2014; 75: 728–738.

10.

Frischer

Weigand

Guo

et al.

Clinical and pathological insights into the dynamic nature of the white matter multiple sclerosis plaque. Ann Neurol 2015; 78: 710–721.

11.

Meltzer

Zubieta

Links

Brakeman

Stumpf

Frost

. MR-based correction of brain PET measurements for heterogeneous gray matter radioactivity distribution. J Cereb Blood Flow Metab 1996; 16: 650–658.

12.

Boussion

Hatt

Lamare

Le Rest

Visvikis

. Contrast enhancement in emission tomography by way of synergistic PET/CT image combination. Comput Meth Prog Biomed 2008; 90: 191–201.

13.

Shidahara

Tsoumpas

Hammers

et al.

Functional and structural synergy for resolution recovery and partial volume correction in brain PET. Neuroimage 2009; 44: 340–348.

14.

Leahy

Cherry

et al.

High-resolution 3D Bayesian image reconstruction using the microPET small-animal scanner. Phys Med Biol 1998; 43: 1001.

15.

Rapisarda

Bettinardi

Thielemans

et al.

Image-based point spread function implementation in a fully 3D OSEM reconstruction algorithm for PET. Phys Med Biol 2010; 55: 4131.

16.

McGinnity

Shidahara

Feldmann

et al.

Quantification of opioid receptor availability following spontaneous epileptic seizures: Correction of [11 C] diprenorphine PET data for the partial-volume effect. Neuroimage 2013; 79: 72–80.

17.

Kim

Shidahara

Tsoumpas

et al.

Partial volume correction using structural–functional synergistic resolution recovery: comparison with geometric transfer matrix method. J Cereb Blood Flow Metab 2013; 33: 914–920.

18.

Shidahara

Tsoumpas

McGinnity

et al.

Wavelet-based resolution recovery using an anatomical prior provides quantitative recovery for human population phantom PET [11C] raclopride data. Phys Med Biol 2012; 57: 3107.

19.

Grecchi

O'Doherty

Veronese

et al.

Multimodal partial volume correction–application to 18F-fluoride PET/CT bone metastases studies. J Nucl Med 2015; 56: 1408–1414.

20.

Veronese

Bodini

García-Lorenzo

et al.

Quantification of [ 11C] PIB PET for imaging myelin in the human brain: a test–retest reproducibility study in high-resolution research tomography. J Cereb Blood Flow Metab 2015; 35: 1771–1782.

21.

Turkheimer

Aston

JAD

Banati

Riddell

Cunningham

. A linear wavelet filter for parametric imaging with dynamic PET. IEEE Trans Med Imag 2003; 22: 289–301.

22.

The dual-tree complex wavelet transform: a new technique for shift invariance and directional filters. In: IEEE digital signal processing workshop, Citeseer, August 1998.

23.

Mallat

. A theory for multiresolution signal decomposition: the wavelet representation. IEEE Trans Pattern Anal Machine Intell 1989; 11: 674–693.

24.

Donoho

Johnstone

Kerkyacharian

et al.

Wavelet shrinkage: asymptopia?

J R Stat Soc Ser B 1995; 57: 301–369.

25.

Turkheimer

Brett

Visvikis

et al.

Multiresolution analysis of emission tomography images in the wavelet domain. J Cereb Blood Flow Metab 1999; 19: 1189–1208.

26.

Wienhard

Schmand

Casey

et al.

The ECAT HRRT: performance and first clinical application of the new high resolution research tomograph. IEEE Trans Nucl Sci 2002; 49: 104–110.

27.

Carson

Barker

Jeih-San

et al.

Design of a motion-compensation OSEM list-mode algorithm for resolution-recovery reconstruction for the HRRT. Nucl Sci Symp Conf Record 2003; 5: 3281–3285.

28.

Stute

Comtat

. Practical considerations for image-based PSF and blobs reconstruction in PET. Phys Med Biol 2013; 58: 3849.

29.

Montgomery

Thielemans

Mehta

et al.

Correction of head movement on PET studies: comparison of methods. J Nucl Med 2006; 47: 1936–1944.

30.

Ikoma

Edison

Ramlackhansingh

et al.

Reference region automatic extraction in dynamic [ 11C] PIB. J Cereb Blood Flow Metab 2013; 33: 1725–1731.

31.

Price

Klunk

Lopresti

et al.

Kinetic modeling of amyloid binding in humans using PET imaging and Pittsburgh Compound-B. J Cereb Blood Flow Metab 2005; 25: 1528–1547.

32.

Hawrylycz

Lein

Guillozet-Bongaarts

et al.

An anatomically comprehensive atlas of the adult human brain transcriptome. Nature 2012; 489: 391–399.

33.

Rizzo

Veronese

Heckemann

et al.

The predictive power of brain mRNA mappings for in vivo protein density: a positron emission tomography correlation study. J Cereb Blood Flow Metab 2014; 34: 827–835.

34.

Rizzo

Veronese

Expert

et al.

MENGA: a new comprehensive tool for the integration of neuroimaging data and the allen human brain transcriptome atlas. PloS One 2016; 11: e0148744.

35.

Jahn

Tenzer

Werner

. Myelin proteomics: molecular anatomy of an insulating sheath. Mol Neurobiol 2009; 40: 55–72.

36.

Hutton

Thomas

Erlandsson

et al.

What approach to brain partial volume correction is best for PET/MRI?

Nucl Inst Meth Phys Res Sec A 2013; 702: 29–33.

37.

Dousset

Grossman

Ramer

et al.

Experimental allergic encephalomyelitis and multiple sclerosis: lesion characterization with magnetization transfer imaging. Radiology 1992; 182: 483–491.

38.

Van Waesberghe

Kamphorst

De Groot

CJA

et al.

Axonal loss in multiple sclerosis lesions: magnetic resonance imaging insights into substrates of disability. Ann Neurol 1999; 46: 747–754.

39.

Kitzler

Zeineh

et al.

Deficient MWF mapping in multiple sclerosis using 3D whole-brain multi-component relaxation MRI. Neuroimage 2012; 59: 2670–2677.

40.

Kolind

Deoni

. Rapid three-dimensional multicomponent relaxation imaging of the cervical spinal cord. Magn Reson Med 2011; 65: 551–556.

41.

Beniac

Luckevich

Czarnota

et al.

Three-dimensional structure of myelin basic protein I. Reconstruction via angular reconstitution of randomly oriented single particles. J Biol Chem 1997; 272: 4261–4268.

42.

Ridsdale

Beniac

Tompkins

et al.

Three-dimensional structure of myelin basic protein II. Molecular modeling and considerations of predicted structures in multiple sclerosis. J Biol Chem 1997; 272: 4269–4275.

43.

Bajaj

LaPlante

Cotero

Fish

Bjerke

Siclovan

et al.

Identification of the protein target of myelin-binding ligands by immunohistochemistry and biochemical analyses. J Histochem Cytochem 2013; 61: 19–30.

44.

Biancalana

Koide

. Molecular mechanism of Thioflavin-T binding to amyloid fibrils. Biochim Biophys Acta 2010; 1804: 1405–1412.

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.00 MB

0.60 MB

Multimodal partial volume correction: Application to [ 11 C]PIB PET/MRI myelin imaging in multiple sclerosis

Abstract

Keywords

Introduction

Materials and methods

The MM-RR method

The MM-RR pipeline

Wavelet decomposition MM-RR algorithm

Scaling factors

MM-RR algorithm

Phantom experiments

Dataset

Positron emission tomography and magnetic resonance imaging

Image analysis and evaluation

MM-RR PET myelin mapping: correlation with myelin associated protein mRNA expressions

Results

Phantom data

Correlation with myelin associated protein mRNA expressions

Application to [11C]PIB PET myelin imaging in patients with MS

Application to [18F]FDOPA PET imaging

Discussion

Conclusion

Footnotes

Funding

Acknowledgments

Declaration of conflicting interests

Authors’ contributions

Supplementary material

References

Supplementary Material

Application to [¹¹C]PIB PET myelin imaging in patients with MS

Application to [¹⁸F]FDOPA PET imaging