Compared to Intensity Modulated Radiation Therapy (IMRT), Stereotactic Body Radiation Therapy (SBRT) is characterized by potent ablative doses and highly conformal dose distributions delivered in a few fractions with a short overall treatment time  . The low number of fractions and the high dose gradient require high-precision treatment delivery techniques to restrict treatment volumes and more accurate PTV margins to avoid radiation damage to normal tissues    .
In order to reduce toxicity of normal tissues, the planning target volume (PTV) should be minimized by keeping margins for setup and inter-fractional position errors as low as possible. Hence, calculation of the setup errors and PTV margins plays an important role in Stereotactic Body Radiation Therapy (SBRT), especially in lung tumor treatment.
Van Herk et al. and Remeijer et al. had introduced various methods to calculate systematic and random errors in radiation therapy   but the methods should be determined by the categories of radiation therapy, such as IMRT, VMAT, SBRT and SRS. Gordon  et al. suggested effective systematic and random errors to correct the residual errors in isocenter and correlations between average isocenter position and individual random errors for radiotherapy with finite fractions. Yukinori  et al. recently proposed variance component analysis, reducing the effect of random component in systematic errors, especially in hypofractionated radiotherapy, to correct the systematic errors introduced by Remeijer.
Variance component analysis has been in use since the 1910s. Searle  et al. brieﬂy summarized the history of this method and described several methods for estimation of variance components in their books. The purpose of this work is to compare standard deviation (SD) of systematic errors and random errors using variance component analysis raised by Yukinori with the methods proposed by Van Herk and Remeijer and to testify the effectiveness of variance component analysis in PTV margins and dose calculation.
In this study, we focused on the verification of variance component analysis raised by Yukinori in Lung Stereotactic Body Radiation Therapy (SBRT) with immobilization devices  . Although previous studies have investigated interfractional uncertainties calculation with analysis of components of variance during lung stereotactic ablative radiotherapy (SABR)  , our study is the first that has rigorously evaluated variance component analysis in Lung SBRT and utilized PTV margins acquired in dose calculation to verify the DVH (dose-volume histogram) variation with three methods (Yukinori method, Van Herk method, Remeijer method).
2. Materials and Methods
2.1. Patient Selection
Twenty patients (12 men, 8 women) receiving lung SBRT between February 2013 and August 2016 at Zhongnan hospital were retrospectively analyzed in this study. Patients were excluded if previous chemo- or thoracic radiotherapy had been administered. In these cases, sixteen targets were primary early-stage lung tumors, and four were metastases of other solid tumors. All patients were treated with prescription dose (50 Gy) in five fractions every other day. In all, the scale of the treatment involved 100 fractions. Patients’ ages ranged from 47 to 70 years with average age of 52.5 years.
2.2. CT Simulation and Treatment Planning
Patients were fixed using R624-SCF immobilization devices  equipped with vacuum bag. As the therapist helps patients to adjust body to the best position, the vacuum bag was deformed to desired shape for patients. The breathing control panel was adjusted to restrict the motion of tumor accordingly. Planning CT images were acquired using a Siemens CT scanner (Somaton Sensation Open, Siemens AG, Erlangen, Germany) with 3 mm CT slice. For each patient, SBRT plans were generated using Eclipse Treatment Planning System (Varian Medical Systems, Palo Alto, CA) with the progressive resolution optimizer 3 (PRO3, ver.11.0.31, Varian Medical Systems, Palo Alto, CA). All patients were planned using nine fields with gantry angles of 10˚, 30˚, 40˚, 75˚, 110˚, 145˚, 310˚, 330˚, and 340˚, using 6 MV SBRT photon beam. For organs at risk (OARs), the following constraints were used: 1) Lung: V20 < 10%; 2) Spinal Cord: Max Dose < 18 Gy; 3) Trachea: Max Dose < 30 Gy. For consistency between plans on the same patient with different margins, plans were all produced by the same medical physicist. All plans were approved by the treating physicians.
2.3. CBCT Acquisition and Registration
All patients were treated on the iX linear accelerator (Varian Medical Systems, Palo Alto, CA). Prior to daily treatment, kV-CBCT was used to verify the position of target area. Analyzed shifts were applied for setup correction and CBCT was repeated. If a large shift (>3 mm) occurred in any direction, the treatment would not start until an additional CBCT was employed for verification after localization. The criterion for image registration between the CT simulation image and CBCT image was soft tissue registration. Due to the size of lung tumor and organ motion in the lung, the rigid registration is not ultilized here with the assistance of immobilization equipment. The setup errors of each patient in three directions were obtained from the variation between CT and CBCT. Figure 1 displays the soft tissue registration in two directions. Alignment was considered optimal when the posterior border of the lung tumor and the anterior lung tumor were in coincidence.
2.4. Analysis of the Setup Variations for the Systematic and Random Errors
1) Yukinori method
Figure 1. The illustrations of soft tissue registrations (Tumor) between CT in simulation and CBCT in treatment. Red curves represent GTV and blue curves represent PTV margins. (a) Cross section of lung; (b) Coronal plane of lung (the red area in the right is the fake GTV which can be ignored).
Yukinori introduced variance component analysis with the calculation methods of standard deviation (SD) of systematic and random errors. Table 1 listed formulas and specific calculation details.
2) Van Herk method
Van Herk proposed the method of calculating standard deviation of systematic and random errors. The standard deviation (SD) of systematic errors were represented by standard deviation (SD) of mean interfractional errors per patient. The standard deviation (SD) of random errors was represented by standard deviation (SD) of root mean square of interfractional errors per patient. Table 2 listed the details of this calculation process. If we consider the rectification of degrees of freedom, Van Herk method is the same to Remeijer method.
3) Remeijer method
The method of calculating systematic and random errors introduced by Remeijer was:
Note: a: number of patients; n: number of measurements per patient; xij: a measured value in fraction j for patient i; M: overall mean of the dataset; mi: individual mean in patient i; σ: standard deviation of random error; ∑: standard deviation of systematic error
2.5. Calculation of PTV Margin from CTV
To discover more general principle in PTV differences, two margin recipes of Stroom and Van Herk were used to calculate the PTV margin form the CTV. ∑ is the standard deviation of total systematic errors and σ is the standard deviation of total random errors. In this work, we consider the interfractional error is the main source of the systematic and random errors.
Table 1. Anova table of Yukinori method.
Note: DF = degree of freedom; SS = sum of squares; MS = mean square; EMS = expected mean square; a = number of patients; n = number of measurements per patient; xij = a measured value in fraction j for patient i; M = overall mean of the dataset; mi = individual mean in patient i; σ2 = random error variance; ∑2 = systematic error variance. Systematic error (∑): Within 95% confidence interval, ∑2 = (MSA − MSE)/n, Random error (σ): Within 95% confidence interval, σ2 = MSE.
Table 2. Calculating table of systematic errors and random errors with Van Herk method.
Note: Quote from VanHerk  .
So Stroom  et al. assumed a 95% dose to 99% of the CTV on average based on tests of realistic plans:
Stroom formula = 2 ∑ + 0.7 σ (3)
Van Herk  et al. assumed a Monte Carlo-based test of 1% TCP loss due to geometric errors for the prostate. The formula is defined as follow:
Van Herk formula = 2.5 ∑ + 0.7 σ − 3mm (4)
2.6. Dose Calculation on Treatment Planning System
In the procedure of dose calculation in Lung Stereotactic Body Radiation therapy, two groups (group 1 and group 2) were divided by Equation (3), Equation (4), respectively. In each group, the dose analysis and comparison among Yukinori, Remeijer and Van Herk were all achieved. However, all PTV margins obtained had to been adjusted before dose recalculation due to the numerical precision of Varian Eclipse treatment planning system. Thus, in group 1, the actual PTV margins in three directions were Yukinori (6 mm, 6 mm, 3 mm), Van Herk (6 mm, 6 mm, 4 mm), Remeijer (6 mm, 6 mm, 4 mm), and in group 2, the PTV margins were Yukinori (4 mm, 4 mm, 1mm), Van Herk (4 mm, 4 mm, 1 mm), Remeijer (5 mm, 4 mm, 1 mm). In group 1 On account of the same numerical precision of Van Herk and Remeijer, we only compare Yukinori method with Remeijer method. In group 2 Yukinori method share the same precision with Van Herk so the comparison between Yukinori method and Remeijer is valid. Four treatment plans (plan441, plan541, plan663, plan664) were generated in the Eclipse system with the progressive resolution optimizer 3. Dose distributions were calculated using the anisotropic analytic algorithm (AAA, ver.11.0.31, Varian Medical Systems, Palo Alto, CA) with a dose calculation grid of 2 mm. R programming language (R Development Core Team 2010) was used for all dose parameters analysis.
3.1. Calculation of Interfractional Errors with Yukinori Method
20 patients’ setup errors were all analyzed with IBM SPSS 19.0. The distributions of interfractional errors including mean errors and standard deviation in lateral (X), longitudinal (Y), vertical (Z) dimensions are shown in Figure 2. It should be noted that the distributions of the errors in three directions do follow the normal distribution.
Table 3 lists the systematic (∑) and random (σ) errors calculated by Yukinori method in X, Y, Z directions. For the standard errors, the lateral, longitudinal, vertical had mean values of 1.11 mm, 0.28 mm, −0.65 mm, respectively, and standard deviations of 2.11 mm, 2.02 mm, 1.06 mm, respectively. For the random errors, the mean values were 1.82 mm, 1.66 mm, 1.25 mm, respectively, and standard deviation were 2.39 mm, 2.14 mm and 1.61 mm, respectively. The standard deviation (SD) of systematic and random errors in three directions were all larger than Se An Oh  et al., analyzing the motion errors of brain tumors with ExacTrac 6D Image Guide System. Although the immobilization devices were used to restrict the location of patients and reduce the motion of lung, the effect of breathing still cannot be ignored.
Figure 2. Distribution of interfractional errors of 20 patients in three dimensions.
Table 3. Systematic errors (∑) and random errors (σ) using Yukinori variance method in three dimensions.
However, variance component analysis can only provide the whole discrepancy of errors in three directions instead of pairwise comparison. In this work, the Dunnett test were used to verify the difference between certain two directions. Table 4 shows pairwise comparison between the interfractional set-up errors in three directions using variance component analysis with Dunnett tests. The variation of errors between X and Y direction is not significant (P > 0.05) while the variation between X and Z, Y and Z are significant (P < 0.05). This indicates the use of immobilization devices indeed reduces breath motion uncertainty in longitudinal direction to some extents.
3.2. Comparison between Yukinori Method and Van Herk Method
Table 5 illustrates the comparisons of systematic and random errors between Yukinori method and Van Herk method. For Yukinori method, standard deviation (SD) of systematic errors in three directions are 2.11 mm, 2.02 mm, 1.06 mm, respectively, and standard deviation (SD) of random errors are 2.39 mm, 2.14 mm, 1.61 mm, respectively. For Van Herk method, standard deviation of systematic errors in three directions are 2.30 mm, 2.19 mm, 1.25 mm, respectively, and standard deviation of random errors are 2.14 mm, 1.93 mm, 1.43 mm, respectively. Compared to the Van Herk method, the systematic errors in Yukinori method were smaller because the random component in systematic errors had been eliminated. As a result of the rectification of degrees of freedom, random errors are bigger in Yukinori method.
Table 6 shows the variation of PTV margins between Yukinori method and Van Herk. To obtain a comprehensive comparison, this study used two CTV-PTV margin formula introduced by Van Herk and Stroom respectively. With Equation (3), PTV margins calculated by Yukinori method in X, Y, Z directions were 5.89 mm, 5.54 mm and 3.25 mm, respectively and PTV margins by Van Herk were 6.10 mm, 5.73 mm and 3.51 mm, respectively. With Equation (4), PTV margins calculated by Yukinori method in three directions were 3.95 mm, 3.55 mm and 0.78 mm, respectively, and PTV margins by Van Herk were 4.25 mm, 3.83 mm and 1.13 mm, respectively. In three directions, PTV margins using Yukinori method were generally a bit smaller than ones using Van Herk method. The rate of change of the two methods in three directions were 3.44%, 3.32% and 7.41% with Equation (3) and 7.06%, 7.31% and 30.97% with Equation (4).
Table 4. Pairwise comparison of interfractional errors in three dimensions.
Table 5. Comparison of standard deviation of systematic errors (∑) and random errors (σ) between Yukinori and Van Herk in three directions.
Table 6. Comparison of PTV margins between Yukinori and Van Herk in three directions.
Note: (3), (4) is Equation (3), Equation (4) respectively; Rate of change = (Van Herk − Yukinori)/Van Herk.
3.3. Comparison between Yukinori Method and Remeijer Method
The application of Remeijer method on errors in Lung SBRT was also achieved and comparisons with Yukinori method are shown in Table 7. The results indicates Yukinori method is indeed a revision of systematic errors for Remeijer method and the systematic errors in Remeijer method are bigger than Yukinori method.
The results in Table 8 show the comparison between these two approaches applied in PTV margins with R624-SCF SBRT immobilization device. To obtain a comprehensive comparison, this study uses two CTV-PTV margin formula introduced by Van Herk and Stroom respectively. With Equation (3), PTV margins calculated by Yukinori method in X, Y and Z directions were 5.89 mm, 5.54 mm, 3.25 mm, respectively, and PTV margins by Remeijer were 6.39 mm, 5.98 mm, 3.69 mm, respectively. With Equation (4), PTV margins calculated by
Table 7. Standard deviation of systematic errors (∑) and random errors (σ) of Yukinori (ANOVA) and Remeijer method.
Table 8. Comparison of PTV margins between Yukinori and Remeijer in three directions.
Yukinori method in three directions were 3.95 mm, 3.55 mm, 0.78 mm, respectively, and PTV margins by Remeijer were 4.57 mm, 4.10 mm, 1.33 mm, respectively. The rate of change of PTV margins with Equation (4) in vertical direction is about 41.35%, the largest among three directions while with Equation (3), the rate of change of PTV margins is also largest in vertical direction, and the reduction rate is nearly 11.92%. Generally, Yukinori method indeed reduced the PTV margins no matter which formula was selected.
3.4. Dose Analysis
Table 9 and Table 10 respectively list the statistical PTV parameters in group 1 and group 2. The mean volume of PTV441 and PTV 541 are 16.79 ± 8.35 cm3 and 17.66 ± 8.69 cm3, which shows significant difference (P < 0.001) between them. The mean volume of PTV 663 and PTV664 are 23.78 ± 3.19 cm3 and 25.19 ± 4.33 cm3, which shows significant difference (P < 0.001) between them. This analysis indicates a meaningful decrease of PTV volume using Yukinori method compared to Van Herk method and Remeijer Method. The boxplot (Figure 3) also displays PTV volume distribution with diverse margins.
During Lung SBRT procedure, due to the influence of setup uncertainty, breathing, heart beat and organ geometric transformation, location deviation of tumor and normal structure, it tend to cause the insufficient dose distribution in target area and excessive dose distribution in normal tissues, which aggravated radiation damage. To improve tumor local control, PTV margins should be calculated as accurate as possible for higher dose gradient, especially in Lung SBRT.
According to AAPM-TG83 report, lung dose parameters (mean dose, V5, V20) are recorded. Table 11 and Table 12 lists the Lung dose parameters in group 1 and group 2. In plan441 of group 1, mean dose, V5 and V20 were 4.75
Figure 3. The boxplot of PTV volumes with four different margins (PTV441, PTV541, PTV663, PTV664).
Table 9. Statistical PTV volume parameters in group 1.
Table 10. Statistical PTV volume parameters in group 2.
Table 11. Lung dose comparison with two PTV margins (Plan441, Plan541).
Table 12. Lung dose comparison with two PTV margins (Plan663, Plan664).
Gy, 23.22% and 6.71% respectively and in plan541, mean dose, V5 and V20 were 4.85 Gy, 23.57% and 6.92% respectively. There is no significant difference (P > 0.05) between these two plans. In plan663 of group 2, mean dose, V5 and V20 are 5.02 Gy, 24.17% and 7.23% respectively and in plan664, mean dose, V5 and V20 are 5.07 Gy, 24.70% and 7.32% respectively. The difference of two plans is also not significant (P > 0.05). Figure 3 illustrates distribution of lung dose parameters in SBRT.
PTV parameters D98, D50 and D2 are also recorded according to AAPM-TG83 report. Table 13 and Table 14 list the PTV parameters with different margins. In plan441 of group 1, D98, D50 and D2 were 49.53 Gy, 51.71 Gy and 53.61 Gy respectively and in plan541, D98, D50 and D2 were 49.48 Gy, 51.77 Gy and 53.70 Gy respectively but the difference is not significant (P > 0.05) between these two plans. In plan 663 of group 2, D98, D50 and D2 were 49.30 Gy, 52.45 Gy and 54.90 Gy respectively and in plan 664, D98, D50 and D2 were 49.09 Gy, 52.69 Gy and 55.00 Gy respectively but the difference is also not significant (P > 0.05) between these two plans. Figure 4 illustrates distribution of PTV dose indices in SBRT.
Figure 5 displays the Dose-Volume Histogram difference between two PTV margins (441 and 541). According to the statistical results, the difference between two PTV margins is not significant (P > 0.05). Figure 6 displays the Dose-Volume Histogram difference between two PTV margins (663 and 664) and there is also no significant difference (P > 0.05) between two PTV margins. The comparison of isodose curves of different PTV margins in group 1 and group 2 is respectively shown in Figure 7 and Figure 8.
Variance component analysis (ANOVA) in radiotherapy introduced by Yukinori is not popular mainly because the validity and effectiveness still require to be verified in large amount of experiments systematically. The purpose of this paper is to verify effectiveness of this methodology on the Lung SBRT, which belongs to hypofractionated radiotherapy referred in the previous paper. We studied the variance component analysis in radiotherapy to explore its effect on PTV margins and cumulative dose distribution in lung SBRT compared with other two conventional methods raised by Van Herk and Remeijer.
However, recent work has demonstrated the limitations of the Van Herk formula for setup errors calculation. Based on the method raised by Van Herk, Gordon  et al. reported convolution method to modify the systematic and
Table 13. Tumor dose comparison with two PTV margins (Plan441, Plan541).
Table 14. Tumor dose comparison with two PTV margins (Plan663, Plan664).
Figure 4. Boxplot of Tumor dose indices with differing PTV margins, (a) records D2 of four different margins, (b) records D50 of four different margins (c) records D98 of four different margins.
random errors, especially for SBRT. This method introduced effective systematic errors, the purpose of which is to add an offset to the isocenter, and effective
Figure 5. DVH difference of target volume between two PTV margins corresponding to two methods (Yukinori and Remeijer), (b) DVH difference of organ at risk (Lung) between two PTV margins corresponding to two methods (Yukinori and Remeijer).
Figure 6. Dose Volume Histogram (DVH) between PTV663 and PTV 664. (a) DVH difference of target volume between two PTV margins corresponding to two methods (Yukinori and Remeijer), (b) DVH difference of organ at risk (Lung) between two PTV margins corresponding to two methods (Yukinori and Remeijer).
Figure 7. The isodose curves of Lung SBRT. *(a) PTV441; (b) PTV541; (c) PTV441 and PTV541.
Figure 8. The isodose curves of Lung SBRT. *(a) PTV663; (b) PTV664; (c) PTV663 and PTV664.
random errors to fully explain correlation between average deviation for isocenter location and each fraction of treatment, and Equation (5) gave the detailed formula:
Yukinori also argued that in Remeijer method, the systematic errors would be overestimated if the number of fractions is limited. In hypofractionated schedules   , like SBRT, SABR, which have been recently applied to cancer in several organs, the overestimation of systematic errors should be carefully considered. Stereotactic Body Radiation Therapy, on account of finite fractions, high level dose per fraction, satisfied the request of hypofractionated radiotherapy so it is meaningful to apply Yukinori method to Lung SBRT, reducing PTV margins compared with Remeijer. According to the explanation given by Yukinori, due to the restriction of treatment fraction and the limitation of sample size, there exists influence of random term in the Remeijer method formula, which caused the overestimation of systematic errors. Thus, Yukinori modified Remeijer formula, listing below:
, , , are standard deviation of systematic and random errors in Yukinori and Remeijer method respectively.
Considering the rectification of degrees of freedom  , we only compared Yukinori method and Remeijer method in the late dose distribution calculation, which indicates no significant difference between them. Van Herk defined that the standard deviation of systematic and random errors are called uncertainty in some situations. Peulen  et al. reported that 8 mm CTV-PTV margin is suitable for mid-ventilation Lung SBRT without immobilization devices. In our work, PTV margin is almost 4 mm, smaller than Peulen, with the immobilization devices. Bregmans  et al. suggested the optimal PTV margin for clinical treatment is 5 mm, which is consistent with our study. On account of the lack of intrafractional errors, the PTV margins will to some extents be smaller than the expected value but it won’t influence the effectiveness of the comparison of PTV margins and Dose Volume Histogram.
Actually, in Lung Stereotactic Body Radiation Therapy (SBRT), the accuracy request for setup errors and PTV margins is more strict so decreased PTV margins calculated with Yukinori method indeed may cut down the dose in normal tissue compared to the Van Herk and Remeijer Method. However, the findings on dose distribution in Lung and PTV shows no significant difference (P > 0.05), which indicates the similarity of the Yukinori method and Remeijer method in dose calculation. One possible reason is that the comparison of PTV margins seems to be insufficient without statistical meaning and it still need more samples to do further verification. Also, the statistical method like Bootstrap with no parameters  can be utilized to obtain the variation between PTV margins. It’s possible that due to the clinical feature of lung radiotherapy with large amount of air cavity, impeding the accumulation of dose, tiny PTV margin reduction won’t cause distinct optimization of dose distribution in TPS according to the Khan’s the physics of radiation therapy  .
In this work, we ignored the intrafrational errors like random motion of tumor and patients on account of the lack of 4D CT simulation and 4D CBCT. However, in the fractional treatment, the image registration method is the soft tissue method (tumor), which can to some extent fully extract a comprehensive effect of setup errors and error of tumor random motion. Some papers have proposed that the CTV-PTV margin formula is not linear because the normal tissues around tumor are inhomogeneous in the actual treatment while the normal tissues in lung around tumor are nearly homogeneous. This research only concentrated on the verification of the new methodology of variance component analysis and we merely selected the most popular CTV-PTV formula raised by Stroom. Although, Van Herk’s CTV-PTV margin formula was mainly used in prostate tumor, the ultilization of it provides a reference to better understand the effectiveness of variance component analysis in Lung SBRT.
In future work, multi-factors variance component analysis will be applied to the calculation of the interfractional setup uncertainties and PTV margins in Stereotactic Body Radiation Therapy (SBRT), to find out whether there exists other factors influencing cumulative dose distribution due to the setup uncertainties and PTV margins. Precision of PTV margins may be dependent on the variation of Treatment Planning System. The examination on brain, liver or other organs should be implemented as well to obtain a more general conclusion. Except for these observations, in this work, we did not find strong relationship between variance component analysis and dose cumulative distribution but expanding sample size is probably a good choice to get better results.
This study evaluated the setup uncertainties of 20 patients treated with Lung SBRT and applied variance component analysis proposed by Yukinori to calculate CTV-PTV margins and corresponding dose distribution. The Yukinori method showed a distinguished effect on the reduction of PTV volumes compared with other two methods (Van Herk and Remeijer). However, with the reduction of PTV margins, there is no evidence showing an outstanding improvement and optimization of dose distribution on PTV and OARs using Yukinori method. Thus, the methodology of variance component analysis, as a novel approach to Stereotactic Body Radiation Therapy, still requires more clinical verification.