An interior point iterative maximum-likelihood reconstruction algorithm incorporating upper and lower bounds with application to SPECT transmission imaging.

The algorithm we consider here is a block-iterative (or ordered subset) version of the interior point algorithm for transmission reconstruction. Our algorithm is an interior point method because each vector of the iterative sequence [x(k)], k = 0, 1, 2, ... satisfies the constraints a(j) < x(j)k < b(j), j = 1, ..., J. Because it is a block-iterative algorithm that reconstructs the transmission attenuation map and places constraints above and below the pixel values of the reconstructed image, we call it the BITAB method. Computer simulations using the three-dimensional mathematical cardiac and torso phantom, reveal that the BITAB algorithm in conjunction with reasonably selected prior upper and lower bounds has the potential to improve the accuracy of the reconstructed attenuation coefficients from truncated fan beam transmission projections. By suitably selecting the bounds, it is possible to restrict the over estimation of coefficients outside the fully sampled region, that results from reconstructing truncated fan beam projections with iterative transmission algorithms such as the maximum-likelihood gradient type algorithm.

Improved image quality and computation reduction in 4-D reconstruction of cardiac-gated SPECT images.

Spatiotemporal reconstruction of cardiac-gated SPECT images permits us to obtain valuable information related to cardiac function. However, the task of reconstructing this four-dimensional (4-D) data set is computation intensive. Typically, these studies are reconstructed frame-by-frame: a nonoptimal approach because temporal correlations in the signal are not accounted for. In this work, we show that the compression and signal decorrelation properties of the Karhunen-Loève (KL) transform may be used to greatly simplify the spatiotemporal reconstruction problem. The gated projections are first KL transformed in the temporal direction. This results in a sequence of KL-transformed projection images for which the signal components are uncorrelated along the time axis. As a result, the 4-D reconstruction task is simplified to a series of three-dimensional (3-D) reconstructions in the KL domain. The reconstructed KL components are subsequently inverse KL transformed to obtain the entire spatiotemporal reconstruction set. Our simulation and clinical results indicate that KL processing provides image sequences that are less noisy than are conventional frame-by-frame reconstructions. Additionally, by discarding high-order KL components that are dominated by noise, we can achieve savings in computation time because fewer reconstructions are needed in comparison to conventional frame-by-frame reconstructions.

Noise characterization of block-iterative reconstruction algorithms: I. Theory.

Researchers have shown increasing interest in block-iterative image reconstruction algorithms due to the computational and modeling advantages they provide. Although their convergence properties have been well documented, little is known about how they behave in the presence of noise. In this work, we fully characterize the ensemble statistical properties of the rescaled block-iterative expectation-maximization (RBI-EM) reconstruction algorithm and the rescaled block-iterative simultaneous multiplicative algebraic reconstruction technique (RBI-SMART). Also included in the analysis are the special cases of RBI-EM, maximum-likelihood EM (ML-EM) and ordered-subset EM (OS-EM), and the special case of RBI-SMART, SMART. A theoretical formulation strategy similar to that previously outlined for ML-EM is followed for the RBI methods. The theoretical formulations in this paper rely on one approximation, namely, that the noise in the reconstructed image is small compared to the mean image. In a second paper, the approximation will be justified through Monte Carlo simulations covering a range of noise levels, iteration points, and subset orderings. The ensemble statistical parameters could then be used to evaluate objective measures of image quality.

Reducing the influence of the partial volume effect on SPECT activity quantitation with 3D modelling of spatial resolution in iterative reconstruction.

Quantitative parameters such as the maximum and total counts in a volume are influenced by the partial volume effect. The magnitude of this effect varies with the non-stationary and anisotropic spatial resolution in SPECT slices. The objective of this investigation was to determine whether iterative reconstruction which includes modelling of the three-dimensional (3D) spatial resolution of SPECT imaging can reduce the impact of the partial volume effect on the quantitation of activity compared with filtered backprojection (FBP) techniques which include low-pass, and linear restoration filtering using the frequency distance relationship (FDR). The iterative reconstruction algorithms investigated were maximum-likelihood expectation-maximization (MLEM), MLEM with ordered subset acceleration (ML-OS), and MLEM with acceleration by the rescaled-block-iterative technique (ML-RBI). The SIMIND Monte Carlo code was used to simulate small hot spherical objects in an elliptical cylinder with and without uniform background activity as imaged by a low-energy ultra-high-resolution (LEUHR) collimator. Centre count ratios (CCRs) and total count ratios (TCRs) were determined as the observed counts over true counts. CCRs were unstable while TCRs had a bias of approximately 10% for all iterative techniques. The variance in the TCRs for ML-OS and ML-RBI was clearly elevated over that of MLEM, with ML-RBI having the smaller elevation. TCRs obtained with FDR-Wiener filtering had a larger bias (approximately 30%) than any of the iterative reconstruction methods but near stationarity is also reached. Butterworth filtered results varied by 9.7% from the centre to the edge. The addition of background has an influence on the convergence rate and noise properties of iterative techniques.

