Review and recommendations on deformable image registration uncertainties for radiotherapy applications

As described in section 2.1, classical algorithms approach the problem of image alignment through a process of iterative optimization. These algorithms search for a global minimum of the solution space, but the choice of similarity metric remains problematic. DSMs aim to improve classical iterative image registration by improving the similarity term. This approach is particularly useful in multi-modal imaging where it has been shown to outperform mutual information (Wu et al 2013, Simonovsky et al 2016). Improvements in difficult monomodal registration problems, low contrast regions and large transformations, have been reported in the literature (Zhao and Jia 2015).

Direct DVF DL algorithms use historic DVFs or artificial DVFs as training data to determine registrations. The optimization phase happens in the training phase, where model parameters are determined. The vast majority of DL models aim for a direct regression of DVF transforms in a supervised approach. Variation between models is primarily a result of algorithm design and methodology.

Reviews (Boveiri et al 2020) cover a range of algorithm architectures. DL architectures include staked auto-encoders (SAEs) (Wang et al 2017, Krebs et al 2018), bayesian frameworks (Deshpande and Bhatt 2019, Khawaled and Freiman 2020, 2022a), implicit neural representations (Wolterink et al 2022) and convolution neural networks (CNNs) (Cao et al 2018, Ferrante et al 2018, Hu et al 2018, Balakrishnan et al 2018, 2019, Kim et al 2019, Kuang and Schmah 2019, Liu et al 2019, Jian et al 2022, Wolterink et al 2022, Xi et al 2022, Liang et al 2023). CNNs have been researched for direct DVF regression, with reported improvements in DVF when coupled with spatial transformer networks (Jaderberg et al 2015). CNN architecture use encoder-decoder networks, rather than a fully connected layer. Such approaches are currently implemented in well-cited solutions (VoxelMorph (Balakrishnan et al 2018, 2019) and U-NET (Liang et al 2023)). Despite the growth of multimodal foundational models in image creation, these reviews do not find application in image registration.

In general, DL training is divided between supervised and unsupervised learning methods (Chen et al 2021). For supervised registration methods, ground truth is either a DVF or a segmentation. The DVF may be created by a conventional DIR algorithm or from synthetic deformations, and the segmentations may be created by manual contouring or other methods. Unsupervised registration methods are further split into training by similarity metrics or generative adversarial networks (GANs) (Mahapatra et al 2018, Elmahdy et al 2019). If similarity metrics are used no ground truth is needed for the learning process but, as in traditional image registration, these models are limited by the same issues as similarity metrics in classical DIR optimization. If GAN is used, a discriminator judges if the warped moving image can be discriminated from the fixed image. When the warped image cannot be distinguished from the fixed image, the registration is deemed to be optimal (Goodfellow et al 2014). GANs show promise for multi-modality DIR problems as they do not require image similarity terms.

One advantage of DL algorithms is improvements in multi-modal registration, which is challenging for classical similarity metrics. Additionally, DL-based algorithms are more computationally efficient (Rohé et al 2017, Cao et al 2018, Balakrishnan et al 2018, 2019).

Non-rigid variations in patient anatomy, such as weight gain or loss, neck flexion and tumour changes can be poorly mapped by rigid and affine registrations. DIR can improve the locally accurate alignment of anatomy (Hill et al 2001). While regularisation is useful to reduce the likelihood of physically unrealistic deformations, the magnitude of anatomical changes may exceed those allowed by an algorithm's settings. This can result in large DIR errors in areas near significant shape changes, particularly in low contrast image regions (Kashani et al 2008) or due to forced anatomical changes such as between external beam RT and brachytherapy (Vásquez Osorio et al 2015) (figure 3(a)).

Zoom In Zoom Out Reset image size

Anatomical changes can be elastic, where the surrounding tissue follows the change and occupies the previous space (e.g. movement, position changes, or displacement) or inelastic, where the surrounding tissue stays in place (e.g. tissue growth, regression or emptying/filling cavities) (Amugongo et al 2022) figures 3(b) and (c). Modelling these changes is challenging (Sonke and Belderbos 2010, Mencarelli et al 2014, Sonke et al 2019). Certain implementations of regularisation can result in significant registration inaccuracies in sites in which naturally sliding boundaries occur, such as a rib bone and its adjacent lung (Sonke et al 2019). Some solutions were proposed to incorporate missing tissue during the DIR (Nithiananthan et al 2012, Vishnevskiy et al 2017, Eiben et al 2018).

