Development of Volumetric Independent Dose Calculation System for Verification of the Treatment Plan in Image-Guided Adaptive Brachytherapy

Purpose: This study aimed to develop a volumetric independent dose calculation (vIDC) system for verification of the treatment plan in image-guided adaptive brachytherapy (IGABT) and to evaluate the feasibility of the vIDC in clinical practice with simulated cases. Methods: The vIDC is based on the formalism of TG-43. Four simulated cases of cervical cancer were selected to retrospectively evaluate the dose distributions in IGABT. Some reference point doses, such as points A and B and rectal points, were calculated by vIDC using absolute coordinate. The 3D dose volume was also calculated to acquire dose-volume histograms (DVHs) with grid resolutions of 1.0 × 1.0 (G1.0), 2.5 × 2.5 (G2.5), and 0.5 × 0.5 mm2 (G0.5). Dosimetric parameters such as D90% and D2cc doses covering 90% of the high-risk critical target volume (HR-CTV) and 2 cc of the organs at risk (OARs) were obtained from DVHs. D90% also converted to equivalent dose in 2-Gy fractions (EQD2) to produce the same radiobiological effect as external beam radiotherapy. In addition, D90% was obtained in two types with or without the applicator volume to confirm the effect of the applicator itself. Validation of the vIDC was also performed using gamma evaluation by comparison with Monte Carlo simulation. Results: The average percentage difference of point doses was <2.28%. The DVHs for the HR-CTV and OARs showed no significant differences between the vIDC and the treatment planning system (TPS). Without considering the applicator volume, the D90% of the HR-CTV calculated by the vIDC decreases with a decreasing calculated dose-grid size (32.4, 5.65, and −2.20 cGy in G2.5, G1.0, and G0.5, respectively). The overall D90% is higher when considering the applicator volume. The converted D90% by EQD2 ranged from −1.29 to 1.00%. The D2cc of the OARs showed that the averaged dose deviation is <10 cGy regardless of the dose-grid size. Based on gamma analysis, the passing rate was 98.81% for 3%/3-mm criteria. Conclusion: The vIDC was developed as an independent dose verification system for verification of the treatment plan in IGABT. We confirmed that the vIDC is suitable for second-check dose validation of the TPS under various conditions.


INTRODUCTION
Image-guided adaptive brachytherapy (IGABT) based on magnetic resonance images (MRIs) has been introduced as a new standard technique to improve the treatment outcome in cervical cancer (1)(2)(3)(4)(5)(6). IGABT delivers a high dose with a small number of fractions after external beam radiotherapy (EBRT). The treatment plan of IGABT is optimized by using dosimetric parameters to meet dose constraints for the organs at risk (OARs) and high-risk clinical target volume (HR-CTV). The dose constraints are normalized to the equivalent dose in 2-Gy fractions (EQD2) to produce the same radiobiological effect as in EBRT. Furthermore, IGABT based on MRIs is more advantageous in the delineation of a region of interest (ROI) such as the HR-CTV and OARs than is high-dose-rate brachytherapy (BT HDR ) based on computed tomography (CT) images (7,8). Because compensation of IGABT for treatment outcomes is difficult, it is more important to generate the correct treatment plan and to deliver accurate doses than with other conventional treatments (9). Thus, at each fraction, IGABT requires re-optimization of the treatment plan and re-defining of the ROIs. Furthermore, the re-optimized treatment plan should be validated for safe and accurate delivery (10).
In general, verification of treatment plan in brachytherapy is performed by comparing the point doses calculated using the treatment planning system (TPS) and an independent dose calculation (IDC) system (11,12). In the previous studies, various IDC systems were reported (13)(14)(15)(16). Formalism of the American Association of Physicists in Medicine Task Group 43 (AAPM TG-43) is usually used to calculate the point doses via templates or were spreadsheet based in these systems. However, these pointdose comparisons alone were not sufficient for the verification of IGABT.
For better validation in IGABT, dosimetric parameters should be calculated using an IDC system. Some commercialized IDC systems, such as BrachyCheck (ROS, California, USA) and DIAMOND for Brachytherapy (PTW, Freiburg, Germany), were recently introduced. BrachyCheck can only calculate a dosevolume histogram (DVH). In addition, DIAMOND calculates point dose and 3D dose distributions, but in order to perform 3D analysis, the dose file must be exported to other analysis software. To solve inconveniences in these commercialized IDC systems, Xianliang et al. reported dose verification software (DVS) that can calculate the point dose and the dosimetric parameters (17). The DVS also utilizes the gamma evaluation, but there are some limitations. In the DVS, the calculation method for the direction of the source is inaccurate in the ring applicator because they defined the source vectors by connecting the existing source dwell-position and next source dwell-position. Moreover, it can use only fixed dose-grid size (1.0 × 1.0 mm 2 ). Therefore, we aimed to develop a volumetric IDC (vIDC) tool suitable for IGABT verification. In the vIDC, the method of determining the direction of the source will be improved, and various dosegrid sizes will be available for selection. We will also add some functions to the vIDC. Firstly, the dosimetric parameters of the HR-CTV will allow the calculation of two types (HR-CTV W or HR-CTV W/O ) with or without applicator volume, and secondly, the calculated 3D dose volumes will be able to be convert to EQD2 to account for the radiobiological effect caused by EBRT.

