Comparison of an in‐house hybrid DIR method to NiftyReg on CBCT and CT images for head and neck cancer

Abstract An in‐house hybrid deformable image registration (DIR) method, which combines free‐form deformation (FFD) and the viscous fluid registration method, is proposed. Its results on the planning computed tomography (CT) and the day 1 treatment cone‐beam CT (CBCT) image from 68 head and neck cancer patients are compared with the results of NiftyReg, which uses B‐spline FFD alone. Several similarity metrics, the target registration error (TRE) of annotated points, as well as the Dice similarity coefficient (DSC) and Hausdorff distance (HD) of the propagated organs at risk are employed to analyze their registration accuracy. According to quantitative analysis on mutual information, normalized cross‐correlation, and the absolute pixel value differences, the results of the proposed DIR are more similar to the CBCT images than the NiftyReg results. Smaller TRE of the annotated points is observed in the proposed method, and the overall mean TRE for the proposed method and NiftyReg was 2.34 and 2.98 mm, respectively (p < 0.001). The mean DSC in the larynx, spinal cord, oral cavity, mandible, and parotid given by the proposed method ranged from 0.78 to 0.91, significantly higher than the NiftyReg results (ranging from 0.77 to 0.90), and the HD was significantly lower compared to NiftyReg. Furthermore, the proposed method did not suffer from unrealistic deformations as the NiftyReg did in the visual evaluation. Meanwhile, the execution time of the proposed method was much higher than NiftyReg (96.98 ± 11.88 s vs. 4.60 ± 0.49 s). In conclusion, the in‐house hybrid method gave better accuracy and more stable performance than NiftyReg.


INTRODUCTION
In a typical adaptive radiotherapy (ART) process, a cone-beam computed tomography (CBCT) scan is usually performed to obtain three-dimensional patient information before treatment. 1 Since rigid registration is not capable of handling the anatomical variations caused by patient motion, organ filling, tumor growth, or Deformable image registration (DIR) has the ability to depict these local deformations and to estimate the geometric displacement of each image pixel. 3 It can be used for image fusion, 4 monitoring the deformation, 5 and dose mapping. 6 In the current study of ART, one common way is to deform the planning CT (pCT) image to a CBCT image, and then recalculate the dose based on the deformed image. Therefore, it is necessary to verify the performance of the DIR algorithms. Currently, there are the geometric, intensity, and deep learning approaches are used for DIR. [7][8][9] As an intensity approach, free-form deformation (FFD) assumes that the local deformation is manipulated by a set of control points. 10 The local deformation in each position is calculated as a weighted sum of the displacements of surrounding control points. 11 The weightings are related to the distances from the pixels to the surrounding control points. 12 FFD is efficient because it optimizes the deformation in control points instead of in each pixel, so it has been used in commercial software such as RayStation, MIM Maestro, and Velocity. [13][14][15][16] RayStation employs an anatomically constrained FFD algorithm, 15 in which a shape-based regularization term and a new penalty term are used to keep the deformation reasonable. Velocity employs a B-Spline FFD algorithm. 16 Although FFD is widely used, it does not guarantee the preservation of topology. 17 Different from FFD, the viscous fluid method assumes the image as a viscous fluid whose deformation is governed by the Navier-Stokes equation. 18,19 It performs better in topology preservation when the displacements are large. However, it costs more in computation time and great care should be taken in selecting the discretization solver for the partial differential equations. 20 As an open-source registration toolkit developed by the Centre of Medical Image Computing at University College London, NiftyReg is an FFD-based method and has been adopted by many research institutions to do medical image registration. 6,21,22 It has been confirmed that NiftyReg is able to map identical structures between the CBCT and CT image in most cases. 21-23 However, it is not able to accurately reproduce large or complex anatomical changes, even though it is possible to improve the NiftyReg performance through further optimization of the registration parameters. 21 There are many parameter options provided in NiftyReg and different combinations of parameter settings will lead to changes in the registration results. Therefore, it is possible to improve the NiftyReg results by optimizing the parameter settings.
To make use of the advantages of different DIR methods, many hybrid methods are proposed. [24][25][26] For example, the landmark technique was introduced in the FFD methods to solve the different transformations in bony structure and soft tissue. 24 In this study, an in-house hybrid method is proposed. It was first tried to use the control points in FFD registration instead of landmark points to accelerate the viscous fluid registration process. The proposed method is applied to transform the pCT images to the CBCT images from 68 patients with head and neck cancer. To verify whether in-house hybrid methods could improve the registration accuracy, multiple different evaluation methods (points, contours, metrics) are employed to analyze the registration accuracy, and the results of our in-house hybrid method are compared with those of NiftyReg. It is hoped that this method can be applied to the field of ART in the future.

