A novel bone suppression algorithm in intensity-based 2D/3D image registration for real-time tumor motion monitoring: Development and phantom-based validation.

BACKGROUND: Real-time tumor motion monitoring (TMM) is a crucial process for intra-fractional respiration management in lung cancer radiotherapy. Since the tumor can be partly or fully located behind the ribs, the TMM is challenging.
PURPOSE: The aim of this work was to develop a bone suppression (BS) algorithm designed for real-time 2D/3D marker-less TMM to increase the visibility of the tumor when overlapping with bony structures and consequently to improve the accuracy of TMM.
METHOD: A BS method was implemented in the in-house developed software for ultrafast intensity-based 2D/3D tumor registration (Fast Image-based Registration [FIRE]). The method operates on both, digitally reconstructed radiograph (DRR) and intra-fractional X-ray images. The bony structures are derived from computed tomography data by thresholding during ray-casting, and the resulting bone DRR is subtracted from intra-fractional X-ray images to obtain a soft-tissue-only image for subsequent tumor registration. The accuracy of TMM utilizing BS was evaluated within a retrospective phantom study with nine different 3D-printed tumor phantoms placed in the in-house developed Advanced Radiation DOSimetry (ARDOS) breathing phantom. A 24 mm craniocaudal tumor motion, including rib eclipses, was simulated, and X-ray images were acquired on the Elekta Versa HD Linac in the lateral and posterior-anterior directions. An error assessment for BS images was evaluated with respect to the ground truth tumor position.
RESULTS: A total error (root mean square error) of 0.87 ± 0.23 mm and 1.03 ± 0.26 mm was found for posterior-anterior and lateral imaging; the mean time for BS was 8.03 ± 1.54 ms. Without utilizing BS, TMM failed in all X-ray images since the registration algorithm focused on the rib position due to the predominant intensity of this tissue within DRR and X-ray images.
CONCLUSION: The BS algorithm developed and implemented improved the accuracy, robustness, and stability of real-time TMM in lung cancer in a phantom study, even in the case of rib interlude where normal tumor registration fails.

INTRODUCTION

Intra‐fractional tumor motion due to breathing can negatively influence the precision of dose delivered to a patient and can lead to cold spots in the tumor as well as hot spots in healthy tissue and thus reduce the effectiveness of the treatment while increasing toxicity to the patients. , , , , Real‐time tumor motion management can reduce this uncertainty by accounting for the motion of the target and therefore improve the accuracy of the delivered dose. ,

There are two approaches to tackle intra‐fractional tumor motion: gating where the beam is only turned on when the target is in the desired location , and tracking where the beam is continuously realigned to the target position. , , A prerequisite for both approaches is a fast and accurate registration between the reference and real tumor position over the time of irradiation. , , , ,

A comprehensive review regarding available tumor motion monitoring (TMM) systems was recently published by Bertholet et al. and the status of its clinical implementation by Anastasi et al. The registration accuracy of TMM, especially in lung tumors, is often compromised by a partial or complete overlap of the tumor with ribs due to the high image intensity of the rib segments compared to the tumor region. , Bone suppression (BS) algorithms integrated into image processing can significantly improve tumor visibility, as reported in several studies. The cancer detection rate, motion tracking error, and computer‐aided localization of lung nodules could be improved with some commercially available softwares, such as OnGuard, SoftView, ClearRead + Detect, , , , or the Samsung BS Software. Block et al. and Roeske et al. tested RapidTrack software (Varian Medical Systems, Palo Alto, CA, USA) in a phantom and patient study and reported a decrease in tracking range error when BS images were used for image registration. Nodule detection visibility could also be improved with a combination of machine learning and pattern recognition algorithms. , , , However, none of the systems was able to provide real‐time bone‐suppressed images from 3D computed tomography (CT) data for TMM.

The aim of this work was to develop a novel method for extracting bone tissue from digitally reconstructed radiographs (DRRs) generated from 3D CT data and the intra‐fractional X‐ray images in real time and to validate its accuracy for 2D/3D image registration with a sophisticated in‐house developed anthropomorphic breathing phantom. , As there is no commercial tumor tracking or radiotherapy system that allows us to conduct research with intrinsic registration algorithms, the in‐house developed software Fast Image‐based Registration (FIRE) , , was further enhanced to integrate a BS algorithm for intensity‐based (marker‐less) tumor registration. FIRE was designed as an open‐source project for Linux and Windows environments to be used by scientists, researchers, and software developers. The main part of FIRE is a 2D/3D registration algorithm, which involves the generation of DRRs from patient CT data.

MATERIALS AND METHODS

FIRE software with a BS algorithm