Volumetric Independent Dose Calculation System
In brachytherapy, it is important to accurately calculate the source position and orientation to improve the accuracy of dose calculations for tandem-and-ring applicators. The vIDC can calculate the appropriate direction of sources in the ring applicator, as shown in Figure 1. To calculate the source direction for the ring applicator, the vIDC obtains various plane vectors of the ring by using three sources' positions. Then, the vIDC selects an appropriate plane vector with analysis of the root mean square error (RMSE) between each source's position and the calculated plane vector. The ring equation of the sphere is obtained using the source positions and the selected plane vector (step 1). The vector direction of each source is calculated by using the tangent between the obtained ring equation and source position (step 2).
The dose calculation for the vIDC was based on the AAPM TG-43U1 formalism. The formalism is classified as a point and line source model, according to the method of source modeling. The line source model was used in the vIDC, and the equation of the model was as follows: where r is the distance from the center of an active source to a point of interest and θ is the polar angle between the source longitudinal axis and the line that connects the center of the active source and the point of interest. The r 0 and the θ 0 denote the reference distance (1 cm) and angle (90 • ), respectively. The air kerma strength (unit: cGy·cm 2 ·h −1 , 1 U = 1 cGy·cm 2 ·h −1 ), dose-rate constant (unit: cGy·h −1 ·U −1 ), geometry factor, radialdose function, and anisotropy function are represented as S K , , G(r, θ ), g(r), and F(r, θ ), respectively. Each value of S K and is defined in the TPS, while g(r) and F(r, θ ) are provided by the manufacturer of the source. The default grid size of the vIDC was 1.0 × 1.0 mm 2 grid (G 1.0 ). It is also possible to select various grid sizes, such as 2.5 × 2.5 (G 2.5 ) and 0.5 × 0.5 (G 0.5 ) mm 2 grids. The calculated dosimetric parameters are converted to EQD2 parameters for evaluating the radiobiological effect of the treatment plan in conjunction with EBRT. EQD2 is calculated as follows: A biologically effective dose (BED) is obtained for calculating EQD2. The BED is used for isoeffective dose calculations, which means that the true biological dose delivered is measured. Each n and d is presented as the number of fractions and dose per fraction, respectively. The α/β ratio means the ratio of "intrinsic radiosensitivity" to "repair capability" of a specified tissue, and it was presented in this paper as the tumor (α/β ratio = 10) and OARs (α/β ratio = 3).

Treatment Planning of Image-Guided Adaptive Brachytherapy Simulated Cases
In order to examine the efficiency of the functions and the clinical feasibility for the vIDC, we retrospectively evaluated with simulated cases treated using IGABT. Four simulated cases with squamous cell carcinoma of the cervix were treated with IGABT using an 192 Ir HDR source and adapted a tandem-andring technique. These cases were approved by the Institutional Review Board (IRB). Treatment plans were generated with the Oncentra Brachy (Elekta, Stockholm, Sweden) TPS, and the 3D dose volume was calculated with a grid size of 1.0 × 1.0 mm 2 . As in EBRT, the dosimetric parameters are used as criteria of the treatment plan ( Table 1). The criteria are defined as the doses covering 90% of the HR-CTV and 2-cc volumes of each OAR (D 90% and D 2cc , respectively). All plans were re-optimized with appropriate dose criteria according to the variation in anatomy and applicator location before treatment.
The point doses at points A and B and the rectum were also calculated to verify the accuracy of dose calculation using the vIDC by comparison with results calculated by the TPS. Point A is defined as 2.0 cm superior to the lateral vaginal fornix and 2.0 cm lateral to the cervical canal. Point B indicates a position 2.0 cm superior to the ring surface and 5.0 cm lateral to the midline. Point A represents dose limits that are delivered to the uterine cervix, whereas point B is used to evaluate the lateral spread of the effective doses, such as nearby doses in the pelvic wall and obturator node. The rectal point dose represents a delivered dose in the rectum. The reference point dose of the TPS was compared with that of the vIDC by using the percentage difference (%diff).

