Abstract
Head motion occurring during brain PET studies leads to image blurring and to bias in measured local quantities. The objective of this work was to implement a correction method for PET data acquired with the mMR synchronous PET/MR scanner. Methods: A list-mode–based motion-correction approach has been designed. The developed rebinner chronologically reads the recorded events from the Siemens list-mode file, applies the estimated geometric transformations, and frames the detected counts into sinograms. The rigid-body motion parameters were estimated from an initial dynamic reconstruction of the PET data. We then optimized the correction for 11C-Pittsburgh compound B (11C-PIB) scans using simulated and actual data with well-controlled motion. Results: An efficient list-mode–based motion correction approach has been implemented, fully optimized, and validated using simulated and actual PET data. The average spatial resolution loss induced by inaccuracies in motion parameter estimates and by the rebinning process was estimated to correspond to a 1-mm increase in full width at half maximum with motion parameters estimated directly from the PET data with a temporal frequency of 20 s. The results show that the rebinner can be safely applied to the 11C-PIB scans, allowing almost complete removal of motion-induced artifacts. The application of the correction method to a large cohort of 11C-PIB scans led to the following observations: first, that more than 21% of the scans were affected by motion greater than 10 mm (39% for subjects with Mini–Mental State Examination scores below 20), and second, that the correction led to quantitative changes in Alzheimer-specific cortical regions of up to 30%. Conclusion: The rebinner allows accurate motion correction at a cost of minimal resolution reduction. Application of the correction to a large cohort of 11C-PIB scans confirmed the necessity of systematically correcting for motion to obtain quantitative results.
The accuracy of PET measurements relies on the patient’s ability to stay still during the data acquisition. Head motion during brain PET studies leads to image blurring, making accurate localization of structures more difficult, as well as leading to quantitative bias in measured local quantities. In dynamic studies, patient motion alters the time–activity curves, measured at each voxel or in regions of interest, and introduces errors in the parameter estimates derived from kinetic modeling. Methods to correct for head motion can be classified into 2 broad categories depending on whether the correction occurs during the list-mode rebinning/reconstruction step or after the reconstruction (1). Most of the postreconstruction approaches rely on image registration algorithms to align each reconstructed time frame with a target, thus neglecting intraframe motion, which is an obvious limitation. In addition, and with the exception of a few implementations (2), these methods do not correct for spatial misalignment between emission data and attenuation data, leading to corrupted reconstructed volumes. Contrary to image registration–based correction, event-by-event correction accounts for intraframe motion and for attenuation–emission mismatches and allows the user to freely define the final framing of the reconstruction image. To the best of our knowledge, the idea to correct for event mispositioning during the rebinning step was first proposed by Menke et al. (3). However, it was only a few years later that Bülher et al. (4) clearly identified possible sources of artifacts that needed to be accounted for during the rebinning for accurate event-by-event correction. Most of the few papers published since then demonstrated the superiority of the event-by-event correction approach over image-based correction methods (5–7).
In this paper, we propose a novel implementation of a list-mode–based correction approach for PET data acquired with the mMR synchronous PET/MR scanner (Siemens Healthcare GmbH) (8). Using realistic Monte Carlo simulations and actual data with well-controlled motion, we assessed the performance of the approach and optimized the whole correction protocol, including derivation of the rigid motion parameters from the PET data only. Finally, we applied the motion correction to a large cohort of 11C-Pittsburgh compound B (PIB) scans from which novel information on the impact of motion was derived.
MATERIALS AND METHODS
The developed rebinner chronologically reads the recorded events from the Siemens list-mode file, applies the estimated geometric transformations, and frames the detected counts into static or dynamic sinograms that can be reconstructed with standard software. Motion correction parameters are stored in a text file, where each line provides the time t at which a motion occurred and the corresponding 6 rigid parameters. In this work, the rigid-body motion parameters were estimated from an initial dynamic reconstruction of the PET data. However, motion parameters derived from simultaneously acquired MR data or using an external tracking device can also be used. The correction method was implemented under Linux Ubuntu 16.04. The source codes are freely available on request. A docker version is also available and automates correction and reconstruction on a deported e7-tools reconstruction machine running Microsoft Windows.
Estimation of the Motion Parameters
To measure the motion parameters, the PET data are first reconstructed using short time-frames and without attenuation or scatter correction for increased registration accuracy (2). The reconstruction is performed using standard ordinary Poisson ordered-subsets expectation maximization (9) with 3 iterations and 21 subsets. By default, a zoom of 2 is used, leading to volume dimensions of 172 × 172 × 127 with a voxel size of 2.09 × 2.09 × 2.03 mm. The 6 rigid transformation parameters that optimize the cross-correlation similarity criterion between each data frame and a reference volume created from a selection of frames is computed. Volume registration and other data manipulation are performed using tools developed at the McConnell Brain Imaging Centre (10). For improved accuracy, the dynamic volume is smoothed by a 3-dimensional (3D) gaussian kernel before the registration. The minimal frame duration that still allows collection of enough counts for accurate motion estimation depends on the counting statistics of the scan, the level of spatial smoothing, and other characteristics such as the tracer pharmacokinetic. The relation between frame duration, smoothing strength, and motion parameter accuracy is investigated in this work.
Rebinner Implementation Details
The rebinner reads the events from the list-mode file and applies the corresponding motion correction parameters to the equation of the corresponding line of responses (LORs). The intersection of this line with the actual scanner ring allows identification of the correct crystal pair and sinogram index. We addressed the LOR discretization issue (4) by spatially oversampling the scanner crystals, allowing for each LOR extra subcrystal combinations to which the simple line transformation is then applied. This approach enables determination of the proportion of events detected in LOR i to be reassigned to LOR j. In the following sections, OS1 means no oversampling (simple LOR approach) whereas OSn (n > 1) means n × n subdivisions of each crystal (tangential × radial), leading for each LOR to n2 direct and n2 – (n mod 2) oblique subcrystal combinations. To account for normalization differences between original and correct LORs, detected counts in each LOR are first normalized by their respective normalization factors before the LOR reassignment. Supplemental Figure 1 illustrates this rebinning process involving crystal oversampling and normalization (supplemental materials are available at http://jnm.snmjournals.org). Event losses caused by the gaps are compensated for by being filled with the neighbor counts before the correction is performed. To account for events leaving the scanner field of view because of motion, after assigning the counts to the correct sinogram bins, the rebinner multiplies the counts by the ratio of the frame duration and the time during which the corresponding LOR was not falling outside the field of view because of motion. Contrary to some correction approaches (5), our implementation does not take advantage of the extra events that are detected because of motion and simply discards them.
Validation Framework
Rebinning Accuracy
We assessed the accuracy of the rebinner operating under different modes using simulated 18F-FDG brain PET data with well-controlled motion. Each scan was a 600-s list-mode acquisition generated with the PET-SORTEO Monte Carlo simulator modeling the geometry and physical characteristics of the Siemens mMR PET/MR system (11) and during which an acute basic rigid motion was applied at t = 200 s. Each motion consisted of either a translation or a rotation with respect to one of the axes (from −10 to +10 mm with a step of 1 mm for the translations and from −10 to 10° with a step of 1° for the rotations). A total of 120 simulated scans with motion were generated this way using the same activity levels, numeric emission phantom, and attenuation map. In addition, 2 reference scans were generated without motion using the same activity distribution and attenuation map, and hence, both reference scans differed only by the noise. Each list-mode scan was then rebinned into a single static 3D sinogram with a span of 1 and a maximum ring difference of 60 without and with motion correction using the true motion parameters. Different combinations of parameters related to the motion correction were tested: crystal oversampling with n varying from 1 (no oversampling) to 4, accounting or not for the difference in normalization factors and compensating or not for data losses. Finally, each sinogram was reconstructed into a 256 × 256 × 127 voxel volume (voxel size, 1.4 × 1.4 × 2.03 mm) with 3D ordinary Poisson ordered-subsets expectation maximization using 3 iterations and 21 subsets and with all corrections applied (9). The normalized root-mean-square error (nRMSE) computed between each reconstructed volume V and the reference volume R and expressed as a percentage of the root sum square of the reference volume [nRMSE(V,R) = 100 ⋅ RMSE(V,R)/root sum square(R)] was the metric used to assess the accuracy of the correction. To remove the differences between the corrected and reference scan that were caused by the noise process only, the nRMSE computed between the 2 reference volumes was subtracted from the nRMSE to obtain a final metric called nRMS0.
Accuracy of Motion Estimates
The accuracy and precision with which the motion parameters can be estimated from an initial dynamic reconstruction were assessed using an actual 1,800-s 11C-PIB list-mode scan selected from the 11C-PIB study described in the following section. We generated first the corresponding static sinogram (no compression, span = 1, prompt and delay events stored in 2 different matrices), thus removing the temporal dimension. This sinogram can be seen as the one obtained from a motionless brain, with the activity distribution being the time-averaged distribution of the original scan. From this static sinogram, motion-free 1,800-s list-mode data were created with varying counting statistics using a nonparametric bootstrap approach which consisted of sampling with replacement of the prompt and delay events from the original scan (12). Four motion-free list-mode files were generated this way, with counting statistics corresponding to 10%, 20%, 50%, and 100% (designated by LM10%, LM20%, LM50%, and LM100%) of the counts from the original scan. Each of these list-mode files was then rebinned into 90 frames of 20 s each (designated by LM10%_20s, LM20%_20s, LM50%_20s, and LM100%_20s). The full-statistics list-mode data (LM100%) were used to generate 3 additional dynamic scans of 30 frames of 60 s (LM100%_60s), 12 frames of 150 s (LM100%_150s), and 6 frames of 300 s (LM100%_300s). Table 1 reports for each dynamic scan the number of frames, the frame duration, the total number of prompts, and the number of net trues in the first and last frames. Finally, using the same 4 list-mode files, we generated 7 additional dynamic sinograms with the same counting statistics, except that rigid motion was included during the rebinning process using our rebinner program. The parameters and time at which the motion was included are reported in Table 2. We deliberately did not impose motion that would have caused significant data loss (Tz [translation along Z], Rx [rotation about X], and Ry [rotation about Y]). Those losses would have induced artifacts in the reconstructed images, impacting the registration process. Also, motion times were chosen such that the motion would always occur between 2 frames for all tested framing parameters. Motion parameter files obtained with a different 3D gaussian smoothing kernel (from 0 to 28 mm in full width at half maximum [FWHM]) applied to each frame before the registration process were then computed for each of the 7 motion-free and 7 motion-corrupted sinograms using our correction program. For each frame, the residual rigid motion parameters were computed as the difference between the estimated and true motion parameters (the latter being identity for the motion-free sinograms). The mean absolute Euclidean distance induced across time by the residual motion to a point located at a radial distance of 7 cm was computed and used as a metric to assess the accuracy of the estimated motion parameters. In addition, in order to estimate the loss of resolution caused by both the rebinning process and the inaccuracies in motion parameter estimates, we simulated the acquisitions of a point source, originally located in the center of the field of view and at a radial distance of 7 cm along the x-axis, and during which the same motion parameters as above were applied (Table 2). The reference spatial resolution for these 2 locations (0 and 7 cm) was measured similarly but without including any motion during the simulation. Each simulated list-mode file was rebinned using motion parameters estimated from the 11C-PIB study (LM10%_20s to LM100%_600s scans with and without motion) and reconstructed using filtered backprojection into a 344 × 344 × 127 element volume with a voxel size of 0.21 × 0.21 × 2.03 mm. The FWHM value was then measured from each reconstructed image following the National Electrical Manufacturers Association procedure (13).
Description of Motion-Free and Motion-Corrupted Sinograms
Motion Parameters Used for Generation of Motion-Corrupted Sinograms
Application to an Actual Hoffman Phantom Scan
A Hoffman brain phantom was filled with a solution of 18F-FDG and placed at the center of the field of view of the scanner. A 60-min acquisition was performed during which the phantom was manually moved on a few occasions after 35 min of motionless acquisition. No attempt was made to control the nature or magnitude of the motion. The period that contained motion (between 30 and 60 min) was rebinned into a single static frame with and without motion correction. Motion parameters were estimated using an initial reconstruction, and the resulting volume was smoothed with a 3D gaussian kernel of 16 mm (FWHM). Different framing parameters were tested: 45 frames of 40 s, 60 frames of 30 s, 72 frames of 25 s, 90 frames of 20 s, 120 frames of 15 s, 180 frames of 10 s, and 360 frames of 5 s. During the period from 4.47 to 30 min, when the phantom was motionless, the same number of disintegrations theoretically occurred as during the period from 30 to 60 min; the period from 4.47 to 30 min was therefore used to generate static scans with and without (reference scan) motion correction and using the same framing parameters as above. Each static scan was reconstructed into a 344 × 344 × 127 element volume with a voxel size of 0.83 × 0.83 × 2.03 mm using 3D ordinary Poisson ordered-subsets expectation maximization with all corrections applied and using 3 iterations and 21 subsets. The image-quality improvement (or the image-quality degradation) resulting from the correction process was assessed by computing, each time, the nRMSE between images obtained from data that originally contained motion (or from data that were originally motion-free) and the reference image. In the case of the motion-corrupted data, an estimate of the nRMSE due to differences that are caused by the noise only was subtracted from the measured nRMSE (nRMS0). Supplemental Figure 2 provides details about the method used to estimate this value. Motion-free images were obtained using the same portion of the original list mode as for the reference image, and therefore this adjustment was not required.
Application to an Actual 11C-PIB Study
In total, 119 patients (Mini–Mental State Examination [MMSE] score, 22 ± 5.7; age, 77 ± 6 y; Alzheimer disease, 16%; no cognitive impairment, 15%; mild cognitive impairment, 57%; vascular dementia, 12%) underwent a 30-min brain PET scan 40 min after injection of 370 (±10%) MBq of 11C-PIB. This human study was approved by the Domain Specific Review Board of the National Healthcare Group, and written informed consent was obtained from all patients. Each list-mode datum was rebinned into a single static frame during which motion correction was applied. Motion parameter estimates were measured from an initial reconstruction of the PET data using 90 frames of 20 s and using the first 6 frames (120 s) to create the target. The resulting volumes were smoothed with a 3D gaussian kernel of 16 mm (FWHM). Motion correction was performed using the normalization, data loss compensation, and an oversampling factor of 3. Selection of these correction parameters was based on results obtained from the previous experiments. Motion- and non–motion-corrected scans were reconstructed using 3D ordinary Poisson ordered-subsets expectation maximization with all corrections applied and using 3 iterations and 21 subsets. SUV ratio volumes were generated using the cerebellum gray matter as the reference region. Regional SUV ratios for 12 Alzheimer-specific regions were measured from the normalized volumes using the subject’s parcellated structural MR images.
RESULTS
Rebinner Accuracy
The nRMS0 results reported in Figure 1 clearly indicate that data loss compensation (see the blue curves for rotations about x and y), normalization (see rotation around z; however, blue and green curves are superimposed), and a minimum oversampling factor of 3 (see translations) are required for an accurate correction. An example of correction is shown in Figure 2, demonstrating visually the impact of a 10-mm motion along the z-axis and the contrast recovery achieved with the correction. These results characterize the accuracy of the correction only when the exact motion parameters are used and do not account for inaccuracies in motion parameter estimates.
nRMS0 error computed between corrected and reference volumes as function of magnitude of motion and for different correction parameters: simple line approach (oversampling of 1 = OS1), OS1 including normalization (OS1+norm), OS1+norm with data loss compensation (OS1+norm+compens), and different levels of crystal oversampling (2, 3, and 4). nRMS0 obtained when no correction was applied is also shown for comparison. For clarity, only results obtained with positive motion are shown.
Uncorrected and corrected simulated 18F-FDG images reconstructed from scan with 10-mm translation along z-axis. Normalization, data loss compensation, and oversampling factor of 3 were used for correction. Simulated reference volume is also shown for comparison.
Accuracy of Motion Estimates
Accuracies in motion estimates reported in Figure 3A show that the higher the number of counts per frame the better the motion estimates. Registration led to better results with motion-free sinograms than when motion was present. Using the full statistics scan, a framing of 20 s and a smoothing of 16 mm led to registration errors of 0.41 and 0.85 mm without and with the presence of motion. Halving the injected dose or using a temporal sampling of 10 s (LM50%_20s) increased the errors to 0.52 and 0.99 mm for the motion-free and motion-corrupted data, respectively. Results shown in Figure 3A correspond to inaccuracies caused by parameter estimates only, whereas Figures 3B–3D show the spatial resolution degradation as measured with the simulated point source, that is, caused both by inaccuracies in motion parameter estimates and by the rebinning process. An initial framing of 20 s of the original scan (LM100%_20 s) led to an average 1-mm increment in FWHM as compared with the motionless point source measurements (reference).
(A) Mean absolute displacement across time of point located at radial distance of 7 cm in central transverse plane caused by inaccuracies in motion parameter estimates as function of gaussian smoothing. Results obtained from motion-free sinograms are shown as solid lines, whereas results obtained from motion-corrupted data are shown as dashed lines. (B–D) Spatial resolution degradation resulting from whole correction process (inaccuracies in motion parameter estimates and resolution loss due to rebinning) as function of counting statistics. “Motion” and “No motion” refer to rebinning of point source using parameters measured from motion-corrupted and motion-free 11C-PIB data, respectively. Reference spatial resolutions measured at center (0 cm) and at radial offset of 7 cm are also shown for comparison.
Hoffman Study
Figure 4 shows the differences (nRMSE) between the corrected and reference images, when the correction was applied to motion-corrupted data and motion-free data. The correction applied to the motion-corrupted data led to a decrease of the nRMSE from 16.7% (no correction) to less than 3%. In addition, when applied to motion-free data, the correction process resulted in an nRMSE of less than 4%. Estimated motion from different framing configurations is given in Supplemental Figure 3. Supplemental Figure 4 shows reconstructed images before and after the correction.
nRMSE computed between corrected and reference images, when applied to “Motion – Corrected” motion-corrupted data and “NoMotion – Corrected” motion-free data. nRMSE computed between “Motion – Uncorrected” corrupted scan and reference scan is also shown for comparison.
Actual Studies
An example of the correction on an actual 11C-PIB scan involving a large level of motion is shown in Figure 5, illustrating the improvement in contrast obtained with the correction. Figure 6 summarizes the motion observed during the correction of the 11C-PIB scans. Scans of subjects with a low MMSE score contained a higher level of motion than scans of subjects with a normal MMSE score (Fig. 6A), hinting that motion may unequally affect the different groups.
(A) Uncorrected and corrected 30-min 11C-PIB study with high level of motion. (B) Corresponding displacement of point located 7 cm radially as function of time. Displacement curve obtained from scan with low level of motion is also shown for comparison.
Motion measured from 119 11C-PIB scans. (A) Mean head displacement over time as function of subject MMSE score. (B) Proportion of scans whose mean and maximum displacements computed over scan duration were above distance given in x-axis. (C) Mean and envelope displacement curves across time. (D) Proportion of scans with displacement greater than 2, 5, and 10 mm as function of time. In all graphs, displacement was computed for point in central plane at radial distance of 7 cm on x-axis.
More than 29% and 9% of the scans with an MMSE score below 20 underwent a mean absolute displacement across time of 5 and 10 mm, respectively, during the collection of the data. Within the same group, in more than 56% and 39% of the scans, a displacement of at least 5 and 10 mm, respectively, was observed at least once (Fig. 6B). As shown in Figure 6D, after 10 min of acquisition, 43%, 16%, and 8% of the scans suffered from displacements above 2, 5, and 10 mm, respectively. Figure 7 shows that the correction induced changes above 5% in more than 20% of the regional SUV ratio measurements.
Change in SUV ratio (SUVR) after correction for different cortical ROIs.
DISCUSSION
In the current work, we proposed a novel event-by-event motion correction approach that is dedicated to the Siemens mMR machine. Our proposed implementation accounts for the sources of artifacts that have been previously identified (4), namely LOR discretization, differences in LOR normalization factors, and data loss. Our results using the simulated 18F-FDG scans with well-controlled motion confirmed the importance of addressing these artifacts for an accurate correction. Each of these sources of error acts differently and with varying importance, depending on the nature of the motion. From this study, we found that an oversampling factor of 3 was sufficient, as no significant improvement was obtained with a higher crystal sampling factor. We also tried to find the optimal trade-off between the temporal frequency of the estimate (impacting the counting statistics within each frame) and the accuracy of the motion estimate. We found that, in the case of our 11C-PIB protocol, accurate motion parameters could be obtained every 20 s, with mean errors of 0.41 mm when no motion was present and 0.85 mm with the presence of motion. This optimization is specific to the spatial distribution and counting statistics of the 11C-PIB data obtained with our protocol. Optimal registration and framing parameters should be reevaluated for any other acquisition protocol and tracers. The corresponding resolution loss induced by motion estimate inaccuracies and by the rebinning process was, on average, below 1 mm and systematically below 1.5 mm (increase in FWHM). To the best of our knowledge, this is the first time that results on resolution loss induced by event-by-event rigid motion correction have been reported. The Hoffman study mainly confirmed, using an actual phantom scan, the previous findings. Application of the correction method to the 119 actual 11C-PIB scans led to the following observations: first, that more than 21% of the scans were affected by motion greater than 10 mm (39% for subjects with an MMSE score below 20), and second, that the correction of motion led to changes in Alzheimer-specific cortical regions of up to 30%, proving to be a great source of variability. Figure 7 indicates that the motion-induced biases cancel out (mean changes per region close to 0). However, it is not certain whether this will hold true if the results are broken down by groups.
CONCLUSION
A list-mode rebinner for the Siemens mMR scanner including rigid-body motion-correction capability was developed and fully validated. The rebinner allows accurate motion correction at a cost of a minimal reduction of resolution caused by inaccuracies in motion estimates and by the rebinning process. Application of the correction to a large cohort of 11C-PIB scans confirmed the necessity of correcting for motion to obtain quantitative results.
DISCLOSURE
This work was supported by the French national ‘invest for the future’ programs (LILI Lyon Integrated Life Imaging: hybrid MR-PET ANR-11-EQPX-0026), the hospital University Institut CESAME (Brain and Mental Health ANR-10-IBHU-0003), and by France Life Imaging (FLI, ANR-11-INBS-0006). The 11C-PIB study was funded by NUHS-CG CIRC Seed Funding and by the SoM Aspiration Fund (World Class) Category. No other potential conflict of interest relevant to this article was reported.
Footnotes
Published online Apr. 13, 2018.
- © 2018 by the Society of Nuclear Medicine and Molecular Imaging.
REFERENCES
- Received for publication November 30, 2017.
- Accepted for publication April 4, 2018.