Clinical data
The pCT and the first treatment CBCT data from 68 patients with head and neck cancer were collected for DIR performance evaluation. All the pCT was acquired with a voxel size of 1.27 mm × 1.27 mm × 3 mm by Siemens SOMATOM Definition AS. The treatment CBCT was reconstructed with a resolution of 0.51 mm × 0.51 mm × 1.99 mm by Varian TrueBeam.

Points and OARs delineation
To assess the warp field error near the bone area, five feature points (denoted as P1-P5) located in different cervical vertebrae are manually labeled in both the CBCT and pCT images of each patient in the same way. P1 and P2 are located on the right and left side of the same vertebra near the oral cavity, P3 and P4 are in the same vertebra near the larynx, and P5 is located in the spinous near the shoulder. Each point is delineated in the axial slice with a small circle. The points with the same index in CBCT and pCT correspond to the same physical position. The pixel position is calculated from the centroid of the small circle. To assess the deformation accuracy in other anatomical structures, five organs at risk (OARs) including the larynx, spinal cord, oral cavity, mandible, and parotid are automatically delineated by MedMind software (MedMind Technology) in both the CBCT and CT images and finally modified by a radiation oncologist.

Data preparation
For each patient, the CBCT data are regarded as the reference image, which would not be warped. Before DIR processing, the input pCT images were processed by the rigid registration. The input pCT for DIR was of the same size and same resolution as the treatment CBCT images. The rigid registration is done by MedMind software.

NiftyReg registration
NiftyReg contains tools for global and deformable image registration. The global registration includes rigid F I G U R E 1 Flowchart of in-house hybrid method registration and affine registration. The deformable registration implementation in NiftyReg uses FFD. Since the deformation is generated from the displacement of control points, its computational consumption can be easily adjusted by setting different numbers of control points. The GPU version of NiftyReg can be obtained from Github. 27 It is developed in C++ codes and CUDA. It uses multi-scale registration (three-level) based on mutual information (MI) optimization. The maximum optimization iteration is set as 30 for each level. The other parameters are set as the default.

In-house hybrid method
Due to the usage of interpolation techniques, NiftyReg has a limited ability to depict fine deformation. On the contrary, viscous fluid registration can provide dense displacement for each voxel but it is time consuming.
To combine the advantages of the two methods, an inhouse hybrid method is proposed. As shown in Figure 1, it applies FFD with sparsely distributed control points to provide a better initial estimation for viscous fluid flow. Then, the viscous fluid flow is applied to obtain the fine warp field estimation. It also uses the multi-scale strategy to implement warp field optimization. Our in-house hybrid method uses MI as the objective function. The fluid velocity is not calculated from successive over-relaxation but simplified to a convolution of force field. 19 The new displacement field is obtained from the current displacement and an increment caused by fluid velocity. Gaussian filtering is further applied to the displacement field to keep the displacement continuous. The pCT image is resampled according to the updated warp field and the optimization steps into the next iteration. The optimization process mentioned above will not stop until the optimization iteration exceeds the maximum optimization iteration or the cost function cannot be further minimized.

Evaluation and statistical tools
After DIR processing, the warping fields define the geometrical displacement between homologous points in the pCT and treatment CBCT images.The warping fields can be used to deform the pCT images to present identical anatomical structures to the CBCT images and to propagate the annotated points and OAR contours defined on pCT to the treatment CBCT.

Intensity-based metric
Similar to the metrics MI, the normalized crosscorrelation (NCC) and the absolute image value differences (which could directly show the registration errors in bony structure) 28,29 are indicators for the general registration performance in the whole overlap area. They are used to quantitatively evaluate the registration error in all anatomical structures.

