A direct ROI quantification method for inherent PVE correction: accuracy assessment in striatal SPECT measurements

Abstract

Purpose

The clinical potential of striatal imaging with dopamine transporter (DAT) SPECT tracers is hampered by the limited capability to recover activity concentration ratios due to partial volume effects (PVE). We evaluated the accuracy of a least squares method that allows retrieval of activity in regions of interest directly from projections (LS-ROI).

Methods

An Alderson striatal phantom was filled with striatal to background ratios of 6:1, 9:1 and 28:1; the striatal and background ROIs were drawn on a coregistered X-ray CT of the phantom. The activity ratios of these ROIs were derived both with the LS-ROI method and with conventional SPECT EM reconstruction (EM-SPECT). Moreover, the two methods were compared in seven patients with motor symptoms who were examined with N-3-fluoropropyl-2-β-carboxymethoxy-3-β-(4-iodophenyl) (FP-CIT) SPECT, calculating the binding potential (BP).

Results

In the phantom study, the activity ratios obtained with EM-SPECT were 3.5, 5.3 and 17.0, respectively, whereas the LS-ROI method resulted in ratios of 6.2, 9.0 and 27.3, respectively. With the LS-ROI method, the BP in the seven patients was approximately 60% higher than with EM-SPECT; a linear correlation between the LS-ROI and the EM estimates was found (r = 0.98, p = 0.03).

Conclusion

The LS-ROI PVE correction capability is mainly due to the fact that the ill-conditioning of the LS-ROI approach is lower than that of the EM-SPECT one. The LS-ROI seems to be feasible and accurate in the examination of the dopaminergic system. This approach can be fruitful in monitoring of disease progression and in clinical trials of dopaminergic drugs.

Appendix

In the LS-ROI method, the generalised inverse is considered for the determination of a set of ROI values from raw tomographic data. The generalised inverse gives the least squares solution to a linear system of equations.

Actually, the tomographic problem can be described using a system of linear equations:

$$ p_{{km}} = {\sum\limits_{ij} {F^{{km}}_{{ij}} } }Y_{{ij}} $$

(4)

where p _km is the measured projection at the m-th angle (m=1, .., M), bin k (k=1,.., K); Y _ij is the number of emitted photons in the object pixel (i, j) and the elements $F^{{km}}_{{ij}}$ describe the acquisition process as well as the geometrical system response (resolution, attenuation and scatter). By using a matrix notation, if we call Y the [IJ×1] object voxel concentration vector, p the [KM×1] projection data set and, finally, F the [KM×IJ] projection matrix, the previous equation can be rewritten as:

$$ FY = p. $$

(5)

Least squares methods seek the minimisation of the functional

$$ \chi ^{2} {\left( Y \right)} = {{\sum\limits_{km} {{\left( {{\sum\limits_{ij} {F^{{km}}_{{ij}} Y_{{ij}} - p_{{km}} } }} \right)}^{2} } }} \mathord{\left/ {\vphantom {{{\sum\limits_{km} {{\left( {{\sum\limits_{ij} {F^{{km}}_{{ij}} Y_{{ij}} - p_{{km}} } }} \right)}^{2} } }} {\sigma ^{2}_{{km,}} }}} \right. \kern-\nulldelimiterspace} {\sigma ^{2}_{{km,}} } $$

(6)

where σ _km is the uncertainty with which p _km is measured. The solution is given by:

$$ Y = {\left( {F^{T} \Phi _{p} ^{{ - 1}} F} \right)}^{{ - 1}} F^{T} \Phi _{p} ^{{ - 1}} p, $$

(7)

where $\Phi ^{{ - 1}}_{p}$ is the inverse of the projection covariance matrix. Φ _p is a diagonal matrix whose diagonal elements are the variances of the projection data or, if an unweighted fit is performed, Φ _p =I.

Equation 7 gives the standard least squares solution to Eq. 4. Application of Eq. 7 to the reconstruction of the tomographic images would be appealing because the solution could be obtained directly, without any iterative procedure, allowing the incorporation of a model of non-stationary factors. Unfortunately, Eq. 7 is extremely ill-conditioned, thus threatening the chance of finding an acceptable physical solution.

Let us suppose that our image space can be divided into a small number N _R of ROIs with constant content, such that Y _ij =X ^α for all pixels (i, j) belonging to the region α (R _α), where α=1, ..,N _R . Equation 4, therefore, becomes:

$$ p_{{km}} = {\sum\limits_{\alpha = 1}^{N_{R} } {{\left( {{\sum\limits_{(i,j) \in R_{\alpha } } {F^{{km}}_{{ij}} } }Y_{{ij}} } \right)}} } = {\sum\limits_{\alpha = 1}^{N_{R} } {X^{\alpha } } }{\sum\limits_{(i,j) \in R_{\alpha } } }F^{{km}}_{{ij}} = {\sum\limits_{\alpha = 1}^{N_{R} } {G^{{km}}_{\alpha } } }X^{\alpha } $$

(8)

where $G^{{km}}_{\alpha }$ are called “sinograms” of the regions and represent the projection along the ray at bin k, angle m, of an object which is equal to 1 for voxels belonging to region α and 0 elsewhere.

If we call X the [N _R ×1] region concentration vector and G the [KM×N _R ] region projection matrix, the (unweighted) least squares solution of Eq. 8 is:

$$ X = {\left[ {G^{T} G} \right]}^{{ - 1}} G^{T} p. $$

(9)

Since we have introduced a lot of a priori information about the solution (which is constant over the ROIs), the ill-conditioning of the “ROI” problem (Eq. 9) is strongly reduced with respect to the original “voxel” problem (Eq. 7). Therefore, the solution can be computed directly from projections by using Eq. 9, once regions have been segmented in some way and the sinograms $G^{{km}}_{\alpha }$ for each region have been obtained from the weighting factors $F^{{km}}_{{ij}}$. Therefore, region activity evaluation is not performed with reconstructed images.

The covariance matrix of the solution can also be computed and is given by:

$$ \Phi _{x} = {\left[ {G^{T} G} \right]}^{{ - 1}} G^{T} \Phi _{p} G{\left[ {G^{T} G} \right]}^{{ - 1}} . $$

(10)