- Email: [email protected]
Dose warping performance in deformable image registration in lung
Physica Medica 37 (2017) 16–23
Contents lists available at ScienceDirect
Physica Medica journal homepage: http://www.physicamedica.com
Original paper
Dose warping performance in deformable image registration in lung Shunsuke Moriya a, Hidenobu Tachibana b,⇑, Nozomi Kitamura c, Amit Sawant d, Masanori Sato a a
Radiological Sciences, Graduate Division of Health Sciences, Komazawa University, Tokyo 1548525, Japan Particle Therapy Division, Research Center for Innovative Oncology, National Cancer Center, Chiba 2778577, Japan c Department of Radiation Oncology, Cancer Institute Hospital of the Japanese Foundation of Cancer Research, Tokyo 1358550, Japan d Department of Radiation Oncology, University of Texas Southwestern Medical Center, 5801 Forest Park Rd., Dallas, TX 75390-9183, USA b
a r t i c l e
i n f o
Article history: Received 14 July 2016 Received in Revised form 13 February 2017 Accepted 20 March 2017
Keywords: Deformable image registration Generalized equivalent uniform dose Four-dimensional computed tomography Lung cancer
a b s t r a c t Purpose: It is unclear that spatial accuracy can reflect the impact of deformed dose distribution. In this study, we used dosimetric parameters to compare an in-house deformable image registration (DIR) system using NiftyReg, with two commercially available systems, MIM Maestro (MIM) and Velocity AI (Velocity). Methods: For 19 non-small-cell lung cancer patients, the peak inspiration (0%)-4DCT images were deformed to the peak expiration (50%)-4DCT images using each of the three DIR systems, which included computation of the deformation vector fields (DVF). The 0%-gross tumor volume (GTV) and the 0%-dose distribution were also then deformed using the DVFs. The agreement in the dose distributions for the GTVs was evaluated using generalized equivalent uniform dose (gEUD), mean dose (Dmean), and threedimensional (3D) gamma index (criteria: 3 mm/3%). Additionally, a Dice similarity coefficient (DSC) was used to measure the similarity of the GTV volumes. Results: Dmean and gEUD demonstrated good agreement between the original and deformed dose distributions (differences were generally less than 3%) in 17 of the patients. In two other patients, the Velocity system resulted in differences in gEUD of 50.1% and 29.7% and in Dmean of 11.8% and 4.78%. The gamma index comparison showed statistically significant differences for the in-house DIR vs. MIM, and MIM vs. Velocity. Conclusions: The finely tuned in-house DIR system could achieve similar spatial and dose accuracy to the commercial systems. Care must be taken, as we found errors of more than 5% for Dmean and 30% for gEUD, even with a commercially available DIR tool. Ó 2017 Associazione Italiana di Fisica Medica. Published by Elsevier Ltd. All rights reserved.
1. Introduction The locations and volumes of organs inside a patient’s body can change significantly during irradiation and between treatments [1,2]. Anatomical changes over an entire course of treatment can compromise the value of the initial treatment plan and therefore the treatment outcomes [3]. Adaptive radiotherapy (ART) has therefore been proposed to overcome this challenge; it allows the plan quality to be maintained by modification or reoptimization of the treatment plan according to changes in a patient’s anatomy [4,5]. Inter- and intra-fraction organ motion is especially significant in lung radiotherapy, and ART can provide a more accurate final dose distribution from the information on the accumulated dose distribution. This is obtained using ⇑ Corresponding author. E-mail addresses: [email protected] (S. Moriya), [email protected] (H. Tachibana), [email protected] (N. Kitamura), Amit.Sawant@ UTSouthwestern.edu (A. Sawant), [email protected] (M. Sato).
fraction-by-fraction dose distributions deformed according to a reference transformation between the planning CT and treatment CT. Deformable image registration (DIR) is therefore an essential tool for ART. The use of four-dimensional computed tomography (4DCT) has also facilitated treatment planning by allowing the respiratory motion of the target and critical organs to be incorporated into the analysis. The tumor volume and shape can be estimated more accurately in a phased 4DCT image than they can from a maximum or mean intensity projection formed from the 4DCT scan. A typical planning also usually considers only a static tumor, even though the tumor will move during the irradiation. An accumulated 4D dose calculation is therefore useful for achieving a more realistic estimation of the dose delivered to the tumor in the lung and the surrounding organs [6–9]. Information on the delivered 4D dose accumulation can help the physician evaluate the treatment and decide when and how a re-plan should be performed. Accurate dose warping is therefore one of the most important factors in the ART methodology.
http://dx.doi.org/10.1016/j.ejmp.2017.03.016 1120-1797/Ó 2017 Associazione Italiana di Fisica Medica. Published by Elsevier Ltd. All rights reserved.
17
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Many investigators have reported on methods for evaluation of the accuracy of DIR, including the use of the Dice similarity coefficient (DSC), landmark comparisons, mean slice-wise Hausdorff distance to agreement, volume differences, receiver operating characteristics, and target registration error [10–14]. However, few studies have directly evaluated dose deformation. Yeo et al. evaluated the accuracy of dose warping using DEFGEL, which is a deformable three-dimensional-dosimetry gel phantom [15–17]. Dosimetry using DEFGEL can provide direct measurements of dose, taking into consideration deformable media. However, the dosimetry information available from gel phantoms is limited, as they lack geometric complexity and are unable to provide assessments for individual patient cases. Additionally, dosimetry with deformable gel phantoms requires specialized and expensive equipment. In addition to dose monitoring in the lung, DIR may be applied in various other fields, such as radiotherapy and nuclear medicine [18,19]. There are several publicly available DIR systems, which may be suitable for contributing to clinical practice and the research environment, provided that the system is at least comparable to commercially available DIR systems. The purpose of this study was twofold; the first part was to build an in-house DIR system with a graphical user interface and image registration using the publically available NiftyReg (freedownloadable software library package, NiftyReg), while the second part of the study compares the in-house NiftyReg DIR system with two commercially available systems, MIM Maestro (MIM Software Inc., Cleveland, OH, USA) and Velocity AI (Varian Medical Systems, Palo Alto, CA, USA), with the evaluations being made using mean dose (Dmean), generalized equivalent uniform dose (gEUD), and three-dimensional (3D) gamma index evaluation. These dosimetric parameters are common in clinical practice and are easy to use, with a previous study having applied the generalized equivalent uniform dose (gEUD) as an index to quantify plan quality [20]. Dmean and gEUD were investigated as indices of the accuracy of dose warping using the DSC representing the spatial accuracy of image registration under the DIR. 2. Materials and methods 2.1. 4D-CT scan acquisition and patients All 4DCT scans were acquired on a 4-slice clinical scanner (Brilliance CT Big Bore, Philips Healthcare, Andover, MA). The respiratory cycle signal was monitored using an abdominal bellows pressure belt system (Philips Medical Systems). Other imaging parameters (e.g. tube voltage [kV] and tube current exposure time product value [mAs]) were set according to site-specific standard imaging protocols. CT data were reconstructed with a field of view of 50 cm on a 512 512 grid with a slice thickness of 2.0 mm. The longitudinal scan length was determined on a scout view. Nineteen patients treated with lung stereotactic body radiation therapy underwent 4DCT scans. The peak inspiration (0%) and peak expiration phases (50%) of the 4DCT image datasets were used in this study. Tumor displacement was defined as the vector sum from the centroid of the 0%gross tumor volume (GTV) to the centroid of the 50%-GTV. The tumor displacement varied from 0.30–1.55 cm across cases. The volumes of the gross tumor volume (GTV) ranged from 0.89– 99.0 cm3 for the 0% phase and 1.06–99.0 cm3 for the 50% phase (Table 1). The mean percentage difference between the GTVs at the 0% and 50% phases was 0.98 ± 0.38% (Maximum: 28.2%). 2.2. Planning and dose warping Treatment planning (structure contouring, treatment field setting, and dose calculation) was performed on the 0%-4DCT image
Table 1 Tumor displacement, GTV obtained at the peak inspiration phase (0%) and the peak expiration phase (50%), and the percentage difference of the volume change in the 19 patients evaluated. Case no.
Tumor displacement [cm]
GTV at 0% (inspiration) phase [cm3]
GTV at 50% (expiration) phase [cm3]
Percentage difference of volume change [%]
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
0.62 1.55 1.51 0.86 1.37 1.22 1.36 0.30 1.54 0.82 1.12 0.96 0.38 0.91 1.07 0.61 1.01 0.77 0.62
2.46 60.81 29.46 9.30 5.35 9.13 25.30 1.58 13.50 0.89 16.32 14.46 99.00 4.58 4.37 4.25 92.43 8.32 14.07
2.50 60.55 29.28 8.19 6.86 9.57 28.13 1.19 13.98 1.06 16.07 16.06 99.84 4.40 3.85 4.03 85.84 8.33 11.62
1.6 0.4 0.6 11.9 28.2 4.8 11.2 24.7 3.6 19.1 1.5 11.1 0.8 3.9 11.9 5.2 7.1 0.1 17.4
dataset using Eclipse treatment planning software (Version 10.0, Varian Medical Systems, Palo Alto, CA, USA). Radiation oncologists were asked to draw contours defining the GTVs on the 50%-4DCT image dataset. A medical physicist with five-years of experience delineated the GTVs on the 0%-4DCT image datasets, under the observation of a radiation oncologist who also checked them all. The mean DSC value representing the inter-observer variability for GTV volume on the 0%-4DCT images was 0.85 ± 0.05 (Range: 0.73–0.92), a value that is comparable to the inter-observer delineation variation found in a previous study [21]. The clinical target volume (CTV) was set as equivalent to the GTV, and the planning target volume (PTV) was created by adding a 5 mm margin to the CTV in all directions. The number of treatment fields ranged from 9 to 13. The aperture shapes of the multileaf collimator were adjusted in the beam’s eye view to cover the PTV with an additional 5 mm margin, to take into consideration the inaccuracy of dose calculation in the border region between the tumor and lung. All plans were generated with the prescription dose (54 Gy/3 fr) covering 95% of the PTV. A dose calculation grid size of 2.5 mm was used for the anisotropic analytical algorithm (AAA) calculations, and a dose volume histogram (DVH) of the GTV was subsequently computed. The 0%-4DCT image dataset was then deformed to the 50%-4DCT image dataset using the three DIR systems, thereby facilitating the computation of the deformation vector field (DVF) for each system. The dose distribution in the 0%4DCT was also deformed according to each of the three DVFs obtained. The deformed dose distributions were then superimposed onto the 50%-4DCT image datasets with the corresponding contoured structures, and the DVHs of the GTVs were computed. 2.3. Dose warping algorithm 2.3.1. In-house deformable image registration system using NiftyReg The in-house DIR system included a graphical user interface and image registration using NiftyReg. The image registration process consisted of three major steps. In the first step, a rigid registration was performed using a block-matching algorithm for the whole body of the patient [22]. The second step entailed a deformation inside the lung using a fast free-form deformation (FFD) algorithm, with the third step involving a larger-scale deformation inside the body. The second step focused on the tumor and lung, which have
18
S. Moriya et al. / Physica Medica 37 (2017) 16–23
the greatest motion, while the third step focused on whole organs, as deformation occurs not only in the tumor and lung, but also in the surrounding organs, such as the heart and liver. Thus, two deformation processes were performed inside the lung. The FFD algorithm consists of a deformation model, an objective function, and an optimization scheme [23]. The deformation model is based on a lattice of control points overlaid onto the reference image. The control point spacing for this lattice was set to a width of 5 voxels. The positional information for the deformation model is calculated using a cubic B-spline interpolation. The objective function uses normalized mutual information (NMI) as a metric and bending energy (BE) as a penalty term. NMI is a similarity measure based on entropy, with the computation being performed through a joint-histogram. Its optimization leads to a maximization of the amount of information shared by both images. The BE is computed from the cubic B-spline deformation model and favors a smooth warping between both images. The optimization scheme was performed to achieve the highest degree of coincidence between the two images according to the results of the objective function. For this optimization, we set the number of iterations to 300 and the weight of the BE penalty term to 0.005%. The DIR was performed from the floating image (0%-4DCT) to the reference image (50%-4DCT), and the resulting displacement vector fields were then used to deform the dose distribution. 2.3.2. MIM Maestro In MIM Maestro, a rigid registration is initially applied, which is followed by a non-rigid registration [24]. In the present study, the rigid registration was defined automatically using the whole body as a starting point for the region of interest (ROI) for deformation. Non-rigid registration was performed using the intensity-based FFD algorithm [24,25]. The deformation model is based on a coarse-to-fine multi-resolution approach, and uses a grid of control points on the static image that is used to search for the best corresponding location in the target volume. For this stage, the user does not need to define whether the voxels should move more or less, as this is automatically determined by the initial rigid registration from which the deformation then proceeds. This intensity matching is constrained, meaning that the voxels are warped according to areas of similarity in and around the region that is being deformed. The image matching metric minimizes intensity differences between the two images. The optimization strategy uses a custom-modified gradient descent algorithm. The regularization minimizes the effects of noise and incorrect correspondence, while still providing a large degree of freedom for each control point, to ensure proper matching with the target volume. 2.3.3. Velocity AI The Velocity AI system can perform registration by means of a rigid and a non-rigid registration using a B-spline algorithm [25]. The B-spline algorithm only defines the deformation on a sparse lattice of nodes overlaid on the image, with the displacement at any voxel being obtained by interpolation from the closest lattice nodes. The B-spline algorithm is based on the Mattes formulation and consists of mutual information and some proprietary method that the vendor does not divulge. The limited-memory Broyden– Fletcher–Goldfarb–Shanno optimizer was used to find the optimal node value, with maximums of 100 iterations and 20 corrections being used as termination conditions for the optimization algorithm [26]. It should be noted that information regarding the parameter settings in this software is not available from the vendor; the user cannot change the parameter settings. However, users can set an ROI for the deformation. In the current study, the in-house and MIM Maestro systems used the whole body for deformation, and therefore, the whole body was set as the ROI for deformation with the Velocity AI system.
2.4. Evaluation of the three deformable image registration systems 2.4.1. Spatial accuracy of deformable image registration Before evaluating dose warping using the dosimetric parameters, the DSC was used to evaluate the spatial accuracy of the DIR. The DSC was calculated using the structure at 50% and the deformed structure at 0%. The DSC is defined as
DSC ¼
2jA \ Bj : jAj þ jBj
ð1Þ
The values of the DSC range from 0 to 1, and are identical to 1 if the A and B volumes are equal with a complete intersection [25]. 2.5. Accuracy of dose warping The accuracy of dose warping was evaluated using Dmean, gEUD, and the 3D gamma index. Additionally, correlation coefficients (R) were calculated between the DSC and the difference of the Dmean, and between the DSC and the difference of the gEUD. The Dmean,org was calculated from the DVH in the 0% plan, and the Dmean,def was calculated from the DVH of the deformed 0%dose distribution and 50% structures. The difference between the Dmean,org and the Dmean,def was measured. The gEUD was converted from multiple points of the DVH curve computed by the treatment planning system as follows:
gEUD ¼ ð
X 1=a V i Dai Þ
ð2Þ
i
where a is a model parameter specific to the normal structure or tumor of interest, and Vi represents the i-th partial volume receiving the dose Di in Gy [27,28]. In the present study, a value of 15 was used for a. The gEUDdef of the GTV from the deformed 0%dose distribution was compared with the gEUDorg from the original dose distribution. The gEUDorg and gEUDdef were also calculated, as were the Dmean,org and the Dmean,def, and the difference between the gEUDorg and the gEUDdef was measured. If the difference was zero, the accuracy of the dose warping could be considered as being very good. Conversely, the greater this difference was, the lower the accuracy of the dose warping. A 3D gamma index comparison was performed to assess the relative performance of the three DIR systems in terms of dosimetric and geometric accuracy. The passing rates were calculated with 3 mm/3% testing criteria and a 10% dose threshold. 3. Results 3.1. Analysis using spatial and dose parameters The DSCs for the three DIR systems are shown in Table 2. A threshold value for a good agreement was set at 0.6 [29]. Of the 19 cases, No. 5 (0.35), No. 9 (0.58), and No. 15 (0.53) resulted in values below the threshold when using Velocity AI. The differences in the Dmean and gEUD values in the three DIR systems are shown in Table 3. The differences in the Dmean and gEUD values in 16 cases showed good agreement, with the spatial agreement being within approximately 3.0%. Fig. 1 shows a typical example from cases with good DSC values and showing good agreement in Dmean and gEUD values, in situations involving: (a) a large respiratory movement; (b) a small target volume; and (c) areas of low image intensity contrast. In the cases showing DSC values less than the threshold value (No. 5, No. 9, and No. 15), as shown in Fig. 2, the differences in Dmean were 11.8%, 4.78% and 1.56% respectively. Additionally, the differences in the gEUD in the three cases were 50.1%, 29.7% and 1.87% respectively. The DSC value for case No. 15 was 0.53;
S. Moriya et al. / Physica Medica 37 (2017) 16–23 Table 2 Comparison of the three deformable image registration programs using the Dice similarity coefficient. Case no.
In-house using the NiftyReg
MIM Maestro
Velocity AI
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 Mean SD
0.89 0.89 0.89 0.88 0.84 0.85 0.82 0.73 0.82 0.70 0.90 0.86 0.93 0.90 0.89 0.86 0.92 0.91 0.90 0.86 0.06
0.86 0.89 0.88 0.87 0.80 0.84 0.91 0.79 0.83 0.81 0.87 0.87 0.93 0.86 0.92 0.92 0.92 0.92 0.88 0.87 0.04
0.80 0.84 0.73 0.82 0.35 0.79 0.82 0.79 0.58 0.65 0.69 0.82 0.93 0.77 0.53 0.84 0.89 0.83 0.80 0.75 0.14
19
however, the dose differences of the Dmean and gEUD were less than 2%. The differences in gEUD were significantly larger than the differences in Dmean. Fig. 3(a) and (b) show the relationships between DSC values and differences in the Dmean, and between DSC values and differences in the gEUD. For the in-house DIR system, the correlation coefficients between DSCs and differences in Dmean (R = 0.76), and between DSCs and differences in gEUD (R = 0.66) show a moderate linear relationship, as they also do for the MIM Maestro (R = 0.76 for Dmean and R = 0.66 for gEUD) and Velocity AI systems (R = 0.85 for Dmean and R = 0.80 for gEUD). From the results of the slopes of the fitted lines, gEUD is demonstrated as being more sensitive than Dmean for the indication of dose differences. 3.2. Situation analysis using dose parameters As shown in Table 3, in almost all cases, the three DIR systems showed a good agreement for gEUD and Dmean, with the accuracy of dose warping being good where the target motion was large, the target volume was small, and the tumor boundary was blurred, as shown in Fig. 1. The Dmean and gEUD for these examples showed a good agreement.
Fig. 1. Three examples of dose warping obtained using the three deformable image registration programs: (a) with a large respiratory movement; (b) with a small target volume; and (c) with an area of low image intensity contrast with good spatial and dosimetric agreement.
20
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Fig. 2. Example of DVH histograms and dose distributions of cases No. 5 (a), No. 9 (b) and No. 15 (c) with less spatial and dosimetric agreements.
We found a low correlation between the accuracy of dose warping and the magnitude of respiratory motion (R = 0.2 for the inhouse DIR system, R = 0.4 for the MIM Maestro, and R = 0.4 for the Velocity AI). Additionally, we also found a low correlation for volume size with the in-house DIR system with NiftyReg (R = 0.2), MIM Maestro (R = 0.4), and Velocity AI systems (R = 0.4). Seventeen of the 19 patients (with the exception of No. 5 and No. 9) in the current study were further divided into two groups, with group A consisting of patients where the tumor boundary was clear, and group B consisting of patients where the tumor boundary was unclear. The differences in gEUD for group A were 0.41% ± 0.41%, 0.34% ± 0.33%, and 1.10% ± 0.78% for the in-house, MIM Maestro, and Velocity AI systems respectively. The differences in the gEUD for group B were 0.64% ± 1.01%, 0.39% ± 0.42%, and 1.09% ± 0.89% for the in-house, MIM Maestro, and Velocity AI systems respectively. The statistical significance of the differences between group A and B was investigated using a t-test at the 1% level, and no significant difference was found between group A and B (p = 0.5, 0.8, and 1.0 for the in-house, MIM Maestro and Velocity AI systems respectively). Thus, there was no relationship
between the accuracy of dose warping and clarity of the tumor boundary. As shown above, almost all cases demonstrated good agreement for Dmean and gEUD; however, the differences in the gEUD in cases No. 5 (50%) and No. 9 (30%) were substantially larger with the Velocity AI than with the other two systems. These two cases shared a similar anatomical situation with their boundary being unclear as the tumors were attached to the chest wall (Fig. 2a and b). Case No. 15 showed poor coverage of the small tumor volume with the Velocity AI in comparison with the other two systems. Fig. 4 summarizes the 3D gamma index comparison for the three DIR systems, showing the path ratios in the deformed dose distributions for the in-house vs. MIM (NM), the in-house vs Velocity AI (NV), and MIM vs. Velocity AI (MV) systems. The path ratios of the NM, NV, and MV were 90.1 ± 7.3%, 87.4 ± 6.1%, and 84.6 ± 5.5% respectively. The in-house and MIM systems showed a better agreement for Dmean and gEUD, resulting in better agreement for the 3D gamma index evaluation, which showed a statistically significant difference between NM and MV.
21
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Fig. 3. Relationship between (a) DSC values and differences of Dmean, and between (b) DSC values and differences of gEUD.
Table 3 Comparison of the three deformable image registration programs using Dmean and gEUD. Case no.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
In-house using the NiftyReg
MIM Maestro
Velocity AI
Dmean [%]
gEUD [%]
Dmean [%]
gEUD [%]
Dmean [%]
gEUD [%]
0.53 0.46 0.32 0.59 0.71 0.20 0.36 0.59 1.00 2.42 0.03 0.44 0.09 0.15 0.12 0.04 0.20 0.02 0.25
0.40 0.51 0.49 0.60 0.90 0.24 0.14 0.65 2.90 3.09 0.03 1.34 0.02 0.24 0.15 0.19 0.22 0.00 0.47
0.08 0.17 0.27 0.49 0.66 0.19 0.24 0.31 1.05 1.06 0.19 0.41 0.14 0.45 0.02 0.11 0.13 0.00 0.26
0.14 0.12 0.54 0.60 0.98 0.17 0.38 0.43 2.75 1.33 0.03 0.89 0.08 0.76 0.05 0.20 0.15 0.02 0.30
0.75 0.67 1.51 0.26 11.76 0.55 0.11 0.50 4.78 1.55 1.23 0.70 0.03 1.47 1.56 0.59 0.17 0.90 0.01
0.71 1.05 2.38 0.44 50.05 0.61 0.15 0.62 29.68 2.02 2.30 1.29 0.03 2.32 1.87 0.96 0.60 1.19 0.10
22
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Fig. 4. Comparison in three-dimensional gamma index evaluation for the three DIR programs (Criteria: 3 mm/3%)
4. Discussion Dmean and gEUD can provide a dosimetric evaluation of the coverage of the GTV in comparisons between the reference dose distribution and the deformed dose distribution. The in-house DIR system using NiftyReg showed either a similar or higher geometrical and dosimetrical accuracy than the MIM and Velocity AI systems. In the 3D gamma index comparison, which showed a relative evaluation between different DIR systems, the in-house DIR system demonstrated a good agreement with MIM. Thus, the in-house system using NiftyReg could be considered to be well tuned to the lung tumor registration. Compared with the in-house and MIM Maestro systems, Velocity AI showed stronger linear relationships between DSC and Dmean and gEUD, as the accuracy of the Velocity AI system was distributed over a wider spatial and dose range. For the in-house and MIM Maestro systems, the DSC values were above 0.7 in all cases, with the differences in Dmean and gEUD being less than approximately 3.0%. However, the DSC values for Velocity AI were less than 0.7 in three cases (Minimum: 0.35), and the differences in Dmean and gEUD for Velocity AI were more than 3% in three cases (Maximum: 11.8% for Dmean and 50.1% for gEUD). If the respiratory tumor motion is large, the DIR system needs to perform a large deformation to match the two images. The magnitude of the deformation is a factor affecting the accuracy of the warped dose distribution [16], and a large deformation of the target may reduce the accuracy of this dose deformation. However, in the current study, only a low correlation was found between the accuracy of dose warping and the magnitude of the respiratory motion. Changes in the accuracy of dose deformation were more pronounced when the tumor volume was small, although we only found a low correlation between dose deformation and volume. This is because for larger tumor volumes, small geometric variations in the deformation have a smaller effect on the accuracy of dose deformation. In other words, the evaluation of dose warping for a large-dose coverage is relatively insensitive, which means that portions of large target volumes with a large variation are neglected more easily than those of small targets. In general, unexpected deformation tends to appear in low image contrast regions, as intensity metrics such as FFD or Bspline are less sensitive to displacement changes in low-contrast
regions [30–32]. However, we found no relationship between accuracy and the definition of tumor boundary in the present study. As demonstrated above, the three DIR systems used in our study could satisfactorily deform the dose distribution, and these DIR systems may represent an improvement over previously reported DIR systems. Compared to previous studies [16,30–32], the results from the three DIR systems examined in this study demonstrate improved accuracy, probably because of upgrades and modifications to the DIR algorithms. This is demonstrated by the fact that the three DIR systems did not show a lower correlation between the accuracy of dose warping and respiratory motion, tumor volume, and an unclear tumor boundary. However, we did encounter situations where the dose warping was poor. In two of our patients the Velocity AI system resulted in large discrepancies in dose deformation when motion was high and the tumor was small, or when the tumor boundary was unclear. The three DIR systems employed in the current study used the same B-spline deformation model. However, the optimization processes driving the deformation differed. The in-house DIR system used a structure-based optimization that operated in two limited regions, such as the lung and the whole body. With MIM Maestro, the region for optimization was automatically selected after a rigid registration of the whole body. Conversely, the Velocity AI system used the whole body to perform the optimization, and therefore, the presence of a large ROI may adversely affect the deformation optimization, with a loss of focus on the tumor being a possible factor behind the reduced accuracy of the dose warping. Thus, it may be difficult to predict when such situations will occur, and sometimes a poor result may be encountered. In terms of the ART and 4D dose calculation, dose warping is an essential requirement, and the accuracy of dose warping can affect the accuracy of the ART and 4D dose calculation. DSC provided a suitable index for measuring the accuracy of DIR, which is the degree of similarity in the anatomy. However, the DSC value is defined as the magnitude of the spatial accuracy. We found a strong relationship between the spatial accuracy factor of the DSC and the two factors Dmean and gEUD, which indicates that these two factors could be used for assessing the accuracy of dose warping and could show the magnitude of dose distribution discrepancies. A deformable gel phantom has been investigated in terms of the accuracy of gel dosimetry. In particular, the dosimetry
S. Moriya et al. / Physica Medica 37 (2017) 16–23
can only provide limited situations where the media in the phantom is equivalent. In the thoracic region, there are several situations which affect the accuracy of dose warping for ART and four-dimensional dose calculation, including a large respiratory movement, a small target volume, and an unclear target boundary [16,30–32]. Thus, our method should provide comprehensive evaluation of dose warping using Dmean and gEUD, and could also be adjusted to clinical situations and to show the magnitude of the impact of dose warping. 5. Conclusions The accuracy of dose warping performed using our specifically tuned in-house DIR system based on publicly available software was compared to two commercially available DIR systems using Dmean, gEUD, and 3D gamma index. Any of these three systems could potentially be used with high accuracy in most patients with lung cancer. The in-house DIR system with better tuning should be able to achieve higher spatial and dose accuracy than the commercially available systems. However, there may be errors of more than 5% and 30% for Dmean and gEUD respectively, even when commercially available DIR tools are used. Even if the same deformation model is used in different systems, the optimization algorithm for the deformation may be different, and different dose warping results may consequently be obtained. Care should be taken regarding not only the deformation model, but also the optimization process, especially in relation to clinical practice and where tumor respiratory motion is large, the tumor volume is small, and the tumor boundary is unclear. Additionally, dose warping accuracy is proportional to spatial accuracy; however, the magnitude of the accuracy of dose warping cannot always be predicted by the magnitude of the DSC value. Dmean and gEUD can demonstrate the accuracy and magnitude of dose warping, with gEUD being a more sensitive parameter than Dmean for the evaluation of target dose coverage from low to high-dose regions. These dosimetric parameters may be used to show the impact of dose warping. Acknowledgments This work was supported by JSPS KAKENHI Grant Number 26713022. References [1] Graf R, Boehmer D, Nadobny J, Budach V, Wust P. Appropriate patient instructions can reduce prostate motion. Radiat Oncol 2012;7:125. [2] Michalski A, Atyeo J, Cox J, Rinks M. Inter- and intra-fraction motion during radiation therapy to the whole breast in the supine position: a systematic review. J Med Imaging Radiat Oncol 2012;56:499–509. [3] van de Bunt L, van der Heide UA, Ketelaars M, de Kort GA, Jurgenliemk-Schulz IM. Conventional, conformal, and intensity-modulated radiation therapy treatment planning of external beam radiotherapy for cervical cancer: the impact of tumor regression. Int J Radiat Oncol Biol Phys 2006;64:189–96. [4] Gomez DR, Chang JY. Adaptive radiation for lung cancer. J Oncol 2011;2011. [5] Men C, Jia X, Jiang SB. GPU-based ultra-fast direct aperture optimization for online adaptive radiation therapy. Phys Med Biol 2010;55:4309–19. [6] Li X, Yang Y, Li T, Fallon K, Heron DE, Huq MS. Dosimetric effect of respiratory motion on volumetric-modulated arc therapy-based lung SBRT treatment delivered by TrueBeam machine with flattening filter-free beam. J Appl Clin Med Phys 2013;14:4370.
23
[7] Chan MK, Kwong DL, Ng SC, Tong AS, Tam EK. Experimental evaluations of the accuracy of 3D and 4D planning in robotic tracking stereotactic body radiotherapy for lung cancers. Med Phys 2013;40:041712. [8] Jung SH, Yoon SM, Park SH, Cho B, Park JW, Jung J, et al. Four-dimensional dose evaluation using deformable image registration in radiotherapy for liver cancer. Med Phys 2013;40:011706. [9] Manescu P, Ladjal H, Azencot J, Beuve M, Testa E, Shariat B. Four-dimensional radiotherapeutic dose calculation using biomechanical respiratory motion description. Int J Comput Assist Radiol Surg 2014;9:449–57. [10] Hardcastle N, van Elmpt W, De Ruysscher D, Bzdusek K, Tome WA. Accuracy of deformable image registration for contour propagation in adaptive lung radiotherapy. Radiat Oncol 2013;8:243. [11] Latifi K, Huang TC, Feygelman V, Budzevich MM, Moros EG, Dilling TJ, et al. Effects of quantum noise in 4D-CT on deformable image registration and derived ventilation data. Phys Med Biol 2013;58:7661–72. [12] Akino Y, Yoshioka Y, Fukuda S, Maruoka S, Takahashi Y, Yagi M, et al. Estimation of rectal dose using daily megavoltage cone-beam computed tomography and deformable image registration. Int J Radiat Oncol Biol Phys 2013;87:602–8. [13] Kim J, Kumar S, Liu C, Zhong H, Pradhan D, Shah M, et al. A novel approach for establishing benchmark CBCT/CT deformable image registrations in prostate cancer radiotherapy. Phys Med Biol 2013;58:8077–97. [14] Peroni M, Spadea MF, Riboldi M, Falcone S, Vaccaro C, Sharp GC, et al. Validation of automatic contour propagation for 4D treatment planning using multiple metrics. Technol Cancer Res Treat 2013;12:501–10. [15] Yeo UJ, Taylor ML, Dunn L, Kron T, Smith RL, Franich RD. A novel methodology for 3D deformable dosimetry. Med Phys 2012;39:2203–13. [16] Yeo UJ, Taylor ML, Supple JR, Smith RL, Dunn L, Kron T, et al. Is it sensible to ‘‘deform” dose? 3D experimental validation of dose-warping. Med Phys 2012;39:5065–72. [17] Yeo UJ, Supple JR, Taylor ML, Smith R, Kron T, Franich RD. Performance of 12 DIR algorithms in low-contrast regions for mass and density conserving deformation. Med Phys 2013;40:101701. [18] Yamamoto T, Kabus S, Lorenz C, Mittra E, Hong JC, Chung M, et al. Pulmonary ventilation imaging based on 4-dimensional computed tomography: comparison with pulmonary function tests and SPECT ventilation images. Int J Radiat Oncol Biol Phys 2014;90:414–22. [19] Suga K. Technical and analytical advances in pulmonary ventilation SPECT with xenon-133 gas and Tc-99m-Technegas. Ann Nucl Med 2002;16:303–10. [20] Louie AV, Rodrigues G, Olsthoorn J, Palma D, Yu E, Yaremko B, et al. Interobserver and intra-observer reliability for lung cancer target volume delineation in the 4D-CT era. Radiother Oncol 2010;95:166–71. [21] Modat M, McClelland J, Ourselin S. Lung registration using the NiftyReg package. Medical Image Analysis for the Clinic-A Grand Challenge, 2010. p. 33–42. [22] Modat M, Ridgway GR, Taylor ZA, Lehmann M, Barnes J, Hawkes DJ, et al. Fast free-form deformation using graphics processing units. Comput Methods Programs Biomed 2010;98:278–84. [23] Jon Piper, Aaron Nelson M, Harper J. Deformable Image Registration in MIM MaestroTM Evaluation and Description. MIM Software Inc 252800 Science Park Drive, Suite 180. 2013:1–5. [24] La Macchia M, Fellin F, Amichetti M, Cianchetti M, Gianolini S, Paola V, et al. Systematic evaluation of three different commercial software solutions for automatic segmentation for adaptive therapy in head-and-neck, prostate and pleural cancer. Radiat Oncol 2012;7:160. [25] Kadoya N, Fujita Y, Katsuta Y, Dobashi S, Takeda K, Kishi K, et al. Evaluation of various deformable image registration algorithms for thoracic images. J Radiat Res 2014;55:175–82. [26] Gay HA, Niemierko A. A free program for calculating EUD-based NTCP and TCP in external beam radiotherapy. Phys Med 2007;23:115–25. [27] Li XA, Alber M, Deasy JO, Jackson A, Jee KWK, Marks LB, et al. The use and QA of biologically related models for treatment planning: Short report of the TG-166 of the therapy physics committee of the AAPM(a)). Med Phys 2012;39:1386–409. [28] Akbarzadeh A, Gutierrez D, Baskin A, Ay MR, Ahmadian A, Riahi Alam N, et al. Evaluation of whole-body MR to CT deformable image registration. J Appl Clin Med Phys 2013;14:4163. [29] Zhong H, Kim J, Li H, Nurushev T, Movsas B, Chetty IJ. A finite element method to correct deformable image registration errors in low-contrast regions. Phys Med Biol 2012;57:3499–515. [30] Gu X, Dong B, Wang J, Yordy J, Mell L, Jia X, et al. A contour-guided deformable image registration algorithm for adaptive radiotherapy. Phys Med Biol 2013;58:1889–901. [31] Zhong H, Kim J, Chetty IJ. Analysis of deformable image registration accuracy using computational modeling. Med Phys 2010;37:970–9.
Contents lists available at ScienceDirect
Physica Medica journal homepage: http://www.physicamedica.com
Original paper
Dose warping performance in deformable image registration in lung Shunsuke Moriya a, Hidenobu Tachibana b,⇑, Nozomi Kitamura c, Amit Sawant d, Masanori Sato a a
Radiological Sciences, Graduate Division of Health Sciences, Komazawa University, Tokyo 1548525, Japan Particle Therapy Division, Research Center for Innovative Oncology, National Cancer Center, Chiba 2778577, Japan c Department of Radiation Oncology, Cancer Institute Hospital of the Japanese Foundation of Cancer Research, Tokyo 1358550, Japan d Department of Radiation Oncology, University of Texas Southwestern Medical Center, 5801 Forest Park Rd., Dallas, TX 75390-9183, USA b
a r t i c l e
i n f o
Article history: Received 14 July 2016 Received in Revised form 13 February 2017 Accepted 20 March 2017
Keywords: Deformable image registration Generalized equivalent uniform dose Four-dimensional computed tomography Lung cancer
a b s t r a c t Purpose: It is unclear that spatial accuracy can reflect the impact of deformed dose distribution. In this study, we used dosimetric parameters to compare an in-house deformable image registration (DIR) system using NiftyReg, with two commercially available systems, MIM Maestro (MIM) and Velocity AI (Velocity). Methods: For 19 non-small-cell lung cancer patients, the peak inspiration (0%)-4DCT images were deformed to the peak expiration (50%)-4DCT images using each of the three DIR systems, which included computation of the deformation vector fields (DVF). The 0%-gross tumor volume (GTV) and the 0%-dose distribution were also then deformed using the DVFs. The agreement in the dose distributions for the GTVs was evaluated using generalized equivalent uniform dose (gEUD), mean dose (Dmean), and threedimensional (3D) gamma index (criteria: 3 mm/3%). Additionally, a Dice similarity coefficient (DSC) was used to measure the similarity of the GTV volumes. Results: Dmean and gEUD demonstrated good agreement between the original and deformed dose distributions (differences were generally less than 3%) in 17 of the patients. In two other patients, the Velocity system resulted in differences in gEUD of 50.1% and 29.7% and in Dmean of 11.8% and 4.78%. The gamma index comparison showed statistically significant differences for the in-house DIR vs. MIM, and MIM vs. Velocity. Conclusions: The finely tuned in-house DIR system could achieve similar spatial and dose accuracy to the commercial systems. Care must be taken, as we found errors of more than 5% for Dmean and 30% for gEUD, even with a commercially available DIR tool. Ó 2017 Associazione Italiana di Fisica Medica. Published by Elsevier Ltd. All rights reserved.
1. Introduction The locations and volumes of organs inside a patient’s body can change significantly during irradiation and between treatments [1,2]. Anatomical changes over an entire course of treatment can compromise the value of the initial treatment plan and therefore the treatment outcomes [3]. Adaptive radiotherapy (ART) has therefore been proposed to overcome this challenge; it allows the plan quality to be maintained by modification or reoptimization of the treatment plan according to changes in a patient’s anatomy [4,5]. Inter- and intra-fraction organ motion is especially significant in lung radiotherapy, and ART can provide a more accurate final dose distribution from the information on the accumulated dose distribution. This is obtained using ⇑ Corresponding author. E-mail addresses: [email protected] (S. Moriya), [email protected] (H. Tachibana), [email protected] (N. Kitamura), Amit.Sawant@ UTSouthwestern.edu (A. Sawant), [email protected] (M. Sato).
fraction-by-fraction dose distributions deformed according to a reference transformation between the planning CT and treatment CT. Deformable image registration (DIR) is therefore an essential tool for ART. The use of four-dimensional computed tomography (4DCT) has also facilitated treatment planning by allowing the respiratory motion of the target and critical organs to be incorporated into the analysis. The tumor volume and shape can be estimated more accurately in a phased 4DCT image than they can from a maximum or mean intensity projection formed from the 4DCT scan. A typical planning also usually considers only a static tumor, even though the tumor will move during the irradiation. An accumulated 4D dose calculation is therefore useful for achieving a more realistic estimation of the dose delivered to the tumor in the lung and the surrounding organs [6–9]. Information on the delivered 4D dose accumulation can help the physician evaluate the treatment and decide when and how a re-plan should be performed. Accurate dose warping is therefore one of the most important factors in the ART methodology.
http://dx.doi.org/10.1016/j.ejmp.2017.03.016 1120-1797/Ó 2017 Associazione Italiana di Fisica Medica. Published by Elsevier Ltd. All rights reserved.
17
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Many investigators have reported on methods for evaluation of the accuracy of DIR, including the use of the Dice similarity coefficient (DSC), landmark comparisons, mean slice-wise Hausdorff distance to agreement, volume differences, receiver operating characteristics, and target registration error [10–14]. However, few studies have directly evaluated dose deformation. Yeo et al. evaluated the accuracy of dose warping using DEFGEL, which is a deformable three-dimensional-dosimetry gel phantom [15–17]. Dosimetry using DEFGEL can provide direct measurements of dose, taking into consideration deformable media. However, the dosimetry information available from gel phantoms is limited, as they lack geometric complexity and are unable to provide assessments for individual patient cases. Additionally, dosimetry with deformable gel phantoms requires specialized and expensive equipment. In addition to dose monitoring in the lung, DIR may be applied in various other fields, such as radiotherapy and nuclear medicine [18,19]. There are several publicly available DIR systems, which may be suitable for contributing to clinical practice and the research environment, provided that the system is at least comparable to commercially available DIR systems. The purpose of this study was twofold; the first part was to build an in-house DIR system with a graphical user interface and image registration using the publically available NiftyReg (freedownloadable software library package, NiftyReg), while the second part of the study compares the in-house NiftyReg DIR system with two commercially available systems, MIM Maestro (MIM Software Inc., Cleveland, OH, USA) and Velocity AI (Varian Medical Systems, Palo Alto, CA, USA), with the evaluations being made using mean dose (Dmean), generalized equivalent uniform dose (gEUD), and three-dimensional (3D) gamma index evaluation. These dosimetric parameters are common in clinical practice and are easy to use, with a previous study having applied the generalized equivalent uniform dose (gEUD) as an index to quantify plan quality [20]. Dmean and gEUD were investigated as indices of the accuracy of dose warping using the DSC representing the spatial accuracy of image registration under the DIR. 2. Materials and methods 2.1. 4D-CT scan acquisition and patients All 4DCT scans were acquired on a 4-slice clinical scanner (Brilliance CT Big Bore, Philips Healthcare, Andover, MA). The respiratory cycle signal was monitored using an abdominal bellows pressure belt system (Philips Medical Systems). Other imaging parameters (e.g. tube voltage [kV] and tube current exposure time product value [mAs]) were set according to site-specific standard imaging protocols. CT data were reconstructed with a field of view of 50 cm on a 512 512 grid with a slice thickness of 2.0 mm. The longitudinal scan length was determined on a scout view. Nineteen patients treated with lung stereotactic body radiation therapy underwent 4DCT scans. The peak inspiration (0%) and peak expiration phases (50%) of the 4DCT image datasets were used in this study. Tumor displacement was defined as the vector sum from the centroid of the 0%gross tumor volume (GTV) to the centroid of the 50%-GTV. The tumor displacement varied from 0.30–1.55 cm across cases. The volumes of the gross tumor volume (GTV) ranged from 0.89– 99.0 cm3 for the 0% phase and 1.06–99.0 cm3 for the 50% phase (Table 1). The mean percentage difference between the GTVs at the 0% and 50% phases was 0.98 ± 0.38% (Maximum: 28.2%). 2.2. Planning and dose warping Treatment planning (structure contouring, treatment field setting, and dose calculation) was performed on the 0%-4DCT image
Table 1 Tumor displacement, GTV obtained at the peak inspiration phase (0%) and the peak expiration phase (50%), and the percentage difference of the volume change in the 19 patients evaluated. Case no.
Tumor displacement [cm]
GTV at 0% (inspiration) phase [cm3]
GTV at 50% (expiration) phase [cm3]
Percentage difference of volume change [%]
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
0.62 1.55 1.51 0.86 1.37 1.22 1.36 0.30 1.54 0.82 1.12 0.96 0.38 0.91 1.07 0.61 1.01 0.77 0.62
2.46 60.81 29.46 9.30 5.35 9.13 25.30 1.58 13.50 0.89 16.32 14.46 99.00 4.58 4.37 4.25 92.43 8.32 14.07
2.50 60.55 29.28 8.19 6.86 9.57 28.13 1.19 13.98 1.06 16.07 16.06 99.84 4.40 3.85 4.03 85.84 8.33 11.62
1.6 0.4 0.6 11.9 28.2 4.8 11.2 24.7 3.6 19.1 1.5 11.1 0.8 3.9 11.9 5.2 7.1 0.1 17.4
dataset using Eclipse treatment planning software (Version 10.0, Varian Medical Systems, Palo Alto, CA, USA). Radiation oncologists were asked to draw contours defining the GTVs on the 50%-4DCT image dataset. A medical physicist with five-years of experience delineated the GTVs on the 0%-4DCT image datasets, under the observation of a radiation oncologist who also checked them all. The mean DSC value representing the inter-observer variability for GTV volume on the 0%-4DCT images was 0.85 ± 0.05 (Range: 0.73–0.92), a value that is comparable to the inter-observer delineation variation found in a previous study [21]. The clinical target volume (CTV) was set as equivalent to the GTV, and the planning target volume (PTV) was created by adding a 5 mm margin to the CTV in all directions. The number of treatment fields ranged from 9 to 13. The aperture shapes of the multileaf collimator were adjusted in the beam’s eye view to cover the PTV with an additional 5 mm margin, to take into consideration the inaccuracy of dose calculation in the border region between the tumor and lung. All plans were generated with the prescription dose (54 Gy/3 fr) covering 95% of the PTV. A dose calculation grid size of 2.5 mm was used for the anisotropic analytical algorithm (AAA) calculations, and a dose volume histogram (DVH) of the GTV was subsequently computed. The 0%-4DCT image dataset was then deformed to the 50%-4DCT image dataset using the three DIR systems, thereby facilitating the computation of the deformation vector field (DVF) for each system. The dose distribution in the 0%4DCT was also deformed according to each of the three DVFs obtained. The deformed dose distributions were then superimposed onto the 50%-4DCT image datasets with the corresponding contoured structures, and the DVHs of the GTVs were computed. 2.3. Dose warping algorithm 2.3.1. In-house deformable image registration system using NiftyReg The in-house DIR system included a graphical user interface and image registration using NiftyReg. The image registration process consisted of three major steps. In the first step, a rigid registration was performed using a block-matching algorithm for the whole body of the patient [22]. The second step entailed a deformation inside the lung using a fast free-form deformation (FFD) algorithm, with the third step involving a larger-scale deformation inside the body. The second step focused on the tumor and lung, which have
18
S. Moriya et al. / Physica Medica 37 (2017) 16–23
the greatest motion, while the third step focused on whole organs, as deformation occurs not only in the tumor and lung, but also in the surrounding organs, such as the heart and liver. Thus, two deformation processes were performed inside the lung. The FFD algorithm consists of a deformation model, an objective function, and an optimization scheme [23]. The deformation model is based on a lattice of control points overlaid onto the reference image. The control point spacing for this lattice was set to a width of 5 voxels. The positional information for the deformation model is calculated using a cubic B-spline interpolation. The objective function uses normalized mutual information (NMI) as a metric and bending energy (BE) as a penalty term. NMI is a similarity measure based on entropy, with the computation being performed through a joint-histogram. Its optimization leads to a maximization of the amount of information shared by both images. The BE is computed from the cubic B-spline deformation model and favors a smooth warping between both images. The optimization scheme was performed to achieve the highest degree of coincidence between the two images according to the results of the objective function. For this optimization, we set the number of iterations to 300 and the weight of the BE penalty term to 0.005%. The DIR was performed from the floating image (0%-4DCT) to the reference image (50%-4DCT), and the resulting displacement vector fields were then used to deform the dose distribution. 2.3.2. MIM Maestro In MIM Maestro, a rigid registration is initially applied, which is followed by a non-rigid registration [24]. In the present study, the rigid registration was defined automatically using the whole body as a starting point for the region of interest (ROI) for deformation. Non-rigid registration was performed using the intensity-based FFD algorithm [24,25]. The deformation model is based on a coarse-to-fine multi-resolution approach, and uses a grid of control points on the static image that is used to search for the best corresponding location in the target volume. For this stage, the user does not need to define whether the voxels should move more or less, as this is automatically determined by the initial rigid registration from which the deformation then proceeds. This intensity matching is constrained, meaning that the voxels are warped according to areas of similarity in and around the region that is being deformed. The image matching metric minimizes intensity differences between the two images. The optimization strategy uses a custom-modified gradient descent algorithm. The regularization minimizes the effects of noise and incorrect correspondence, while still providing a large degree of freedom for each control point, to ensure proper matching with the target volume. 2.3.3. Velocity AI The Velocity AI system can perform registration by means of a rigid and a non-rigid registration using a B-spline algorithm [25]. The B-spline algorithm only defines the deformation on a sparse lattice of nodes overlaid on the image, with the displacement at any voxel being obtained by interpolation from the closest lattice nodes. The B-spline algorithm is based on the Mattes formulation and consists of mutual information and some proprietary method that the vendor does not divulge. The limited-memory Broyden– Fletcher–Goldfarb–Shanno optimizer was used to find the optimal node value, with maximums of 100 iterations and 20 corrections being used as termination conditions for the optimization algorithm [26]. It should be noted that information regarding the parameter settings in this software is not available from the vendor; the user cannot change the parameter settings. However, users can set an ROI for the deformation. In the current study, the in-house and MIM Maestro systems used the whole body for deformation, and therefore, the whole body was set as the ROI for deformation with the Velocity AI system.
2.4. Evaluation of the three deformable image registration systems 2.4.1. Spatial accuracy of deformable image registration Before evaluating dose warping using the dosimetric parameters, the DSC was used to evaluate the spatial accuracy of the DIR. The DSC was calculated using the structure at 50% and the deformed structure at 0%. The DSC is defined as
DSC ¼
2jA \ Bj : jAj þ jBj
ð1Þ
The values of the DSC range from 0 to 1, and are identical to 1 if the A and B volumes are equal with a complete intersection [25]. 2.5. Accuracy of dose warping The accuracy of dose warping was evaluated using Dmean, gEUD, and the 3D gamma index. Additionally, correlation coefficients (R) were calculated between the DSC and the difference of the Dmean, and between the DSC and the difference of the gEUD. The Dmean,org was calculated from the DVH in the 0% plan, and the Dmean,def was calculated from the DVH of the deformed 0%dose distribution and 50% structures. The difference between the Dmean,org and the Dmean,def was measured. The gEUD was converted from multiple points of the DVH curve computed by the treatment planning system as follows:
gEUD ¼ ð
X 1=a V i Dai Þ
ð2Þ
i
where a is a model parameter specific to the normal structure or tumor of interest, and Vi represents the i-th partial volume receiving the dose Di in Gy [27,28]. In the present study, a value of 15 was used for a. The gEUDdef of the GTV from the deformed 0%dose distribution was compared with the gEUDorg from the original dose distribution. The gEUDorg and gEUDdef were also calculated, as were the Dmean,org and the Dmean,def, and the difference between the gEUDorg and the gEUDdef was measured. If the difference was zero, the accuracy of the dose warping could be considered as being very good. Conversely, the greater this difference was, the lower the accuracy of the dose warping. A 3D gamma index comparison was performed to assess the relative performance of the three DIR systems in terms of dosimetric and geometric accuracy. The passing rates were calculated with 3 mm/3% testing criteria and a 10% dose threshold. 3. Results 3.1. Analysis using spatial and dose parameters The DSCs for the three DIR systems are shown in Table 2. A threshold value for a good agreement was set at 0.6 [29]. Of the 19 cases, No. 5 (0.35), No. 9 (0.58), and No. 15 (0.53) resulted in values below the threshold when using Velocity AI. The differences in the Dmean and gEUD values in the three DIR systems are shown in Table 3. The differences in the Dmean and gEUD values in 16 cases showed good agreement, with the spatial agreement being within approximately 3.0%. Fig. 1 shows a typical example from cases with good DSC values and showing good agreement in Dmean and gEUD values, in situations involving: (a) a large respiratory movement; (b) a small target volume; and (c) areas of low image intensity contrast. In the cases showing DSC values less than the threshold value (No. 5, No. 9, and No. 15), as shown in Fig. 2, the differences in Dmean were 11.8%, 4.78% and 1.56% respectively. Additionally, the differences in the gEUD in the three cases were 50.1%, 29.7% and 1.87% respectively. The DSC value for case No. 15 was 0.53;
S. Moriya et al. / Physica Medica 37 (2017) 16–23 Table 2 Comparison of the three deformable image registration programs using the Dice similarity coefficient. Case no.
In-house using the NiftyReg
MIM Maestro
Velocity AI
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 Mean SD
0.89 0.89 0.89 0.88 0.84 0.85 0.82 0.73 0.82 0.70 0.90 0.86 0.93 0.90 0.89 0.86 0.92 0.91 0.90 0.86 0.06
0.86 0.89 0.88 0.87 0.80 0.84 0.91 0.79 0.83 0.81 0.87 0.87 0.93 0.86 0.92 0.92 0.92 0.92 0.88 0.87 0.04
0.80 0.84 0.73 0.82 0.35 0.79 0.82 0.79 0.58 0.65 0.69 0.82 0.93 0.77 0.53 0.84 0.89 0.83 0.80 0.75 0.14
19
however, the dose differences of the Dmean and gEUD were less than 2%. The differences in gEUD were significantly larger than the differences in Dmean. Fig. 3(a) and (b) show the relationships between DSC values and differences in the Dmean, and between DSC values and differences in the gEUD. For the in-house DIR system, the correlation coefficients between DSCs and differences in Dmean (R = 0.76), and between DSCs and differences in gEUD (R = 0.66) show a moderate linear relationship, as they also do for the MIM Maestro (R = 0.76 for Dmean and R = 0.66 for gEUD) and Velocity AI systems (R = 0.85 for Dmean and R = 0.80 for gEUD). From the results of the slopes of the fitted lines, gEUD is demonstrated as being more sensitive than Dmean for the indication of dose differences. 3.2. Situation analysis using dose parameters As shown in Table 3, in almost all cases, the three DIR systems showed a good agreement for gEUD and Dmean, with the accuracy of dose warping being good where the target motion was large, the target volume was small, and the tumor boundary was blurred, as shown in Fig. 1. The Dmean and gEUD for these examples showed a good agreement.
Fig. 1. Three examples of dose warping obtained using the three deformable image registration programs: (a) with a large respiratory movement; (b) with a small target volume; and (c) with an area of low image intensity contrast with good spatial and dosimetric agreement.
20
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Fig. 2. Example of DVH histograms and dose distributions of cases No. 5 (a), No. 9 (b) and No. 15 (c) with less spatial and dosimetric agreements.
We found a low correlation between the accuracy of dose warping and the magnitude of respiratory motion (R = 0.2 for the inhouse DIR system, R = 0.4 for the MIM Maestro, and R = 0.4 for the Velocity AI). Additionally, we also found a low correlation for volume size with the in-house DIR system with NiftyReg (R = 0.2), MIM Maestro (R = 0.4), and Velocity AI systems (R = 0.4). Seventeen of the 19 patients (with the exception of No. 5 and No. 9) in the current study were further divided into two groups, with group A consisting of patients where the tumor boundary was clear, and group B consisting of patients where the tumor boundary was unclear. The differences in gEUD for group A were 0.41% ± 0.41%, 0.34% ± 0.33%, and 1.10% ± 0.78% for the in-house, MIM Maestro, and Velocity AI systems respectively. The differences in the gEUD for group B were 0.64% ± 1.01%, 0.39% ± 0.42%, and 1.09% ± 0.89% for the in-house, MIM Maestro, and Velocity AI systems respectively. The statistical significance of the differences between group A and B was investigated using a t-test at the 1% level, and no significant difference was found between group A and B (p = 0.5, 0.8, and 1.0 for the in-house, MIM Maestro and Velocity AI systems respectively). Thus, there was no relationship
between the accuracy of dose warping and clarity of the tumor boundary. As shown above, almost all cases demonstrated good agreement for Dmean and gEUD; however, the differences in the gEUD in cases No. 5 (50%) and No. 9 (30%) were substantially larger with the Velocity AI than with the other two systems. These two cases shared a similar anatomical situation with their boundary being unclear as the tumors were attached to the chest wall (Fig. 2a and b). Case No. 15 showed poor coverage of the small tumor volume with the Velocity AI in comparison with the other two systems. Fig. 4 summarizes the 3D gamma index comparison for the three DIR systems, showing the path ratios in the deformed dose distributions for the in-house vs. MIM (NM), the in-house vs Velocity AI (NV), and MIM vs. Velocity AI (MV) systems. The path ratios of the NM, NV, and MV were 90.1 ± 7.3%, 87.4 ± 6.1%, and 84.6 ± 5.5% respectively. The in-house and MIM systems showed a better agreement for Dmean and gEUD, resulting in better agreement for the 3D gamma index evaluation, which showed a statistically significant difference between NM and MV.
21
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Fig. 3. Relationship between (a) DSC values and differences of Dmean, and between (b) DSC values and differences of gEUD.
Table 3 Comparison of the three deformable image registration programs using Dmean and gEUD. Case no.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
In-house using the NiftyReg
MIM Maestro
Velocity AI
Dmean [%]
gEUD [%]
Dmean [%]
gEUD [%]
Dmean [%]
gEUD [%]
0.53 0.46 0.32 0.59 0.71 0.20 0.36 0.59 1.00 2.42 0.03 0.44 0.09 0.15 0.12 0.04 0.20 0.02 0.25
0.40 0.51 0.49 0.60 0.90 0.24 0.14 0.65 2.90 3.09 0.03 1.34 0.02 0.24 0.15 0.19 0.22 0.00 0.47
0.08 0.17 0.27 0.49 0.66 0.19 0.24 0.31 1.05 1.06 0.19 0.41 0.14 0.45 0.02 0.11 0.13 0.00 0.26
0.14 0.12 0.54 0.60 0.98 0.17 0.38 0.43 2.75 1.33 0.03 0.89 0.08 0.76 0.05 0.20 0.15 0.02 0.30
0.75 0.67 1.51 0.26 11.76 0.55 0.11 0.50 4.78 1.55 1.23 0.70 0.03 1.47 1.56 0.59 0.17 0.90 0.01
0.71 1.05 2.38 0.44 50.05 0.61 0.15 0.62 29.68 2.02 2.30 1.29 0.03 2.32 1.87 0.96 0.60 1.19 0.10
22
S. Moriya et al. / Physica Medica 37 (2017) 16–23
Fig. 4. Comparison in three-dimensional gamma index evaluation for the three DIR programs (Criteria: 3 mm/3%)
4. Discussion Dmean and gEUD can provide a dosimetric evaluation of the coverage of the GTV in comparisons between the reference dose distribution and the deformed dose distribution. The in-house DIR system using NiftyReg showed either a similar or higher geometrical and dosimetrical accuracy than the MIM and Velocity AI systems. In the 3D gamma index comparison, which showed a relative evaluation between different DIR systems, the in-house DIR system demonstrated a good agreement with MIM. Thus, the in-house system using NiftyReg could be considered to be well tuned to the lung tumor registration. Compared with the in-house and MIM Maestro systems, Velocity AI showed stronger linear relationships between DSC and Dmean and gEUD, as the accuracy of the Velocity AI system was distributed over a wider spatial and dose range. For the in-house and MIM Maestro systems, the DSC values were above 0.7 in all cases, with the differences in Dmean and gEUD being less than approximately 3.0%. However, the DSC values for Velocity AI were less than 0.7 in three cases (Minimum: 0.35), and the differences in Dmean and gEUD for Velocity AI were more than 3% in three cases (Maximum: 11.8% for Dmean and 50.1% for gEUD). If the respiratory tumor motion is large, the DIR system needs to perform a large deformation to match the two images. The magnitude of the deformation is a factor affecting the accuracy of the warped dose distribution [16], and a large deformation of the target may reduce the accuracy of this dose deformation. However, in the current study, only a low correlation was found between the accuracy of dose warping and the magnitude of the respiratory motion. Changes in the accuracy of dose deformation were more pronounced when the tumor volume was small, although we only found a low correlation between dose deformation and volume. This is because for larger tumor volumes, small geometric variations in the deformation have a smaller effect on the accuracy of dose deformation. In other words, the evaluation of dose warping for a large-dose coverage is relatively insensitive, which means that portions of large target volumes with a large variation are neglected more easily than those of small targets. In general, unexpected deformation tends to appear in low image contrast regions, as intensity metrics such as FFD or Bspline are less sensitive to displacement changes in low-contrast
regions [30–32]. However, we found no relationship between accuracy and the definition of tumor boundary in the present study. As demonstrated above, the three DIR systems used in our study could satisfactorily deform the dose distribution, and these DIR systems may represent an improvement over previously reported DIR systems. Compared to previous studies [16,30–32], the results from the three DIR systems examined in this study demonstrate improved accuracy, probably because of upgrades and modifications to the DIR algorithms. This is demonstrated by the fact that the three DIR systems did not show a lower correlation between the accuracy of dose warping and respiratory motion, tumor volume, and an unclear tumor boundary. However, we did encounter situations where the dose warping was poor. In two of our patients the Velocity AI system resulted in large discrepancies in dose deformation when motion was high and the tumor was small, or when the tumor boundary was unclear. The three DIR systems employed in the current study used the same B-spline deformation model. However, the optimization processes driving the deformation differed. The in-house DIR system used a structure-based optimization that operated in two limited regions, such as the lung and the whole body. With MIM Maestro, the region for optimization was automatically selected after a rigid registration of the whole body. Conversely, the Velocity AI system used the whole body to perform the optimization, and therefore, the presence of a large ROI may adversely affect the deformation optimization, with a loss of focus on the tumor being a possible factor behind the reduced accuracy of the dose warping. Thus, it may be difficult to predict when such situations will occur, and sometimes a poor result may be encountered. In terms of the ART and 4D dose calculation, dose warping is an essential requirement, and the accuracy of dose warping can affect the accuracy of the ART and 4D dose calculation. DSC provided a suitable index for measuring the accuracy of DIR, which is the degree of similarity in the anatomy. However, the DSC value is defined as the magnitude of the spatial accuracy. We found a strong relationship between the spatial accuracy factor of the DSC and the two factors Dmean and gEUD, which indicates that these two factors could be used for assessing the accuracy of dose warping and could show the magnitude of dose distribution discrepancies. A deformable gel phantom has been investigated in terms of the accuracy of gel dosimetry. In particular, the dosimetry
S. Moriya et al. / Physica Medica 37 (2017) 16–23
can only provide limited situations where the media in the phantom is equivalent. In the thoracic region, there are several situations which affect the accuracy of dose warping for ART and four-dimensional dose calculation, including a large respiratory movement, a small target volume, and an unclear target boundary [16,30–32]. Thus, our method should provide comprehensive evaluation of dose warping using Dmean and gEUD, and could also be adjusted to clinical situations and to show the magnitude of the impact of dose warping. 5. Conclusions The accuracy of dose warping performed using our specifically tuned in-house DIR system based on publicly available software was compared to two commercially available DIR systems using Dmean, gEUD, and 3D gamma index. Any of these three systems could potentially be used with high accuracy in most patients with lung cancer. The in-house DIR system with better tuning should be able to achieve higher spatial and dose accuracy than the commercially available systems. However, there may be errors of more than 5% and 30% for Dmean and gEUD respectively, even when commercially available DIR tools are used. Even if the same deformation model is used in different systems, the optimization algorithm for the deformation may be different, and different dose warping results may consequently be obtained. Care should be taken regarding not only the deformation model, but also the optimization process, especially in relation to clinical practice and where tumor respiratory motion is large, the tumor volume is small, and the tumor boundary is unclear. Additionally, dose warping accuracy is proportional to spatial accuracy; however, the magnitude of the accuracy of dose warping cannot always be predicted by the magnitude of the DSC value. Dmean and gEUD can demonstrate the accuracy and magnitude of dose warping, with gEUD being a more sensitive parameter than Dmean for the evaluation of target dose coverage from low to high-dose regions. These dosimetric parameters may be used to show the impact of dose warping. Acknowledgments This work was supported by JSPS KAKENHI Grant Number 26713022. References [1] Graf R, Boehmer D, Nadobny J, Budach V, Wust P. Appropriate patient instructions can reduce prostate motion. Radiat Oncol 2012;7:125. [2] Michalski A, Atyeo J, Cox J, Rinks M. Inter- and intra-fraction motion during radiation therapy to the whole breast in the supine position: a systematic review. J Med Imaging Radiat Oncol 2012;56:499–509. [3] van de Bunt L, van der Heide UA, Ketelaars M, de Kort GA, Jurgenliemk-Schulz IM. Conventional, conformal, and intensity-modulated radiation therapy treatment planning of external beam radiotherapy for cervical cancer: the impact of tumor regression. Int J Radiat Oncol Biol Phys 2006;64:189–96. [4] Gomez DR, Chang JY. Adaptive radiation for lung cancer. J Oncol 2011;2011. [5] Men C, Jia X, Jiang SB. GPU-based ultra-fast direct aperture optimization for online adaptive radiation therapy. Phys Med Biol 2010;55:4309–19. [6] Li X, Yang Y, Li T, Fallon K, Heron DE, Huq MS. Dosimetric effect of respiratory motion on volumetric-modulated arc therapy-based lung SBRT treatment delivered by TrueBeam machine with flattening filter-free beam. J Appl Clin Med Phys 2013;14:4370.
23
[7] Chan MK, Kwong DL, Ng SC, Tong AS, Tam EK. Experimental evaluations of the accuracy of 3D and 4D planning in robotic tracking stereotactic body radiotherapy for lung cancers. Med Phys 2013;40:041712. [8] Jung SH, Yoon SM, Park SH, Cho B, Park JW, Jung J, et al. Four-dimensional dose evaluation using deformable image registration in radiotherapy for liver cancer. Med Phys 2013;40:011706. [9] Manescu P, Ladjal H, Azencot J, Beuve M, Testa E, Shariat B. Four-dimensional radiotherapeutic dose calculation using biomechanical respiratory motion description. Int J Comput Assist Radiol Surg 2014;9:449–57. [10] Hardcastle N, van Elmpt W, De Ruysscher D, Bzdusek K, Tome WA. Accuracy of deformable image registration for contour propagation in adaptive lung radiotherapy. Radiat Oncol 2013;8:243. [11] Latifi K, Huang TC, Feygelman V, Budzevich MM, Moros EG, Dilling TJ, et al. Effects of quantum noise in 4D-CT on deformable image registration and derived ventilation data. Phys Med Biol 2013;58:7661–72. [12] Akino Y, Yoshioka Y, Fukuda S, Maruoka S, Takahashi Y, Yagi M, et al. Estimation of rectal dose using daily megavoltage cone-beam computed tomography and deformable image registration. Int J Radiat Oncol Biol Phys 2013;87:602–8. [13] Kim J, Kumar S, Liu C, Zhong H, Pradhan D, Shah M, et al. A novel approach for establishing benchmark CBCT/CT deformable image registrations in prostate cancer radiotherapy. Phys Med Biol 2013;58:8077–97. [14] Peroni M, Spadea MF, Riboldi M, Falcone S, Vaccaro C, Sharp GC, et al. Validation of automatic contour propagation for 4D treatment planning using multiple metrics. Technol Cancer Res Treat 2013;12:501–10. [15] Yeo UJ, Taylor ML, Dunn L, Kron T, Smith RL, Franich RD. A novel methodology for 3D deformable dosimetry. Med Phys 2012;39:2203–13. [16] Yeo UJ, Taylor ML, Supple JR, Smith RL, Dunn L, Kron T, et al. Is it sensible to ‘‘deform” dose? 3D experimental validation of dose-warping. Med Phys 2012;39:5065–72. [17] Yeo UJ, Supple JR, Taylor ML, Smith R, Kron T, Franich RD. Performance of 12 DIR algorithms in low-contrast regions for mass and density conserving deformation. Med Phys 2013;40:101701. [18] Yamamoto T, Kabus S, Lorenz C, Mittra E, Hong JC, Chung M, et al. Pulmonary ventilation imaging based on 4-dimensional computed tomography: comparison with pulmonary function tests and SPECT ventilation images. Int J Radiat Oncol Biol Phys 2014;90:414–22. [19] Suga K. Technical and analytical advances in pulmonary ventilation SPECT with xenon-133 gas and Tc-99m-Technegas. Ann Nucl Med 2002;16:303–10. [20] Louie AV, Rodrigues G, Olsthoorn J, Palma D, Yu E, Yaremko B, et al. Interobserver and intra-observer reliability for lung cancer target volume delineation in the 4D-CT era. Radiother Oncol 2010;95:166–71. [21] Modat M, McClelland J, Ourselin S. Lung registration using the NiftyReg package. Medical Image Analysis for the Clinic-A Grand Challenge, 2010. p. 33–42. [22] Modat M, Ridgway GR, Taylor ZA, Lehmann M, Barnes J, Hawkes DJ, et al. Fast free-form deformation using graphics processing units. Comput Methods Programs Biomed 2010;98:278–84. [23] Jon Piper, Aaron Nelson M, Harper J. Deformable Image Registration in MIM MaestroTM Evaluation and Description. MIM Software Inc 252800 Science Park Drive, Suite 180. 2013:1–5. [24] La Macchia M, Fellin F, Amichetti M, Cianchetti M, Gianolini S, Paola V, et al. Systematic evaluation of three different commercial software solutions for automatic segmentation for adaptive therapy in head-and-neck, prostate and pleural cancer. Radiat Oncol 2012;7:160. [25] Kadoya N, Fujita Y, Katsuta Y, Dobashi S, Takeda K, Kishi K, et al. Evaluation of various deformable image registration algorithms for thoracic images. J Radiat Res 2014;55:175–82. [26] Gay HA, Niemierko A. A free program for calculating EUD-based NTCP and TCP in external beam radiotherapy. Phys Med 2007;23:115–25. [27] Li XA, Alber M, Deasy JO, Jackson A, Jee KWK, Marks LB, et al. The use and QA of biologically related models for treatment planning: Short report of the TG-166 of the therapy physics committee of the AAPM(a)). Med Phys 2012;39:1386–409. [28] Akbarzadeh A, Gutierrez D, Baskin A, Ay MR, Ahmadian A, Riahi Alam N, et al. Evaluation of whole-body MR to CT deformable image registration. J Appl Clin Med Phys 2013;14:4163. [29] Zhong H, Kim J, Li H, Nurushev T, Movsas B, Chetty IJ. A finite element method to correct deformable image registration errors in low-contrast regions. Phys Med Biol 2012;57:3499–515. [30] Gu X, Dong B, Wang J, Yordy J, Mell L, Jia X, et al. A contour-guided deformable image registration algorithm for adaptive radiotherapy. Phys Med Biol 2013;58:1889–901. [31] Zhong H, Kim J, Chetty IJ. Analysis of deformable image registration accuracy using computational modeling. Med Phys 2010;37:970–9.