Warp-based metric
The target registration error (TRE) of the manually annotated points in CBCT and pCT images is used to quantitatively evaluate the accuracy of the warp fields. The Dice similarity coefficient (DSC) between the delineated OAR volume in CBCT and that in warped pCT is a quantitative indicator for the agreement between the deformed pCT OAR volume and the delineated CBCT ones. The DSC is defined as 22 : where A is the binary mask of an OAR in CBCT and B denotes the warped binary mask of the same OAR in pCT. Each OAR is evaluated separately. The Hausdorff distance (HD) is defined as the maximum of the closest distance between two volumes. The closest distance is computed for each vertex of the two volumes.
The HD is very sensitive to outliers, since the most mismatched point is the sole determining criterion of the distance.
A paired t-test was used for statistical analyses of the two methods, and the p-value ≤ 0.05 was considered statistically significant. Figure 2 shows the image value difference between the CBCT images and the warped pCT images given by the two DIR methods. Indicated by the color bar, a larger image difference is shown in a brighter color. It is evident that more bright areas appear in the NiftyReg results, marked by the red circle in Figure 2j-i, whereas the inhouse hybrid method gives darker areas in corresponding green circles. In addition, NiftyReg presents large image discrepancy in the purple circle. In the shoulder area marked by orange circles, large image value differences appear in both methods. As shown in Table 1, the MI and NCC of the proposed DIR are slightly better than those of NiftyReg. The absolute image differences of the proposed method are smaller than the NiftyReg result (p < 0.001). Five feature points located in different cervical vertebrae are manually labeled in both the CBCT data and CT data ( Figure 3). Table 2 lists the mean TRE of the five different positions from P1 to P5. The overall mean TRE in the in-house hybrid method is lower than the NiftyReg result (2.34 ± 1.02 mm vs. 2.98 ± 1.21 mm). For P1 and P2, NiftyReg gives lower TRE than the in-house hybrid method, but the difference was not statistically significant. Figure 4 shows the visual differences between the delineated OARs on the CBCT images and the warped OARs from pCT. In the axial results (a), (d), and (g), the NiftyReg result (shown in the green area) is larger than the blue ground truth. The in-house hybrid method gives the red result closer to the ground truth. In Figure 4c, unrealistic shape deformation (green part) appears in the NiftyReg result on the bottom of the larynx. For the other OARs such as the oral cavity,mandible,spinal cord, and parotid, the two DIR methods give similar results. Qualitative statistics of DSC and HD in different OARs are presented in Figure 5. As shown in Figure 5, the in-house hybrid algorithm performs better. The average DSC of the studied OAR results given by the in-house hybrid method ranged from 0.77 to 0.91, higher than the NiftyReg results (ranging from 0.76 to 0.90), and the HD result of the in-house hybrid method ranged from 1.71 to 3.80 mm, lower than the NiftyReg results, which ranged from 2.00 to 5.19 (p < 0.05).

RESULTS
In addition, it is necessary to compare the execution time of the DIR algorithms. This helps to choose the appropriate DIR algorithm in different RT applications. The average execution time of NiftyReg to register one set of pCT data is about 5 s (4.60 ± 0.49 s), which is much less than the 97 s (96.98 ± 11.88 s) given by the in-house hybrid method.

DISCUSSION
In this study, an in-house hybrid method, combining FFD and a viscous fluid registration method, and NiftyReg were applied to transform pCT to treatment CBCT images. To quantitatively assess their registration accuracy, three intensity-based similarity metrics were employed to analyze the registration accuracy in the overall overlap. The TRE of annotated points, DSC, and HD of deformed OARs were used to indicate the warp error in different structures. Indicated by the analysis of the different similarity metrics, the in-house hybrid method gave warped pCT results that were more similar to the CBCT image. From the TREs of the five annotated points, NiftyReg gave more variable performance than the in-house hybrid method. According to the DSC and HD of different OARs, the in-house hybrid method performs better than NiftyReg.
As an intensity DIR,FFD uses a set of control points to generate the warp field for each voxel. It gives no consideration to the deformation difference between bony and soft tissue. Therefore, it has a limited ability to handle the deformations of complex structures, such as rigid bony structures surrounded by soft tissues with large deformations. On the contrary, viscous fluid registration provides a dense warp field and it has the potential to depict local warp fields in detail. Enlightened by the hybrid DIR F I G U R E 2 Axial, coronal, and sagittal images between registered computed tomography (CT) and cone-beam computed tomography (CBCT) image. (a-c) CBCT, (d-f) registered CT, (g-i) difference of in-house hybrid method, and (j-l) difference of NiftyReg methods, 24-26 the in-house hybrid method uses FFD registration with sparse control points to provide the initial warp field for the subsequent viscous fluid registration to accelerate the dense warp field optimization process. Due to the neck position changes in the different examinations, FFD is a good choice to correct the large offset in both the head bony structure and the neck bony structure at the same time. That is why FFD shows good results in head and neck applications. 1,14,30 The image value difference is a direct similarity indicator between the warped pCT and CBCT images. The in-house hybrid method gives a mean image difference around 60.70, lower than the 67.25 in the NiftyReg results. As shown in the intensity-based metrics, the in-house hybrid method gives a higher mean MI and NCC. Seen from both the statistical analysis and the visual comparison, the in-house hybrid method gives more stable results than NiftyReg. Because no unrealistic deformation of NiftyReg in the HN part was reported in previous studies, the current parameter choices in NiftyReg need to be optimized to reduce the deformation errors. Similar to the analysis based on landmarks, 31,32 the TREs of annotated points are used as an indicator for registration error. As shown in Table 2, the TREs of the NiftyReg method results in this study ranged from 2.04 to 4.41 mm, which are similar to the results in other studies. 14,30 The 4 mm TRE was reported for mono modality image registration, 14 and the TRE calculated from 14 landmarks was 2.44 mm. 30 In this study, the TREs of the in-house hybrid method (2.13 mm in P1, 2.17 mm in P2,2.27 mm in P3,2.16 mm in P4,2.96 mm in P5) are within the voxel diagonal of the pCT voxel diagonal 3.49 mm (1.27 mm × 1.27 mm × 3 mm). According to AAPM TG-132, the target registration error tolerance should be no more than the maximum volume dimension. [13] Therefore, it is possible to use the inhouse hybrid method in medical registration. Figure 5 shows the mean DSC and HD for different OARs in HN data. The in-house hybrid method presents a better mean DSC than NiftyReg in the larynx and spinal cord. A lower HD from the in-house hybrid method is found in most of the analyzed OARs. These DSC statistics are in agreement with results in previous studies. 23,30,33,34 Bastien showed that the mean DSC for the parotid gland of 15 patients was 0.75. 30