The anatomical changes caused by natural patient motion, such as respiration, muscle contraction, and blood flow can lead to image artifacts (Nehmeh and Erdi 2008, Zhang et al 2012, Spin-Neto and Wenzel 2016, Giganti et al 2022). For example, motion artifacts during the image acquisition can result in implausible anatomy (Yamamoto et al 2008, Persson et al 2010) and implants such as prostheses in the imaging area can lead to streaking or voids (Ritter et al 2009, Fontenele et al 2018, Lee et al 2021). As these artifacts disrupt the true image intensity gradients of the patient tissue several papers have demonstrated decreased intensity-based DIR quality in their presence (Serban et al 2008, Sonke and Belderbos 2010, Fusella et al 2016) (figure 4).

Zoom In Zoom Out Reset image size

Sensitivity of DIR algorithms to image noise, resolution (Constable and Henkelman 1991, Verdun et al 2015, Zhao et al 2016, Sarrut et al 2017), field of view (Barber et al 2020) and image contrast (Mencarelli et al 2014, Barber et al 2020, Dowling and O'Connor 2020) has been demonstrated in the literature. However, other studies find that the effect of image noise has only minor effects on DIR results for computed tomography (CT) to CT registrations (Nesteruk et al 2022).

Research on the implementation of iterative image reconstruction algorithms has shown reduced noise and improved image quality for both CT and cone-beam CT (CBCT) (Held et al 2016, Giacometti et al 2019, Jarema and Aland 2019, Greffier et al 2020, Loi et al 2020), which may allow for improved quality intensity-based DIR.

Multimodal DIR offers considerable clinical benefit in contour propagation (Söhn et al 2008, Vásquez Osorio et al 2012, Barber et al 2020, Zachiu et al 2020). However, multimodal DIR remains challenging, and similarity metrics must be selected with care.

For example, magnetic resonance imaging (MRI) to CT registration in the lung is difficult because of low contrast and resolution in MRI (Yang et al 2015) and in the prostate, lack of a clear boundary of the prostate gland in CT may lead to failures in MR-CT DIR (Zhong et al 2015). In the HN, limited soft tissue contrast and dental artifacts in CT images compared to MR influence the DIR uncertainty (Nix et al 2017, Kiser et al 2019). Additionally, gradient contrast artifacts in MRI may impair the DIR quality between different image modalities (figure 5) (Vásquez Osorio et al 2012). McKenzie at el. found monomodal registration from synthetic CT (generated from the MRI) to CT to be more accurate than the multimodal registration from the original MRI to CT for large deformations of HN patients (McKenzie et al 2020). Of course, the synthetic CT generation also faces uncertainties, for example the resulting Hounsfield units (HUs) differ between CT and synthetic CT. Boulanger et al report a mean absolute error of 76 HU in head and liver, and 42 HU in the pelvic area in average over multiple methods generating synthetic CTs (Boulanger et al 2021). Geometric differences between structures can also appear (Palmér et al 2021).

Zoom In Zoom Out Reset image size

Geometric measures are not independent and self-correlate. Loi et al found a linear relationship between mean distance to agreement (MDA) and dice similarity coefficient (DSC) (Loi et al 2018). Also, a correlation between distance discordance metric (DDM) and Harmonic energy (HE) has been found (Kierkels et al 2018). Reporting multiple measures is still useful despite being redundant. For example, the DSC limitations can be critically analysed in conjunction with other metrics, such as MDA for different structures and volumes (Jena et al 2010, Brock et al 2017, Loi et al 2018). Combining different geometrical metrics can improve the understanding of the overall quality of the DIR for a specific application.

Different implementations and specific ways to use the same measure can lead to vastly different results. For example variations of up to 50% in DSC, 50% in Hausdorff distance (HD) and 200% in MDA were found between the same structure sets, evaluated by different institutions (Gooding et al 2022). Comparing results from different studies and centres should therefore be taken with care. The correlation between geometric and dosimetric measures was found to be low (Hvid et al 2016, Pukala et al 2016, Poel et al 2021, Nash et al 2022, Kamath et al 2023).

The TG-132 report and other recommendations suggest imaging the patient in the treatment position whenever possible to minimise the magnitude of the required deformation during registration (Brock et al 2017, Barber et al 2020).

DIR has been used to derive motion-corrected images from 4D CT scans (Wolthaus et al 2008) with average landmark-position differences of 0.5 mm for all directions in the tumour region. DIR is also used to reconstruct time-resolved 4D MRI (Nie et al 2020), with reported centre of mass differences of 2.9±0.6 mm. We expect the geometrical uncertainties of propagated images to be similar to those of structure propagation, considering both utilise the same input data.

With MRI linac or CBCT-based online adaptation becoming more commonly available, the interest in deforming images between fractions for dose calculation and optimization is increasing (Kraus et al 2017, Tenhunen et al 2018, Irmak et al 2020, Byrne et al 2021). In these workflows, the calculated dose distribution is unlikely to be accurate considering the spatial uncertainties in the deformed CT, especially in areas with large density changes.

To correct for density changes that are not represented by the deformed image such as moving air in the gastro-intestinal organs, the density in these areas is often overwritten with the density of air or water (van Timmeren et al 2020). Research investigating the impacts of these overwrites on photon RT has found these impacts to be not clinically relevant (Pham et al 2022), except for very large air cavities (Thapa et al 2019). For protons, these density corrections are likely more relevant.

To avoid the use of DIR and manual density overwrites, direct dose calculation on the MRI or CBCT images has been investigated. The generation of synthetic CT images from MRI is reviewed elsewhere (Owrangi et al 2018, Hoffmann et al 2020, Boulanger et al 2021). Methods of scatter correction to make CBCT usable for dose calculation are widely explored (Kurz et al 2016, Giacometti et al 2019, Jarema and Aland 2019, Lalonde et al 2020, Trapp et al 2022).

Follow-up images after intervention can be registered to a planning CT to understand the relation and location of local failure such as recurrence or necrosis with a planned dose distribution and planning structures (Chang et al 2018, Kamal et al 2020, Abdel-Aty et al 2022). In these cases, dramatic changes are observed caused by the time between images, surgical intervention, or other medical issues. Systematic studies quantifying the impact of DIR uncertainties for intervention follow-up are rare and further work needs to be done to quantify and account for them in patterns of failure.

For treatment planning, structures are commonly defined on the planning CT. Structure definition can be challenging on CT due to low contrast compared to other imaging modalities such as MRI. Including multiple imaging modalities for contouring can lead to a reduction in inter-observer variability (Caldwell et al 2001, Farina et al 2017, Hall et al 2018). Though it is common to merge CT with PET, MRI or other images for contouring, the effect of deformable registration errors is not well investigated. (Barber et al 2020) therefore suggested using rigid registration wherever possible.

DIR is also used in atlas-based auto-segmentation, which is increasingly used in clinics to assist contouring. In this case, DIR is applied between images from different patients (Vrtovec et al 2020). Research showing a time benefit in using altas-based contours also show the necessity of manual corrections (Gooding et al 2013, Cardenas et al 2019, Welgemoed et al 2023). To our knowledge, there are no systematic studies on the impact of DIR implementation and DIR uncertainty for atlas-based segmentation. Studies do however investigate the impact of atlas selection (Schipaanboord et al 2019) or institution-specific implementation (Gooding et al 2013). As an alternative to atlas-based auto-segmentation, DL-based auto-segmentation was developed. Different auto-segmentation methods are reviewed elsewhere (Cardenas et al 2019, Schipaanboord et al 2019, Vrtovec et al 2020, Harrison et al 2022). The details of auto segmentation methods are out of the scope of this study.

DIR is used to map contours between different breathing phases, or intrafraction changes in patient treatment positions, to reduce the time needed for contouring. In clinical practice, the propagated structures are visually checked and corrected if necessary (Gaede et al 2011, Peroni et al 2013, Liu et al 2016, Ma et al 2017, Willigenburg et al 2022).

Structure propagation can speed up recontouring for repeated imaging of a patient (Sonke et al 2019). This can be used for evaluation of recalculated doses or adaptive planning, but it is especially important for online adaptive workflows. The clinical availability of regular or daily imaging, such as scheduled repeated CT or daily CBCT, and the implementation of online adaptive workflows has led to multiple studies on the quality of deformed structures for adaptive planning.

Table 2 summarises recent studies investigating geometric DIR uncertainties for structure propagation for different anatomical sites and imaging modalities. The structures have been evaluated using geometrical measures introduced in table 1. Commercial and open access algorithms show similar performance (Scaggion et al 2020b). The majority of studies incorporate reference structures of a single expert physician. Variations between structures defined by different physicians are observed in many studies, and these inter-expert structure variations are currently de-facto the clinically accepted variability. Research comparing DIR propagated structure uncertainties to physician-to-physician uncertainties, has demonstrated results approaching inter-expert contour variation (Riegel et al 2016, Woerner et al 2017, Rigaud et al 2019, Nash et al 2022).

Currently, there is no consensus on the use of DIR propagated structures for plan adaptation in the literature. Some authors conclude that propagated structures can be used for reoptimization and/or dose evaluation (Beasley et al 2016, Hvid et al 2016, Qiao et al 2019, Nenoff et al 2021b, Nash et al 2022), while others found that manual corrections are still necessary (Li et al 2017, Christiansen et al 2020, 2021). Generally, the literature agrees that a visual inspection of the DIR propagated structures remains necessary for dose evaluation and optimization. Furthermore, it has been observed that for most organs at risk (OARs) geometric uncertainties correlated only weakly to dosimetric errors (Hvid et al 2016, Pukala et al 2016, Nash et al 2022).

There are a small number of studies evaluating the dosimetric effect of uncertainties in propagated structures for dosimetric evaluation or plan optimization during adaptive RT (table 3). For pancreas stereotactic body radiotherapy (SBRT), physician-drawn structures were compared to propagated structures by MIM and Precision DIR algorithms (Magallon-Baro et al 2022). They compared uncorrected propagated structures with physician-drawn structures in 0.5, 1 and 3 cm distance rings from the target. They found that replanning with uncorrected propagated structures improves the target coverage and OAR sparing compared to no adaptation. For the majority of fractions, manual correction of propagated structures could be avoided or be limited to the region closest to the target. Ray at al. evaluated the use of automatic deformed CTVs compared to physician defined CTVs and proposed a framework to determine PTV margins based on automatic deformed CTVs for adaptive planning (Ray et al 2020). Nash et al showed that even large geometrical structure differences rarely had a statistically significant impact on OAR dose-volume-histograms (DVH) parameters and concluded that DIR propagated structures are suitable for dose evaluation (Nash et al 2022). Similar conclusions were found by Hvid et al (2016).

Also for proton therapy the dosimetric impact of using propagated structures for proton dose evaluation and optimization has been investigated: Qiao et al and Elmahdy et al investigated prostate structures propagated from CT to CT with the open source DIRs in Elastix (Elmahdy et al 2019, Qiao et al 2019). They gave an extensive geometrical evaluation (included in table 2) and dosimetric evaluation (included in table 3) that showed that DIR propagated structures can be used for optimization in online-adaptive intensity-modulated proton therapy (IMPT). Similar conclusions were found for lung cancer patients by Nenoff et al, showing that daily IMPT optimization on CT based on propagated, uncorrected structures was better than no adaptation (Nenoff et al 2021b). Daily manual recontouring on each CT gives a small additional benefit for some patients and OARs. They also investigated if including the inter-algorithm variation between structures propagated with DIR in the adaptive IMPT optimization could improve the adapted plan against structure uncertainties (Nenoff et al 2022). They found that adaptation on propagated, uncorrected structures showed a benefit over no adaptation for MRI-to-MRI registrations for pancreas and liver patients and CT-to-CT registrations for HN patients. Only for the HN patients including structure propagation uncertainties in the optimization significantly improved the adapted plan. Recently, Smolders et al compared the effect of different auto-segmentation methods, among those DIR based structure propagation, for the dosimetric quality of online adaptive proton therapy plans. They found the dosimetric influence of using automatic contours for the optimization to be small for OARs and larger for targets, with DIR propagated structures performing best for both OARs and targets (Smolders et al 2023a).

A commonly proposed use of dose accumulation is for 4D treatment planning or the dose reconstruction of the dose in a moving area. Both, 4D optimisation (Graeff et al 2013, Engwall et al 2018, Spautz et al 2023) and 4D dose evaluation (Zhang et al 2019, Meijers et al 2020) require DIR and therefore show DIR-related uncertainties. Both are mostly applied in anatomical areas affected by breathing motion, registering all phases of a 4D CT or 4D MRI scan into a reference phase or average image (Rosu and Hugo 2012, Engwall et al 2018, Meijers et al 2020). Most 4D dose optimisation and dose calculation studies do not investigate DIR-facilitated dosimetric uncertainties (table 3). Those who do, report different metrics between different studies, to quantify these uncertainties. For example, (Ribeiro et al 2018) found differences in the target V95% of up to 11.34% for 4D dose accumulation of liver cancer proton therapy. In contrast, (Sarudis et al 2019) found only dose deviations of ≤3.0% between different visually acceptable DIRs in 4D lung volumetric modulated arc therapy (VMAT) dose accumulations. Mogadas et al tested five open-source registration algorithms on lung and liver SBRT, using the delta D95% of the target using 4D dose reconstruction compared to the static plan. For lung metastases, accumulated dose distributions were similar regardless of the DIR algorithm. In contrast, for liver metastases, accumulated dose distributions strongly varied, due to large DIR uncertainties in low contrast regions (Mogadas et al 2018).

Another DIR-facilitated application has been to map doses re-calculated (or re-optimized) on 3D images from different fractions on the planning CT, to get an estimation of the total delivered treatment dose (Chetty and Rosu-Bubulac 2019, Ziegler et al 2019, Nenoff et al 2020). This technique has been extended to evaluate the validity of treatment plans with reduced margins (Wu et al 2009, van Kranen et al 2016, Lowther et al 2020b, van der Bijl et al 2022). The reporting of the uncertainties is very application-dependent and not standardised, which makes direct comparisons challenging (table 3). Examples of the magnitude of dose mapping uncertainties for interfractional changes include (Nenoff et al 2020), reporting differences caused by DIR uncertainties of 8.7% in the CTV of accumulated proton doses and (Wang et al 2018) reporting a maximum dose variation of 3.3 Gy for hollow organs in the abdomen for interfraction dose mapping. More recently, (Huesa-Berral et al 2022) reported a dosimetric uncertainty between fractions below 2 Gy in tumour and OAR in lung SBRT. This study also concluded that inter-fraction variations dominated and that dose accumulation for these patients should prioritise day-to-day changes over respiratory motion.

There have been several propositions to also predict uncertainties on the geometrical and the dosimetric level. The inter-algorithm variability was proposed to be used for geometric as well as dosimetric DIR uncertainties (Nenoff et al 2020, Amstutz et al 2021b). Probabilistic unsupervised DL methods have also been proposed to predict the variance of DVFs in interfraction datasets (Gong et al 2022, Smolders et al 2022b, 2022a).

Dose mapping and accumulation have been used in work on treatment method combinations and patient re-irradiation. Application of the technology presents the possibility of greater outcome modelling in combined methodologies, and long term outcomes in re-irradiation. Research regarding the combination of external beam RT and brachytherapy was done for cervical cancer patients (Vásquez Osorio et al 2015, Swamidas et al 2020, Zeng et al 2020). Van Heerden did not find clinically relevant improvements when using DIR for dose accumulation compared to adding uniform external beam RT doses or overlapping high dose volumes (van Heerden et al 2017).

In recent years, improved survival has led to an increase in the numbers of re-irradiations (Nieder et al 2013, Andratschke et al 2022) with particular focus made on cancers of the brain (glioma), lung, HN, abdomen, pelvis and spine (Abusaris et al 2012, Mantel et al 2013, De Ruysscher et al 2014, Nieder et al 2016). Dose from previous treatments can be deformed to the current anatomy to evaluate potential dose overlap (Meijneke et al 2013, Nix et al 2022). Thereby, being used to define safe dose tolerances in those previously treated regions (Embring et al 2021, Andratschke et al 2022, Brooks et al 2022, Nix et al 2022). In addition, DIR-facilitated dose warping can be used to correlate places of local failure with previously planned and/or delivered dose distributions (Boman et al 2017, McVicar et al 2018, Skjøtskift et al 2018, Embring et al 2021, Nix et al 2022). Registration algorithms are challenged by dramatic anatomical changes caused by the time between treatments, often months or years, not to mention sequels of treatments such as fibrosis resulting from radiation or surgery (Nix et al 2022, Vasquez Osorio et al 2023b). Systematic studies about the DIR uncertainties in the re-irradiation setting are rare, but some reports indicate that DIR uncertainty increases with the magnitude of anatomical changes, in particular for lung radiographic changes after SBRT (Mahon et al 2020). DIR uncertainty is only one of multiple uncertainty factors which makes the definition of organ constraints for re-irradiation challenging. The lack of standardised toxicity scoring or cumulative DVHs over multiple treatments, partially influenced by DIR uncertainty remain reasons why the recovery of organs over time is not well quantified. The calculation of biologically effective dose can improve the understanding of normal tissue responses over time (Brooks et al 2022, Nix et al 2022) and allow a better estimation of safe dose constraints during re-irradiation.

Currently, tumour control probability (TCP) and normal tissue complication probability (NTCP) models are built on planned doses. They are however designated to correlate to delivered doses which can differ from the planned dose. Dose accumulation of reconstructed doses on repeated images, requiring DIR in most anatomical areas, is the closest surrogate to the delivered dose that is available. The impact of DIR uncertainty on the accumulated doses directly affects the outcome calculation (Nenoff et al 2021a, Smolders et al 2023b). Deformation-free methods (Niemierko 1997, Niebuhr et al 2021) have their own (not well quantified) uncertainties. Niebuhr et al found larger differences when assuming a registration error of 3 mm, compared to changing alpha-beta values for prostate RT. (Niebuhr et al 2021) There is more research needed to fully understand and quantify the impact of DIR uncertainty for outcome calculation.

Conventional outcome modelling simplifies the planned dose distribution to a single value, often using DVH statistics. Voxel-based analysis techniques that maintain the spatial distribution of doses have been used to explore local correlations between dose and treatment outcomes. Voxel-based analysis (figure 6) relies on DIR to 'spatially normalise' dose distributions into a common reference anatomy (Palma et al 2020, Shortall et al 2021). In summary, DIR is first performed between the planning CTs of each patient and an arbitrarily selected reference CT scan. The DIR result is then used to map the dose distributions to the reference anatomy allowing the local dose to be correlated with the studied outcome. The region is evaluated with statistical modelling, often quantifying the improvement in model discrimination when the dose to the identified region is included in a multivariable predictive model (including other demographic and clinical variables). The region is then used to generate hypotheses which are then tested and validated in external cohorts aiming at generating dose constraints to ultimately improve treatment outcomes.

Zoom In Zoom Out Reset image size

With voxel-based techniques, doses to anatomical subregions have been linked to outcomes, such as the dose to the base of the heart to overall survival in lung RT (McWilliam et al 2017, Green et al 2020), the inferior–anterior hemi-anorectum dose to rectal bleeding in prostate RT (Dréan et al 2016) and the cricopharyngeus muscle, cervical oesophagus and the base of the brainstem dose to dysphagia in HN RT (Monti et al 2017).

Several measures to evaluate DIR uncertainty for voxel-based analysis have been proposed (Palma et al 2020, Shortall et al 2021, McWilliam et al 2023). Quantified DIR uncertainties are often incorporated in the analysis by treating them as random errors and blurring the mapped dose distributions (McWilliam et al 2017, Beasley et al 2018, Green et al 2020, Vasquez Osorio et al 2023a). Therefore, DIR uncertainties can result in a decrease of significance for small radiosensitive regions and local changes in their shapes.