This repository code to reproduce the results in [1] which proposes the use of the Stochastic Primal-Dual Hybrid Gradient algorithm (SPDHG) [2] for PET imaging with non-smooth priors (see also [3]). The results are on real clinical PET data in span 1 from the Siemens Biograph mMR based at UCLH Macmillan Cancer Centre. SPDHG is a direct generalization of the popular Primal-Dual Hybrid Gradient algorithm (PDHG) also known as the Chambolle-Pock algorithm [4, 5].
[1] M. J. Ehrhardt, P. J. Markiewicz, and C.-B. Schönlieb, “Faster PET reconstruction with non-smooth priors by randomization and preconditioning,” Phys. Med. Biol., 2019. 10.1088/1361-6560/ab3d07
Example results for various priors on FDG data are shown below.
You need to install ODL and NiPET [6] to run the examples.
We recommend to use conda to manage the correct versions of the packages. If you want to use conda, you could start like
conda create --name pet python=2.7
source activate pet
Then get some packages that nipet NiPET need.
pip2 install numpy dicom
Download the PET data set available here, e.g.
wget https://zenodo.org/record/1472951/files/amyloidPET_FBP_TP0.zip path_to_data
Clone the GIT repository and install (use path_to_data/umap when asked for the attenuation maps)
git clone https://github.com/pjmark/NIPET.git path_to_nipet
cd path_to_nipet
pip2 install --no-binary :all: --verbose .
(most recent) ODL
Clone the GIT repository and install
git clone https://github.com/odlgroup/odl.git path_to_odl
cd path_to_odl
pip2 install -e .
The code might also be compatible with older versions of ODL which are slightly easier to install:
pip2 install odl
Assuming you have installed the dependencies as below, you can run various examples as outlined below.
The numerical examples in [1] can be reproduced with
- Maximum Likelihood (ML) example
- ML example which breaks OSEM
- Maximum A-Posteriori (MAP) example with (isotropic) total variation (TV) prior
- MAP example with anisotropic TV prior
- MAP example with directional TV prior (uses MRI)
- MAP example with total generalized variation (TGV) prior
and the figures be recreated with
[1] M. J. Ehrhardt, P. J. Markiewicz, and C.-B. Schönlieb, “Faster PET reconstruction with non-smooth priors by randomization and preconditioning,” Phys. Med. Biol., 2019. 10.1088/1361-6560/ab3d07
[2] A. Chambolle, M. J. Ehrhardt, P. Richtárik, and C.-B. Schönlieb, “Stochastic Primal-Dual Hybrid Gradient Algorithm with Arbitrary Sampling and Imaging Applications,” SIAM J. Optim., vol. 28, no. 4, pp. 2783–2808, 2018.
[3] M. J. Ehrhardt, P. J. Markiewicz, P. Richtárik, J. Schott, A. Chambolle, and C.-B. Schönlieb, “Faster PET Reconstruction with a Stochastic Primal-Dual Hybrid Gradient Method,” in Proceedings of SPIE, 2017, vol. 10394, pp. 1–12.
[4] A. Chambolle and T. Pock, “A First-Order Primal-Dual Algorithm for Convex Problems with Applications to Imaging,” J. Math. Imaging Vis., vol. 40, no. 1, pp. 120–145, Dec. 2011.
[5] T. Pock and A. Chambolle, “Diagonal Preconditioning for First Order Primal-Dual Algorithms in Convex Optimization,” in Proceedings of the IEEE International Conference on Computer Vision, 2011, pp. 1762–1769.
[6] P. J. Markiewicz et al., “NiftyPET: a High-throughput Software Platform for High Quantitative Accuracy and Precision PET Imaging and Analysis,” Neuroinformatics, vol. 16, no. 1, pp. 95–115, 2018.