F I G U R E 4
Examples of contour propagation using two deformable image registration (DIR) methods. The blue region is the delineation by cone-beam computed tomography (CBCT), the red is warped organs at risk (OARs) generated by in-house hybrid method, and the green is generated by NiftyReg. (a-c) The brain stem, larynx, and spinal cord results; (d-f) the oral cavity; (g-i) mandible, parotid for right parotid, 0.79 for left parotid, 0.86 for larynx). 33 Nobnop applied three different DIR methods to register MV CT images and kV CT images in the HN area, 34 and the DSC in the right parotid, left parotid, and spinal cord ranged from 0.54 to 0.80. Li tested many different types of DIR methods on CBCT and CT data from HN patients and the mean DSC of the FFD method fell in the range of 0.7 to 0.9. 23 According to the DSC tolerance for clinical application (from 0.8 to 0.9) in AAPM TG-132, 13 both the in-house hybrid method and NiftyReg can be used for some OARs. The HD of the two methods in this study ranges from 1.76 to 5.16 mm, which is in accordance with the 5 mm reported in the previous study. 13

F I G U R E 5
The Dice similarity coefficient (DSC) and Hausdorff distance (HD) of organ at risk (OAR) contours given by in-house hybrid method and NiftyReg The prominent drawback of the in-house hybrid method is that it consumes too much time, 21 times longer than NiftyReg. This is mainly because the NiftyReg used in this study was a GPU version, and was tested with the extra support of NVIDIA GeForce GTX 1660. Meanwhile, the in-house hybrid method was developed on a CPU and lacks parallel computing, and the testing environment is Intel i5-9400F CPU (2.9 GHz) with 16GB RAM. The in-house hybrid method could be optimized into a CUDA version and the time consumption could potentially be reduced to 30 ms in our estimation. Since the viscous fluid registration process in the in-house hybrid method needs to optimize the local deformation in each pixel position, its computation takes longer than that of NiftyReg. It should be mentioned that the impact of DIR on dose calculations has not been analyzed in this study. This is the major drawback but it will be studied in further research.

CONCLUSION
In this study, an in-house hybrid DIR method and NiftyReg are evaluated to register pCT images to CBCT images in the head and neck region.The in-house hybrid method gives a better intensity-based similarity, and the registration error in most of the annotated points (P3,P4, P5), and DSC and HD in most of the tested OAR (larynx, spinal cord, oral cavity, mandible, and parotid) are better than NiftyReg. Furthermore, the in-house hybrid method presents more stable performance than NiftyReg, and it has the potential to be used in the field of adaptive radiotherapy.

AUTHOR CONTRIBUTIONS
Chunling Jiang: Data curation, methodology, project administration, and writing-original draft. Yuling Huang: Data curation, formal analysis, and investigation. Shenggou Ding: Investigation and project administration. Xiaochang Gong and Shaobin Wang: Data curation and investigation. Xingxing Yuan: Formal analysis. Jingao Li: Methodology and writing-review. Yun Zhang: Conceptualization,formal analysis,and writing-review and editing.