Script-Based Automatic Intensity-Modulated Radiotherapy Planning with Robust Optimization for Craniospinal Irradiation: A Retrospective Study


 Background

Craniospinal irradiation (CSI) is essential for treating central nervous system malignancies. However, it is very challenging and time-consuming to design a robust treatment plan capable of delivering a uniform dose to the entire treatment volume for standard linear accelerators. This study proposes a script-based automatic planning method with robust optimization for CSI to reduce sensitivity to setup errors and increase planning efficiency.
Methods

Ten CSI patients with planning target volume (PTV) lengths between 49.8 and 85.0 cm were retrospectively studied. IMRT plans with robust optimization were generated by the automatic planning script. The plans’ sensitivities to positional inaccuracy were evaluated by deliberately shifting the middle beamsets ±5 mm longitudinally. Non-robust optimization plans were also created and compared to the robust ones. The planning time was recorded to evaluate the efficiency of the automatic CSI planning.
Results

For the ten patients enrolled in the study, homogeneous dose coverage were achieved for both robust and non-robust plans with average conformity index value of 0.827 (±0.028) and 0.841 (±0.023) (p = 0.047), homogeneity index values of 0.076 (±0.014) and 0.078 (±0.027) (p = 0.114), respectively. In robust plans, the variations of D1% and D99% were between -0.17–0.99% and -1.18–0.53% for cranial PTV, and -1.02–7.2% and -4.78–0.31% for spinal PTV, while in non-robust plans, the variations were much more significant with corresponding values of -0.53–8.31% and -5.17–0.57% for cranial PTV, and 0.53–40.14% and -46.06% to -2.01% for spinal PTV. The total planning time of a robust optimization plan was about 9-13 minutes by the automatic planning script.
Conclusions

The script-based automatic CSI planning with robust optimization could efficiently create a uniform dose coverage for the entire target volume regardless of PTV length, with the advantage of being insensitive to setup errors and standardizing the planning process.


Background
Craniospinal irradiation (CSI) is an essential component in the curative treatment of central nervous system malignancies with the cerebrospinal uid transmission, including most medulloblastomas, some ependymomas, some germinomas, and central nervous system lymphomas [1][2][3]. Delivering a uniform dose to the entire planning target volume (PTV) using standard linear accelerators (Linac) is very challenging because of long and irregular target volumes and beam divergence. Traditionally, CSI treatment plans using a three-dimensional conformal radiation therapy (3D-CRT) technique include two opposing brain elds and one or two posterior spinal elds. In the 3D-CRT technique, it is di cult to deliver a homogeneous dose, particularly at the eld junction areas. Hence, isocenters are usually shifted by 1 or 2 cm weekly to feather excessive hot or cold spots. High-precision conformal radiotherapy, such as intensity-modulated radiation therapy (IMRT), volumetric modulated arc therapy (VMAT), helical tomotherapy, and proton therapy, have been explored in CSI to improve dose homogeneity and reduce doses to organs at risk (OARs) [4][5][6][7][8][9][10][11][12][13]. Comparing to 3D-CRT, IMRT and VMAT have another advantage of being less susceptible to patient setup errors due to intentionally overlapped adjacent elds [5,13]. Many approaches, such as jagged-junction IMRT [14], gradient-optimized [15], auto-feathering [16], and robustoptimized [17] have also been applied in optimization to further improve the robustness to setup errors.
For the jagged-junction method, adjacent elds were intentionally overlapped, and eld edges were staggered in a certain step size. The robustness against setup errors depended on the step size. The gradient-optimized method has been recognized as very robust against setup errors by producing a slow, linear, and complementary dose gradient at the eld junctions. However, this method is time-consuming due to additional contour segmentation and complex plan optimization and is inclined to generate a stepshape dose distribution. The auto-feathering algorithm integrated into Eclipse version 15.6 treatment planning systems (TPS) (Varian Medical System, USA) showed comparable robustness to the gradientoptimized method. However, it is not available in other TPS. Recently, robust optimization with multiple independent isocenters has been available in RayStaion treatment planning systems (TPS), version 9A, (RaySearch laboratory, Sweden). The method implemented in RayStation is minimax optimization, which tries to nd the best solution in the worst-case scenario [18]. CSI can particularly bene t from such robust optimization due to eld matching uncertainty can be handled during optimization. In IMRT and VMAT planning, another challenge for CSI is to determine the number and location of isocenters for the long PTV and the overlap length between two adjacent elds. In this study, we will describe a script-based automatic planning method with robust optimization for CSI, which can effectively generate a homogenous dose coverage for PTV and is insensitive to eld matching errors so as to simplify and standardize the CSI planning process.

