| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | RSS | TABLE OF CONTENTS |
|
|
|
||||||||
|
||||||||||
Clinical Investigations |
1 Division of Nuclear Medicine, Department of Radiology, Harvard Medical School and Brigham and Women’s Hospital, Boston, Massachusetts
2 Ernest Orlando Lawrence Berkeley National Laboratories, Berkeley, California
 |
ABSTRACT |
|---|
 TOP ABSTRACT INTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
Key Words: quantitative 82Rb cardiac PET • factor analysis of dynamic sequences • compartment analysis
|
INTRODUCTION |
|---|
|
TOPABSTRACTINTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
15 min), allowing rest and stress imaging under virtually identical conditions and decreasing the total time required to scan each patient. In addition to the poor spatial resolution due to the relatively long positron range of the 82Rb positron emitter, the major challenges of parametric 82Rb imaging are estimation of an accurate input function from noisy data without arterial blood sampling and the absolute quantitation of myocardial perfusion. One approach, proposed by Lin et al. (8,9), to deal with noise is to filter the reconstructed volume before region-of-interest quantitation using wavelet transforms. Compartmental modeling for 82Rb cardiac PET has also been investigated (3,8,10), and 2 models have been proposed by Gould (10): a model that reduces the 2 compartments to one unknown, flow, which is calculated from myocardial uptake and arterial input function, and a compartmental model that fits the model equations to observed myocardial and blood time–activity curves (TACs) to estimate the 2 unknowns. The physiologic model of 82Rb kinetics in the myocardium has also been compared with reduced-order models by Coxson et al. (11). In this study, we used a generalized form of the least-squares factor analysis of dynamic sequences (FADS) to obtain a robust estimate of left and right ventricular (LV+RV) input functions automatically, without the need to draw volumes of interest (VOIs), and developed a compartment model based on orthogonal grouping to estimate on a voxel-by-voxel basis, rather than within selected VOIs, extraction of 82Rb in the myocardium. We validated our approach by Monte Carlo simulation studies and demonstrated clinical feasibility by patient studies.
|
MATERIALS AND METHODS |
|---|
|
TOPABSTRACTINTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
Furthermore, although these techniques increase the range of situations in which unique FADS solutions are achieved, they do not ensure a unique solution in all cases. We have previously developed a technique that is more general than previously reported approaches and can be used in a variety of applications, with a wide range of ejection fractions and blood TACs (18). In this work, we address the nonuniqueness problem by penalizing its main effect, spatial overlap between factor images. In our approach, FADS is performed first by any factor analysis method such as the apex-seeking (19) or least-squares (17) approaches. After FADS, a second step is used to minimize the overlap between factor images and, hence, increase the probability of a unique solution. This second step is described and tested in Monte Carlo simulations of realistic dynamic cardiac 82Rb PET studies, and the feasibility of the approach is demonstrated in patients.
Theory.
A dynamic sequence of medical images can be represented by an N x M matrix A. N is the number of voxels in one dynamic image and M is the number of time frames. The factor model of the dynamic data assumes that the data matrix can be represented by the following equation:
 | (Eq. 1) |
 | (Eq. 2) |
Nonnegativity constraints were applied by adding another penalty term, fn(C,F), to the objective function (Eq. 3). fn(C,F) is a quadratic function that heavily penalizes negative values of C and F and yields a positive solution (17):
 | (Eq. 3) |
 | (Eq. 4) |
Therefore, the total objective function that was used in the first step was:
 | (Eq. 5) |
Because the factor model described by Equation 1 is not mathematically unique, the result of the optimization described by Equation 5 cannot be used for quantitative analysis. To address this problem, we generalized our previous approach (18) to allow implementation as a postprocessing method following any factor model technique. In this approach, GFADS, we modify the factors and factor images obtained in the first step by minimizing the nonnegative term funi(R) that penalizes the overlap between images of factor coefficients while keeping CF constant:
 | (Eq. 6) |
The result of the minimization is C' = CR, where R is the rotation matrix determined during the optimization, so that the results obtained in the first step hold during optimization in the second step:
 | (Eq. 7) |