The open‐source software FIRE , , was developed at the Medical University of Vienna for research in the field of TMM in real time. It was based on a ray‐casting algorithm that was executed on a graphics processing unit (GPU) to render DRRs from patient CT data. A discrete line integral (summation) of Hounsfield units (HUs) was applied with a constant step size along the ray path to calculate the intensities of the detector pixels. The following merit functions for 2D/3D tumor registration were implemented in FIRE and can be selected via the graphical user interface (GUI) depending on the specific application:

Stochastic rank correlation (SRC) for intra‐modal registration of images with monotonous contrast dependency (e.g., different X‐ray energies).

Mutual information (MI) , for intermodal image registration.

Normalized cross‐correlation (NCC) for intra‐modal image registration only with linear contrast dependency.

For 2D/3D registration, a region of interest (ROI) spanning the area of the tumor motion was selected via the GUI of FIRE for the DRR and intra‐fractional X‐ray image. Thus, the following registration focused on tumor motion and neglected influences from neighboring image segments. A screenshot of the GUI interface is provided in Supporting Information (Figure S1).

The derivative‐free NEWUOA optimizer algorithm based on Powell's method was utilized to find the closest match between DRR and X‐ray images (Figure 1) by generating the translation and rotation vector for the volume until the best similarity was reached. The initial optimization radius (InitRhoStart) was set to 10.0, the function tolerance to 1 × 10−3, and the maximum number of 300 iterations per registration were defined as stop criteria. Several iterations were required until the optimum alignment (registration result) of the volume was found. For the assessment of the tumor motion trajectory, due to computational speed, the procedure was performed in the two most dominant degrees of freedom (DoF) regarding tumor motion (craniocaudal [CC]/lateral [LAT] or CC/posterior–anterior [PA] translation) for each intra‐fractional X‐ray image of the series acquired for the evaluation of the BS algorithm as described later in Section 3.2.

FIGURE 1

Workflow of tumor motion monitoring (TMM) and the intensity‐based 2D/3D registration of Fast Image‐based Registration (FIRE). Computed tomography (CT) data were loaded into the graphics processing unit (GPU) memory before the ray‐casting of a digitally reconstructed radiography (DRR) was performed (a). A merit function (e.g., normalized cross‐correlation [NCC], stochastic rank correlation [SRC], or mutual information [MI]) was selected via the graphical user interface (GUI) for similarity measurement between the region of interest (ROI) of the DRR and the intra‐fractional X‐ray image (b). The resulting value was used by the optimizer algorithm to generate a volume transformation vector T that contains desired values for translation and rotation of the volume (c). A ray‐casting was performed, and the entire process (a–c) was repeated until the best match between the ROI of the DRR and the X‐ray image was reached

All registrations were performed on an HP Zbook 15G3 laptop computer (four core CPU, 3.4 GHz CPU‐clock, and 16 GB RAM) with an integrated Nvidia Quadro M2000 graphics card (768 cores, 1.1 GHz clock, and 4 GB RAM).

To improve the accuracy and robustness of the TMM with FIRE, a BS method was implemented. For BS, an HU threshold similar to a patient rib value was used to indicate bone tissue during ray‐casting. For voxels above the threshold, a weighting factor was applied for the line integral to control suppression intensity. Since the voxel filter threshold and the weighting factor were patient specific, a slider for each parameter was implemented in the GUI of FIRE to allow optimized manual adjustment via direct visual feedback on the screen. Therefore, residual rib segments (high weighting) or dark segments (low weighting) in the BS DRR were avoided.

In the first step of BS, a bone DRR (DRRBone) was rendered, and a weighting factor to adjust the suppression strength was applied. Next, the image was subtracted from the original DRR (DRROrig), which contained all tissue compartments, to obtain the bone‐suppressed DRR (DRRBS), as shown in Equation (1).

Afterwards, the intensity transfer function (ITF) was applied on DRROrig and DRRBS to convert both images into the X‐ray imaging domain. The subsequent subtraction led to the simulated bone tissue X‐ray image XRAYBone (Equation 2).

Finally, XRAYBone was subtracted from the intra‐fractional X‐ray image (XRAYOrig), and the result represented the BS X‐ray image (XRAYBS) (Equation 3).

Optimization and evaluation of the BS algorithm

For optimization and evaluation of the BS algorithm, a tumor phantom (16.3 mm × 16.0 mm × 10.7 mm, 1.043 cm3) was 3D printed from polymethyl methacrylate material (24 ± 14 HUs). The printed tumor was placed inside the lung of the in‐house developed Advanced Radiation DOSimetry (ARDOS) breathing phantom. A set of CT scans was acquired on a Philips Brilliance Big Bore CT scanner (Philips Medical Systems, Inc., Cleveland, OH, USA) using the following settings: 120 kVp, 50 mAs, and slice thickness of 0.75 mm. Additionally, X‐ray images of a CC tumor motion of 50 mm with a 2 mm step size were acquired on the MedPhoton ImagingRing (MedPhoton GmbH, Salzburg, Austria) at PA and LAT imaging angles (120 kVp, 25 mA, and 20 ms). All CT and X‐ray image acquisitions were performed with and without the ribcage inside of the ARDOS phantom to obtain a ground truth scenario for BS.