Monte Carlo Simulation for Validation of Volumetric Independent Dose Calculation System
Validation of vIDC was performed with Monte Carlo (MC) simulations using GATE (Geant4 Application for Tomographic Emission, Version 8.1). The grid size in MC was 1.0 × 1.0 mm 2 , and the 192 Ir source (HDR 192 Ir mHDR-v2) has been modeled, as described by Granero et al. (18). A virtual treatment plan that included a tandem-and-ring technique was established within the simulated water phantom. The plan was imported into the TPS and the vIDC, and the respective calculation engine was used to obtain the dose distribution in the axis plane of the ring applicator's center. With gamma evaluation, dose distributions calculated by the TPS and the vIDC were compared with the MC results.  Table 2 shows the percentage differences in the calculated point doses between TPS and vIDC at the reference points to evaluate the accuracy of the developed vIDC. All average percentage differences of vIDC were <-2.28% at reference points of the same absolute coordinates, compared with the clinical TPS.

Comparison of the Dose-Volume Histogram and Dosimetric Parameters
The vIDC facilitated comprehensive 3D dose evaluations for clinical use, even when the G 1.0 was used. Figure 3 shows averaged DVHs from the calculated 3D dose volumes by the vIDC and TPS using the same dose-grid size of 1.0 mm (G 1.0 ). As shown in Figure 3, no dose differences were observed for the HR-CTV W/O , bladder, and sigmoid in the DVHs. However, a dose difference with a high-dose range exceeding 1,000 cGy occurred in HR-CTV W . In addition, a dose difference of <10 cGy was observed at the same volume in the intermediate dose range for the rectum. Table 3 shows the dosimetric parameters between the vIDC and TPS for the HR-CTV with or without the applicator volume for various dose-grid sizes. The volumes of the HR-CTV  calculated by the vIDC were decreased by decreasing the dosegrid size, as shown in Table 3. The target volume calculated with G 2.5 was larger in the vIDC than TPS. The minimum absolute percentage difference for the target volume was 1.43% in HR-CTV W/O calculated with G 1.0 . The percentage difference for the dosimetric parameter also decreases by decreasing the dosegrid size. For the G 2.5 , the differences of D 100% and D 90% were more than −6%. The differences decreased to −3.42 and −0.44% when G 0.5 was used. The standard deviations for both evaluated dosimetric parameters were within 40 cGy for D 100% and 23 cGy for D 90% . Overall, the standard deviations were similar for the vIDC and TPS regardless of grid sizes. In addition, the D 100% values of vIDC were the same in HR-CTV W and HR-CTV W/O . However, all standard deviations of D 90% used the same grid size were higher in HR-CTV W/O than HR-CTV W. The differences of D 90% were also greater in HR-CTV W/O . Table 4 indicates the dosimetric parameters and volumes calculated with the TPS and vIDC for all OARs such as the bladder, rectum, and sigmoid under the condition of various dose-grid sizes. All OARs were acceptable for the D 2cc dose criterion shown in Table 1. For the calculated volumes for all OARs, the smallest difference (<1.82%) compared with that in vIDC and TPS was observed in the G 1.0 used in TPS. The volume differences for all OARs in the other grid were more than 4% for G 0.5 and 8% for G 2.5 . Compared with the dosimetric parameters of TPS using G 1.0 , all evaluated parameters of vIDC for the bladder and sigmoid have dose differences within ∼-3% for G 0.5 and G 1.0 , but there was a difference of more than −4% for G 2.5 . On the contrary, there was a dose difference for the rectum of <-3% for G 2.5 and more than −4% for G 0.5 and G 1.0 .

Comparison of the Equivalent Dose in 2 Gy per Fraction
In order to evaluate the radiobiological effect, the EQD2 of dosimetric parameters for the HR-CTV and OARs was analyzed in each faction. Table 5 shows the EQD2 for the HR-CTV and OARs in each fraction for one simulated case. For HR-CTV W/O , EQD2 differences in D 90% calculated using the vIDC and TPS ranged from −1.00 to 1.29%. For each fraction, the differences of EQD2 for dosimetric parameters were higher in OARs than HR-CTV W/O . The maximum difference of EQD2 was −9.15% in the sigmoid and was observed in D 0.1cc of the first fraction.

