Improving contrast between gray and white matter of Logan graphical analysis' parametric images in positron emission tomography through least-squares cubic regression and principal component analysis

LGA is derived from a set of differential equations describing the in vivo behavior of an administered radiotracer. LGA [2] is described by the following equation:

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\displaystyle \int }_{0}^{t}C(u){du}}{C(t)}={\text{}}{DVR}\\ \quad \cdot \,\left(\displaystyle \frac{{\displaystyle \int }_{0}^{t}{C}^{{\rm{R}}}(u){du}+\tfrac{1}{{k}_{2}^{{\rm{R}}}}{C}^{{\rm{R}}}(t)}{C(t)}\right)+{int},\end{array}\end{eqnarray} \tag{ 1 }$$

where C(t) and C^R(t) denote the radioactivity in the target and reference tissue, respectively. ${k}_{2}^{R}$ denotes the dissociation rate constant in the reference region. The target tissue, with the target receptors, is the region of interest, whereas the reference region is chosen to have a negligible concentration of target receptors. For ¹¹ C-PiB PET data used in this study, the cerebellum gray matter region was used as the reference region [13].

The y-intercept term (int) in (1) becomes constant after a time t^*. Then, a linear relationship between $\tfrac{{\int }_{0}^{t}C(u){du}}{C(t)}$ and $\tfrac{{\int }_{0}^{t}{C}^{R}(u){du}+\tfrac{1}{{k}_{2}^{R}}{C}^{R}(t)}{C(t)}$ can be obtained from (1). Applying linear regression to (1), the DVR can be estimated for time T > t^*.

The limitation of OLS-based LGA estimates is that OLS only accounts for noise in the response variable, but both the predictor and response variables of the LGA in (1) contain noise due to the noisy term C(t) in their denominators. As a linear regression method, LSC was developed to account for noise in both the predictor and response variables [14]. Therefore, if applied to LGA, LSC can account for noise in both the predictor and response variables by minimizing the weighted squared residuals in both the variables. In addition, LSC also accounts for noise correlation in the variables. Accounting for noise correlation is also important as the noise in the LGA predictor and response variables are highly correlated [5].

Mathematical details of LSC as a regression method can be found in [14]. For specific details of LSC in LGA applications, see [10].

PCA is a dimensionality reduction technique used mostly for feature extraction and noise filtering. PCA transforms variables into new sets of variables that are linear functions of the original variables [15]. The new sets of variables are uncorrelated and are defined by sets of orthogonal basis vectors and principal components that optimally describe the variance in the data. The principal components and corresponding variance are obtained by solving for the eigenvalues and eigenvectors of the covariance matrix of the data. The eigenvectors are the principal components and eigenvalues the corresponding variances.

The original variables of the noise filtering applications used in this study are recovered using lower-dimensional projections, resulting in denoised estimates of the original variables.

The number of principal components retained to estimate the denoised data in this study were set such that the specified percentage variance could be realized. pecifically, we evaluated the estimates of LSC–PCA and OLS–PCA methods obtained from the data denoised at 9 percentage variance levels of 91% to 99%. This is because in real practical situations, the variance level at which the best trade-off between denoising and information loss is not known. It is therefore necessary to investigate the two methods over a range of percentage variance levels.

A set of kinetic parameters—the exchange rate constants K₁, k₂₋₋₄ and non-displaceable distribution volume (V_ND )—obtained from other studies, such as [11] and [16], was used to make the first simulation data set, for the CFN radiotracer.

• •  
  11 values of DVR were set in the range of [0.0, 4.0], covering a range of [0.0, 0.35] for k₃.
• •  
  $[{K}_{1}\,\,{k}_{4}]=[0.1835\,{\text{mL cm}}^{-1}\,{{\rm{\min }}}^{-1}\,\,\,\,0.115\,{{\rm{\min }}}^{-1}]$.
• •  
  [k₂ k₃] = [K₁/V_ND (DVR − 1) · k₄].
• •  
  V_ND = 1.59 mLcm⁻¹.