Hence, the problem of nonuniqueness of the FADS solution is addressed in our GFADS approach by minimizing the major effect of nonuniqueness, overlap between the factor images. This is accomplished by optimizing the rotation matrix R by a minimization of funi under nonnegativity constraints imposed by fn(CR,R–1F). The total objective function that was used in the second step (funi + bfn) was minimized by a simplex algorithm (20), and convergence was obtained after 100–500 iterations. The value of the penalty parameter, b, 0.1 in this work, was not found to be critical for the convergence to a solution when chosen in the range [0.05, 0.5].
VOI Analysis
LV, RV, and myocardial TACs were also estimated using a VOI-based approach. The initial VOIs used to define the activity distributions in the LV, RV, or myocardium for the Monte Carlo simulations were also used for estimating the TACs. These volumes were first reduced by one eighth (50% in each dimension) to mimic usual clinical practice to minimize spillover among LV, RV, and myocardial tissues. For the patient studies, VOIs were drawn by a skilled operator over several slices of the 3-dimensional (3D) volume of the LV, RV, or myocardium obtained by summing the first 24 dynamic frames, and counts within these volumes were estimated at each frame to generate the corresponding VOI TAC.
Compartment Analysis
Because 82Rb is only partially extracted by the myocardium, a 2-compartment model is required for accurate estimation of myocardial extraction fraction (10). The 2 compartments of the model are the "free rubidium space" (blood perfusing the myocardium plus interstitial space) and the "trapped rubidium space" (muscle of the myocardium). The main parameters of the model are the kinetic transport constants k1 (mL/min/g) and k2 (min–1), which denote the extraction (forward) and egress (backward) rates of transport between the metabolically trapped space (myocardium) and the freely diffusible space (blood pool), respectively.
Kinetic Model.
The TAC in each voxel was modeled as a combination of 3 contributions: the contribution from myocardial tissue, modeled using a 2-compartment model, and contributions from LV and RV cavities, modeled as fractions of measured LV and RV functions. In cardiac imaging it is necessary to explicitly model the contribution from the RV to obtain stable and robust voxel parameter estimates (10). This model can be expressed as:
 | (Eq. 8) |
Orthogonal Grouping.
The estimation of parametric images on a voxel-by-voxel basis yields very noisy images due to high levels of noise in TACs derived from single voxels. Using a filtering approach reduces noise but can result in a degradation of spatial resolution of the parametric images. In this work we used a new method, orthogonal grouping, which consisted of a 3D grouping of voxels with similar time course into several groups of size G. The average TACs within the groups were used for calculation of the kinetic parameters. This approach leads to reduced levels of noise with minimal effects on spatial resolution. The orthogonal grouping consisted of 2 steps:
, where = 1... 3 that described the location of a given dynamic voxel TACi(t) in the 3D data subspace S, and where P denotes the projection operator onto the th base vector. . First, Euclidian distances dij = 
P(TACi) – P(TACj) were calculated between each 2 dynamic voxels. The dynamic voxel for which the sum of distances to all other voxels was highest was determined. Next, the G – 1 closest dynamic voxels were found. These G dynamic voxels (initial 1 and G – 1) constituted the first group. The procedure was repeated until all dynamic voxels were grouped. G was chosen iteratively to be the smallest value that ensures convergence for all groups (G was 30 for the Monte Carlo simulations and 100 for the patent studies). Finally, the estimated kinetic parameters were remapped back to the original volume.
Monte Carlo Simulations of Dynamic 82Rb PET Studies
Kinetic Modeling of Dynamic 82Rb Studies.
LV and RV TACs were measured in 10 human rest and dipyridamole stress dynamic cardiac 82Rb PET studies using small VOIs positioned over the LV and RV cavities. The patients used to generate the LV and RV input functions were imaged at different time points in the life of the 82Rb generator. There were 6 subjects without evidence of obstructive coronary artery disease (CAD) as determined by the absence of regional perfusion defects on the dipyridamole stress images (summed stress score [SSS] = 0). In addition, there were 4 patients with known obstructive CAD and extensive perfusion abnormalities on the dipyridamole stress images (SSS > 18). This yielded LV and RV TACs that varied greatly in shape and represented a wide range of physiologic conditions.
Because 82Rb is a partially extracted tracer, the myocardial activity detected in a small region is the sum of the activity trapped in the cell and the free circulating activity. When modeling the myocardial uptake, we assumed that the myocardial tissue contains approximately 25% of arterial blood (10). The LV TAC, I(t), was used as input function for a 2-compartment kinetic model to generate a pure myocardial extraction fraction (extravascular space) TAC:
 | (Eq. 9) |
