|Home | About | Journals | Submit | Contact Us | Français|
We study the noise characteristics of an image reconstruction algorithm that incorporates a model of the non-stationary detector blurring (DB) for a mouse-imaging positron emission tomography (PET) scanner. The algorithm uses ordered subsets expectation maximization (OSEM) image reconstruction, which is used to suppress statistical noise. Including the non-stationary detector blurring in the reconstruction process (OSEM(DB)) has been shown to increase contrast in images reconstructed from measured data acquired on the fully-3D MiCES PET scanner developed at the University of Washington. As an extension, this study uses simulation studies with a fully-3D acquisition mode and our proposed FORE+OSEM(DB) reconstruction process to evaluate the volumetric contrast versus noise trade-offs of this approach. Multiple realizations were simulated to estimate the true noise properties of the algorithm. The results show that incorporation of detector blurring (FORE+OSEM(DB)) into the reconstruction process improves the contrast/noise trade-offs compared to FORE+OSEM in a radially dependent manner. Adding post reconstruction 3D Gaussian smoothing to FORE+OSEM and FORE+OSEM(DB) reduces the contrast versus noise advantages of FORE+OSEM(DB).
The Micro Crystal Element Scanner (MiCES)  is a small animal positron emission tomography (PET) scanner dedicated for mouse-imaging that is under development at the University of Washington. We have previously proposed a pragmatic approach that uses a fully-3D (high sensitivity) data acquisition mode followed by Fourier rebinning (FORE) to reduce the fully-3D data set to a computationally simpler data set of a stacked set of 2D sinograms. The contiguous image planes are then reconstructed using ordered subsets expectation maximization (OSEM) and a model of non-stationary detector blurring (DB) to form a volumetric image . The details of the FORE+OSEM(DB) approach and relation to the more accurate fully-3D maximum a posteriori (MAP) regularized methods [3,4] in the context of mouse-imaging PET scanners have been previously discussed . In brief, the use of iterative methods allows modeling of physical effects (e.g., statistical noise, detector blurring, attenuation, etc.). Including these effects in the image reconstruction, in principle, lead to a more accurate and precise reconstructed estimate of the true activity distribution.
Detector blurring, i.e. detector resolution loss, arises from several sources including the statistics of the scintillation process where the high energy annihilation photon is converted to a large number of low energy photons, each with a certain probability of detection. There is also Compton scatter in the scintillation crystal as well as the finite size of the scintillation crystal. These effects, however, change with position in the scanner as illustrated by Fig. 1. For lines of response (LORs) with increasing radial distance from the center of the scanner, there is an increasing broadening and skewing of the intrinsic detector blurring due to increased probability of interaction, changes in optical photon optics, and differential changes in the likelihood of secondary interactions with respect to angle. These effects have been discussed in more detail elsewhere [2,5-7,10].
Several studies have shown improvements in reconstructed image contrast by including a model of the non-stationary detector blurring in the reconstruction algorithms for both clinical and preclinical scanners [2-5, 7-12].
Qi and Leahy investigated the noise properties of maximum a posteriori (MAP) reconstruction method, compensating physical effects such as photon acollinearity and spatially variant blurring effects in the detector by incorporating those effects into the system model [3-5]. They also compared the performance of MAP algorithm with filtered backprojection by both simulation and experiments . While they showed that MAP improved resolution over FBP by both simulation and in-vivo study in , noise characteristics through different radial positions were not fully investigated.
Less well understood are the trade-offs in image quality resulting from changes in the noise properties resulting from including additional a priori information in the reconstruction process . One of the studies has been conducted by Alessio et. al . They explored a 3D model of detector PSFs which varied in both radial and axial locations, incorporated the radially variant 2D blurring functions into OSEM reconstruction, and evaluated the performance with multiple figures of merit by simulations as well as measured phantoms. They proved that the inclusion of detector PSFs removes spatial variance caused by different radial and axial positions and improves quantitative accuracy with a reduction of tumor bias over OSEM.
While they used 3D and 2D models of detector PSFs in simulation and reconstruction respectively for a clinical PET scanner, in this study we investigate how modeling non-stationary 2D detector blurring followed by reconstruction with 1D PSFs affects the noise structure in reconstructed images for a small animal PET scanner. As an extension of our previous study , we investigate the characteristics [5, 9] of the OSEM(DB) method by simulation studies. We use a fully 3D acquisition mode followed by Fourier rebinning and reconstruction to evaluate the volumetric contrast versus noise trade-offs of including the non-stationary detector blurring in the system model. To reduce computational burden, we utilized analytical simulation methods that incorporate radially variant and axially invariant 2D detector response functions that had been estimated by Monte-Carlo simulations [7,15]. To model the effect of non-stationary detector blurring in the OSEM algorithm, we have added a factorized system matrix, as proposed by Qi and Leahy [3,4], to the ASPIRE reconstruction library [2,16] using the following model, Y = PDBPgeomPprangeX where X is the N x 1 image vector (e.g. typically N = 1282), Prange is an NxN image space operator describing the positron range blurring, Pgeom is an MxN operator describing the geometrical component of the probability of detection of an event from image voxel Xi in data element Yj, PDB is an MxM data space operator describing the detector bluring and Y is a M x 1 vector of the data measurements. In this study, radially varying 1D detector blurring PSFs were utilized to account for the largest component of detector blurring degradation. Image reconstruction was followed by 3D Gaussian smoothing for regularization. Multiple realizations (i.e., 50) were simulated to estimate the true noise properties of the algorithm.
The full MiCES scanner consists of 4 rings (12 cm inner diameter) of modules, with each ring comprised of 18 detector modules. The scanner utilizes a total of 72 photo multiplier tubes (PMTs) (1 per module × 4 rings), each coupled to a 22 × 22 array of 0.8 × 0.8 × 10 mm discrete mixed lutetium silicate (MLS) crystals. There is a 0.1 mm inter-crystal gap between adjacent crystals. In modeling the scanner for this study, the target detector modules were simplified as shown in Fig. 2. In the simplified model, each ring was divided into 396 equally spaced discrete crystals (22 crystals × 18 blocks) along the circumference of the ring. In the axial direction, 4 rings were split into 88 crystals (22 crystals × 4 blocks). Keeping the same detector ring diameter, the crystal cross-sections were approximately 1×1 mm2. No gap was considered between adjacent crystals and the length of the crystals was 10 mm.
Before investigating the noise characteristics of the FORE+OSEM(DB) method, we determined the behavior of FORE+OSEM and FORE+OSEM(DB) with increasing iteration to determine reasonable stopping points for each method. We used two test phantoms that contained cold and hot spheres with a 6 mm diameter in a cylindrical background of 50 mm diameter and 88 mm length as shown in Fig. 3(a) and (b). The activity ratio of the hot sphere to the background object was 2 to 1. For the cold sphere, the ratio was 0.5 to 1.
The analytical simulation tool (ASIM)  was used for this study. The phantom objects were analytically forward-projected into a set of fully 3D sinograms for the combinations of different detector rings, where each sinogram was 198 by 198 (distance bins by angle bins). The sinogram dimensions were selected to reflect the specifications of MiCES system. The slice thickness and distance bins were 0.5 and 0.3 mm respectively. The oblique projection planes were convolved with a 2D positron range point spread function (PSF). The empirical exponential functions by Derenzo were used for implementing the blurring effects of 18F positron range . Then, those projection planes were blurred again with a radially-varying and axially constant 2D detector PSFs estimated with the SimSET photon-tracking simulation tool [19,20]. Fig. 4 shows 2D PSFs at the center, 10, and 20 mm off-centered position for detector blurring. Poisson noise was then added to each of the 3D oblique sinogram bins. The projections were scaled so that the rebinned 2D sinograms have 100M and 10M total counts. Each oblique sinogram bin was Poisson distributed around its mean value leading to realistic noise levels. The noise from random and scattered events was not included in these simulations since our initial estimates are that the scatter fraction will approximately 12% for the MiCES scanner, which we have confirmed with measurements from commercial system with a similar geometry. For the MiCES scanner a typical FDG injection of ~200μCi will lead to a randoms fraction of less than 3%. Based on these measures we are assuming that the Poisson assumption is still valid.
The noise-applied oblique sinograms were resorted by Fourier rebinning into a 2D sinogram stack. Then the 2D sinograms were reconstructed using OSEM and OSEM(DB) with a total of up to 450 updates (9 subsets × 50 iterations). To model the effect of spatially varying detector blurring (DB) in the FORE+OSEM(DB) algorithm, we have generated 1D PSFs using SimSET as shown in Fig. 5 and utilized a dynamic linked library for the factorized system matrix in the ASPIRE reconstruction package . The final images had 128×128×175 voxels with 0.5 mm axial thickness and 0.47 mm pixel distance in a transverse direction.
After each iteration (9 subsets × 1 iteration), regions of interest (ROIs) the same size as the hot and cold spheres were used to calculate mean values of the spheres. A background ROI was also used with the same size as the contrast spheres for the mean of the background. These data were then used to determine the number of iterations used for further experiments.
The test phantom for this study consisted of three small spheres (0.004 cc or 2 mm diameter) in a 25 mm diameter by 88 mm long cylinder, reflecting the entire field of view (FOV) of the MiCES system. The sphere-to-background ratio was set to 4:1 and the sphere centers had a 0, 10, 20 mm offset from the center of the phantom in the transaxial direction. In the same manner as the section II.B, 50 realizations of the sinogram were simulated by ASIM for the case with the PSF for 18F positron range, detector blurring, and Poisson noise level (100M total counts). The 50 realizations were reconstructed with FORE+OSEM and FORE+OSEM(DB) using the parameters determined above. The resulting images were compared with and without post-reconstruction smoothing by a 3D Gaussian filter with 0.47, 0.94, 1.41, 1.88, 2.34, 2.82, and 3.29 mm FWHM. For comparison, the same sinogram data was also reconstructed with standard filtered-backprojection (FORE+FBP) using a Hanning filter with a cutoff at 1.0 , 0.8, 0.6, 0.4, and 0.2 of the Nyquist frequency.
From these images, the noise properties of FORE+OSEM, FORE+OSEM(DB), and FORE+FBP were investigated by estimating the ensemble mean and variance of the ROI values of the background values and of the contrast spheres. In addition images of the variance for each voxel across the 50 realizations were calculated. For quantitative evaluation, bias of ROI mean values were calculated by eq. (1) and eq. (2)
where M is the number of realizations, T is the true value for the ROI, fi,j is a pixel value at index i in the reconstructed image ROI k, and N is the number of pixels in the ROI k.
The results of the contrast versus iteration study are shown in Fig. 6. For both hot (Fig. 6(a)) and cold (Fig. 6(b)) test contrast spheres at both median position(Fig. 3(a)) and FOV boundary(Fig. 3(b)), the values calculated with FORE+OSEM becomes stable with smaller number of iterations than needed for FORE+OSEM(DB). Since the estimated values of the spheres become stable around 72 iterations (FORE+OSEM) and 90 iterations (FORE+OSEM(DB)) as shown in Fig. 6, the number of updates for further study were chosen to be 72 for OSEM (9 subsets and 8 iterations) and 90 for FORE+OSEM(DB) (9 subsets and 10 iterations).
The 50 realizations of the test object were reconstructed with the seven smoothing parameters described above.
For an anecdotal comparison at one set of smoothing parameters (Hanning 0.6 for FORE+FBP and 3D post filtering for FORE+OSEM and FORE+OSEM(DB) with 0.94 mm FWHM), Fig. 7 shows transverse sections through reconstructed realizations at the two count levels. Even though the ideal comparison should be performed among reconstructed images with matched resolution or contrast, this is difficult due to several factors such as the spatially varying noise and resolution, and the lack of a precise match in noise or contrast for any set of parameters. Therefore, Fig. 7 shows a triplet of images with the reconstruction parameters that result in the closest variance of ~ 3 for the point at a 20 mm radial offset with 100 M events (Fig. 11(c) left).
Figs. 8 and and99 show transverse and coronal sections through the mean and variance images of the 50 realizations. Fig. 10 (a) and (b) show profiles through the mean and variance images in 100M and 10M total counts respectively. The smoothing parameters utilized in Fig. 7 were also applied in the images in Figs. 8--1010.
For the spheres at all three radial distances, the bias decreased at the expense of increased noise. In addition the OSEM(DB) showed a relative further reduction in bias for a given level of variance at lower levels of smoothing (i.e. increased variance) compared to FBP and OSEM. This relative advantage of OSEM(DB) increased with radial distance.
Detection bluring (DB) in a fully 3D PET scanner is described by a seven-dimensional function; there are three dimensions describing the location of a point source and there are four dimensions, thus four directions of blurring, in the data space . Modeling the detector blurring in the reconstruction has the potential to improve the accuracy of the estimated images, but with consequent changes in image noise and noise correlations. In this study we investigated the noise characteristics of the FORE+OSEM(DB) reconstruction method for the MiCES mouse imaging PET scanner. We used a simplified 2D mode acquisition and reconstruction process to evaluate volumetric contrast versus noise trade-offs of including the non-stationary detector blurring in the system model. A more realistic model would include all seven dimensions of the detector blurring, but as has been shown the radial part of the radially-dependent blurring is the most significant component . Thus, while the actual impact of modeling the complete detector blurring in the reconstruction will differ from what we find, the general trends should be correct.
We first note the change in convergence properties of FORE+OSEM(DB) as compared to FORE+OSEM as shown in Fig. 6. There is an improvement in the level of contrast reached for both hot and cold objects, although it takes more iterations, or update (iterations x subsets) to reach the higher contrast levels. From these results we chose an increased number of updates of OSEM(DB) compared to OSEM for the subsequent analyses.
From the aggregate images of mean and variance (Figs. 8--10)10) some features are evident: The variance is related to the object value for FORE+OSEM as expected  which can be seen by comparing the mean and variance images for FORE+OSEM. In addition there is a small ring artifact at the edge of the background cylinder for the FORE+OSEM(DB) image. The display levels in Fig. 8 (right) emphasize the artifact, but the true magnitude of the effect can be appreciated by comparison with the profile through the mean image in Fig. 10 (left). This effect has been noted by others, and we currently do not have a complete explanation. We hypothesize that FORE+OSEM(DB) is in essence performing an iterative deconvolution of a point spread function which will cause high frequency image features (boundaries of regions) to hyper-resolve. We used the term “hyper-resolve” because FORE+OSEM(DB) leads to a small ring artifact (Fig. 10 (left)) at the boundary edges which were exaggerated over the true values.
In addition, for small objects (i.e. less than the edge effect size), the ringing(or Gibbs-type) artifact may lead to an overestimation and thus artifactually increased contrast recovery. In patients there are no edges but there may be artifactually increased contrast recovery for small hot spots. Some rational for this is given by Snyder et. al..
We purposely choose to simulate data with a different data blurring operation than the one we model in the reconstruction in order to mimic the reality that the reconstruction system model will always contain approximations.
The variance images in Fig. 9 shows that the FORE+OSEM(DB) algorithm produced images with noise that decreased slightly as radial distance increased, while FORE+OSEM images were not sensitive to the radial position. This effect is likely caused by the radially-varying detector PSFs incorporated in FORE+OSEM(DB) method, such that the second order statistics have a radially variant convergence rate and also converge to values that vary radially. The result of the specific shape of the PSF (Fig. 5) is that, per update, the algorithm to converge more slowly at the edge of the field of view, leading to less noise at the edge than at the center. Meanwhile, the results in Fig. 6 show that the convergence speed for hot and cold objects at increased radius is similar for the same objects at the median radial position. We note that the background value, however, takes longer to converge. This could be an effect of the ringing artifact discussed earlier.
The mean and variance images do not capture the effects of noise correlations within an image. To better understand the noise across realizations we looked at the mean and variance of an ROI placed over the small spheres. Our results (Fig. 11) agree with earlier findings [e.g. 5,8-11], that show that the incorporation of detector blurring into the system matrix (FORE+OSEM(DB)) improves contrast, or in our analysis reduces bias, compared to FORE+OSEM. In addition we find that incorporation of detector blurring also induces a decreased bias at all radial locations at lower levels of smoothing (i.e. increased variance) compared to FORE+FBP and FORE+OSEM. This relative advantage of FORE+OSEM(DB) increased with radial distance. For higher levels of smoothing, the three algorithms produced roughly equivalent trade-offs in noise versus bias. These results are for a geometry appropriate for a mouse-imaging PET scanner, but share scale-invariant properties with whole-body clinical PET scanners [5,8,9].
With a more accurate model (i.e. 7D) of the detector blurring incorporated in the system model used in image reconstruction, we would expect a more accurate image, i.e. with less bias. The overall gain however is still likely to be reduced when image smoothing is increased, as the relative gains in improved contrast will be overshadowed by the applied regularization.
The overall benefit of this approach then will depend on the specific task, whether it is estimation, detection, localization, classification, or some combination of these metrics. Given the particular metric, incorporation of detector blurring in the system model used in image reconstruction will favor metrics that have better performance with images that are less smoothed or regularized. In the case that these metrics favor more regularized (smoothed) images there would not seem to be an advantage to FORE+OSEM(DB) compared to FORE+OSEM or FORE+FBP.
This work was supported by National Institutes of Health under grants R01-CA74135, R01-CA86892, R01-EB0217 and R01-CA115870. We acknowledge helpful discussions with Drs Jeffrey Fessler and Richard Leahy.