Introduction
Computed tomography (CT) imaging is influenced by significant image artefacts induced by patient motion particularly in radiation therapy.Reference Keall, Mageras and Balter1 Patients are usually simulated with immobilisation devices, such as helical computed tomography (HCT) or axial computed tomography (ACT), which are used afterwards to outline the tumour and critical structures in treatment planning. CT image quality determines the accuracy of the volumes of these structures, and the accuracy of the CT numbers is invalidated by motion artefacts that affect the accuracy of the dose calculation algorithms, considering tissue heterogeneity. Patient motion leads generally to spread out of the CT number distributions of the CT images from the different imaging modalities.Reference Ali, Jackson, Alsbou and Ahmad2 Different approaches are used to manage patient motion at different stages of treatment with radiation therapy.Reference Keall, Mageras and Balter1 Breath holding and four-dimensional CT (4D-CT) techniques are used during simulation to reduce image artefacts due to patient motion.Reference Bernatowicz, Keall, Mishra, Knopf, Lomax and Kipritidis3 In 4D-CT, the projections are acquired at different breathing phases and several image sets are reconstructed with reduced motion artefacts.Reference Vedam, Keall, Kini, Mostafavi, Shukla and Mohan4 The 4D-CT images are used to generate an internal target volume (ITV) using selected respiratory phases that include a margin around the clinical target volume (CTV) that accounts for patient motion. Beam gating is used in dose delivery where dose is delivered to the patient within the breathing window coinciding with the phases that are used for the ITV by tracking an infrared marker.Reference Giraud, Morvan and Claude5 Other systems track a fiducial marker that is implanted surgically inside the tumour for localisation of the tumour and beam gating such as the CyberKnife system.Reference Gibbs and Loo6 Other studies investigated the use of the multi-leaf collimators to track tumour motion and beam gating.Reference Fast, Kamerling, Ziegenhein, Menten, Bedford, Nill and Oelfke7
In recent years, the applications of deformable image registration (DIR) have been increasing in radiation therapy.Reference Yang, Brame and El Naqa8 This involves correction of variations in patient anatomy due to variations in tumour volume over the course of treatment with radiation therapy.Reference Yan, Vicini, Wong and Martinez9 The changes in the shape and filling of organs at risk such as stomach, bladder or rectum affect dose sparing of these structures. Motion management particularly breathing motion is one of the main challenges to perform adaptive radiation therapy.Reference Yan, Vicini, Wong and Martinez9, Reference Wu, Li, Wu and Yin10 The goal of this study is to investigate quantitatively the displacement vector fields (DVFs) obtained from different DIR algorithms in HCT, ACT and cone beam CT (CBCT) of a mobile phantom and their correlation with motion parameters using controlled cyclic motion patterns. Four DIR algorithms including DemonsReference Thirion11–Reference Pennec, Cachier, Ayache, Taylor and Colchester13, Fast-DemonsReference Wang, Dong and O’Daniel12, Reference Yang, Li, Low, Deasy and El Naqa14, Horn–SchunckReference Horn and Schunck15, Reference Meinhardt-Llopis, Perez and Kondermann16 and Lucas–KanadeReference Bruhn, Weickert and Schnörr17, Reference Lucas and Kanade18 from the DIRART software were selected to register CT images of this mobile phantom.Reference Yang, Brame and El Naqa8 This study investigates the use of DVF to extract respiratory motion parameters and internal motion margins based on actual tumour motion.
Materials and Methods
Phantom setup and imaging
A thorax phantom mounted on a mobile platform (Standard Imaging, Inc., Middleton, WI, USA) was imaged using HCT, ACT and CBCT imaging techniques. The thorax phantom included three targets that are manufactured from water-equivalent material and embedded in low-density foam to simulate lung lesions. The targets had well-known sizes and lengths (10, 20 and 40 mm) along the direction of motion (y axis) as shown in Figure 1. The phantom was moved with controlled motion patterns using a range of motion amplitudes (0–20 mm) and frequencies (0·125–0·5 Hz). A thorax imaging technique was used to acquire the HCT images with a CT simulator (GE Discovery-CT-590RT, General Electric Healthcare, Milwaukee, WI, USA) with the following parameters: pitch of 1·375, 2·5 mm slice thickness, 120 kVp and 440 mA. Similar scanning parameters were used in ACT imaging. The CBCT images were acquired with an on-board imager mounted on a TrueBeam STx radiotherapy machine (Varian Medical Systems, Inc., Palo Alto, CA, USA). The imaging parameters used for CBCT images were 2 mm slice thickness, 125 kVp and 264 mA.
(a) Experimental setup of the thorax phantom with motion platform. (b) Coronal view of the thorax phantom with the three mobile targets embedded in the lung tissue.
Deformable image registration algorithms
In this study, four DIR algorithms from the DIRART8 research software developed at the University of Washington–Saint Louis were used: DemonsReference Thirion11, Reference Pennec, Cachier, Ayache, Taylor and Colchester13, Fast-DemonsReference Wang, Dong and O’Daniel12, Reference Yang, Li, Low, Deasy and El Naqa14, Horn–SchunckReference Horn and Schunck15, Reference Meinhardt-Llopis, Perez and Kondermann16 and Lucas–Kanade.Reference Bruhn, Weickert and Schnörr17, Reference Lucas and Kanade18 The DIRART software provides different algorithms that are coded with MATLAB (MathWorks, Inc., Natick, MA, USA) and other software tools to perform adaptive radiation therapy.10 The CT images of the mobile phantom were registered to the corresponding CT images of the stationary phantom, which was considered as the reference image set. These selected algorithms perform autonomous intensity-based DIR. The Horn–Schunck algorithm solves the optical flow equation where the variations in the intensity of an image are represented by a distribution of velocities from the movement of different voxels in the image.Reference Meinhardt-Llopis, Perez and Kondermann16 It uses a global method where a constraint optical flow equation is convolved with a Gaussian low-pass filter. In contrast, the Lucas–Kanade algorithm uses small motions of the different voxels in the images and the image intensity remained constant, where a voxel shifts like its neighbours.Reference Bruhn, Weickert and Schnörr17, Reference Lucas and Kanade18 The Demons algorithm uses a diffusion model involving the physics principles of the Maxwell’s demons in fluids.Reference Thirion11, Reference Pennec, Cachier, Ayache, Taylor and Colchester13 The Fast-Demons algorithm extends the Demons by considering Newton’s third law of motion in the diffusing model that allows both images to diffuse in each other in contrast to the demons forces that diffuse only the moving image through the static image.Reference Wang, Dong and O’Daniel12, Reference Yang, Li, Low, Deasy and El Naqa14
Results
CBCT images
Figure 2 shows that the maximal DVF increased with increasing motion amplitude of the mobile targets, which represented the displacements needed to shift the superior side of mobile target to match with the superior side of the stationary target. In addition, the minimal DVFs decreased with increasing motion amplitude that represented the displacements required to shift the inferior side of mobile targets to match the stationary target. The DVFs obtained from the different DIR algorithms displayed similar patterns; however, different displacements were calculated by the various algorithms. Large DVFs for image voxels were obtained only at the interference of the water-equivalent mobile targets with the low-density material surrounding them in the thorax phantom. Otherwise, small DVFs were obtained for image voxels surrounded by homogenous material.
CT number profiles of CT-numbers along the direction of motion for small and medium targets in (a); and the large target in (b) for the indicated motion amplitude in CBCT images. The DVF profiles along the direction of motion for the small and medium targets in the first column and large target in the second column calculated by the Demons (c, d), Fast Demons (e, f), Horn–Schunck (g, h) and Lucas–Kanade (i, j) deformable image registration algorithms for CBCT images. Abbreviations: CT, computed tomography; CBCT, cone beam computed tomography; DVF, displacement vector field.
Figures (3a and 3b) show that the maximal DVF increased linearly with motion amplitude for the small and medium targets in the first column and the large target in the second column. The minimal DVF decreased linearly with motion amplitude as shown in Figures (3c and 3d). The minimal DVF calculated by the Fast-Demons algorithm did not decrease linearly in a consistent way with increasing motion amplitude particularly with the small and medium targets with high-gradient variation in the CT number level and stronger artefacts induced by motion as shown in Figure 3a. Table 1 shows significant correlation between the maximal and minimal DVFs and the motion amplitudes with the Pearson correlation coefficient (<0·05) in CBCT for most of the DIR algorithms. However, the p value of the Pearson coefficient for minimal DVFs was >0·05 only for the Fast-Demons algorithms, which generally performed poorer than other algorithms. The optical flow-based algorithms such as Horn–Schunck and Lucas–Kanade produced consistent DVF distributions where minimal and maximal DVF changed linearly with the motion amplitude. However, the DVFs from the Demons and Fast-Demons algorithms, which uses attraction forces, were less consistent with their correlation with the motion amplitudes compared to the optical flow algorithms.
Maximal DVF versus motion amplitude for the small and meduim targets (a) and the large target (b) imaged with CBCT. Minimal DVF versus motion amplitude for the small and meduim targets (c) and the large target (d) from CBCT images. Abbreviations: CBCT, cone beam computed tomography; DVF, displacement vector field.
p Values from the Pearson correlation coefficeint of the maximal and minimal DVF with the motion amplitude for the different mobile targets: DVF12 for the small and medium targets, and the large target DVF3 for the Demons, Fast-Demons, Horn–Schunck and Lucas–Kanade DIR algorithms using CBCT, HCT and ACT
Abbreviations: DVF, displacement vector field; ACT, axial computed tomography; HCT, helical computed tomography; CBCT, cone beam computed tomography.
HCT images
Figure 4 shows that maximal DVFs calculated by the different DIR algorithms increased linearly with increasing motion amplitude in HCT images, while the minimal DVFs decreased linearly with increasing motion amplitude. These maximal and minimal DVFs represented the largest shifts needed to match both sides of a mobile target with that of the stationary target. The maximal and mininal DVFs in HCT were generaly smaller than the corresponding DVFs obtained from CBCT because of less image artefacts induced by motion. Most of the DIR algorithms calculated nearly similar maximal and minimal DVFs from the deformation of the HCT images for a certain motion amplitude. Similar DVFs were calculated for the small, meduim and large targets by the different DIR algorithms in HCT in contrast to CBCT where different DVFs were obtained from the different algorithms. The Pearson corrleation coefficients were <0·05 for maximal and minimal DVFs obtained for all the mobile targets using different DIR algorithms in HCT as listed in Table 1. However, the p values in HCT were larger than the corresponding values in CBCT. The Horn–Schunck and Lucas–Kanade correlations were superior to the Demons and Fast-Demons correlations in HCT.
Maximal DVF versus motion amplitude for the small and meduim targets (a); and the large target (b) imaged with HCT. Minimal DVF versus motion amplitude for the small and meduim targets (c); and the large target (d) from HCT images. Abbreviations: DVF, displacement vector field; HCT, helical computed tomography.
ACT images
Figure 5 shows that the maximal and minimal DVFs obtained from the different DIR algorithms increased or decreased gradually with motion amplitude; however, they did not produce a linear relationship as seen in HCT and CBCT. In ACT images, the profiles did not show a regular pattern of increasing DVFs with motion amplitude. The DVFs calculated by the different DIR algorithms were generally smaller than those calculated in HCT and CBCT. The image artefacts induced by motion were the strongest in ACT compared with other modalities; and thus, the different DIR algorithms did not deform well the CT images of the mobile targets to reproduce the size, shape and CT values in the ACT images of the stationary targets. In ACT, the Pearson correlation coefficients of the maximal and minimal DVFs were >0·05 for certain DIR algorithms as listed in Table 1. The p values in ACT were generally larger than the corresponding values in CBCT and HCT images, indicating poor correlation.
(a) Maximal DVF versus motion amplitude for the small and meduim targets imaged with ACT, (b) maximal DVF versus motion amplitude for the large target imaged with ACT, (c) minimal DVF versus motion amplitude for the small and meduim targets in (c) and large target in (d) from the different DIR algoritms: Demons, Fast-Demons, Horn–Schunk (HS) and Lucas–Kanade (LK). Abbreviations: DVF, displacement vector field; ACT, axial computed tomography.
Discussion
In CBCT and HCT, the maximal and minimal DVFs calculated by the different DIR algorithms changed linearly with motion amplitude; however, in ACT, the corresponding DVF did not correlate well with the motion amplitude as represented by the Pearson correlation algorithm. This was due to the patterns of motion artefacts induced in the different CT imaging modalities. In CBCT, image artefacts induced by motion led mainly to elongation of the mobile targets along the direction of motion. In HCT, the motion artefacts were mainly elongation, shrinkage or spatial shifting in the image of the mobile targets that were handled easily by the different DIR algorithms. In ACT, the motion artefacts were not regular as in CBCT and HCT where the image of the mobile targets split and the CT number distributions spread-out irregularly. Most of the DIR algorithms were not able to reproduce the shape and volume of the stationary targets in ACT. The maximal and minimal DVFs can be used to determine quantitatively the margins to be added to the CTV to obtain the ITV from CBCT and HCT images, which is given as follows:
$${\rm{ITV}} = {\rm{CTV}} + F{\,^*}{\rm{DV}}{{\rm{F}}_{\max }}\left( {x,y,z} \right)$$
where DVFmax is the maximal DVF, F is a correction factor that can be obtained from the correlation curves of the DVF with motion amplitude as shown in Figures 3 and 4 associated with each DIR algorithm for CBCT and HCT, respectively. This approach can be used as an alternative to the ITV obtained from 4D-CT imaging where the images are sorted in different respiratory phases using an external marker.
The image artefacts induced by cyclic motion were mainly dependent on motion amplitude in CBCT. The dependence on frequency and phase varied with the different imaging modalities, considering the time needed to acquire the projections during CT scanning. The scanning time in CBCT covered nearly 6–12 respiratory cycles, where the image artefacts were mainly dependent on motion amplitude. HCT took the shortest scanning time that caused the least image artefacts induced by motion. The individual targets were scanned in a fraction of a respiratory cycle in HCT with the continuous motion of the imaging couch. Thus, the image artefacts in HCT depended on the motion frequency and phase. In ACT, scanning the mobile targets takes a longer time compared with scanning by HCT where partial volumes were imaged as the couch was indexed at different motion phases. Thus, the motion artefacts in ACT were complicated with their dependence on the different motion parameters that included amplitude, frequency and phase. The correlation between the DVF and motion amplitude particularly provides a useful parameter that can be used to manage patient motion or correct image artefacts in HCT and CBCT images based on internal tumour motion instead of external marker motion attached to the patient chest wall using 4D-CT. The motion amplitude of a certain target is given as follows:
$$A = C\,^*DV{F_{\max }}$$
where A is the motion amplitude, DVFmax is the maximal DVF along the direction of motion, C is a calibration factor that can be obtained from the correlation curve of DVF with motion amplitude for a certain DIR algorithm. More information about motion frequency and phase are needed to reconstruct the internal motion trajectory of the tumour that might be extracted from the DVF distributions in future studies.
Conclusions
This study demonstrated that the DVF correlated well with the motion amplitude of the mobile phantom in CBCT and HCT images. However, in ACT, the DVF did not correlate well with motion amplitudes. The correlation of the DVF with motion amplitude provides valuable information about unknown motion parameters of the mobile organs in real patients particularly using CBCT and HCT images. The maximal or minimal DVF can be used to extract motion parameters such as amplitude, which can be used for the management of patient motion as an alternative technique of the regular 4D-CT and beam gating based on the motion of an external marker. The DVF can be used to extract the internal margins needed for the ITV based on actual internal tumour motion that can be used in treatment planning.
Acknowledgements
None.