Myocardial extraction TACs were generated for rest studies with k1 values of 0.4, 0.7, and 1.2 mL/min/g and for stress studies with k1 values of 0.7, 1.4, and 2.1 mL/min/g. Because 82Rb is an analog of potassium that is taken up rapidly and retained for a long time within myocardial cells (10), we used very low k2 values of 0.01 and 0.02 min–1 in rest and stress studies. The selected values of k1 and k2 spanned the range of values reported in the literature in rest and stress studies (10). The lowest two k1 values were associated with k2 = 0.01 min–1, whereas the highest k1 value was associated with k2 = 0.02 min–1. Typical myocardial TACs obtained using our kinetic model were compared with those measured using small myocardial VOIs in patients, and several were found to match patient TACs quite well. Since 3 pairs of (k1, k2) values were used with 10 pairs of (LV, RV) input functions, a set of 30 unique realistic simulated dynamic studies was obtained. The TACs for the liver and soft tissue were assumed to be delayed constant functions, identical to each other. Sample LV and RV input functions, as well as 3 myocardial TACs, are shown in Figure 1 (solid lines).
|
We modeled a 2D PET scanner similar to the GE DST Discovery PET/CT scanner (GE Healthcare) with a uniform 3-cm-thick bismuth germanate (BGO) detector (energy resolution, 20% at 511 keV; axial extent, 15 cm; ring diameter, 88 cm). A stack of annular septa (septal thickness, 0.1 cm) was explicitly modeled in the Monte Carlo simulation in the axial direction, yielding twenty-eight 3.3-mm-thick slices (voxel size, 3.3 x 3.3 x 3.3 mm3). Events reaching the detector were binned into 128 projections using only the direct planes. Nine orders of Compton scatter were modeled. Positron range was also modeled using the analytic approach developed by Palmer and Bronwell (24), which uses ß-decay energy spectra and empiric range formulae.
Random coincidences, not presently modeled in SimSET, were also modeled in the simulation software that we used in this work. Random coincidences were estimated from singles using:
 | (Eq. 10) |
Nine billion photons were generated when simulating each structure (LV, RV, myocardium, soft tissues). This yielded essentially noise-free sinograms that were used as the basis for generating noisy sinograms. Scatter, known from simulation, as well as random coincidences, estimated using the singles-rate formula (Eq. 10), was subtracted. The resulting sinograms were precorrected for attenuation by multiplying each ray-sum by the inverse of the corresponding attenuation factor, and reconstructed using ordered-subsets expectation maximization (OSEM, 5 iterations, 4 subsets) (26). Stationary resolution was assumed in the forward projector and no postfiltering was performed. This yielded 30 studies, each consisting of twenty-four 5-s frames (128 x 128 x 28 voxels each).
Dynamic 82Rb Patient Studies
Thirteen patients were included in this analysis. Five patients had evidence of obstructive CAD, and 8 had no evidence of CAD (Table 1). Rest and stress dynamic protocols were followed on 8 subjects, whereas rest dynamic data were available on 4 subjects and stress dynamic data were available on 1 subject. Patients undergoing the rest/stress protocol were injected with a bolus of 2.22 ± 0.19 GBq (60 ± 5 mCi) of 82Rb (Bracco Diagnostics) in 14 ± 6 mL of saline and imaged dynamically for 6 min at rest and, 10 min later, during dipyridamole stress (intravenous infusion of 0.14 mg/kg/min for 4 min). The dynamic imaging protocol consisted of twenty-four 5-s frames followed by eight 30-s frames. The initial fast frames were used to capture the rapid washin and washout of 82Rb from the LV and RV cavities.
|
Data Analysis
Estimation performance was assessed for the Monte Carlo–simulated data by computing known true values. The LV and RV TACs estimated with GFADS or by VOI analysis were compared with the LV and RV input functions used to generate the dynamic studies. Likewise, the myocardial TACs estimated with GFADS and VOI were compared with the true myocardial TAC. For a given Monte Carlo simulation, the total GFADS (or VOI) estimation error 
GFADS of the LV, RV, or myocardial TAC was computed over the 24 time frames, and the average error over the 30 Monte Carlo studies was expressed as:
 | (Eq. 11) |
