Bayesian Tomographic Reconstruction Using Riemannian MCMC Stefano Pedemonte 12 , Ciprian Catana 1 , Koen Van Leemput 123 1 Athinoula A. Martinos Center for Biomedical Imaging, MGH/Harvard, MA, USA 2 Department of Computer Science and Department of Neuroscience and Biomedical Engineering, Aalto University, FI 3 Department of Applied Mathematics and Computer Science, Technical University of Denmark, DK Abstract. This paper describes the use of Monte Carlo sampling for tomographic image reconstruction. We describe an efficient sampling strategy, based on the Riemannian Manifold Markov Chain Monte Carlo algorithm, that exploits the peculiar structure of tomographic data, en- abling efficient sampling of the high-dimensional probability densities that arise in tomographic imaging. Experiments with positron emission tomography (PET) show that the method enables the quantification of the uncertainty associated with tomographic acquisitions and allows the use of arbitrary risk functions in the reconstruction process. 1 Introduction A tomographic imaging device produces an indirect measurement of a spatially- dependent quantity. In the case of positron emission tomography, the interaction of photons with the array of detection crystals provides information about the rate of nuclear decay. In the early days, the tomographic reconstruction prob- lem was addressed with a mathematical formulation arising in the context of integral geometry: the Radon transform. Such idealized formulation provides an exact inversion formula under the assumptions of noiseless measurements and infinitesimally small detection elements. Starting with Shepp and Vardi [1], the formulation of the reconstruction problem in the probabilistic framework has enabled the development of algorithms that account for the uncertainty asso- ciated to the measurements and employ accurate models of the characteristics of the imaging devices. While in the context of tomographic imaging the re- search community has focused on the development of algorithms to compute a best guess of the unknown parameters, under the maximum likelihood (ML) and maximum a posteriori (MAP) point estimation criteria, in the context of low-level vision, starting with Geman and Geman [2], there has been a grow- ing interest in the fully Bayesian approach to probabilistic reasoning, based on Markov Chain Monte Carlo (MCMC) sampling. In contrast to point-estimation criteria, which aim to maximize a density function, the fully Bayesian approach aims to characterize the entire posterior probability density induced by the imag- ing experiment. This enables the quantification of the uncertainty, as reported