The data was simulated for 60 min, in 27-frame scans (6 × 10 s, 3 × 30 s, $5\,\times \,1\,{\rm{\min }}$, $5\,\times \,2.5\,{\rm{\min }}$, and $8\,\times \,5\,{\rm{\min }}$). Two-tissue compartment PET data were simulated using a clinically measured plasma time-activity curve. The reference region was formed with a minimal noise level, and DVR = 1. As per the physiological assumptions of the reference region, the delivery (K₁) and clearance (k₂) rates were set identical to those of the target tissues.

For ¹¹ C-PiB simulations, the kinetic parameters, K₁, k₂₋₋₄ and V_ND , obtained from [17], were used to simulate the two-tissue compartmental model PET data.

• •  
  11 values of DVR were set in the range of [0.0, 4.0], covering a range of [0.0, 0.045] for k₃.
• •  
  $[{K}_{1}\,\,{k}_{2}]=[0.262\,{\text{mL cm}}^{-1}\,{{\rm{\min }}}^{-1}\,\,\,\,0.121\,{{\rm{\min }}}^{-1}]$.
• •  
  [k₃ k₄] = [(DVR − 1) · k₄ 0.015].
• •  
  V_ND = K₁/k₂.

For contrast comparison in simulation data, a second set of hypothetical white matter regions were simulated. For ¹¹ C-PiB Aβ-negative images, white matter regions can have high retention of the radiotracer due to nonspecific retention [18, 19], and slow accumulation and slower clearance [18]. Therefore, the white matter regions were simulated with the values of K₁ (0.190 mLcm⁻¹) and k₂ (0.100 ${{\rm{\min }}}^{-1}$) lower than those of the gray matters, and DVR values high than those of the gray matters. Thus, to simulate white matter with high DVR values than the garay matters, we added the 11 unit vector, [0.4: 0.055: 0.95], to the 11 unit vector of the gray matter DVRs. This gave the lowest white matter DVR a value of 1.6727 (the lowest for gray matter is 1.2727), and the highest of 4.95 (the highest for gray matter is 4). The ratio of these gray and white matter DVR values are within the ranges reported in a previous study [20].

The scan period and number of frames for the ¹¹ C-PiB simulation were set to be similar to those of the clinical ¹¹ C-PiB PET images used in this study. The scan period was set at 70 min with 25-frame scans (6 × 10 s, 3 × 20 s, $2\,\times \,1\,{\rm{\min }}$, $2\,\times \,3\,{\rm{\min }}$, and $12\,\times \,5\,{\rm{\min }}$). The reference region was formed as for CFN; with minimal noise level; DVR = 1; K₁ and k₂ identical to those of the target tissues.

For the two radiotracers, a noise-free tTAC was formed for each DVR value. Statistical noise was added to the noise-free tTACs to form noisy tTACs. For each noise-free tTAC, 1024 noisy tTACs were simulated as a slice of 32 × 32 voxels. Four methods, LSC–PCA, OLS–PCA, LSC, and conventional OLS-based LGA, were then used to estimate the DVR values from the noisy tTACs, and the results were compared. The true value for the dissociation rate of the reference region (${k}_{2}^{R}$) used to construct the simulation data was used for all methods.

Following the common approach [3, 7, 11], the noise model was based on counting statistics via the half-life of the radionuclide and scanning time, such that the noise increases in the later time frames, which have lower count rates. The noise was assumed to have a zero-mean Gaussian distribution with variance proportional to the true curves. This approach is well-described in [7]. The noise was scaled to within the magnitude range observed in the actual voxel-based PET data.

The simulated data of the two radiotracers, CFN and ¹¹ C-PiB, were denoised individually by PCA. The gray and white matter of the simulation data for ¹¹ C-PiB were denoised together. For each data set (1024 simulated tTACS for each of the 11 true DVRs), DVR estimation from the noisy tTACs were performed thrice at each percentage variance level (i.e., 27 times) for the PCA based methods. The three estimates at each variance level were differed by a slight change of <5% in the magnitude of noise level. The magnitude of the noise level refers to the percentage ratio of the means of "standard deviation in the noisy tTACs to the mean of the noise-free tTAC".