The SD of the average error was also calculated over the 30 Monte Carlo studies after averaging over the time frames. Likewise, the average errors associated with myocardial tissue extraction (k1) and egress (k2), as well as RV and LV contributions (fvi, rvi), were computed over the 30 Monte Carlo simulations. The rate constants k1, k2, as well as fvi, and rvi, were averaged over each tissue compartment using the known simulated activity distribution as a mask, and errors were calculated using the known parameters’values used to generate the dynamic studies. Means and SDs of the errors were calculated over the 30 studies.
|
RESULTS |
|---|
|
TOPABSTRACTINTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
|
Patient clinical data at rest and stress, as well as the SSS and CFR, are given in Table 1. In all patient studies,GFADS yielded LV and RV input functions, myocardial TACs, and corresponding factor images of the expected form (Fig. 3). Furthermore, the 2-compartment kinetic model yielded parametric maps of myocardial tissue extraction and egress as well as spillover fractions from the blood pool into the myocardial tissue (Fig. 4). Note that these images were obtained despite the short acquisition time of the dynamic frames (5 s each) and without postreconstruction filtering; this implies that the new approach is robust to noise. As expected, the k2 values measured in patients using the LV and RV input functions estimated with GFADS were very low and ranged between 0.012 and 0.028 min–1. The corresponding k1 values ranged between 0.52 mL/min/g in patients with suspected CAD and 1.03 mL/min/g in patients without suspected CAD at rest. The corresponding range of k1 values at stress was 0.74–1.84 mL/min/g. Furthermore, k1 values estimated in anterior, inferior, lateral, and septal walls did not differ significantly (NS, paired t test) in patients with normal myocardial blood perfusion (SSS = 0). The flow values (k1) estimated using the LV and RV input functions estimated by VOI analysis were systematically greater (8.3%–27.6%) than those obtained using GFADS LV and RV input functions (P < 0.05). k2 values associated with VOIs were slightly greater (5%–10%) than those with GFADS; however, fv (0.39 ± 0.07) and rv (0.31 ± 0.02) parameters were significantly greater (13.4%–31.7%) with VOIs as compared with GFADS (P < 0.05). The CFR was also calculated as the ratio of flow at peak stress and flow at rest (k1 stress/k1 rest). CFR values ranged from 1.58 to 2.3 in patients with no prior known CAD and from 1.39 to 1.51 in patients with known prior CAD. k1 and k2 values were not strongly correlated at rest or stress (r2 < 0.69).
|
|
|
DISCUSSION |
|---|
|
TOPABSTRACTINTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
In the Monte Carlo simulations, positron range was modeled using an analytic approach (24). Modeling positron range is important in PET Monte Carlo simulations of 82Rb, because the full width at half maximum of the positron range is 2.6 mm as compared with 1 mm for 18F, so the effect of positron range on spatial resolution and image quality is greater for 82Rb than for 18F. Likewise, modeling random coincidences is of paramount importance in the case of 82Rb cardiac PET because of the high counting rates after injection of 2.22 GBq (60 mCi) of 82Rb, and the torso attenuation. Dead time is another important factor that affects singles and coincidences differently; although this effect was not modeled in this study, we are currently working on including it in the Monte Carlo simulation.
In the Monte Carlo simulations, patient-derived TACs were used to generate spatially uniform LV and RV compartments, and GFADS yielded accurate estimates of the original simulated TACs with very little overlap between LV, RV, and myocardium. In patient studies, one factor might not be sufficient to describe the dynamic activity distribution in the LV or RV when the activity does not vary uniformly within these structures. In this case, more factors would be needed to fully describe the dynamic ventricular activities. We are presently investigating the impact of patient hemodynamics on the shape of LV and RV input functions. The kinetic approach presented in this work was designed to model such situations; the grouping algorithm does not necessarily result in spatially contiguous clusters as the distance metric is not a spatial distance vector. This allows different locations within a given structure to have different kinetic parameters, a desirable property when nonuniform dynamic activity distributions, such as those in diseased patients, are to be modeled. However, in the presence of respiratory motion, the lack of spatial constraints in the orthogonal grouping can lead to similar time curves that represent different tracer kinetics. One potential solution to this problem is to use spatial constraints in the orthogonal analysis. Note that differences in tissue count recovery due to partial-volume effect or cardiac motion do not affect the estimation of k1 and k2 by orthogonal grouping, as shown earlier, but only the estimation of fv and rv.
In all patient studies, GFADS yielded LV and RV input functions of the expected form and the 2-compartment kinetic model yielded parametric maps of myocardial tissue extraction and egress as well as spillover fractions from the blood pool into the myocardial tissue. The flows estimated with our approach were well within the range reported previously using an approach based on wavelet transforms, and validated in H215O (9), and were consistent with the clinical status and catheterization results in the limited number of patients considered. CFR was systematically higher in subjects with no prior known CAD (2.1 ± 0.3) than in patients with known prior CAD (1.4 ± 0.1). The estimated CFR values in subjects without obstructive CAD (Table 1) are somewhat lower than those usually found in the literature (2.5). This may be related to the fact that these subjects had multiple coronary risk factors (e.g., hypertension, dyslipidemia, diabetes) and were not healthy volunteers. Furthermore, as flow increases, myocardial extraction may no longer be proportional to flow (as in the low-flow regime) and peak flow may be underestimated in dipyridamole studies. This flow dependency of the first-pass extraction fraction can be modeled in our kinetic modeling approach after estimating the LV and RV input functions with GFADS if it is known a priori and, hence, would not significantly affect the accuracy of the flow estimates. It can also be included as a postprocessing step.
In subjects with known CAD who underwent dipyridamole stress, k1 was 30%–50% lower when the VOI drawn on the parametric image corresponded to a region irrigated by a stenotic vessel than when the VOI was drawn over a large myocardial area. Therefore, our approach, based on factor and compartment analyses, allowed discrimination between normal and diseased myocardial blood flow in different subjects, as well as in the same subject. Although GFADS taken alone does not allow the discrimination between different myocardial regions in the same patient (because only one factor is currently used to model the myocardial TAC), the combination of the generalized factor analysis and compartment analysis approaches allows this aim to be achieved by using GFADS to estimate the LV and RV input functions and then by estimating k1 values on a voxel-by-voxel basis.
|
CONCLUSION |
|---|
|
TOPABSTRACTINTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
|
ACKNOWLEDGMENTS |
|---|
|
FOOTNOTES |
|---|
For correspondence or reprints contact: Georges El Fakhri, PhD, Division of Nuclear Medicine, Department of Radiology, Brigham and Women’s Hospital, 75 Francis St., Boston, MA 02115.
E-mail: elfakhri{at}bwh.harvard.edu
|
REFERENCES |
|---|
|
TOPABSTRACTINTRODUCTIONMATERIALS AND METHODSRESULTSDISCUSSIONCONCLUSIONREFERENCES |
|---|
This article has been cited by other articles:
|
|
 |