Weighting factor, ITF, and contrast

For the evaluation of the weighting factor, the pixel intensity of a line profile was used to evaluate an optimal BS result. Values in the range of 0.6–1.0 were tested to find the optimum setting.

To perform BS on intra‐fractional X‐ray images, the BS DRR needed to be converted into the X‐ray imaging domain. Therefore, a nonlinear ITF of the selected ROI in the DRR and X‐ray image was sampled. Since this function depended on multiple variables (mean pixel intensity, contrast, and imaging angle), resampling was performed automatically before each BS in FIRE.

To evaluate the contrast improvements of a BS X‐ray image in comparison to the original X‐ray images, the contrast between ribs and lung, tumor and lung, and ribs and tumor was measured. These contrast measurements were performed for 36 combinations of bone/lung ROIs (Figure 2: blue/yellow rectangles), six combinations of tumor/lung ROIs (Figure 2: green/yellow rectangles), and six combinations of bone/tumor ROIs (Figure 2: blue/green rectangles). The formula for evaluating the contrast C of the mean intensity of an ROI (I) with respect to the mean intensity of the background (I) is given in Equation (4), and the contrast improvement CI is represented by the contrast ratio of the BS and original X‐ray image (Equation 5). The definitions of contrast and contrast improvement were used according to the work of Menten et al.

FIGURE 2

Location of regions of interest (ROIs) in the original X‐ray image used for contrast evaluation (lung, bone, and tumor segments are presented in yellow, blue, and green colors, respectively). Six ROIs of lung tissue (yellow rectangles) were defined at the corners and the vertical center of the ROI encompassing the whole motion range (red rectangle). Additionally, six ROIs of rib tissue (blue rectangles) were defined within the rib segments and one ROI within the tumor. The positions were chosen to investigate the influence of bone suppression (BS) on the bone/lung, bone/tumor, and tumor/lung contrast. The entire ROI (red rectangle) was used for fitting the intensity transfer function (ITF) to establish an intensity transformation between the digitally reconstructed radiograph (DRR) and the X‐ray imaging domain

Robustness and performance of TMM utilizing BS

For the evaluation of registration errors, robustness and performance of TMM utilizing BS, the influence of three different merit functions, that is, NCC, SRC, and MI, was investigated. , Error evaluation of the registration with ribcage against registration without ribcage was performed by analysis of the mean absolute error (MAE), root mean square error (RMSE), minimum and maximum error. For the evaluation of the real‐time capability regarding BS, the mean time for 100 consecutive DRR renderings was measured for different ROI sizes (full size, rib segment, and tumor motion segment) and ray‐step increments (2.5, 1.0, and 0.1 mm). A ray‐step increment is the resolution of the CT volume along a virtual X‐ray path used for DRR generation. High resolution is required to reduce artifacts in the DRR, especially in tissue with high HU values; however, it requires more computational time. Thus, we introduced a fourth option of ray‐step increments, which we called “dynamic bone oversampling” (DBOS). This method automatically adjusts the ray‐step increment from 1.0 mm in non‐bone tissue to 0.1 mm within bone. Therefore, a high‐quality DRR used for bone subtraction from the intra‐fractional X‐ray image was generated by a minimum amount of additional runtime.

Validation of TMM utilizing BS in a phantom study

TMM including BS was validated in a phantom study with nine different tumor inserts focusing on accuracy, robustness, and runtime. Therefore, the best settings identified in the evaluation of the BS algorithm were used for further validation of the TMM.

Tumor inserts of different sizes were extracted from patient CT data and 3D printed, resulting in the following inserts:

small tumor insert (9.3 mm × 10.1 mm × 14.4 mm, volume: 0.665 cm3);

medium tumor insert (22.6 mm × 24.0 mm × 18.8 mm, volume: 5.345 cm3);

large tumor insert (29.1 mm × 36.2 mm × 26.0 mm, volume: 11.611 cm3).

For each tumor size, three replicas with different HU ranges (low, medium, high) were printed to allow the evaluation of TMM with respect to the X‐ray attenuation. The medium HU range was defined in a similar radiological range as real patient tumors, whereas the low and high HU ranges were used to include HU range variability among different patients in the evaluation. The list of all tumor phantoms, including properties, is given in Table 1.