To have equal set of data for all methods, the LSC and OLS estimates were performed 27 times at noise levels corresponding to those of the PCA methods. The mean and standard deviations were then calculated for 27 estimates for each 1024 simulated tTACs for each true DVR value.

For ¹¹ C-PiB data, in addition to the averages and standard deviations, the DVR parametric estimates were compared numerically in terms of the contrast between the simulated gray and white matter regions. The contrast index was calculated as [21]:

$$\begin{eqnarray}&&{ \mathcal C }=\displaystyle \frac{{\mu }_{W}-{\mu }_{G}}{{\sigma }_{W}+{\sigma }_{G}},\end{eqnarray} \tag{ 2 }$$

where μ_G and μ_W denote the mean DVR of the gray and white matter regions, respectively, and σ_G and σ_W denote the standard deviations of DVR in the gray and white matter, respectively.

A cohort of 12 (11 Aβ-negative and 1 Aβ-positive) subjects was included in this study. The subjects underwent 70 min ¹¹ C-PiB PET scans of 25 frames (6 × 10 s, 3 × 20 s, $2\,\times \,1\,{\rm{\min }}$, $2\,\times \,3\,{\rm{\min }}$, and $12\,\times \,5\,{\rm{\min }}$). The scan procedure was carried out according to the Alzheimer's Disease Neuroimaging Initiative (ADNI) ¹¹ C-PiB standards. The study was approved by the Ethics Committee of Kindai University Hospital, and written informed consent was obtained from all participants.

The obtained PET data were in the form of voxels of size 2.1 × 2.1 × 3.4 mm³ in matrices of 128 × 128 × 47 slices. DVR parametric images were then generated from the PET data using the LSC–PCA and OLS–PCA methods.

As with the simulation data, DVR estimates were performed at the 9 different percentage variance levels. However, the estimates were only performed once at each percentage variance since there is no noise adjustment. Contrasts were then calculated and compared for LSC–PCA and OLS–PCA for these estimates. Contrasts were calculated between the four main gray matter cortices (frontal, temporal, occipital, and parietal) and white matter (corona radiata).

All data analyses were performed in MATLAB (The MathWorks, Inc., Natick, MA, USA).

figure 1 presents the simulated tTACs, mimicking the two radiotracers, 1(a) CFN and 1(b) ¹¹ C-PiB. Noise-free tTACs are represented by the thick curves, whereas the thin curves represent the noisy tTACs. For each noise-free curve, out of the 1024 simulated noisy curves, only one is shown in each of the subfigures in figure 1. The lowermost thick curve denotes the reference region, and the upper 11 curves denote the 11 DVR values.

Zoom In Zoom Out Reset image size

figure 2 shows examples of PCA-denoised tTACs for the two radiotracers, 2(a) CFN and 2(b) ¹¹ C-PiB. The percentage labels of the tTACs, PCA (91%)—PCA (99%), in the legends denote the tTACs denoised with the number of PCs that respectively explain 91%—99% variance in the data. The corresponding noise-free and noisy (unprocessed by PCA) tTACs are also included. The average number of PCs that denoised the data at the specified percentage variance levels are shown in table 1. Table 1 shows that the CFN data required more PCs to attain the specified percentage variance value, compared with ¹¹ C-PiB. The total number of PCs for the CFN and ¹¹ C-PiB data are 27 and 25, respectively, corresponding to their number of time points. For both radiotracers, the tTACs denoised at 98% variance level are showing to be still noisy. At 95% and 91%, noise reduction could be seen in the tTACs. One can consider lowering the number of PCs, but that would equate to lowering the amount of explained variance, and hence information loss.

Zoom In Zoom Out Reset image size