Accelerating the EMML algorithm and related iterative algorithms by rescaled block-iterative methods.

Analysis of convergence of the algebraic reconstruction technique (ART) shows it to be predisposed to converge to a solution faster than simultaneous methods, such as those of the Cimmino-Landweber type, the expectation maximization maximum likelihood method for the Poisson model (EMML), and the simultaneous multiplicative ART (SMART), which use all the data at each step. Although the choice of ordering of the data and of relaxation parameters are important, as Herman and Meyer have shown, they are not the full story. The analogous multiplicative ART (MART), which applies only to systems y=Px in which y>0, P= or >0 and a nonnegative solution is sought, is also sequential (or "row-action"), rather than simultaneous, but does not generally exhibit the same accelerated convergence relative to its simultaneous version, SMART. By dividing each equation by the maximum of the corresponding row of P, we find that this rescaled MART (RMART) does converge faster, when solutions exist, significantly so in cases in which the row maxima are substantially less than one. Such cases arise frequently in tomography and when the columns of P have been normalized to have sum one. Between simultaneous methods, which use all the data at each step, and sequential (or row-action) methods, which use only a single data value at each step, there are the block-iterative (or ordered subset) methods, in which a single block or subset of the data is processed at each step. The ordered subset EM (OSEM) of Hudson et al. is significantly faster than the EMML, but often fails to converge. The "rescaled block-iterative" EMML (RBI-EMML) is an accelerated block-iterative version of EMML that converges, in the consistent case, to a solution, for any choice of subsets; it reduces to OSEM when the restrictive "subset balanced" condition holds. Rescaled block-iterative versions of SMART and MART also exhibit accelerated convergence.

Choice of initial conditions in the ML reconstruction of fan-beam transmission with truncated projection data.

We investigate the effects of initial conditions in the iterative maximum-likelihood (ML) reconstruction of fan-beam transmission projection data with truncation. In an iterative ML reconstruction, the estimate of the transmission reconstructed image in the previous iteration is multiplied by some factors to obtain the current estimate. Normally, a flat initial condition (FIC) or an image with equal positive pixel values is used as initial condition for an ML reconstruction. Usage of FIC has also been perceived as a way of preventing any bias on the reconstruction which may have come from the initial condition. When projection data have truncation, we show that using an FIC in an ML iterative reconstruction can introduce a bias to the reconstruction inside the densely sampled region (DSR), whose projection data have no truncation at any angle. To reduce this bias, we propose to use the largest right singular vector (LRSV) of the system matrix as the initial condition, and demonstrate that the bias can be reduced with the LRSV. When data truncation is reduced, the LRSV approaches the FIC. This result does not contradict to the use of FIC when projection data are not truncated. We also demonstrate that the reconstructed transmission image using LRSV as initial condition provides a more accurate attenuation coefficient distribution than that using FIC. However, the improvement is mostly in the area outside the DSR.

Convergent block-iterative algorithms for image reconstruction from inconsistent data.

It has been shown that convergence to a solution can be significantly accelerated for a number of iterative image reconstruction algorithms, including simultaneous Cimmino-type algorithms, the "expectation maximization" method for maximizing likelihood (EMML) and the simultaneous multiplicative algebraic reconstruction technique (SMART), through the use of rescaled block-iterative (BI) methods. These BI methods involve partitioning the data into disjoint subsets and using only one subset at each step of the iteration. One drawback of these methods is their failure to converge to an approximate solution in the inconsistent case, in which no image consistent with the data exists; they are always observed to produce limit cycles (LCs) of distinct images, through which the algorithm cycles. No one of these images provides a suitable solution, in general. The question that arises then is whether or not these LC vectors retain sufficient information to construct from them a suitable approximate solution; we show that they do. To demonstrate that, we employ a "feedback" technique in which the LC vectors are used to produce a new "data" vector, and the algorithm restarted. Convergence of this nested iterative scheme to an approximate solution is then proven. Preliminary work also suggests that this feedback method may be incorporated in a practical reconstruction method.

Block-iterative methods for image reconstruction from projections.