|
  V. Vaccarino, J. Votaw, T. Faber, E. Veledar, N. V. Murrah, L. R. Jones, J. Zhao, S. Su, J. Goldberg, J. P. Raggi, et al. Major Depression and Coronary Flow Reserve Detected by Positron Emission Tomography Arch Intern Med, October 12, 2009; 169(18): 1668 - 1676. [Abstract] [Full Text] [PDF]
|
|
|
|
 |
|
 G. El Fakhri, A. Kardan, A. Sitek, S. Dorbala, N. Abi-Hatem, Y. Lahoud, A. Fischman, M. Coughlan, T. Yasuda, and M. F. Di Carli Reproducibility and Accuracy of Quantitative Myocardial Blood Flow Assessment with 82Rb PET: Comparison with 13N-Ammonia PET J. Nucl. Med., July 1, 2009; 50(7): 1062 - 1071. [Abstract] [Full Text] [PDF]
|
|
|
|
 |
|
 F. M. Bengel, T. Higuchi, M. S. Javadi, and R. Lautamaki Cardiac positron emission tomography. J. Am. Coll. Cardiol., June 30, 2009; 54(1): 1 - 15. [Abstract] [Full Text] [PDF]
|
|
|
|
 |
|
 M. M. Hajjiri, M. B. Leavitt, H. Zheng, A. E. Spooner, A. J. Fischman, and H. Gewirtz Comparison of Positron Emission Tomography Measurement of Adenosine-Stimulated Absolute Myocardial Blood Flow Versus Relative Myocardial Tracer Content for Physiological Assessment of Coronary Artery Stenosis Severity and Location J. Am. Coll. Cardiol. Img., June 1, 2009; 2(6): 751 - 758. [Abstract] [Full Text] [PDF]
|
|
|
|
|
|
M. F. Di Carli, S. Dorbala, J. Meserve, G. El Fakhri, A. Sitek, and S. C. Moore Clinical Myocardial Perfusion PET/CT J. Nucl. Med., May 1, 2007; 48(5): 783 - 793. [Abstract] [Full Text] [PDF]
|
|
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | RSS | TABLE OF CONTENTS |
| JOURNAL OF NUCLEAR MEDICINE TECHNOLOGY | THE JOURNAL OF NUCLEAR MEDICINE |