TABLE 1

List of 3D‐printed tumor inserts used for the validation of bone suppression (BS) in Fast Image‐based Registration (FIRE)

Phantom‐ID	Size	Material	Mean HU
S‐PLA‐LOW	Small	Polylactic acid	−137 ± 15
S‐PLA‐MEDIUM	Small	Polylactic acid	64 ± 20
S‐NYL‐HIGH	Small	Nylon	86 ± 4
M‐NYL‐LOW	Medium	Nylon	−53 ± 41
M‐NYL‐MEDIUM	Medium	Nylon	55 ± 13
M‐NYL‐HIGH	Medium	Nylon	80 ± 9
L‐PLA‐LOW	Large	Polylactic acid	−31 ± 30
L‐PLA‐MEDIUM	Large	Polylactic acid	10 ± 23
L‐PLA‐HIGH	Large	Polylactic acid	72 ± 36

Note: The Phantom‐ID summarizes information about tumor size, material, and an increased or decreased Hounsfield unit (HU) range (low/high).

Abbreviations: NYL, nylon; PLA, polylactic acid.

For 3D printing, a custom filament printer (Original Prusa i3 MK3S, Prusa Research a.s., Praha, Czech Republic) was used. Printing temperature and speeds were selected according to the filament manufacturer's and 3D printer specifications. Samples were printed from two different filament materials: polylactic acid and nylon 12.

The individual tumors were inserted into the ARDOS breathing phantom without residual air gaps. , A setup can be seen in Figure S2. A linear tumor motion of 24 mm in the CC direction was simulated according to the work of Seppenwoolde et al., which presents an average breathing amplitude of 12 ± 2 mm for the 20 investigated patients in clinical practice. The tumor was moved in discrete steps of 2.0 mm followed by an X‐ray image acquisition after each step to obtain the ground truth tumor position.

For DRR generation, volumetric CT data of the ARDOS breathing phantom including the tumor inserts were acquired on the Siemens Somatom Definition AS CT (Siemens Healthcare AG, Erlangen, Germany) with the following settings: 120 kVp and a slice thickness of 2.0 mm. The voxel size was 0.8 × 0.8 × 2.0 mm3 after reconstruction, and for the purpose of TMM, volumes were resampled to 1.0 × 1.0 × 1.0 mm3. The scans were acquired with and without the ribcage inside of the ARDOS phantom.

For image registration, X‐ray images were acquired (120 kVp, 5 mA) on the Elekta XVI imaging system (Elekta Versa HD Linac). Figure S2 shows the measurement setup in the treatment room. In summary, four image series consisting of 13 X‐ray images (start position and 12 step increments of 2 mm) were acquired for each of the nine tumor phantoms. For error evaluation of the TMM, three different configurations were used:

ARDOS with ribcage and TMM without BS (TMMRibs),

ARDOS with ribcage and TMM with BS (TMMRibs+BS),

ARDOS without ribcage and TMM without BS (TMMnoRibs).

The scenario TMMnoRibs was used as the reference configuration since the missing ribcage represented a perfect BS result that was compared with our BS method in scenario TMMRibs+BS. TMMRibs represented the initial configuration of FIRE with known influence of the ribs on the registration result.

For clinical relevance, registrations were performed in two DoFs: the CC and LAT directions in PA imaging and the CC and PA directions in LAT imaging. The registration error was evaluated in Matlab (Release R2019b; The MathWorks, Inc., Natick, MA, USA) in the CC and LAT directions or the CC and PA directions, corresponding to the imaging angle. Additionally, the RMSE and 2D (Euclidean) error were evaluated.

The mean runtime of the registrations for the scenarios TMMRibs+BS and TMMnoRibs was measured for comparison. Additionally, the mean runtime for the BS process (ITF sampling, intensity adaptation, image subtraction) was measured for the scenario TMMRibs+BS.

To evaluate a potential dependency of the TMM error with respect to tumor size and HUs, 2D registration errors were grouped for each tumor, and a Kruskal–Wallis test was performed to test for significant differences in the distribution of the registration error.

RESULTS

The evaluation results of the optimal weighting factor are shown in Figure 3. From the pixel intensity of the vertical line profile for different values of the weighting factor, the value of 0.8 was found to be optimal. Values below the optimum showed a residual rib intensity in the line profiles, whereas values above introduced artifacts in the form of dark rib segments.

FIGURE 3

Pixel intensities in the digitally reconstructed radiograph (DRR) before (solid blue line) and after bone suppression (BS) (dashed lines) of the vertical line profile across the ribs (blue line, right side) for evaluation of the weighting factor w. Optimized BS without residual ribs or the introduction of artifacts in the form of dark segments was performed in the case of w = 0.8

