Accelerating iterative coordinate descent using a stored system matrix

Abstract

Purpose The computational burden associated with model‐based iterative reconstruction (MBIR) is still a practical limitation. Iterative coordinate descent (ICD) is an optimization approach for MBIR that has sometimes been thought to be incompatible with modern computing architectures, especially graphics processing units (GPUs). The purpose of this work is to accelerate the previously released open‐source FreeCT_ICD to include GPU acceleration and to demonstrate computational performance with ICD that is comparable with simultaneous update approaches. Methods FreeCT_ICD uses a stored system matrix (SSM), which precalculates the forward projector in the form of a sparse matrix and then reconstructs on a rotating coordinate grid to exploit helical symmetry. In our GPU ICD implementation, we shuffle the sinogram memory ordering such that data access in the sinogram coalesce into fewer transactions. We also update NS voxels in the xy‐plane simultaneously to improve occupancy. Conventional ICD updates voxels sequentially (NS = 1). Using NS > 1 eliminates existing convergence guarantees. Convergence behavior in a clinical dataset was therefore studied empirically. Results On a pediatric dataset with sinogram size of 736 × 16 × 13860 reconstructed to a matrix size of 512 × 512 × 128, our code requires about 20 s per iteration on a single GPU compared to 2300 s per iteration for a 6‐core CPU using FreeCT_ICD. After 400 iterations, the proposed and reference codes converge within 2 HU RMS difference (RMSD). Using a wFBP initialization, convergence within 10 HU RMSD is achieved within 4 min. Convergence is similar with NS values between 1 and 256, and NS = 16 was sufficient to achieve maximum performance. Divergence was not observed until NS > 1024. Conclusions With appropriate modifications, ICD may be able to achieve computational performance competitive with simultaneous update algorithms currently used for MBIR.

The Reconstruction of Magnetic Particle Imaging: Current Approaches Based on the System Matrix

“…Besides avoiding the high initial cost for decomposition based on direct solvers, a further advantage of an iterative algorithm like the Kaczmarz method is that it allows the incorporation of physical constraints, such as a real non-negative particle concentration c in the iteration procedure. Usually, the number of Kaczmarz iterations is chosen between 1 and 10 [ 53 ]. It converges rapidly due to the orthogonality of the system matrix, and results in a reasonable image quality [ 54 ].…”

Recent developments of the reconstruction in magnetic particle imaging

“…The Kaczmarz method is a row-action method that operates on rows independently. Convergence has been reported to occur within ten iterations of the Kaczmarz method [ 52 ]. Therefore, the state-of-the-art method for SM-based reconstruction is the Tikhonov regularization solved by the Kaczmarz method [ 53 ].…”

Accelerated 3D image reconstruction with a morphological pyramid and noise-power convergence criterion

“…Reduction of MBIR computational burden is an active area of research that encompasses a variety of methods, including hardware modifications tailored to MBIR (Fang and Mueller 2005, Jia et al 2011, Ravi et al 2019 and algorithmic approaches for convergence acceleration or reduction of per-iteration runtime (Ramani and Fessler 2012, Nien and Fessler 2016, Hsieh et al 2019. Recent work has explored approaches beyond conventional optimization to obtain the minimization of classic MBIR cost functions, such as penalized weighted least squares (PWLS), via combinations of conventional optimization methods and machine learning techniques.…”