Patient Selection and Structure De nition
Ten patients with PTV lengths between 49.8 and 85.0 cm, with media length of 66.5 cm, who received CSI were retrospectively included in this study. All the patients were positioned supine, with arms resting at both sides and immobilized by thermoplastic masks in our regular procedure. Computed tomography (CT) images were acquired using a 16-slice GE Revolution CT scanner (General Electric Medical Systems, USA) of 3-mm thickness and transferred to the RayStation TPS, version 9A. All plans were designed for an Elekta Versa HD TM Linac (Elekta Medical System, UK) equipped with a high-resolution multileaf collimator with a leaf width of 5 mm at the isocenter.
The cranial clinical target volume (CTV) included the whole brain within the inner table of the skull, while the spinal CTV included the entire subarachnoid space that encompasses the extensions along with the nerve roots laterally. The cranial PTV (PTV_C) was created by expanding the cranial CTV by 5 mm, while the spinal PTV (PTV_S) was created by expanding the spinal CTV by 7 mm. The OARs to be delineated included the eyes, lenses, cochlea, parotid, salivary glands, larynx, esophagus, thyroid gland, breasts in females, lungs, heart, kidneys, intestines, and gonads. The actual prescription doses were 23.4 -36 Gy to PTV with 13 -20 fractions for the ten patients. However, the same dose of 23.4 Gy was prescribed to each patient with 13 fractions, and every plan was normalized to 95% of PTV in this retrospective study.
Automatic CSI Planning Work ow An automatic CSI planning script based on Python language was created in RayStation TPS. The work ow was described below and showed in Figure 1.

Isocenter setting and beamsets arrangement
The number of isocenters depends on the length of PTV_S, which is denoted as "L." If L was < 38 cm, two isocenters were needed, wherein the rst one was placed in PTV_C and the second in the center of PTV_S. If L was between 38 cm and 42 cm, two isocenters were still needed; however, the second one was moved 10 cm toward the patient's back. If L was between 42 and 72 cm, three isocenters were needed, wherein the rst one was still placed in PTV_C, and the second and third isocenters were placed in the midline of PTV_S. If L was between 72 and 80 cm, three isocenters were needed, but the second and third were moved 10 cm toward the patient's back. Since L is generally smaller than 80 cm for a majority of patients, cases in which L is larger than 80 cm were not considered. A case with L between 42 and 72 cm was shown in Figure 2 to illustrate the determination of isocenters and overlapping of adjacent elds. The coordinates of the center of PTV_C and PTV_S were extracted and denoted as (x c , y c , z c ) and (x s , y s , z s ), respectively, with the positive x, y, and z coordinates toward the patient's left arm, abdomen, and head, respectively. All the isocenters have the same x and y coordinates; therefore, only a longitudinal shift of couch was performed during the treatment. The coordinates of the three isocenters can be expressed as follows: Iso1 (x, y, z) = (x s , y s , z c ), Iso2 (x, y, z) = (x s , y s , z s + (L-F)/2), Iso3 (x, y, z) = (x s , y s , z s -(L-F)/2), where "L" is the length of PTV_S and can be derived by the GetBoundingBox () built-in function in RayStation. "F" is the achievable largest eld length assuming that the maximum eld length of standard Linac is 40 cm. However, considering beam divergence and penumbra, F was set to 38 cm.
The cranial beamset with the rst isocenter (Iso1) consisted of 6 elds with gantry angles of 180, 260, 300, 60, and 100 degrees and Y2 jaw limit not exceeding the fth cervical vertebra level to avoid shoulders irradiation. Both the upper and lower spinal beamsets with the second isocenter (Iso2) and the third isocenter (Iso3), respectively, consisted of three posterior or near posterior elds with gantry angles of 140, 180, and 220 degrees. The cranial and upper spinal beamsets overlapped at the neck from the rst to the fth cervical vertebra. The upper and lower spinal beamsets overlap at the middle spinal with a length equal to 2×F-L. If the overlapping length was set to be at least 4 cm and F was known to be 38 cm, the maximum value of L was 72 cm. However, when the second and third isocenters were moved 10 cm toward the patient's back, F was magni ed by 1.1 times, and L could reach almost 80 cm.

