## Abstract

### Background & purpose

### Materials & methods

### Results

### Conclusion

## 1. Introduction

## 2. Materials and methods

where $\mathcal{A}$ is the projection transform, $\mathbf{x}$ is the (4D) image to be reconstructed and $\mathbf{p}$ are the log-transformed projections. $\mathcal{W}$ is either the identity transform or (as in this work) a diagonal matrix containing weights [

in the standard 3D isotropic version. Here

*t*is the phase index, indicating that the total TV is the sum of the TV for the individual phases. We investigated three different choices for regularization.

### 2.1 Sparse Frequency Regularization (SFR)

$\mathcal{D}$ denotes a box filter, downsampling the image by a factor of two in the spatial dimensions. ${\mathcal{F}}_{t}$ is the Fourier transform along the temporal dimension. ${\lambda}_{\mathit{TV}}$, ${\lambda}_{\mathit{ATV}}$ and ${\lambda}_{\mathcal{F}}$ are unitless weights, controlling the strength of each regularisation parameter. The total variation term ensures that the images are sparse spatially, while the Fourier term enforces that the change in each pixel can be represented with a few Fourier coefficients. The use of total variation on multiple scales was previously proposed by Huang et al. [

### 2.2 4D total variation

where ${\mathit{TV}}_{t}(\mathbf{x})$ is the temporal total variation, given by

Again, ${\lambda}_{\mathit{TV}}$ and ${\lambda}_{4D}$ are unitless weights, controlling strength of the regularisation.

### 2.3 Prior image constrained compressed sensing (PICCS)

where ${\mathbf{x}}_{t}$ is the 3D image corresponding to phase

*t*of $\mathbf{x}$. The motivation behind PICCS is that most of the image is stationary, so the difference term ${\mathbf{x}}_{t}-{\mathbf{x}}_{\mathbf{p}}$ should mainly contain zero. While the original article used a constrained approach Lauzier et al. [

### 2.4 McKinnon-Bates

### 2.5 RTK ROOSTER

### 2.6 Patient scans

### 2.7 Extraction of breathing curves

### 2.8 Virtual XCAT phantom

### 2.9 Implementation

### 2.10 Selection of subsets

*N*subsets, and then perform a step of an iterative algorithm using only this subset. As the cost of an iteration is usually proportional to the size of the data, this provides a theoretical

*N*times speed up. Recently, ordered subsets have been combined with Nesterov’s method [

### 2.11 Forward and backward projection

*n*equally spaced points through the volume, using trilinear interpolation. For the backprojection, the projections belonging to each phase were loaded into a stacked 2D texture. A GPU thread was used per image voxel and backprojection was performed by finding the corresponding point on each projection and sampling with bilinear interpolation. This approach is similar to what is used in other GPU reconstruction frameworks such as RTK [

### 2.12 Parameters and evaluation

Eriksson D, Bindel D, Shoemaker C. Surrogate Optimization Toolbox (pySOT) [Internet]. Cited [2018 Mar 2]. Available from: github.com/dme65/pySOT.

## 3. Results

### 3.1 XCAT phantom

${\lambda}_{\mathit{TV}}$ | ${\lambda}_{\mathit{ATV}}$ | Temporal | Runtime | ${\mathit{SSIM}}_{\mathrm{min}}$ | |
---|---|---|---|---|---|

TV4D | $0.20$ | – | ${\lambda}_{\mathit{TV}4D}=0.62$ | $45.3\phantom{\rule{0.25em}{0ex}}\mathrm{s}$ | $0.912$ |

PICCS | $0.24$ | – | ${\lambda}_{\mathit{PICCS}}=0.13$ | $49.4\phantom{\rule{0.25em}{0ex}}\mathrm{s}$ | $0.903$ |

SFR | $0.12$ | $0.42$ | ${\lambda}_{\mathit{SFR}}=0.49$ | $56.8\phantom{\rule{0.25em}{0ex}}\mathrm{s}$ | $\mathbf{0}.\mathbf{916}$ |

FDK | – | – | – | – | $0.858$ |

MKB | – | – | – | – | $0.786$ |

ROOSTER (10 iterations) | – | – | – | $140\phantom{\rule{0.25em}{0ex}}\mathrm{s}$ | $0.857$ |

ROOSTER | – | – | – | $1344\phantom{\rule{0.25em}{0ex}}\mathrm{s}$ | $0.915$ |

### 3.2 Patients

## 4. Discussion

## Conflict of interest

## Acknowledgements

## Appendix – Algorithm

Algorithm 1. Momentum accelerated ordered subsets method with augmented Lagrangian multipliers for penalized weighted least squares (PWLS). |
---|

1 |

2 function OS-MOM(f) |

3 ${\mathbf{x}}^{0}\leftarrow \mathbf{0}$ |

4 $\mathbf{g}\leftarrow {\mathcal{A}}^{\u2020}{\mathcal{W}}^{\u2020}\mathcal{WA}\mathbf{1}$ ▷SQS diagonal preconditioner |

5 ${t}^{0}\leftarrow 1$ |

6 for $k\leftarrow 1$ to iterations do |

7 for $l\leftarrow 1$ to subsets do ▷Loop over subsets |

8 $m\leftarrow l+k\times \mathrm{subsets}$ |

9 ${\mathbf{y}}^{m\ast}\leftarrow \mathrm{subsets}\xb7{\mathcal{A}}_{\mathbf{l}}^{\u2020}({\mathbf{f}}_{\mathbf{l}}-{\mathcal{A}}_{\mathbf{l}}{\mathbf{x}}^{m-1})\xf8\phantom{\rule{0.25em}{0ex}}\mathbf{g}$ ▷ $\xf8$: elementwise-division |

10 ${\mathbf{y}}^{m}$ denoise(${\mathbf{y}}^{m\ast}\text{,}\mathbf{g}$) ▷Denoise using algorithm 2 |

11 ${t}^{m}\leftarrow \frac{1}{2}\left(1+\sqrt{1+4{t}^{m-1}\xb7{t}^{m-1}}\right)$ |

12 ${\mathbf{x}}^{m}\leftarrow {\mathbf{y}}^{m}+\frac{{t}^{m-1}-1}{{t}^{m}}\left({\mathbf{y}}^{m}-{\mathbf{y}}^{m-1}\right)$ ▷ Momentum acceleration step |

13 return ${\mathbf{x}}^{\mathrm{final}}$ |

Algorithm 2. Arrow Hurwicz denoising 28 , 38 . |
---|

1 function denoise$({\mathbf{x}}^{0}\text{,}\mathbf{g})$ |

2 ${\tau}^{1}\leftarrow {10}^{-3}$ |

3 $\gamma \leftarrow 0.35$ |

4 ${\sigma}^{1}\leftarrow \frac{1}{16{\tau}^{1}}$ |

5 for $k\leftarrow 1$ to iterations do |

6 ${\mathbf{p}}^{k}\leftarrow {y}^{k-1}+{\sigma}^{k}\mathcal{K}{\overline{\mathbf{x}}}^{k-1}$ ▷ $\mathcal{K}$ calculates the finite difference |

7 ${y}_{i}^{k}\leftarrow \frac{\frac{{\mathbf{p}}_{i}^{k}}{1+{\sigma}^{k}\alpha}}{\mathrm{max}\left(1\text{,}\left|\frac{{\mathbf{p}}_{i}^{k}}{1+{\sigma}^{k}\alpha}\right|\right)}$ |

8 ${\mathbf{x}}^{k\ast}\leftarrow {\mathbf{x}}^{k-1}-{\tau}^{k}{\mathcal{K}}^{\u2020}{\mathbf{y}}^{k}+\frac{{\tau}^{k}}{\lambda}{\mathbf{x}}^{0}\phantom{\rule{0.25em}{0ex}}\xf8\phantom{\rule{0.25em}{0ex}}\mathbf{g}$ |

9 ${\mathbf{x}}^{k}\leftarrow {\mathbf{x}}^{k\ast}\phantom{\rule{0.25em}{0ex}}\xf8\phantom{\rule{0.25em}{0ex}}\mathbf{g}$ |

10 ${\theta}^{k}\leftarrow \frac{1}{\sqrt{1+2\gamma \tau}}$ |

11 ${\tau}^{k+1}\leftarrow {\theta}^{k}{\tau}^{k}$ |

12 ${\sigma}^{k+1}\leftarrow \frac{{\sigma}^{k}}{{\theta}^{k}}$ |

13 return $\mathcal{P}({\mathbf{x}}^{\mathrm{final}})$ ▷ $\mathcal{P}$ sets all negative elements to $0$ |

## Supplementary data

- Supplementary data 1

- Supplementary data 2

## References

- Accuracy and inter-observer variability of 3D versus 4D cone-beam CT based image-guidance in SBRT for lung tumors.
*Radiat Oncol.*2012; 7: 81 - Difference in performance between 3D and 4D CBCT for lung imaging: a dose and image quality analysis.
*J Appl Clin Med Phys.*2016; 17: 97-106 - Iterative image reconstruction for CBCT using edge-preserving prior.
*Med Phys.*2008; 36: 252-260 - Performance studies of four-dimensional cone beam computed tomography.
*Phys Med Biol.*2011; 56: 6709-6721 - Registration-based reconstruction of four-dimensional cone beam computed tomography.
*IEEE Trans Med Imaging.*2013; 32: 2064-2077 - Estimating 4D-CBCT from prior information and extremely limited angle projections using structural PCA and weighted free-form deformation for lung radiotherapy.
*Med Phys.*2017; 44: 1089-1104 - A technique for estimating 4D-CBCT using prior knowledge and limited-angle projections.
*Med Phys.*2013; 40: 121701 - Motion-aware temporal regularization for improved 4D cone-beam computed tomography.
*Phys Med Biol.*2016; 61: 6856-6877 - 5D respiratory motion model based image reconstruction algorithm for 4D cone-beam computed tomography.
*Inverse Prob.*2015; 31: 115007 - High-quality four-dimensional cone-beam CT by deforming prior images.
*Phys Med Biol.*2013; 58: 231-246 - Simultaneous motion estimation and image reconstruction (SMEIR) for 4D cone-beam CT.
*Med Phys.*2013; 40: 101912 - 4D cone-beam CT reconstruction using multi-organ meshes for sliding motion modeling.
*Phys Med Biol.*2016; 61: 996-1020 - Optimal short scan convolution reconstruction for fan beam CT.
*Med Phys.*1982; 9: 254-257 Huang Y, Taubmann O, Huang X, Haase V, Lauritsch G, Maier A. A new scale space total variation algorithm for limited angle tomography. In: CT-Meeting 2016 proceedings (the 4th international meeting on image formation in X-ray computed tomography); 2016. p. 149–52.

- 4D cone beam CT via spatiotemporal tensor framelet.
*Med Phys.*2012; 39: 6943-6946 - Prior image constrained compressed sensing (PICCS): a method to accurately reconstruct dynamic CT images from highly undersampled projection data sets.
*Med Phys.*2008; 35: 660-663 - Practical cone-beam algorithm.
*J Opt Soc Am A.*1984; 1: 612 - Prior image constrained compressed sensing: implementation and performance evaluation.
*Med Phys.*2011; 39: 66-80 - Towards imaging the beating heart usefully with a conventional CT scanner.
*IEEE Trans Bio-med Eng.*1981; 28: 123-127 - The Reconstruction Toolkit (RTK), an open-source cone-beam CT reconstruction toolkit based on the Insight Toolkit (ITK).
*J Phys Conf Ser.*2014; 489 (012079) Mory C, Jacques L. A modified 4D ROOSTER method using the Chambolle–Pock algorithm. In: Proc. 3rd intl. conf. on image formation in X-ray CT; 2014. p. 191–3.

- Clinical use of iterative 4D-cone beam computed tomography reconstructions to investigate respiratory tumor motion in Lung Cancer patients.
*Acta Oncol.*2014; 53: 1107-1113 - Image matching as a diffusion process: an analogy with Maxwell’s demons.
*Med Image Anal.*1998; 2: 243-260 Zijp L, Sonke JJ, van Herk M. Extraction of the respiratory signal from sequential thorax cone-beam X-ray images. In: International conference on the use of computers in radiation therapy; 2004. p. 507–9.

- 4D XCAT phantom for multimodality imaging research.
*Med Phys.*2010; 37: 4902 - Fast reconstruction of low dose proton CT by sinogram interpolation.
*Phys Med Biol.*2016; 61: 5868-5882 - Combining ordered subsets and momentum for accelerated X-ray CT image reconstruction.
*IEEE Trans Med Imaging.*2015; 34: 167-178 - Studies in linear and non-linear programming.in: Stanford Mathematical Studies in the Social Sciences. vol. II. Stanford University Press, Stanford, California1958
- Accelerated fast iterative shrinkage thresholding algorithms for sparsity-regularized cone-beam CT image reconstruction.
*Med Phys.*2016; 43: 1849-1872 - Introductory lectures on convex optimization: a basic course.Springer US, New York2004
- Image quality assessment: from error visibility to structural similarity.
*IEEE Trans Image Process.*2004; 13: 600-612 - Mixture surrogate models based on Dempster-Shafer theory for global optimization problems.
*J Global Optim.*2011; 51: 79-104 Eriksson D, Bindel D, Shoemaker C. Surrogate Optimization Toolbox (pySOT) [Internet]. Cited [2018 Mar 2]. Available from: github.com/dme65/pySOT.

- An investigation of 4D cone-beam CT algorithms for slowly rotating scanners.
*Med Phys.*2010; 37: 5044 - A novel technique for markerless, self-sorted 4D-CBCT: feasibility study.
*Med Phys.*2012; 39: 1442 - Image quality in thoracic 4D cone-beam CT: a sensitivity analysis of respiratory signal, binning method, reconstruction algorithm, and projection angular spacing.
*Med Phys.*2014; 41 (041912) - A hybrid reconstruction algorithm for fast and accurate 4D cone-beam CT imaging.
*Med Phys.*2014; 41 (071903) - A first-order primal-dual algorithm for convex problems with applications to imaging.
*J Math Imaging Vision.*2011; 40: 120-145

## Article info

### Publication history

### Identification

### Copyright

### User license

Creative Commons Attribution – NonCommercial – NoDerivs (CC BY-NC-ND 4.0) |## Permitted

### For non-commercial purposes:

- Read, print & download
- Redistribute or republish the final article
- Text & data mine
- Translate the article (private use only, not for distribution)
- Reuse portions or extracts from the article in other works

## Not Permitted

- Sell or re-use for commercial purposes
- Distribute translations or adaptations of the article

Elsevier's open access license policy