Hierarchical prior for Bayesian deconvolution of radioactive sources with Poisson statistics

The problem of deconvolution is a classical problem in Gamma imaging. Using the conventional statistical model introduced by Shepp and Verdi (1982), many authors have adopted a Bayesian approach to address this problem and to estimate the different numbers of Gamma rays emitted lambda = (lambda(s))(s is an element of T) at each pixel s. As prior on lambda to regularize the picture, one uses a Gauss Markov random Field. However, this approach has one major drawnback: the choice of the regularization parameter beta. If beta is too high, the discontinuities are oversmoothed, and if beta is too low, the picture is not regularized enough. In this paper, we introduce a new hierarchical prior model for lambda, in which beta is not constant over the picture. This hierarchical prior model uses a Markov random field to describe spatial variation of the logarithm of the smoothing parameter log beta = (log beta(s))(s is an element of T) in a second random field which describes the spatial variation in lambda. The coupled Markov random fields are used as prior distributions. Similar ideas have occurred in Aykroyd (1996), but our prior model is quite different. Our new hierarchical prior model is applied for the problem of deconvolution of radioactive sources in Gamma imaging. The estimation of lambda and beta is based on their joint posterior density, following a Bayesian framework. This estimation is performed using a new SAGE EM algorithm (Hero and Fessler, 1995), where the parameters lambda and beta are updated sequentially. Our new prior model is tested on synthetic and real data and compared to the conventional Gauss Markov random field prior model : our algorithm increases significantly the results obtained by using a classical prior model.