Figure 4 demonstrates the variation in pixel intensities (blue markers) surrounding the ITF (red curve) originating from differences in image sharpness, image noise, and beam‐hardening artifacts, which were not considered during ray‐casting because of the strict real‐time requirements. Additionally, a residual registration error contributed to the spread of pixel intensities.

FIGURE 4

Measured intensity transfer function (ITF) (red line) for the conversion of the digitally reconstructed radiograph (DRR) into the intra‐fractional X‐ray image domain required for bone suppression (BS). For curve fitting, pixel intensities of the region of interest (ROI) were sampled from the DRR and intra‐fractional X‐ray image (blue markers)

For the contrast enhancement evaluation of the tumor due to BS, six ROIs for lung tissue (15 × 15 pixels), six ROIs for bone tissue (15 × 15 pixels), and one ROI in the tumor region (11 × 11 pixels) were selected in original and BS X‐ray images. Contrast improvements (Equation 5) regarding the median contrast of the original and BS X‐ray images were found to be 0.986, 0.056, and ‐0.718 for tumor/lung, bone/lung, and bone/tumor, respectively. Thus, the contrast between tumor and lung tissue remained unchanged, whereas the contrast of bone and lung tissue was strongly reduced (Figures 5 and 6).

FIGURE 5

(a) Vertical line profiles across the ribs of the bone suppression (BS) (dashed red line) and original (solid blue line) X‐ray image acquired on the MedPhoton Imaging Ring. The intensity of peak plateaus (blue line), which were introduced from the high attenuating rib tissue, was reduced to a minimum, whereas lung tissue was unaffected (overlapping regions of blue and dashed red line). The original and bone‐suppressed X‐ray images are shown in (b) and (c), respectively.

FIGURE 6

Contrast in original and bone suppression (BS) X‐ray images of different tissue segments (bone, tumor, and lung). The contrast of bone and lung was decreased (middle column), while the contrast of tumor and lung was unaffected (left column)

An overview of the results regarding the robustness of the TMM is presented in Table 2, which shows the mean number of DRR renderings required for the tumor registration for the three tested merit functions (NCC, SRC, MI) in 1DoF/2DoF depending on the ray‐step size and DBOS setting. The results for a step size of 0.1 mm were not evaluated since the high runtime for ray‐casting would not allow for tumor registration in real time. The resulting values focus on the efficiency of the registration algorithm in the form of the mean number of required DRR renderings for a single tumor registration and therefore present a measure that is independent of the performance of the CPU and GPU; for NCC 16.0/30.7 and SRC 15.6/29.7 renderings in 1DoF/2DoF were needed, whereas for MI, TMM failed in every case when BS was applied. In terms of registration error (RMSE, MAE), NCC outperformed SRC in all tested scenarios (Table S2).

TABLE 2

Tumor motion monitoring (TMM) performance and robustness analysis with respect to the ray‐casting step size, presence of ribs, degree of freedom in registration, and different merit functions (normalized cross‐correlation [NCC], stochastic rank correlation [SRC], or mutual information [MI])

			Mean number of DRR renderings CC (1DoF)			Mean number of DRR renderings CC/LAT (2DoF)
Ray‐step size	Ribcage in phantom	BS (DBOS)	NCC	SRC	MI	NCC	SRC	MI
2.5 mm	No	No	16.8	17.0	×	32.2	32.5	×
	Yes	No	×	×	×	×	×	×
	Yes	Yes (none)	×	×	×	×	×	×
	Yes	Yes (0.1)	17.1	17.8	×	33.7	×	×
1.0 mm	No	No	13.9	14.3	16.3	28.4	29.7	33.7
	Yes	No	×	×	×	×	×	×
	Yes	Yes (0.1)	16.0	15.6	×	30.7	29.7	×

Note: A step size of 0.1 mm was not applicable in terms of real‐time ray‐casting; therefore, the corresponding results were not evaluated. All cases where tumor registration failed were marked with “×.”

Abbreviations: BS, bone suppression; CC, craniocaudal; DBOS, dynamic bone oversampling; DoF, degrees of freedom; DRR, digitally reconstructed radiography; LAT, lateral.

For the evaluation of the average ray‐casting runtime of 100 consecutive renderings (including memory transfer), the time was measured for the recommended configuration (ray‐step size 1.0 mm, DBOS 0.1 mm). DRR rendering of the ROI took 2.4 ms for the ribs (410 × 140 pixels) and 0.8 ms for the tumor motion region (34 × 119 pixels). An extended overview regarding runtimes with respect to different ray‐step sizes is presented in Table S1.

