Development of clinical application program for radiotherapy induced cancer risk calculation using Monte Carlo engine in volumetric-modulated arc therapy

Background The purpose of this study is to develop a clinical application program that automatically calculates the effect for secondary cancer risk (SCR) of individual patient. The program was designed based on accurate dose calculations using patient computed tomography (CT) data and Monte Carlo engine. Automated patient-specific evaluation program was configured to calculate SCR. Methods The application program is designed to re-calculate the beam sequence of treatment plan using the Monte Carlo engine and patient CT data, so it is possible to accurately calculate and evaluate scatter and leakage radiation, difficult to calculate in TPS. The Monte Carlo dose calculation system was performed through stoichiometric calibration using patient CT data. The automatic SCR evaluation program in application program created with a MATLAB was set to analyze the results to calculate SCR. The SCR for organ of patient was calculated based on Biological Effects of Ionizing Radiation (BEIR) VII models. The program is designed to sequentially calculate organ equivalent dose (OED), excess absolute risk (EAR), excess relative risk (ERR), and the lifetime attributable risk (LAR) in consideration of 3D dose distribution analysis. In order to confirm the usefulness of the developed clinical application program, the result values from clinical application program were compared with the manual calculation method used in the previous study. Results The OED values calculated in program were calculated to be at most approximately 13.3% higher than results in TPS. The SCR result calculated by the developed clinical application program showed a maximum difference of 1.24% compared to the result of the conventional manual calculation method. And it was confirmed that EAR, ERR and LAR values can be easily calculated by changing the biological parameters. Conclusions We have developed a patient-specific SCR evaluation program that can be used conveniently in the clinic. The program consists of a Monte Carlo dose calculation system for accurate calculation of scatter and leakage radiation and a patient-specific automatic SCR evaluation program using 3D dose distribution. The clinical application program that improved the disadvantages of the existing process can be used as an index for evaluating a patient treatment plan.


