Solution of inverse problems in image processing by wavelet expansion

Joint NDT image restoration and segmentation using Gauss-Markov-Potts prior models and variational Bayesian computation

In this paper, we propose a method to simultaneously restore and to segment piecewise homogeneous images degraded by a known point spread function (PSF) and additive noise. For this purpose, we propose a family of nonhomogeneous Gauss-Markov fields with Potts region labels model for images to be used in a Bayesian estimation framework. The joint posterior law of all the unknowns (the unknown image, its segmentation (hidden variable) and all the hyperparameters) is approximated by a separable probability law via the variational Bayes technique. This approximation gives the possibility to obtain practically implemented joint restoration and segmentation algorithm. We will present some preliminary results and comparison with a MCMC Gibbs sampling based algorithm. We may note that the prior models proposed in this work are particularly appropriate for the images of the scenes or objects that are composed of a finite set of homogeneous materials. This is the case of many images obtained in nondestructive testing (NDT) applications.

A fast multilevel algorithm for wavelet-regularized image restoration

We present a multilevel extension of the popular "thresholded Landweber" algorithm for wavelet-regularized image restoration that yields an order of magnitude speed improvement over the standard fixed-scale implementation. The method is generic and targeted towards large-scale linear inverse problems, such as 3-D deconvolution microscopy. The algorithm is derived within the framework of bound optimization. The key idea is to successively update the coefficients in the various wavelet channels using fixed, subband-adapted iteration parameters (step sizes and threshold levels). The optimization problem is solved efficiently via a proper chaining of basic iteration modules. The higher level description of the algorithm is similar to that of a multigrid solver for PDEs, but there is one fundamental difference: the latter iterates though a sequence of multiresolution versions of the original problem, while, in our case, we cycle through the wavelet subspaces corresponding to the difference between successive approximations. This strategy is motivated by the special structure of the problem and the preconditioning properties of the wavelet representation. We establish that the solution of the restoration problem corresponds to a fixed point of our multilevel optimizer. We also provide experimental evidence that the improvement in convergence rate is essentially determined by the (unconstrained) linear part of the algorithm, irrespective of the type of wavelet. Finally, we illustrate the technique with some image deconvolution examples, including some real 3-D fluorescence microscopy data.

Weight assignment for adaptive image restoration by neural networks

This paper presents a scheme for adaptively training the weights, in terms of varying the regularization parameter, in a neural network for the restoration of digital images. The flexibility of neural-network-based image restoration algorithms easily allow the variation of restoration parameters such as blur statistics and regularization value spatially and temporally within the image. This paper focuses on spatial variation of the regularization parameter.We first show that the previously proposed neural-network method based on gradient descent can only find suboptimal solutions, and then introduce a regional processing approach based on local statistics. A method is presented to vary the regularization parameter spatially. This method is applied to a number of images degraded by various levels of noise, and the results are examined. The method is also applied to an image degraded by spatially variant blur. In all cases, the proposed method provides visually satisfactory results in an efficient way.

A Bayesian wavelet-based multidimensional deconvolution with sub-band emphasis

This paper proposes a new algorithm for wavelet-based multidimensional image deconvolution which employs subband-dependent minimization and the dual-tree complex wavelet transform in an iterative Bayesian framework. In addition, this algorithm employs a new prior instead of the popular l(1) norm, and is thus able to embed a learning scheme during the iteration which helps it to achieve better deconvolution results and faster convergence.

Multigrid tomographic inversion with variable resolution data and image spaces

A multigrid inversion approach that uses variable resolutions of both the data space and the image space is proposed. Since the computational complexity of inverse problems typically increases with a larger number of unknown image pixels and a larger number of measurements, the proposed algorithm further reduces the computation relative to conventional multigrid approaches, which change only the image space resolution at coarse scales. The advantage is particularly important for data-rich applications, where data resolutions may differ for different scales. Applications of the approach to Bayesian reconstruction algorithms in transmission and emission tomography with a generalized Gaussian Markov random field image prior are presented, both with a Poisson noise model and with a quadratic data term. Simulation results indicate that the proposed multigrid approach results in significant improvement in convergence speed compared to the fixed-grid iterative coordinate descent method and a multigrid method with fixed-data resolution.

A segmentation-based regularization term for image deconvolution

This paper proposes a new and original inhomogeneous restoration (deconvolution) model under the Bayesian framework for observed images degraded by space-invariant blur and additive Gaussian noise. In this model, regularization is achieved during the iterative restoration process with a segmentation-based a priori term. This adaptive edge-preserving regularization term applies a local smoothness constraint to pre-estimated constant-valued regions of the target image. These constant-valued regions (the segmentation map) of the target image are obtained from a preliminary Wiener deconvolution estimate. In order to estimate reliable segmentation maps, we have also adopted a Bayesian Markovian framework in which the regularized segmentations are estimated in the maximum a posteriori (MAP) sense with the joint use of local Potts prior and appropriate Gaussian conditional luminance distributions. In order to make these segmentations unsupervised, these likelihood distributions are estimated in the maximum likelihood sense. To compute the MAP estimate associated to the restoration, we use a simple steepest descent procedure resulting in an efficient iterative process converging to a globally optimal restoration. The experiments reported in this paper demonstrate that the discussed method performs competitively and sometimes better than the best existing state-of-the-art methods in benchmark tests.