Plan Optimization
The initial optimization functions are listed in Table 1. PTV dose constraints were set to robustness with independent isocenters position uncertainty of 5 mm in superior and inferior directions. The total scenarios to compute were 27. Two auxiliary structures, Ring_PTV_C and Ring_PTV_S, were involved in the optimization to improve the target conformity. The Ring_PTV_C was derived by rstly expanding PTV_C with large margins of 2 cm in anterior, posterior, left and right directions, and small margins of 0.5 cm in all six directions, then subtracting the smaller expanded PTV_C from the larger expanded PTV_C. The Ring_PTV_S was derived in the same way but with large margins of 4 cm in anterior, posterior, left and right directions, and small margins of 1 cm in anterior, left, right, superior, and inferior, and 2 cm in posterior directions. The larger posterior margin in Ring_PTV_S is due to posterior or near posterior elds for PTV_S. Function values of each optimization objective were computed after one round of optimization. If the function value is out of range of 0.002-0.006 and the organ type is not target, the script will tighten or loosen dose constraints until all the values are within the range. Then optimization and adjustment were repeated until the function values were within the range, or the iteration number reached 10. At the end of the optimization, a very highly weighted max dose constraint was added to the external body, then the last round of optimization was performed. The non-robust plans had the same initial optimization functions except that PTV dose constraints were set to non-robustness.

Setup Error Sensitivity
The sensitivity of the automatic CSI plans to setup error was evaluated by deliberately shifting Iso2 5 mm superiorly and then inferiorly. Doses were recalculated with applied shifts but without altering any planning parameters. A dose comparison was performed for the shifted and un-shifted plans.

Plan evaluation
The dose homogeneity and conformity referred to the PTV were determined using the homogeneity index (HI) [19] and conformity index (CI) [20] through the following equations: HI = (D2-D98) /D50, CI = (V T,P × V T,P ) / (V T × V P ), where D50 is the median dose to the PTV, and D2 and D98 are the near maximum and minimum doses that cover 2% and 98% volume of the PTV, respectively. V T,P is the volume of the target enclosed by the prescription dose. V P is the volume of tissues including the target covered by the prescription dose, and V T is the volume of the target. D1, D2, D98, and D99 of PTV_C and PTV_S, D1 of Lens, and mean dose of lungs, heart, and kidneys were also evaluated, where Di is the dose that covers i% volume of PTV. The total planning time was recorded to evaluate the e ciency of the automatic CSI planning.

Statistical methods
Wilcoxon matched-pairs signed-rank test was performed with p < 0.05 to identify differences between robust and non-robust optimization plans. The computations were performed with the use of the SPSS statistics program (Version 22, IBM SPSS Inc., Chicago, IL).

Results
Dose distribution, Statistics and Planning Time Table 2 summarizes the dose statistics to PTV and OARs for both robust and non-robust plans.
Compared with the non-robust plans, the robust plans had slightly inferior CI value of 0.827±0.028 vs. 0.841±0.023, p = 0.047, and slightly higher lung mean dose of 15.9%±2.0% vs 15.5%±2.2%, p = 0.013. The D1, D2, and D99 of PTV were also slightly inferior to those of non-robust plans (p < 0.05). One representative patient's dose distribution and DVH of both robust and non-robust plans are shown in Figure 3. A longitudinal dose pro le of one patient's cranial beamsets and spinal beamsets of the robust plan was also depicted in Figure 4.
Isocenter setting, auxiliary ring structure derivation, initial planning setup, and optimized function settings were quickly performed by the automatic planning script in about half a minute. During the optimization, the maximum optimization round number was set to ten, with each round taking approximately 1-2 minutes. However, for most of the patients, the round number usually was 3 to 5. The total planning times for each patient were approximately 9-13 minutes. The simulation of the ±5-mm longitudinal error resulted in overdose or underdose mainly in the eld junction areas. The dose differences of D1, D2, D98, and D99 of PTV_C and PTV_S between shifted and non-shifted plans were illustrated by boxplots in Figure 5. In robust plans, the variations of D1 and D99 were between -0.17% to 0.99% and -1.18% to 0.53% for PTV_C, and -1.02% to 7.2% and -4.78% to 0.31% for PTV_S, while in non-robust plans, the variations were much more signi cant with corresponding values of -0.53% to 8.31% and -5.17% to 0.57% for PTV_C, and 0.53% to 40.14% and -46.06% to -2.01% for PTV_S. The dose variations of two lines, which were placed in the neck and spinal elds junction in sagittal pro le for a robust plan were depicted in Figure 6.