The simultaneous MART algorithm (SMART) and the expectation maximization method for likelihood maximization (EMML) are extended to block-iterative versions, BI-SMART and BI-EMML, that converge to a solution in the feasible case, for any choice of subsets. The BI-EMML reduces to the "ordered subset" EMML of Hudson et al. (1992, 1994) when their "subset balanced" property holds.

Noniterative compensation for the distance-dependent detector response and photon attenuation in SPECT imaging.

A filtering approach is described, which accurately compensates for the 2D distance-dependent detector response, as well as for photon attenuation in a uniform attenuating medium. The filtering method is based on the frequency distance principle (FDP) which states that points in the object at a specific source-to-detector distance provide the most significant contribution to specified frequency regions in the discrete Fourier transform (DFT) of the sinogram. By modeling the detector point spread function as a 2D Gaussian function whose width is dependent on the source-to-detector distance, a spatially variant inverse filter can be computed and applied to the 3D DFT of the set of all sinogram slices. To minimize noise amplification the inverse filter is rolled off at high frequencies by using a previously published Wiener filter strategy. Attenuation compensation is performed with Bellini's method. It was observed that the tomographic point response, after distance-dependent filtering with the FDP, was approximately isotropic and varied substantially less with position than that obtained with other correction methods. Furthermore, it was shown that processing with this filtering technique provides reconstructions with minimal degradation in image fidelity.

Iterative image reconstruction algorithms based on cross-entropy minimization.

The related problems of minimizing the functionals F(x)=alphaKL(y,Px)+(1-alpha)KL(p,x) and G(x)=alphaKL(Px,y)+(1-alpha)KL(x,p), respectively, over the set of vectors x=/>0 are considered. KL(a, b) is the cross-entropy (or Kullback-Leibler) distance between two nonnegative vectors a and b. Iterative algorithms for minimizing both functionals using the method of alternating projections are derived. A simultaneous version of the multiplicative algebraic reconstruction technique (MART) algorithm, called SMART, is introduced, and its convergence is proved.

The effect of intrinsic attenuation correction methods on the stationarity of the 3-D modulation transfer function of SPECT.

The application of stationary restoration techniques to SPECT images assumes that the modulation transfer function (MTF) of the imaging system is shift invariant. It was hypothesized that using intrinsic attenuation correction (i.e., methods which explicitly invert the exponential radon transform) would yield a three-dimensional (3-D) MTF which varies less with position within the transverse slices than the combined conjugate view two-dimensional (2-D) MTF varies with depth. Thus the assumption of shift invariance would become less of an approximation for 3-D post- than for 2-D pre-reconstruction restoration filtering. SPECT acquisitions were obtained from point sources located at various positions in three differently shaped, water-filled phantoms. The data were reconstructed with intrinsic attenuation correction, and 3-D MTFs were calculated. Four different intrinsic attenuation correction methods were compared: (1) exponentially weighted backprojection, (2) a modified exponentially weighted backprojection as described by Tanaka et al. [Phys. Med. Biol. 29, 1489-1500 (1984)], (3) a Fourier domain technique as described by Bellini et al. [IEEE Trans. ASSP 27, 213-218 (1979)], and (4) the circular harmonic transform (CHT) method as described by Hawkins et al. [IEEE Trans. Med. Imag. 7, 135-148 (1988)]. The dependence of the 3-D MTF obtained with these methods, on point source location within an attenuator, and on shape of the attenuator, was studied. These 3-D MTFs were compared to: (1) those MTFs obtained with no attenuation correction, and (2) the depth dependence of the arithmetic mean combined conjugate view 2-D MTFs.(ABSTRACT TRUNCATED AT 250 WORDS)

Image-restoration algorithms for a fully connected architecture.

We describe the implementation of a technique for achieving image superresolution using a fully connected network of simple processors operating in an iterative mode. We show that an updating scheme can be specified that ensures convergence for the serial (asynchronous) updating case. With the appropriate hardware, parallel (synchronous) updating becomes of particular interest because of the potential for accelerated convergence; it is this approach that we envisage implementing in optical hardware. For this case also, we present a convergent scheme that can be related to a regularized form of the Gerchberg-Papoulis algorithm.

Limit of continuous and discrete finite-band Gerchberg iterative spectrum extrapolation.

It is known that the Gerchberg method for iteratively extrapolating the spectrum of an object of finite support has a different limit when the infinite frequency band used in the iterations is replaced by a finite one. It is shown that, for both continuous and discrete data, this limit involves an optimal extrapolation within the finite band to minimize the component of the limit that lies beyond the object support.