Background
Recently, owing to the development of radiotherapy technology such as the multi-leaf collimators (MLC), the use of intensity-modulated radiotherapy (IMRT) and volumetric modulated arc therapy (VMAT) is increasing [1][2][3]. Especially, VMAT is a treatment technique with the advantage of short treatment time because it provides complex intensity modulated beams by employing dynamic motion of a MLC with rotating gantry. However, VMAT is more likely to induce head-scattered radiation because it uses a higher monitor unit (MU) than conventional therapies such as conformal radiotherapy (CRT) [2,4,5]. Therefore, secondary cancer risk (SCR) in VMAT can be increased due to low doses to healthy tissues induced by scatter and leakage radiation from gantry head.
Currently, the SCR receives increasing interest and there are many studies related to the SCR [4,[6][7][8][9][10][11][12][13][14]. Researches regarding SCR have been reported that SCR increases as the volume of organs receiving radiation increases and as the distance from the source reduces [11,15]. In Hall et al. 's studies, it was reported that SCR is approximately 1.0-1.8% higher (10 years survival rate) in IMRT than CRT [13,16]. It is also reported that secondary scattered doses at locations approximately at 20-50 cm from the iso-center were approximately 1.0-5.0 mGy per 1 Gy for planning target volume (PTV) in IMRT [17].
The most commonly used methods to analyze the effects of secondary cancers involve the utilization of human phantoms to measure secondary doses using devices such as the thermo-luminescent dosimeter (TLD) and radio-photoluminescent glass dosimeter (RPLD), and the calculation of doses using commercialized treatment planning system (TPS) [4,8,[18][19][20]. Recently, studies have also been conducted to calculate scattered radiation for in-field and out-field, which has a significant impact on SCR using physics-based analysis algorithms [21][22][23][24]. Nevertheless, the existing evaluation methods have limitations. First, the method of evaluating by placing the measuring element inside the phantom cannot measure many points, and only the part with holes made for measurement can be measured. Because organ doses are calculated from values of predetermined points due to the structure of the phantom, there is a possibility that it may be evaluated differently depending on the person measuring it [25]. For example, in the case of the thyroid gland, some researchers measure 3 points to evaluate the dose to the thyroid gland, while others measure 5 points to define the thyroid dose. In addition, if the location of the pre-defined hole is not the center of the organ, the dose from the thyroid gland may vary depending on the position selected. Second, the method of measuring the amount of head-scattered radiation according to the distance from the iso-center has little variation in measurement value [26]. However, it is difficult to directly evaluate the patient organ dose because it is a method of measuring the effect of head-scattered radiation [10]. Third, using the TPS dose distribution, the estimated dose that has a large effect on the SCR is not accurate. The reason for the inaccuracy that occurs in TPS is mostly the effect of scatter and leakage radiation in out-of-field [27][28][29][30][31][32].
The Monte Carlo simulation is another method for SCR analysis, and many studies have been conducted [33][34][35][36]. The organ dose estimation utilizing Monte Carlo simulation has the advantage that the dose can be calculated accurately. In a previous study, we have developed a phantom-based dose calculation system for VMAT using GATE Monte Carlo that directly apply scatter dose, difficult to calculate using TPS [37]. And, the system developed in previous studies has been improved to allow direct calculation of patient dose by adding new modules. We verified that the Monte Carlobased dose calculation system for VMAT calculated the dose distributions for patient-specific treatment plan accurately.
Several studies regarding SCR from stray radiation have reported that it is especially important to predict SCR in order to provide patient's long-term health information in VMAT treatment with high MU [8,9]. However, there was no method of patient-specific SCR evaluation program through accurate dose calculation using Monte Carlo in VMAT treatment. In this study, we developed a clinical application program to apply a SCR evaluation process into the workflow of patient treatment. In addition, the automatic quality assurance (QA) system [37] for VMAT based on GATE Monte Carlo simulation which is useful for accurate dose calculation, has added and automated a module that enables calculation using patient specific CT data for direct evaluation in a system developed through previous research. By combining SCR evaluation program and program that improved the disadvantages of the existing process can be used as an index for evaluating a patient treatment plan.

Keywords: TPS, Monte Carlo, Secondary cancer risk, Clinical application program
Monte Carlo based automatic QA system for VMAT to improve the accuracy, we intend to apply clinical evaluation of SCR that may be caused by radiotherapy as an indicator of patient treatment.

Patient-specific dose calculation using Monte Carlo engine
We utilized a previously developed clinically available automatic VMAT QA system for this study [37]. The system is designed to re-calculate the beam sequence of equipment and treatment plan using the Monte Carlo engine, so it is possible to accurately calculate and evaluate scatter and leakage radiation, difficult to calculate in TPS [38]. The basic program was designed for QA, so the module was additionally configured to enable dose calculation using patient CT data. The CT contrast phantom (Electron Density phantom 062, CIRS, Norfolk, VA) was utilized to calibrate the electron density and CT number.
Typically, clinical programs perform dose calculation by means of electron density correlations based on the Hounsfield unit (HU) and Compton scattering. In this study, more accurate dose calculations were performed through stoichiometric calibration [39][40][41]. The mass density and the chemical composition of the tissue were obtained by the method developed by Schneider et al. as a necessary pre-calculation stage for the Monte Carlobased secondary cancer calculation system [41][42][43]. The mass density and the chemical composition are shown in Table 1. They are divided into 27 human tissues with CT numbers between − 1050 and 4000 and are interpolated between sections.
The patient CT data utilized in the Monte Carlo-based SCR evaluation program is converted to analyze 7.5 or metadata format. The patient's CT data was maintained as the 1.27 × 1.27 × 2.5 mm 3 voxel value of the original raw data for dose calculation. The commercialized TPS performs calculation after reconstructing the pixel values of patient CT data into a preset calculation grid size of 3 mm 3 . In the case of Monte Carlo calculation, after the physical calculation of the dose deposit was completed based on the voxel value of 1.27 × 1.27 × 2.5 mm 3 , it was reconstructed into a pixel size of 2 mm 3 to represent the  dose distribution [44,45]. Information regarding the number of fields, position of the gantry, field weight, and MLC of the patient treatment plan is converted to a macro file and stored in the system [37]. GATE v8.1, officially released in 2018, was employed for system configuration because it has the advantage of easy geometry configuration and fast calculation time with parallel computing using job split. Optimization was improved compared to the previous study, so the dose calculation time was approximately 10 h per patient [46]. The dose distribution after the dose calculation is stored in a 3D dose distribution DICOM format compatible with TPS.

Patient-specific calculation of radiation-induced cancer
The automatic SCR evaluation program created with a simple graphical user interface (GUI) using MATLAB (R2016a, MathWorks Inc., Natick, MA, USA) was set to re-calculate the treatment plan with the Monte Carlo engine and analyze the results to calculate SCR. The results calculated by the Monte Carlo engine are first converted into dose-volume histogram (DVH) in the application program. The DVH is a volume-based graph of dose distribution and is often used as the primary data for dose analysis in clinical and research applications [47][48][49]. The results of organ equivalent dose (OED) are constructed to be calculated based on DVH. The SCR calculation part of program is divided into two categories as shown in Fig. 1: OED and SCR. After selecting the organs and the parameters to be evaluated, the OED was calculated and the results were utilized to sequentially calculate the excess absolute risk (EAR), excess relative risk (ERR), and the lifetime attributable risk (LAR). The organs, gender, baseline cancer risk, age and biological factors for calculating EAR, ERR, and LAR were made up of a selection list to allow easy selection and calculation. The SCR was evaluated using the Biological Effect of Ionizing Radiation (BEIR) VII model [8,9,50]. The OED is a dose calculation method to evaluate the effects of SCR caused by radiotherapy. The OED corresponds to a dose response proportional to the population cancer incidence based on the same gender and age. In other words, it is a concept that has the same risk for SCR when having the same OED, and is configured to be able to calculate using 3D dose distribution. The OED has 3 models (Linear, Bell shape, Plateau model) [26]. In the case of SCR, we used a plateau model for specific OED because the linear type of dose-response as in the low dose region does not appear in the high dose region. This is because there is no effect of cancer risk in dead cells due to cell killing in the high dose region.
The OED was calculated as a plateau model by (1) based on tissue dose distribution [13,[51][52][53]. V is the total volume of the organ for calculating the OED, V i is the corresponding volume, D i is the volume dose, α and δ are dose-response factors of specific organs [54,55]. The EAR is the rate of disease in the exposed population minus the rate of disease in the unexposed population. It means the number of people with secondary cancer caused by exposure to radiation. The ERR is the same as EAR, but is a relative concept. The EAR and ERR of secondary cancer were calculated with (2), the LAR were evaluated based on two calculated EAR and ERR [8].
The biological parameters utilized in (2) were based on the data obtained from the Hiroshima bombing and the Chernobyl accident [56]. D, s, e, and a denote exposure dose, gender, exposure age; e* = (e−30)/10 for e < 30 and 0 for e > 30, and attained age, respectively [50]. The four parameters are the factors that have the highest effect on the results of the SCR. LAR was calculated from EAR and ERR utilizing (3).  The baseline cancer risk was calculated from the annual report of cancer statistics in Korea in 2016 [8,57]. And for SCR, LAR is a single estimate that can be expressed as one value and is calculated as a weighted average value for ERR and EAR. By the recommendation of the BEIR VII Committee, the weight of 0.7 for the ERR and 0.3 for the EAR were calculated on a logarithmic scale [9]. We set the life expectancy to 90 years and calculated the change in SCR according to the exposure age from 30 to 80 years. The automatic SCR evaluation program based on the GUI is designed to achieve the desired results by changing several parameters using 3D dose distribution. In order to evaluate clinical application of the program, we randomly selected five treatment plans using VMAT that were actually treated. For the evaluation of various treatment sites, five treatment plans were assigned for each disease site of brain, lung, abdomen, prostate and spine.

Results
The Fig. 2 shows a schematic diagram of the workflow that our developed clinical application program is used in actual clinical practice. The left side is the existing flow and the right side shows the flow when the developed clinical application program is applied. The five treatment plans selected to evaluate the program were calculated step by step according to the Monte Carlo dose calculation system developed in the study. The patientrelated data is converted in the program to apply to the re-calculation using the Monte Carlo engine, and the previous results were calculated as the value for SCR by applying each biological parameter. The most time-consuming procedure in the system was the dose re-calculation using the Monte Carlo engine, which took about 10 h. The OED calculation results are shown in Table 2. The OEDs calculated using the Monte Carlo engine were compared with OEDs calculated using the dose from the TPS. The OED results using the Monte Carlo based program developed in this study showed a difference of about 0-13.3% compared to the OED results of TPS. The out-of-field organ located away from the treatment region showed a larger difference in OED calculation. Table 3 shows the LAR of each organ according to the age difference exposure from the ages of 30-80 years old based on the attained age of 90 years old and varied between 0.017 and 6.294%. The minimum value of 0.017% appeared in the lens of brain disease site and the maximum value of 6.294% appeared in the colon of abdomen disease site. For example, in the case of LAR on the right lens, based on an exposure age 80 to attained age 90, this means that 17 people per 100,000 population is at risk for secondary cancer. Figure 3 shows the relationship between cumulative baseline cancer risk and LAR with increasing age at exposure. The baseline cancer risk is the rate of cancer incidence by age without the effects of radiation. And it was converted to a cumulative baseline cancer risk to compare with LAR. As with LAR, the cumulative baseline cancer risk was expressed as the sum from exposure age to attained age. By comparing the two values, we tried to easily confirm the effect of secondary cancer caused by radiation. The cumulative baseline cancer risk and LAR expected to be caused by radiation decreased with increasing exposure age. In addition, radiotherapyinduced SCR is relatively high at a young age. The results show that exposure to radiation at a young age is relatively dangerous. Through Fig. 3, we can confirm some things. In radiotherapy of prostate and spine, bladder and liver could be seen that the incidence of secondary cancer caused by radiation was very small than the cumulative baseline cancer risk. This means that compared to other organs, the bladder and liver do not have to consider the effect of secondary cancer. On the other hand, it was found that the SCR of esophagus was significantly higher than cumulative baseline cancer risk in lung radiation therapy. These results suggest that the effects for esophagus should be carefully considered when treating lung site. Of course, these results are only one of the indicators for analyzing the treatment plan, but it does not mean that other organs need not be considered. Depending on the patient, the effects of esophagus may not be significant in lung radiation therapy. So, we think that a SCR evaluation program that can be applied directly to each patient is needed. We should consider the risk of secondary cancer from radiation when evaluating these organs. The LAR result calculated by the conventional manual calculation method showed a maximum difference of 1.24% compared to the result of the developed clinical application program. We confirmed that the clinical application program proceeded without problems in the procedure for evaluating treatment plans in clinical practice.

Discussion
The radiotherapy-induced secondary cancer risk has been recognized and evaluated in different studies [4,[6][7][8][9][10][11][12][13][14]. The developed clinical application program is designed to evaluate SCR caused by radiation through analysis of the patient's treatment plan. However, existing SCR evaluation methods had some the disadvantage in measurement method. The disadvantages include a method of measuring the amount of stray radiation that varies depending on the distance from iso-center, or methods of measuring organ doses by inserting measuring elements such as TLD and RPLD into the human phantom. The method utilizing the point dose measurement device and human phantom has a disadvantage that it is difficult to calculate the accurate organ dose depending on the fixed measurement location [25]. Therefore, dose effect on the volume cannot be accurately evaluated. And the measurement results are standard data using phantoms, actual radiotherapy present high dose difference depending on treatment plan or patient, even with the same location of the disease.
In addition, the method using the TPS has a disadvantage of impossibility to calculate the effects of scatter and leakage radiation, which produces an important effect on SCR. For example, in the case of beam data measurement conducted by commissioning of TPS, 6 MV energy measured a range up to 20.5 cm based on 1.5 cm depth (D max ) and 10 × 10 cm 2 field size. Therefore, it is estimated that no dose is present at the outside of 10.25 cm off the central axis in beam data. However, the actual scatter and leakage radiation exist, and the results of the Monte Carlo simulation indicate that a dose of approximately  . 3 Comparison of LAR and cumulative baseline cancer risk according to exposure age for OAR (treatment site). LAR lifetime attributable risk, OAR organs at risk 1% is calculated even at a location 13.0 cm from the central axis. And at the same position, TPS showed a value of 0.4%. When actually measured, a value over 1% was measured at the same location. For these reasons, in the case of radiotherapy in the same disease site, the OED value of the organ may vary depending on the patient or the treatment plan. In addition, OED showed a greater difference between Monte Carlo and TPS as the distance from the treatment field increased. As the distance from the treatment region increases, the difference in OED can be considered for the following reasons. First, it is the limitation of the application for TPS commissioning as described above. And second is the use of OED based on the plateau model. The plateau model is a biological model that explains that the dose and risk are not proportional to the dose above the threshold, as in low dose. Increasing the distance from the treatment region means a difference in the low dose, and thus has a greater effect.
We believe that a patient-specific SCR evaluation program that can be applied directly to a patient is required because the effect of the dose is different for each patient to predict results with one or two measurements. The difficulty in system configuration was plan re-calculation utilizing the Monte Carlo engine, and required approximately 10 h to evaluate the patient to optimize the effects from scatter and leakage radiation. However, the most time-consuming Monte Carlo dose calculation is constantly being improved. Currently, general clinical application is difficult in workflow due to time limit, so it can function as a prototype program that is applied to special treatment such as stereotactic radio surgery (SRS), stereotactic body radiation therapy (SBRT) or pediatric treatment. We are now using the 88 node high performance computing (HPC) system we designed, but we are planning to upgrade the HPC system so that it can be calculated within 3 h soon. After that, the program is expected to be clinically applicable to the VMAT plan.
In SCR evaluation, it is important to perform a biological evaluation of the effect of the overall OAR volume. There have been many studies on neutrons and scatter radiation from gantry head to evaluate the effects of secondary cancer [58,59]. The calculation method of secondary cancer by scatter radiation from gantry head without evaluating the volume of a specific organ may have high uncertainty. The human phantom and Monte Carlo calculation are required to evaluate the volume of a specific organ. So, we developed a clinical application program that can evaluate organ doses for individual patients. We have constructed and evaluated a clinical patient-specific automatic SCR evaluation program that can accurately calculate scatter and leakage radiation. The dose calculation of the treatment plan for the program test showed that the scattered dose difference between the Monte Carlo calculation system and the TPS was larger as distance increased from the treatment region. To verify the accuracy of the calculated SCR evaluation, no significant difference within 1.24% compared with the manual calculation method was observed when compared with the results using the excel sheet from the previous study [4,9]. The disadvantage of the existing method is that the procedure is complex and must be checked manually to verify the accuracy of the calculation results. In addition, it is necessary to observe the effect of changes in the region to be evaluated and various factors, but there are many parameters to be checked in the sheet because the related equation is complex. And it takes a lot of time because it has to go through a process such as going through another person's verification to confirm the calculated data. For this reason, the existing system had many difficulties in workflow to be applied to each patient in clinical trials.
The automatic patient-specific SCR evaluation program has the advantage of enabling accurate evaluation of scatter and leakage radiation from gantry head for each patient treatment plan and convenient SCR evaluation for various parameter changes.

Conclusions
The effects of radiotherapy-induced secondary cancer risk have been evaluated in several studies [4,[6][7][8][9][10][11][12][13][14], and the development of evaluation tools for treatment plan to help with various clinical decisions has been carried out. To improve the disadvantage of the previous study, the system combines Monte Carlo dose calculation system and SCR evaluation program. The automatic SCR evaluation program can be useful for OED calculation and SCR evaluation for each patient.
By utilizing the developed clinical application program, we plan to conduct a patient-specific SCR evaluation and further study how to establish the means for a SCR database for patients in South Korea. The data will be used as a new evaluation index for SCR in South Korea.