Discussion
We proposed a script-based automatic IMRT planning method for CSI, and this method could e ciently create clinically acceptable plans, regardless of PTV length. With robust optimization, our study showed great results in the plan quality. For the ten enrolled patients, homogeneous doses to PTV were achieved with CI values in a range of 0.78-0.86 and HI values in a range of 0.06-0.09. These results were similar to those reported by Cao et al. [14] in which a jagged-junction IMRT method was applied. They reported good homogeneity and conformity with HI values between 0.07 and 0.09 and CI values between 0.79 and 0.8. In our study, all the robust plans were very robust against elds matching errors that the maximum variations of D1 was not exceeding 7.20% and D99 not exceeding -4.78% for PTV_S. However, the maximum variations of D1 and D99 in non-robust plans can reach 40.14% and -46.06% for PTV_S. A study by Maddalo et al. [16] compared several VMAT strategies, including automatic feathering integrated into TPS, staggered overlap, and gradient optimization, and evaluated their robustness. They found that the automatic feathering algorithm and the gradient-optimized method yielded comparable plan robustness with dose differences of D2 between 5.7%-6.0% and D98 between 6.1%-7.6%. Our robust plans yielded better results with dose differences of D2 and D98 within 1% for PTV_C, D2 between -1.1% to 6.0%, and D98 between -3.63% to 0.63% for PTV_S.
Some studies have investigated the effect of the length of the overlap area on setup error. Theoretically, a smooth, slowly changing gradient is achieved with an increasing length of overlap, which leads to smaller dose variation due to setup errors. This was also con rmed in our study, as shown in Figure 6, the long overlap length in spinal junction induced more slightly dose variation than that of the short overlap length in neck junction. Myers et al. [15] recommended the junction overlap should be longer than at least 3.5 cm to account for major setup deviations. Seppälä et al. [5] suggested that the overlap should be at least 4 cm in their dynamic split eld IMRT plans. In our study, the overlap length depended on the patient's anatomy geometry. At the neck junction, the overlap length ran across the rst to the fth cervical vertebra, with values between 4 to 9 cm for most patients. At the spinal junction, the overlap length was equal to 2×F-L with a value of at least 4 cm. The patient-speci c overlap length in our study has the advantages of maximizing the overlapping length and minimizing the number of isocenters.
It is noted that the dose differences of D1, D2, D98, and D99 of PTV_S between shifted and non-shifted plans were much more signi cant than those values of PTV_C, as shown in Figure 5. This is attributed to the fact that the ±5-mm longitudinal error resulted in overdose or underdose mainly in the eld junction areas, which are both located in PTV_S with the beamsets arrangement in this study.
In this study, the spinal target was treated by three elds with gantry angles of 215, 180, and 145 degrees; thus, the low dose volume irradiated was not increased much when compared to 3D-CRT. Another advantage of these posterior or near posterior beams was that when the length of PTV_S was slightly longer than 38 cm or 72 cm, SSD was extended by moving the isocenter toward the patient's back to enlarge the achievable largest eld length instead of adding another isocenter, which would greatly increase the treatment time. Although these posterior or near posterior beams would probably sacri ce the homogeneity that may generate hot spots, it usually occurs in the patient's back muscles where there are no critical organs.
During the robust optimization, only ±5 mm longitudinal setup errors were considered; this is because the lateral and vertical errors were accounted for by the PTV margin expanding from CTV. Additionally, the amplitude of setup errors in cranial-spinal direction for CSI usually is < 5 mm after careful setup [21,22].
Larger errors can be accounted for by setting larger setup uncertainty during robust optimization.
However, this would enlarge irradiation volume, which would deteriorate the conformity and homogeneity and increase OARs' doses. As can be seen in Table 2, the robust plans had inferior CI values and higher lung doses compared to the non-robust plans, despite the difference were very slight.

Conclusion
The script-based automatic CSI planning with robust optimization could e ciently obtain a homogeneous dose coverage to the entire target volume with a standardized procedure. It greatly simpli es the planning process while maintaining the advantage of being insensitive to setup errors. This study was approved by the West China Hospital ethics committee. We con rm that all methods were carried out in accordance with relevant guidelines and regulations.

Consent for publication
Patients' CT images were anonymized for dosimetry studies, informed consent was not required.

Availability of data and materials
The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.
Xuetao Wang: Methodology, Data Curation; Writing; Changhu Li: Methodology, Data Curation, Writing; Jianling Zhao: Data Curation and preparation of tables; Guangjun Li: Data Curation and preparation of gures. Sen Bai: Experiments guiding and manuscript revising. All authors read and approved the nal manuscript.

Figure 1
The work ow of automatic intensity-modulated planning of craniospinal irradiation   Abbreviations: PTV = planning target volume; Di = the dose that cover i% volume of PTV