figure 3 compares the percentage recovered DVR in the averages of the values estimated from the simulated noisy tTACs using four methods, LSC–PCA, OLS–PCA, LSC, and OLS. In other words, the difference from 100% denotes the percentage bias. For LSC–PCA and OLS–PCA, the averages are for estimates from tTACs denoised at 9 different percentage variance levels, estimated thrice at each level (27 estimates). LSC and OLS averages are from tTACs at 27 different noise levels. As expected, LSC–PCA estimates have smaller error bars than the LSC estimates. This is especially noticeable for ¹¹ C-PiB estimates in figure 3(b), in which some larger variances seen in LSC estimates are well stabilized in LSC–PCA. This demonstrates PCA's effectiveness in reducing the variance in the estimates —the purpose for which it is needed in combination with LSC. In figure 3(a), we see minimal effects in the bias and standard deviations in the estimates for both LSC–PCA compared to LSC, and OLS–PCA compared to OLS. In this case, LSC–PCA maintains its advantage of smaller biases, whereas the bias in OLS–PCA estimates approach those for OLS. The difference in the observations for the two radiotracer data could be due to that CFN data required more number of PCs, thus reintroducing the noise in the data. Radiotracer-focused studies may help determine the ranges of percentage variance that should be retained for various radiotracers. Tables 2 and 3 show the average values of the estimated DVR for the respective CFN and ¹¹ C-PiB radiotracer, and respectively corresponding to figures, 3(a) and figure 3(b).

Zoom In Zoom Out Reset image size

figure 4 compares the contrast between the simulated gray and white matter regions. The comparisons include all calculated DVR values (three times at 9 percentage variance levels for LSC–PCA and OLS–PCA, and at 27 different noise levels for LSC and OLS). Statistical comparisons of these contrast are shown in table 4, using the Friedman test [22] and Post-Hoc multiple comparison. A statistical significant difference is observed between LSC–PCA and all other methods. On the other hand, no statistical significant difference is observed between OLS–PCA and LSC.

Zoom In Zoom Out Reset image size

figure 5 compares the contrast between gray matter cortices and white matter for clinical ¹¹ C-PiB DVR images obtained by LSC–PCA and OLS–PCA, calculated for the 11 Aβ-negative subjects, at 9 different percentage variance levels. The Friedman test achieved an overall p value of 0.008, and a Pot-Hoc p value of 0.187 between LSC-PCA and OLS-PCA.

Zoom In Zoom Out Reset image size

figure 6 shows axial slices of the DVR parametric images obtained using four methods, LSC–PCA, OLS–PCA, LSC, and OLS. PCA-based methods were denoised with the number of PCs that explained 95% variance for these images. The two images in the upper row are from two Aβ-negative subjects, while the bottom row is from the Aβ-positive subject. The means and standard deviations of the number of PCs that denoised the data at the specified percentage variance levels for the 12 clinical ¹¹ C-PiB brain images are shown in table 5. These data shows the expected trend of the number of PCs increasing with the amount explained percentage variance in the data.

Zoom In Zoom Out Reset image size

In figure 6, noise components in the images estimated by the LSC method can be seen appearing as random noisy pixels/regions with high DVR values in the brain region (see the black arrows in figure 6). Images obtained by the LSC–PCA method exhibit a rather fair distribution of the DVR values with no sharp increases in the pixels/regions in the main brain region. However, some bright pixels could be seen on the outskirts of the LSC–PCA images in the middle and bottom row. These could be attributed to the LSC, as this is part of the variance concept which is being addressed by combining LSC and PCA. Further analysis should be carried out to investigate if these bright pixels could be observed within the brain region in LSC–PCA images. The images obtained by OLS–PCA appear similar to those for LSC–PCA. The LSC images thus appear slightly chunky (compared to the LSC–PCA and OLS–PCA images), for these specific settings, which can hinder the contrast between the brain structures. The lowest DVR estimates are observed in OLS-based images, as expected.

figure 7 shows DVR parametric images in which the difference between LSC–PCA and OLS–PCA could be seen at some levels of percentage variance for PCA denoising. In this figure we can see that as the percentage variance to be explained is increased to 97% from 95%, the OLS–PCA estimates get biased. Calculated contrast values for this subject for DVR estimates from the data denoised at 97% variance are shown in table 6, showing lower values for OLS–PCA. Another case for an image showing biased estimates for OLS–PCA is shown in "Additional file" (figure 2A). These results reveals the effect of noise reintroduction in the data, when many PCs are used, on the OLS–PCA estimates.

Zoom In Zoom Out Reset image size