DISCUSSION
The source direction determined by Xianliang et al. (17) was obtained by connecting each source's position. This method is not suitable for obtaining accurate source direction, especially in using ring applicators. Therefore, we initially calculated the sphere's ring equation, and then we calculated the source direction using the tangent of the sphere equation at each source position. We showed the point doses in reference points and the doses in dosimetric parameters calculated by Xianliang's method in Tables 6, 7. For five reference points, the point doses using the vIDC were similar with those obtained with by Xianliang's method. In addition, the doses of dosimetric parameters calculated by Xianliang's method and the vIDC using G 1.0 show significant differences, as shown in Table 7.
Overall, the dose percentage differences of dosimetric parameters are higher in Xianliang's method than vIDC. In particular, the percentage difference in D 90% increased ∼3 times (1.83% for the vIDC and to 5.55% for Xianliang's method) compared with the dose calculated with TPS. The percentage differences also increased with the other dosimetric parameters for OARs. However, the percentage difference compared with TPS was found to be smaller in vIDC than Xianliang's method. Thus, we can conclude that the tangential method used in this study is more accurate than Xianliang's method. Maybe this method can be applied to conventional brachytherapy. Compared with the calculated point dose with TPS, the calculated point doses with vIDC at the right and left point B were the dose percentage difference of more than −1%. However, these differences cannot be defined as the inaccuracy of dose calculation by vIDC, because there was a small difference in the calculated doses between the vIDC and TPS. The average dose differences between vIDC and TPS for all reference points were <0.03 Gy. Therefore, it is indicated that the dose calculations obtained with the vIDC are accurate on the basis of the dose differences.
For the effect of grid size, the percentage difference in doses of all dosimetric parameters for OARs was the smallest in G 1.0 used in TPS. On the other hand, D 100% and D 90% for the target volume had the smallest percentage difference with G 0.5 regardless of whether applicator volume was used. The difference is probably influenced by the dose gradient of the contoured ROI. As the dose-grid size decreases, accurate dose verification is possible, especially in regions with steep dose gradients. However, the dose difference between G 1.0 and G 0.5 is relatively small (<10 cGy) and can be neglected. Moreover, when the dose-grid size becomes two times smaller, the time required for dose calculation increases four times. Therefore, we recommend determining the dose-grid size, taking into account the desired dose calculation accuracy and existing time constraints.
In this study, we noted the dosimetric difference due to the applicator volume. Potter et al. (19) reported that the applicator volume does not affect the dosimetric parameters if the target volume is sufficiently larger than the applicator volume. However, the authors also mentioned that this dosimetric effect of the applicator volume needs to be investigated because it has not been clearly established. Thus, we demonstrated the dosimetric effect with and without applicator volume using the vIDC, and there was almost no the dosimetric effect, as shown in Table 3. The dose differences were <4 cGy in all D 90 . On the basis of our results, we confirmed that the dosimetric effect does not need to be considered when the volume of the HR-CTV is sufficiently large, such as cervical cancer cases.
IGABT is known to improve treatment efficiency by delivering additional dose to local lesion control after EBRT. The IGABT is usually used at high doses (D 90% of HR-CTV > 550 cGy) in low fractions (<5 fractions) using the 192 Ir source with an average energy of 0.38 MeV. In contrast, the EBRT delivers normally the total dose of 45 Gy in 1.8 Gy per fraction delivered using a high-energy photon beam (20,21). Because of the characteristic difference of these delivery doses, it is difficult to accurately determine the same radiobiological effect between IGABT and EBRT for treatment dose. Therefore, the prescribed doses in each treatment were determined by conversion to EQD2 to enable accurate optimization of the delivered doses. In our study, radiobiological evaluation was performed by calculating EQD2, and this evaluation will enable safer and more accurate treatment. Gamma evaluation analyzed with MC simulation results was able to verify the vIDC. The gamma passing rate has no difference, although slightly better in vIDC than TPS. The point dose and dosimetric parameters of the vIDC indicated similarities with the TPS. In view of these results, we confirmed that the vIDC has almost the same performance as TPS. Thus, our vIDC can be used to validate the treatment plan as a second dose check that re-calculates the treatment plan and detects errors in the TPS. In this study, however, the vIDC was only evaluated in the tandem-and-ring cases. We will validate this vIDC in other treatment cases for clinical application in a future study.

CONCLUSION
The vIDC was developed for use as a second dose check to verify the clinical IGABT treatment plans. This study revealed that vIDC has the potential to verify the dose volumes calculated by TPS. In addition, we have been found that the doses in the dosimetric parameters are more affected by the dose grid size, not the applicator volume. Finally, the vIDC can improve the safety of IGABT by verifying the treatment plans with a variety of conditions.

DATA AVAILABILITY STATEMENT
All datasets generated for this study are included in the article/supplementary material.