A fast and stable maximum a posteriori conjugate gradient reconstruction algorithm.
We have derived a maximum a posteriori (MAP) approach for iterative reconstruction based on a weighted least-squares conjugate gradient (WLS-CG) algorithm. The WLS-CG algorithm has been shown to have initial convergence rates up to 10x faster than the maximum-likelihood expectation maximization (ML-EM) algorithm, but WLS-CG suffers from rapidly increasing image noise at higher iteration numbers. In our MAP-CG algorithm, the increasing noise is controlled by a Gibbs smoothing prior, resulting in stable, convergent solutions. Our formulation assumes a Gaussian noise model for the likelihood function. When a linear transformation of the pixel space is performed (the "relaxation" acceleration method), the MAP-CG algorithm obtains a low-noise, stable solution (one that does not change with further iterations) in 10-30 iterations, compared to 100-200 iterations for MAP-EM. Each iteration of MAP-CG requires approximately the same amount of processing time as one iteration of ML-EM or MAP-EM. We show that the use of an initial image estimate obtained from a single iteration of the Chang method helps the algorithm to converge faster when acceleration is not used, but does not help when acceleration is applied. While both the WLS-CG and MAP-CG methods suffer from the potential for obtaining negative pixel values in the iterated image estimates, the use of the Gibbs prior substantially reduces the number of pixels with negative values and restricts them to regions of little or no activity. We use SPECT data from simulated hot-sphere phantoms and from patient studies to demonstrate the advantages of the MAP-CG algorithm. We conclude that the MAP-CG algorithm requires 10%-25% of the processing time of EM techniques, and provides images of comparable or superior quality.