Bayesian wavelet-based image deconvolution: a GEM algorithm exploiting a class of heavy-tailed priors

Image deconvolution is formulated in the wavelet domain under the Bayesian framework. The well-known sparsity of the wavelet coefficients of real-world images is modeled by heavy-tailed priors belonging to the Gaussian scale mixture (GSM) class; i.e., priors given by a linear (finite of infinite) combination of Gaussian densities. This class includes, among others, the generalized Gaussian, the Jeffreys, and the Gaussian mixture priors. Necessary and sufficient conditions are stated under which the prior induced by a thresholding/shrinking denoising rule is a GSM. This result is then used to show that the prior induced by the "nonnegative garrote" thresholding/shrinking rule, herein termed the garrote prior, is a GSM. To compute the maximum a posteriori estimate, we propose a new generalized expectation maximization (GEM) algorithm, where the missing variables are the scale factors of the GSM densities. The maximization step of the underlying expectation maximization algorithm is replaced with a linear stationary second-order iterative method. The result is a GEM algorithm of O(N log N) computational complexity. In a series of benchmark tests, the proposed approach outperforms or performs similarly to state-of-the art methods, demanding comparable (in some cases, much less) computational complexity.

Wavelet-based reconstruction for limited-angle X-ray tomography

The aim of X-ray tomography is to reconstruct an unknown physical body from a collection of projection images. When the projection images are only available from a limited angle of view, the reconstruction problem is a severely ill-posed inverse problem. Statistical inversion allows stable solution of the limited-angle tomography problem by complementing the measurement data by a priori information. In this work, the unknown attenuation distribution inside the body is represented as a wavelet expansion, and a Besov space prior distribution together with positivity constraint is used. The wavelet expansion is thresholded before reconstruction to reduce the dimension of the computational problem. Feasibility of the method is demonstrated by numerical examples using in vitro data from mammography and dental radiology.

An EM algorithm for wavelet-based image restoration

This paper introduces an expectation-maximization (EM) algorithm for image restoration (deconvolution) based on a penalized likelihood formulated in the wavelet domain. Regularization is achieved by promoting a reconstruction with low-complexity, expressed in the wavelet coefficients, taking advantage of the well known sparsity of wavelet representations. Previous works have investigated wavelet-based restoration but, except for certain special cases, the resulting criteria are solved approximately or require demanding optimization methods. The EM algorithm herein proposed combines the efficient image representation offered by the discrete wavelet transform (DWT) with the diagonalization of the convolution operator obtained in the Fourier domain. Thus, it is a general-purpose approach to wavelet-based image restoration with computational complexity comparable to that of standard wavelet denoising schemes or of frequency domain deconvolution methods. The algorithm alternates between an E-step based on the fast Fourier transform (FFT) and a DWT-based M-step, resulting in an efficient iterative process requiring O(N log N) operations per iteration. The convergence behavior of the algorithm is investigated, and it is shown that under mild conditions the algorithm converges to a globally optimal restoration. Moreover, our new approach performs competitively with, in some cases better than, the best existing methods in benchmark tests.

Wavelet based multiresolution expectation maximization image reconstruction algorithm for positron emission tomography

Maximum Likelihood (ML) estimation based Expectation Maximization (EM) [IEEE Trans Med Imag, MI-1 (2) (1982) 113] reconstruction algorithm has shown to provide good quality reconstruction for positron emission tomography (PET). Our previous work [IEEE Trans Med Imag, 7(4) (1988) 273; Proc IEEE EMBS Conf, 20(2/6) (1998) 759] introduced the multigrid (MG) and multiresolution (MR) concept for PET image reconstruction using EM. This work transforms the MGEM and MREM algorithm to a Wavelet based Multiresolution EM (WMREM) algorithm by extending the concept of switching resolutions in both image and data spaces. The MR data space is generated by performing a 2D-wavelet transform on the acquired tube data that is used to reconstruct images at different spatial resolutions. Wavelet transform is used for MR reconstruction as well as adapted in the criterion for switching resolution levels. The advantage of the wavelet transform is that it provides very good frequency and spatial (time) localization and allows the use of these coarse resolution data spaces in the EM estimation process. The MR algorithm recovers low-frequency components of the reconstructed image at coarser resolutions in fewer iterations, reducing the number of iterations required at finer resolution to recover high-frequency components. This paper also presents the design of customized biorthogonal wavelet filters using the lifting method that are used for data decomposition and image reconstruction and compares them to other commonly known wavelets.

Wavelet domain image restoration with adaptive edge-preserving regularization