The mean time for tumor registration (2DoF, NCC as merit function, ray‐step size 1.0 mm, DBOS 0.1 mm) was measured for all three tumor ROIs (ROIsmall: 48 × 64, ROImedium: 64 × 87, ROIlarge: 78 × 120 pixels). The respective registration times were 30.6 (±2.5) ms, 32.9 (±2.9) ms, and 42.8 (±4.0) ms (no ribs, no BS, no DBOS), and 140.0 (±12.5) ms, 144.1 (±10.9) ms, and 162.9 (±14.9) ms (with ribs and BS + DBOS). The number included the time for BS (ITF sampling, ROI conversion, and subtraction), which was 6.9 (±1.1) ms, 7.6 (±1.1) ms, and 9.5 (±1.1) ms for the three ROIs. An overall registration time of 149.0 (±16.3) ms, including 8.0 (±1.5) ms for BS, was found. In summary, a disproportionate dependency between registration time and ROI size was found, and consequently, the ROI size of very large tumors may interfere with real‐time requirements in a worst‐case scenario.

Validation of the TMM using BS in a phantom study

From the evaluation of the robustness and performance described above, we found that TMM utilizing BS showed the best results for a ray‐step size of 1.0 mm, DBOS step size of 0.1 mm, and NCC as the merit function. Therefore, these settings were applied for all subsequent measurements of the phantom study.

Table 3 summarizes the RMSE of the TMM with BS in both imaging scenarios (PA and LAT) for the corresponding directions. The registration succeeded on all X‐ray images, and total errors of 0.87 ± 0.23 mm and 1.03 ± 0.26 mm were found for PA and LAT imaging, respectively. In the reference scenario when ribs were removed from the phantom, total errors of 0.12 ± 0.06 mm and 0.16 ± 0.05 mm were found for PA and LAT imaging, respectively. More details are provided in Table S2.

TABLE 3

Root mean square error (RMSE) of tumor motion monitoring (TMM) performed with Fast Image‐based Registration (FIRE) including bone suppression (BS) with respect to the ground truth tumor motion (24 mm linear tumor motion simulated in the Advanced Radiation DOSimetry [ARDOS] breathing phantom)

RMSE (mm) of TMM (NCC) with ribs and BS
	PA imaging			LAT imaging
Tumor	CC	LAT	2D	CC	PA	2D
L‐PLA‐LOW	0.80	0.56	0.97	0.79	0.20	0.81
L‐PLA‐MEDIUM	0.57	0.34	0.67	0.67	0.32	0.74
L‐PLA‐HIGH	0.60	0.36	0.70	0.77	0.47	0.91
M‐NYL‐LOW	0.72	0.46	0.85	1.60	0.17	1.60
M‐NYL‐MEDIUM	0.46	0.51	0.69	1.21	0.17	1.22
M‐NYL‐HIGH	0.51	0.43	0.67	0.97	0.24	1.00
S‐PLA‐LOW	1.31	0.07	1.31	1.11	0.33	1.16
S‐PLA‐MEDIUM	0.88	0.10	0.89	0.87	0.22	0.96
S‐NYL‐HIGH	0.96	0.51	1.09	0.83	0.40	0.92
Mean RMSE	0.76 ± 0.27	0.37 ± 0.18	0.87 ± 0.22	0.98 ± 0.29	0.28 ± 0.11	1.03 ± 0.26

Note: Registration was performed in 2 degrees of freedom (DoF) in posterior–anterior (PA) and lateral (LAT) images (90° and 0° gantry angle) acquired on the Elekta XVI imaging system for nine different tumors (three different tumor sizes with three different mean Hounsfield units [HUs] each).

Abbreviations: CC, craniocaudal; NCC, normalized cross‐correlation; NYL, nylon; PLA, polylactic acid.

Figure 7 shows the distribution of the absolute error for each 2D/3D tumor registration regarding tumor size, mean HU, and imaging angle (PA and LAT imaging in the left and right columns, respectively). The results are plotted for two phantom configurations, one with ribs and BS (blue markers) and one with the ribs removed from the phantom as a reference (green markers).

FIGURE 7

Distribution of the absolute error of tumor motion monitoring (TMM) with Fast Image‐based Registration (FIRE) for each simulated tumor trajectory within the Advanced Radiation DOSimetry (ARDOS) breathing phantom. The resulting values are shown in two dimensions, craniocaudal (CC)/lateral (LAT) (left column) and CC/posterior–anterior (PA) (right column), for the corresponding imaging angle. The absolute error is given for all tumor phantoms (three different sizes with three different mean Hounsfield units [HUs]) for the motion scenario with ribs (blue markers) and without ribs (green markers)

