Gradient-based iterative methods often converge slowly for tomographic image reconstruction and image restoration problems, but can be accelerated by suitable preconditioners. Diagonal preconditioners offer some improvement in convergence rate, but do not incorporate the structure of the Hessian matrices in imaging problems. Circulant preconditioners can provide remarkable acceleration for inverse problems that are approximately shift-invariant, i.e. for those with approximately block-Toeplitz or block-circulant Hessians. However, in applications with nonuniform noise variance, such as arises from Poisson statistics in emission tomography and in quantum-limited optical imaging, the Hessian of the weighted least-squares objective function is quite shiftvariant, and circulant preconditioners perform poorly. Additional shift-variance is caused by edge-preserving regularization methods based on nonquadratic penalty functions. This paper describes new preconditioners that approximate more accurately the Hessian matrices of shiftvariant imaging problems. Compared to diagonal or circulant preconditioning, the new preconditioners lead to significantly faster convergence rates for the unconstrained conjugate-gradient (CG) iteration. We also propose a new efficient method for the line-search step required by CG methods. Applications to positron emission tomography (PET) illustrate the method.

We present a new class of algorithms for edge-preserving restoration of piecewise-smooth images measured in nonGaussian noise under shift-variant blur. The algorithms are based on minimizing a regularized objective function, and are guaranteed to monotonically decrease the objective function. The algorithms are derived by using a combination of two previously unconnected concepts: A. De Pierro's convexity technique for optimization transfer, and P. Huber's iteration for M-estimation. Convergence to the unique global minimum is guaranteed for strictly convex objective functions. The convergence rate is very fast relative to conventional gradient-based iterations. The proposed algorithms are flexibly parallelizable, and easily accommodate nonnegativity constraints and arbitrary neighborhood structures. Implementation in Matlab is remarkably simple, requiring no cumbersome line searches or tolerance parameters. Keywords: Image restoration, non-Gaussian noise, deconvolution, Bayesian meth...

Fourier-based reprojection methods have the potential to reduce the computation time in iterative tomographic image reconstruction. Interpolation errors are a limitation of Fourier-based reprojection methods. We apply a min-max interpolation method for the nonuniform fast Fourier transform (NUFFT) to minimize the interpolation errors. Numerical results show that the min-max NUFFT approach provides substantially lower approximation errors in tomographic reprojection and backprojection than conventional interpolation methods.

In a wide range of applied problems of 2D and 3D imaging a continuous formulation of the problem places great emphasis on obtaining and manipulating the Fourier transform in Polar coordinates. However, the translation of continuum ideas into practical work with data sampled on a Cartesian grid is problematic. In this article we develop a fast high accuracy Polar FFT. For a given two-dimensional signal of size N × N, the proposed algorithm’s complexity is O(N^2 log N), just like in a Cartesian 2D-FFT. A special feature of our approach is that it involves only 1D equispaced FFT’s and 1D interpolations. A central tool in our method is the pseudo-Polar FFT, an FFT where the evaluation frequencies lie in an oversampled set of nonangularly equispaced points. We describe the concept of pseudo-Polar domain, including fast forward and inverse transforms. For those interested primarily in Polar FFT’s, the pseudo-Polar FFT plays the role of a halfway point—a nearly-Polar system from which conversion to Polar coordinates uses processes relying purely on 1D FFT’s and interpolation operations. We describe the conversion process, and give an error analysis of it. We compare accuracy results obtained by a Cartesian-based unequally-sampled FFT method to ours, both algorithms using a small-support interpolation and no pre-compensating, and show marked advantage to the use of the pseudo-Polar initial grid.

Abstract—A low-complexity recursive procedure is presented for minimum mean squared error (MMSE) estimation in linear regression models. A Gaussian mixture is chosen as the prior on the unknown parameter vector. The algorithm returns both an approximate MMSE estimate of the parameter vector and a set of high posterior probability mixing parameters. Emphasis is given to the case of a sparse parameter vector. Numerical simulations demonstrate estimation performance and illustrate the distinctions between MMSE estimation and MAP model selection. The set of high probability mixing parameters not only provides MAP basis selection, but also yields relative probabilities that reveal potential ambiguity in the sparse model. 1 I.

A low-complexity recursive procedure is presented for model selection and minimum mean squared error (MMSE) estimation in linear regression. Emphasis is given to the case of a sparse parameter vector and fewer observations than unknown parameters. A Gaussian mixture is chosen as the prior on the unknown parameter vector. The algorithm returns both a set of high posterior probability mixing parameters and an approximate MMSE estimate of the parameter vector. Exact ratios of posterior probabilities serve to reveal potential ambiguity among multiple candidate solutions that are ambiguous due to observation noise or correlation among columns in the regressor matrix. Algorithm complexity is linear in the number of unknown coefficients, the number of observations and the number of nonzero coefficients. If hyperparameters are unknown, a maximum likelihood estimate is found by a generalized expectation maximization algorithm. Numerical simulations demonstrate estimation performance and illustrate the distinctions between MMSE estimation and maximum a posteriori probability model selection.

Preconditioning methods can accelerate the convergence of gradient-based iterative methods for tomographic image reconstruction and image restoration. Circulant preconditioners have been used extensively for shiftinvariant problems. Diagonal preconditioners offer some improvement in convergence rate, but do not incorporate the structure of the Hessian matrices in imaging problems. For inverse problems that are approximately shiftinvariant (i.e. approximately block-Toeplitz or blockcirculant Hessians), circulant or Fourier-based preconditioners can provide remarkable acceleration. However, in applications with nonuniform noise variance (such as arises from Poisson statistics in emission tomography and in quantum-limited optical imaging), the Hessian of the (penalized) weighted least-squares objective function is quite shift-variant, and the Fourier preconditioner performs poorly. Additional shift-variance is caused by edge-preserving regularization methods based on nonquadratic penalty ...

Since the data sizes in fully 3D PET imaging are very large, iterative image reconstruction algorithms must converge in very few iterations to be useful. One can improve the convergence rate of the conjugategradient (CG) algorithm by incorporating preconditioning operators that approximate the inverse of the Hessian of the objective function. If the 3D cylindrical PET geometry were not truncated at the ends, then the Hessian of the penalized least-squares objective function would be approximately shift-invariant, i.e. G 0 G would be nearly block-circulant, where G is the system matrix. We propose a Fourier preconditioner based on this shift-invariant approximation to the Hessian. Results show that this preconditioner significantly accelerates the convergence of the CG algorithm with only a small increase in computation. I. Introduction Statistical methods for tomographic image reconstruction from fully 3D PET scans are particularly promising since the number of unknown parameters (...

In this paper, we focus on the restoration of images that have incomplete data in either the image domain or the transformed domain or in both. The transform used can be any orthonormal or tight frame transforms such as orthonormal wavelets, tight framelets, the discrete Fourier transform, the Gabor transform, the discrete cosine transform, and the discrete local cosine transform. We propose an iterative algorithm that can restore the incomplete data in both domains simultaneously. We prove the convergence of the algorithm and derive the optimal properties of its limit. The algorithm generalizes, unifies, and simplifies the inpainting algorithm in image domains given in [8] and the inpainting algorithms in the transformed domains given in [7,16,19]. Finally, applications of the new algorithm to super-resolution image reconstruction with different zooms are presented. 1

Abstract not found