In this paper, we consider a wavelet based edge-preserving regularization scheme for use in linear image restoration problems. Our efforts build on a collection of mathematical results indicating that wavelets are especially useful for representing functions that contain discontinuities (i.e., edges in two dimensions or jumps in one dimension). We interpret the resulting theory in a statistical signal processing framework and obtain a highly flexible framework for adapting the degree of regularization to the local structure of the underlying image. In particular, we are able to adapt quite easily to scale-varying and orientation-varying features in the image while simultaneously retaining the edge preservation properties of the regularizer. We demonstrate a half-quadratic algorithm for obtaining the restorations from observed data.

Adaptive regularized constrained least squares image restoration

In noisy environments, a constrained least-squares (CLS) approach is presented to restore images blurred by a Gaussian impulse response, where instead of choosing a global regularization parameter, each point in the signal has its own associated regularization parameter. These parameters are found by constraining the weighted standard deviation of the wavelet transform coefficients on the finest scale of the inverse signal by a function r which is a local measure of the intensity variations around each point of the blurred and noisy observed signal. Border ringing in the inverse solution is proposed decreased by manipulating its wavelet transform coefficients on the finest scales close to the borders. If the noise in the inverse solution is significant, wavelet transform techniques are also applied to denoise the solution. Examples are given for images, and the results are shown to outperform the optimum constrained least-squares solution using a global regularization parameter, both visually and in the mean squared error sense.

A multiresolution restoration method for cardiac SPECT imaging

In this study we present a multiresolution based method for restoring cardiac SPECT projections. Original projections were decomposed into a set of sub-band frequency images by using analyzing functions localized in both the space and frequency domain. This representation allows a simple denoising and restoration procedure by discarding high-frequency channels and performing inversion only in low frequencies. The method was evaluated in bull's eye reconstructions of a realistic cardiac chest phantom with a custom-made liver insert and 99mTc liver-to-heart activity ratios (LHAR) of 0:1, 1.5:1, 2.5:1, and 3.5:1. The cardiac phantom in free air was used as the reference standard. Reconstructions were performed by filtered backprojection using (1) no correction; (2) restoration without attenuation correction; (3) attenuation correction without restoration; and (4) restoration and attenuation correction. The attenuation correction was carried out with the Chang's method for one iteration. Results were compared with those obtained using an optimized prereconstruction Metz filter. Quantitative analysis was performed by calculating the normalized chi-square measure and mean +/- s.d. of bull's eye counts. In reconstructions with high liver activity (LHAR > 2), attenuation correction without restoration severely distorted the polar maps due to the spill-over of liver activity into the inferior myocardial wall. Both restoration methods when combined with an attenuation correction compensated this artifact and yielded uniform polar maps similar to that of the standard reference. There was no visual or quantitative difference between the performance of Metz filtering and multiresolution restoration. However, the main advantage of the multiresolution method is that it states a more concise and straightforward approach to the restoration problem. Multiresolution based methods does not require information about the object image or optimization processes, such as in conventional nuclear medicine restoration filters. In addition, the method is easy to implement using DFT techniques and potentially can be extended to noniterative spatially shift-invariant restorations in SPECT.

Total least-squares reconstruction with wavelets for optical tomography

In a previous paper [Zhu et al., J. Opt. Soc. Am. A 14, 799 (1997)] an iterative algorithm for obtaining the total least-squares (TLS) solution of a linear system based on the Rayleigh quotient formulation was presented. Here we derive what to our knowledge are the first statistical properties of this solution. It is shown that the Rayleigh-quotient-form TLS (RQF-TLS) estimator is equivalent to the maximum-likelihood estimator when noise terms in both data and operator elements are independent and identically distributed Gaussian. A perturbation analysis of the RQF-TLS solution is derived, and from it the mean square error of the RQF-TLS solution is obtained in closed form, which is valid at small noise levels. We then present a wavelet-based multiresolution scheme for obtaining the TLS solution. This method was employed with a multigrid algorithm to solve the linear perturbation equation encountered in optical tomography. Results from numerical simulations show that this method requires substantially less computation than the previously reported one-grid TLS algorithm. The method also allows one to identify regions of interest quickly from a coarse-level reconstruction and restrict the reconstruction in the following fine resolutions to those regions. Finally, the method is less sensitive to noise than the one-grid TLS and multigrid least-squares algorithms.

A wavelet-based multiresolution regularized least squares reconstruction approach for optical tomography

In this paper, we present a wavelet-based multigrid approach to solve the perturbation equation encountered in optical tomography. With this scheme, the unknown image, the data, as well as the weight matrix are all represented by wavelet expansions, thus yielding a multiresolution representation of the original perturbation equation in the wavelet domain. This transformed equation is then solved using a multigrid scheme, by which an increasing portion of wavelet coefficients of the unknown image are solved in successive approximations. One can also quickly identify regions of interest (ROI's) from a coarse level reconstruction and restrict the reconstruction in the following fine resolutions to those regions. At each resolution level a regularized least squares solution is obtained using the conjugate gradient descent method. This approach has been applied to continuous wave data calculated based on the diffusion approximation of several two-dimensional (2-D) test media. Compared to a previously reported one grid algorithm, the multigrid method requires substantially shorter computation time under the same reconstruction quality criterion.