The Kruskal–Wallis test examining the potential differences in the 2D error distribution caused by different tumor sizes and HU ranges showed no significant differences (chi square = 5.75, p = 0.67, df = 8 and chi square = 4.04, p = 0.85, df = 8) for PA and LAT imaging.

DISCUSSION

Respiration‐induced tumor motion compromises treatment delivery accuracy, especially in the lung. TMM using intra‐fractional X‐ray imaging offers a potential solution. However, in clinical practice, the accuracy and robustness of the TMM are affected by the overlapping of the tumor with bony tissue in the projection images. The proposed BS algorithm implemented in our open‐source software FIRE for 2D/3D image registration software improved the accuracy and robustness of TMM. The method was designed to incorporate volume data of the planning CT. It does not require a double X‐ray exposure as in dual energy subtraction, , which would lead to an increased imaging dose in the patient.

In contrast to existing techniques for image registration, , our proposed method utilizes ray‐casting of the bone structure for subtraction in the intra‐fractional X‐ray images after intensity adaptation. Therefore, the following 2D/3D registrations of TMM will be performed without the interfering contrast of bone tissue. Ray‐casting is based on the planning CT. Thus, it does not require an additional daily cone‐beam CT.

The BS algorithm required manual adjustment of the patient‐specific weighting factor for optimized subtraction. Lower values contributed to residual rib structures, whereas higher values introduced dark artifacts within rib segments. The weighting factor as well as the bone threshold value was adjusted manually by visual inspection of the GUI and therefore optimized to suppress ribs, including bone marrow and cartilage, in a single operation. This approach performed well for the phantom rib structure, and the first promising tests performed on patient data (compare Figures S3 and S4) but will need additional refinements when applied to other bony structures, such as the vertebral column. Since FIRE was designed as a research tool, fine tuning of parameters was limited to manual adjustment in the GUI. However, within a clinical routine, an automatic adjustment of the weighting factor and voxel filter threshold should be implemented to save time.

Three different merit functions were used for the evaluation of TMM accuracy and robustness. With respect to BS, MI performed worse in contrast to NCC and SRC, which showed lower RMSE in 1DoF and 2DoF. Additionally, the results in Table 2 showed successful TMM in the case of MI combined with a ray‐casting step size of 1.0 mm for images without ribcage, whereas registration failed for a step size of 2.5 mm. Steep gradients created by MI result in a higher sensitivity to image artifacts. Thus, a step size of 2.5 mm leads to a noisy merit function, which introduces additional local minima. This led to a decreased optimizer accuracy during registration. In terms of robustness, we recommend utilization of NCC as a merit function of choice for the proposed BS method. However, MI is mainly used for registration of images that do not share a similar intensity function (e.g., CT to magnetic resonance) and is not required for CT to X‐ray registration.

The subtraction of bony tissue required an intermodal image conversion between the DRR and the X‐ray imaging domain. In FIRE, the patient‐ and imaging system‐specific ITF was sampled automatically before each performed BS for the selected ROI regarding the entire tumor motion. To maintain the real‐time capability of the BS algorithm, image sharpness and noise, as well as beam hardening effects, were neglected during ray‐casting. Such differences in image sharpness between DRR and X‐ray images led to residual artifacts at the edges of ribs after subtraction of DRRBone. Therefore, a higher registration error in the CC direction was observed since the ribs of the phantom were in horizontal alignment. In addition, the residual registration error contributed to a spread out of corresponding pixel values of the DRR and the X‐ray. This issue is resolved using a polynomial fit to obtain the ITF used for BS.

To avoid aliasing artifacts in DRRBS, we found that a ray‐step size of 0.1 mm was needed. Such a small step size would compromise the real‐time capability of ray‐casting. Thus, the DBOS algorithm was developed to alter the ray‐step size depending on the tissue type (0.1 mm in bony tissue; 1.0 mm outside). In general, we recommend a ray‐step size equal to the voxel size of the volume data, since in this configuration every voxel along the ray‐path will be sampled at least once. In between the voxel centers, the GPU was configured to calculate interpolated values of neighboring voxels. In this way, aliasing artifacts were reduced, and the required image quality for BS was achieved.

The proposed BS algorithm was validated with the ARDOS breathing phantom. This phantom allows the comparison to the ground truth tumor motion used for registration error evaluation as it can be used without ribcage. Three different tumor volumes were extracted from real patient CT data, and 3D tumor phantoms were printed, mimicking both the geometry and HU in a clinical situation. A limitation was the homogenous tumor material and restriction to isolated tumors without attachment to organs. Additionally, the ARDOS breathing phantom utilized in this work represented a simplified representation of the human torso, and further limitations are discussed in the work of Kostiukhina et al.

Another limitation of the study is the investigation of roughly spherical and isolated tumors moving in the CC direction. Since the merit functions of intensity‐based registration methods are focusing on the strongest image intensities independent of the tumor shape, registration may fail when tumors are attached to organs like heart, chest wall, or vertebrae owing to dominating contrast of organ tissue compared with the tumor tissue.

The main advantage of the validation using our phantom was the removal of the ribcage to obtain a ground truth BS in addition to the ground truth tumor motion. However, the ribcage position was fixed in our phantom setup, whereas rib motion and deformation were present in patients due to breathing. To evaluate the performance of the algorithm in situ, a validation study with real patient data will be conducted as a next step. To prepare our BS method for rib motion, we performed rib registration for every intra‐fractional X‐ray image immediately before each tumor registration to reduce the effect of rib motion. To minimize registration time and maximize robustness, a rigid registration method was implemented in FIRE, which neglects deformation of the ribcage caused by breathing. In a clinical setting, this problem can be mitigated by using 4D CT data correlated to a breathing signal for DRR generation.

In the presence of ribs, TMM using FIRE without BS failed due to the high rib contrast compared to lung tissue, which effected the optimizer to got stuck in local minima. Using TMM with BS succeeded in all test scenarios, independent of the distance of the tumor position between CT and X‐ray. For both phantom configurations (ribs with BS, no ribs), the higher RMSE found in LAT imaging was a consequence of the decreased tumor contrast due to the higher background intensity caused by lung and soft‐tissue compartments of the phantom mimicking the mediastinum.

Evaluation of the registration speed of FIRE including BS confirmed the real‐time capability. In 2DoF, the mean time required for tumor registration was 149 ± 16 ms, which allowed TMM at a rate of 6.7 Hz. Here, the time for BS (ROI rendering, ITF sampling, intensity adaption, and ROI subtraction) only adds 8 ± 2 ms. Registration of the rib structure dominated the runtime due to the large size of the rib ROI compared to the ROI of tumor motion. Furtado et al. reported a total error of 1.6 ± 0.3 mm for gold marker‐based validation of FIRE in 3DoF. Adding BS to FIRE improved the total error by 45.6% to 0.87 ± 0.23 mm for PA and by 35.6% to 1.03 ± 0.26 mm for LAT imaging.

Remmerts de Vries et al. reported successful marker‐less TMM during radiotherapy in 71% of the treatment time on average. In the remaining 29% where TMM failed, one challenging factor was presented by the over‐projection of the spine and the resulting image saturation.

Tanaka et al. focused on the registration error of marker‐less tumor registration and utilized an artificial neural network (ANN) for BS. The authors reported registration errors of 1.9 ± 1.7 mm and 0.4 ± 0.3 mm for the original and BS X‐ray images, respectively. The low error in the latter scenario represents a remarkable result and demonstrates the power of BS with ANN. Unfortunately, these methods tend to suffer from extended training and computation times (e.g., 15 s per image in the case of Tanaka et al.) and are not feasible for real‐time applications in general. Thus, these methods are mainly suitable for image post‐processing.

The ground truth scenario (phantom without ribcage) yielded total errors of 0.12 ± 0.06 mm and 0.16 ± 0.05 mm in PA and LAT, respectively. These results highlight the potential of BS with respect to TMM. In our study, several approximations were introduced to achieve real‐time requirements. Our approach neglected the simulation of beam hardening, sharpness, and noise during DRR rendering. With increasing GPU performance in the future, an implementation of such details seems feasible and promises a potential reduction of the registration error.

The validation of the 2D/3D registration with BS was performed with a phantom study with static images with the tumor in different positions simulating a linear tumor motion of 24 mm in the CC direction, which was selected as a representative scenario for the mean tumor motion found in patients in clinical practice. This demonstrates the feasibility of such an approach, which can be used for both tumor tracking as well as gating. For an upcoming patient study, the acquisition of the breathing phase with a surface scanner will allow us to select the corresponding volume from 4D CT data for ray‐casting. Therefore, the ROI size will be reduced, and tumor registration of a larger tumor motion than 24 mm will be possible.

FIRE was designed to foster research in the field of real‐time TMM, which at the same time enables software modifications and subsequent algorithm testing. Currently, medical device certification is not available. However, due to the promising results, an ethics application is planned to conduct a patient study including 4D CT data and breathing phase detection to optimize and validate the method on patient data in the near future.

CONCLUSION

The BS algorithm developed and implemented in an in‐house developed software for 2D/3D registration improved the accuracy, robustness, and stability of real‐time TMM in lung cancer in a phantom study. To evaluate the potential benefits of this promising method in clinical routine, a patient study will be conducted to confirm the benefit of BS for real‐time TMM.

CONFLICTS OF INTEREST

The authors declare they have no conflicts of interest.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.