Physica Medica 31 (2015) 1e8
Contents lists available at ScienceDirect
Physica Medica journal homepage: http://www.physicamedica.com
Original paper
Normal tissue complication probability models for severe acute radiological lung injury after radiotherapy for lung cancer M. Avanzo a, *, M. Trovo b, C. Furlan b, L. Barresi a, A. Linda c, J. Stancanello d, L. Andreon b, E. Minatel b, M. Bazzocchi c, M.G. Trovo b, E. Capra a a
Medical Physics, CRO Aviano National Cancer Institute, 33081 Aviano, Italy Radiation Oncology Department, CRO Aviano National Cancer Institute, 33081 Aviano, Italy Institute of Diagnostic Radiology, Department of Medical and Biological Sciences, University of Udine, 33100 Udine, Italy d MRI Applications and Workflow, General Electric, 78533 Buc, France b c
a r t i c l e i n f o
a b s t r a c t
Article history: Received 28 February 2014 Received in revised form 7 October 2014 Accepted 8 October 2014 Available online 30 October 2014
Purpose: To derive Normal Tissue Complication Probability (NTCP) models for severe patterns of early radiological radiation-induced lung injury (RRLI) in patients treated with radiotherapy (RT) for lung tumors. Second, derive threshold doses and optimal doses for prediction of RRLI to be used in differential diagnosis of tumor recurrence from RRLI during follow-up. Methods and materials: Lyman-EUD (LEUD), Logit-EUD (LogEUD), relative seriality (RS) and critical volume (CV) NTCP models, with DVH corrected for fraction size, were used to model the presence of severe early RRLI in follow-up CTs. The models parameters, including a/b, were determined by fitting data from forty-five patients treated with IMRT for lung cancer. Models were assessed using Akaike information criterion (AIC) and area under receiver operating characteristic curve (AUC). Threshold doses for risk of RRLI and doses corresponding to the optimal point of the receiver operating characteristic (ROC) curve were determined. Results: The a/bs obtained with different models were 2.7e3.2 Gy. The thresholds and optimal doses curves were EUDs of 3.2e7.8 Gy and 15.2e18.1 Gy with LEUD, LogEUD and RS models, and md of 0.013 and 0.071 with the CV model. NTCP models had AUCs significantly higher than 0.5. Occurrence and severity of RRLI were correlated with patients’ values of EUD and md. Conclusions: The models and dose levels derived can be used in differential diagnosis of tumor recurrence from RRLI in patients treated with RT. Cross validation is needed to prove prediction performance of the model outside the dataset from which it was derived. © 2014 Associazione Italiana di Fisica Medica. Published by Elsevier Ltd. All rights reserved.
Keywords: Lung tumors IMRT Lung injury NTCP Tomotherapy
Introduction Radiological radiation-induced lung injury (RRLI) is the increase in lung tissue density in patients treated with radiotherapy (RT) as seen on follow-up chest x-rays or computed tomography (CT) [1]. It is originated by the replacement of lung tissue with exudates or fibrotic tissues mostly confined into the radiation fields [2]. These changes are categorized as acute before 6 months (pneumonitis), and late after more than 6 months (lung fibrosis) [3e5]. The majority of patients with such radiographic changes remain asymptomatic [3]. Understanding doseeresponse relationship of tissues is essential for optimization and evaluation of RT treatment plans. Several studies indicated that appearance and degree of RRLI are correlated
* Corresponding author. Tel.: þ39 0434659175. E-mail address:
[email protected] (M. Avanzo).
with dosimetric surrogates of the 3D dose distribution such as mean lung dose (MLD) or percent lung volumes receiving more than a certain dose [1,2,5e7]. Models of Normal Tissue Complication Probability (NTCP) correlate the occurring of side-effects after RT with dose-volume data, fractionation scheme and nondosimetric variables [8e12]. For these reasons NTCP models are currently used as objective functions in the inverse planning process and to compare treatment plans [13]. A number of authors have investigated local doseeresponse relationships from changes in CT density in lung cancer RT patients [1,5,7,14] The doseeresponse for RRLI assessed from CT follow-up scans has been described in the past through NTCP modeling for breast cancer patients treated with 3D-conformal radiotherapy (3D-CRT) using the end-point of mild RRLI [10]. However, at best of our knowledge, appearance of severe patterns of early RRLI has never been described trough NTCP models.
http://dx.doi.org/10.1016/j.ejmp.2014.10.006 1120-1797/© 2014 Associazione Italiana di Fisica Medica. Published by Elsevier Ltd. All rights reserved.
2
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
The aim of the present study was to derive threshold dose levels for early severe RRLI to be used as optimization constraints during planning of the treatment and for evaluation of the safety of proposed treatment plans. For this purpose, we modeled the risk of severe patterns of early RRLI using previously established NTCP models with DVHs corrected for fraction size. The best-fit parameters of doseeresponse models, including the fractionation parameter, a/b, were derived from fitting observed incidences of severe patterns of acute RRLI in patients treated with RT for lung tumors. The capability of models to predict severe RRLI was then tested through a series of statistical tests. The second aim was to derive dose levels for classification of patients at high risk of RRLI. In highly conformal RT the differential diagnosis between benign changes on follow-up CT from progression or recurrence is a difficult task because areas of dense consolidation usually develop around the treated tumor [4]. Identification of patients at high risk for RRLI can help the physician in differential diagnosis of RRLI from tumor recurrence or persistence. Methods and materials Patient treatment Forty-five Non-small Cell Lung Cancer (NSCLC) patients were treated with radical RT at our institution. Each patient underwent simulation using a 32-slice CT (Toshiba Aquilion LB, Toshiba Medical Systems Europe, Zoetermeer, the Netherlands). During the planning CT scan, patients were lying on a flat table in a supine and overhead arm position and allowed to breathe freely. No intravenous contrast was used. The lung was segmented on the CT images automatically in the Oncentra Masterplan (Elekta, Stokholm, Sweden) treatment planning system (TPS) by use of a thresholding method at a CT density between 300 and 400 Hounsfield units and edited manually based on visual inspection. The intrathoracic gross tumor was therefore automatically excluded from the lung during the segmentation. The planning target volume (PTV) included the gross tumor volume with an additional margin of 10 mm in the craniocaudal and 5 mm in the axial directions. Patients and RT treatment characteristics are shown in Table 1. Forty-two patients were treated with Helical Tomotherapy (HT) [15] delivered with the Hi-Art II System (Accuray Inc, Sunnyvale, CA). Three patients were treated with IMRT delivered with a Varian CLINAC Trilogy system (Varian Medical Systems, Palo Alto, CA). Intensity modulation was achieved using a sliding window technique. Before treatment, each patient underwent kilovolt conebeam CT imaging for set-up correction. Tomotherapy treatments were planned using the Tomotherapy TPS (Accuray Inc, Sunnyvale, CA), and IMRT and 3D-CRT treatments using the Eclipse (Varian Medical Systems, Palo Alto, CA). In order to derive differential dose volume histograms (DVHs) of all the treatments in the same format, the Tomotherapy treatment plans were exported to Eclipse. DVHs of the ipsilateral lung and PTV were calculated within Eclipse and exported in text format with a dose resolution of 5 cGy. Dose received from imaging procedures was not considered in calculations. All of the Tomotherapy treatment plans were exported to the Eclipse TPS. Differential dose-volume histograms (DVHs) of the ipsilateral lung were calculated with the Eclipse TPS and exported with a dose resolution of 5 cGy. The following dosimetric parameters were recorded for each patient: percent volumes of the lung receiving more than 5 Gy and 20 Gy, mean dose to the lung, average total dose to the PTV, average dose per fraction to PTV.
Patient follow-up The dosimetric and follow-up data were reviewed with the approval of our Institutional Review Board and after obtaining informed consent from each patient. Diagnostic CT scans were acquired during follow-up at 45 days, 3 months and 6 months after the end of RT and were read by a radiologist and a radiation oncologist in consensus. The lung injuries occurring within 6 months from the completion of radiation therapy were defined as acute RRLI and were classified into 5 categories [4,16]: (1) diffuse consolidation; (2) patchy consolidation and ground-glass opacities (GGO); (3) diffuse GGO; (4) patchy GGO; (5) no evidence of increasing density. Diffuse consolidation (Fig. 1a) and patchy consolidation and GGO (Fig. 1b) were defined as patterns of severe RRLI in this study. When more than one CT scan within the first 6 months from RT were available on the same patient, the patient outcome was assumed as the lowest category of RRLI seen on the CT scans. The following patient characteristics that were previously studied as risk factors for RRLI [2,3,14,16] were recorded for each patient: sex, age, smoking status at the time of treatment (smoker vs non-smoker), tumor location (right vs left lung, lower vs upper/ middle), administration of chemotherapy or corticosteroids. The ManneWhitney U test and the chi-square tests were used to assess the correlation between radiographic lung injury and pretreatment clinical or dosimetric variables. Tests of statistical significance were two-sided. Bivariate correlation of dosimetric factors and grade of RRLI were assessed using Spearman's rank correlation. A 5% significance level was considered for the analysis. Models for NTCP In order to account for the variability in dose fractionation, the doses Di in the differential ipsilateral lung DVHs of were corrected into equivalent doses at 2 Gy/fraction, EQD2,i, using the linear quadratic model [17]:
0
1
BDNi B
EQD2;i ¼ Di @
2
þ abC C þ ab A
(1)
where N is the number of fractions in the treatment, i is the subvolume of the organ irradiated with dose Di in the differential DVH, and a/b ratio the parameter of the linear-quadratic model for the lung. To account for the effect of irradiated volume and dose distribution on risk of side effect, we used the equivalent uniform dose (EUD). This is the dose that, when delivered uniformly to the organ or tissue, is assumed to produce the same NTCP as the inhomogeneous dose distribution [18].
EUD ¼
X
!n 1 n
vi ,EQD2;i
(2)
i
where vi is the i-th fractional sub-volume of the organ irradiated with dose EQD2,i and the parameter n describes the volume effect of the irradiated organ or tissue. For dose response modeling, we used previously established NTCP models [8,9,19] whose equations are described in appendix A. Two models based on EUD were chosen to model RRLI: the LymanEUD (LEUD) and the Logit-EUD (LogEUD) [9] models. These models rely on the hypothesis that all volume elements inside an organ have the same importance for the appearance of the side effect. The relative seriality (RS) or Kallman s-model [19] and the population-averaged critical volume (CV) model [8] were also used.
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
3
Table 1 Dosimetric and patient factors for the patients used for the study. Patient n
Age/gender
Tumor location (side/lobe)
PTV size (cc)
Tech-nique
Mean PTV dose
Frac-tions
V20Gy%
V5Gy%
Mean lung dose (Gy)
EUD (Gy)
Smoke/Chemo/ steroids
Time between RT and last CT (m)
RRLI grade
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45
83/F 78/M 52/F 71/M 77/F 74/F 83/F 63/M 72/M 63/M 86/M 78/M 68/M 81/M 89/M 73/M 54/F 86/M 57/M 80/M 74/M 68/F 85/M 80/M 67/M 72/M 84/M 77/M 67/M 68/F 61/F 49/F 82/M 79/M 65/F 73/M 70/M 69/M 70/F 64/F 53/M 82/F 64/M 77/F 74/M
R/U R/U L/U R/U L/L L/U R/L R/U L/U R/L L/L L/L R/U R/L R/L R/U R/U L/L R/U R/U R/U L/L L/L L/U L/L R/U R/U L/U R/U R/L R/U L/U R/L R/L R/U R/U L/U R/U R/U L/L R/U L/L L/L L/U L/L
33 373 379 165 86 447 243 192 62 780 19 40 783 148 291 23 117 155 370 75 577 85 63 67 549 385 85 98 487 46 23 197 438 169 350 58 178 107 126 571 1170 64 351 20 5
Tomo Tomo IMRT Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo IMRT Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo Tomo IMRT
53.0 63.5 58.3 53.0 63.0 60.9 56.3 61.4 51.4 59.6 53.5 52.2 61.0 55.3 51.9 54.2 54.5 58.4 59.9 54.7 58.4 55.6 41.4 55.0 59.8 61.3 53.0 46.9 59.7 54.2 53.5 58.3 57.6 62.4 54.9 53.8 30.8 40.2 56.8 59.0 61.6 54.1 61.6 55.3 53.1
5 25 25 5 8 25 8 25 5 25 5 5 25 8 8 5 5 8 25 8 25 8 5 8 25 25 8 5 25 5 5 25 25 25 25 5 5 5 8 25 25 5 25 5 8
12.0 54.0 54.0 20.0 16.0 33.0 43.0 57.0 15.0 61.0 5.0 9.0 38.0 21.0 27.0 4.0 22.0 33.0 46.0 11.0 48.0 30.0 10.0 15.0 59.0 35.0 19.0 16.0 41.0 19.0 12.0 39.0 47.0 38.0 38.0 16.0 4.0 14.0 24.0 52.0 35.0 17.0 41.0 11.0 1.0
27.0 64.0 61.0 50.0 47.0 67.0 72.0 97.0 47.0 92.0 30.0 28.0 52.0 53.0 62.0 23.0 34.0 69.0 87.0 19.0 64.0 63.0 52.0 43.0 84.0 48.0 38.0 58.0 69.0 52.0 29.0 53.0 64.0 58.0 78.0 42.0 26.0 52.0 59.0 88.0 50.0 42.0 66.0 29.0 14.0
6.0 27.0 30.0 12.0 11.0 17.0 22.0 28.0 10.0 30.0 5.0 6.0 21.0 12.0 15.0 4.0 10.0 17.0 26.0 6.0 24.0 16.0 8.0 9.0 29.0 17.0 10.0 11.0 22.0 12.0 7.0 18.0 23.0 18.0 20.0 9.0 4.0 9.0 14.0 30.0 18.0 10.0 21.0 6.0 2.0
16.0 30.0 34.1 26.5 18.9 17.9 39.0 27.0 20.7 30.4 9.5 12.8 24.4 21.4 25.0 8.3 28.1 29.7 26.9 11.8 27.0 25.6 13.2 16.1 29.5 19.0 17.7 21.4 23.9 27.1 16.3 20.7 24.8 19.1 19.4 20.0 6.5 16.6 24.4 32.3 20.9 23.5 22.7 14.2 2.7
Y/N/N Y/Y/N N/Y/N Y/N/N N/N/N Y/Y/N N/N/N N/Y/N Y/N/N Y/Y/N Y/N/N N/N/N N/N/N Y/N/Y N/N/N Y/N/N Y/Y/Y N/N/N Y/Y/N N/N/N N/Y/N N/N/N N/N/N N/N/N Y/Y/N N/Y/Y N/N/Y N/N/N N/Y/Y N/N/N N/N/N Y/Y/N N/N/N N/Y/Y N/N/N N/N/N N/N/N N/N/N N/N/Y N/Y/N Y/Y/N N/N/N Y/Y/N N/N/N N/N/N
3 1.5 6 6 6 3 3 6 1.5 6 6 6 6 3 3 3 3 3 3 3 6 6 1.5 1.5 1.5 6 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5
5 2 2 5 2 1 1 4 5 1 5 1 4 2 5 5 2 1 1 4 5 1 5 2 1 1 2 2 2 1 1 4 1 2 4 4 5 5 5 2 1 5 5 5 4
In these models, normal tissues are composed of a large number of microscopic functional subunits (FSUs) which are independent. In organs with serial architecture the side-effect occurs if a single FSU is damaged similarly to a metal chain under stress. In organs with parallel architecture the complication occurs if a number of FSUs, or functional reserve, is destroyed. Organs may also have a mixed structure, containing both serial and parallel sub-structures [19]. In the RS model, the organ at risk of complication is considered as an arbitrary combination of serial and parallel series of FSUs. In the population averaged CV model, the organ is damaged when a significant volume, termed critical volume, is damaged from radiation. The inter-patient variability of the critical volume in a group of patients is accounted for in the model by assuming that it is distributed following a log-log-normal distribution [8]. Parameter estimation and data fitting The radiobiological parameters maximizing the logarithm of the maximum likelihood function [9] were determined by fitting the occurrence of early severe RRLI in our patients with the NTCP
models. The maximum likelihood method was used to find the model parameters which yielded the best fit to clinical outcomes. The quantity log-likelihood was maximized by using the NelderMead simplex algorithm [20].
L¼
N Y
NTCPi i ð1 NTCPi Þð1epi Þ ep
(3)
i¼1
where NTCPi and NTCP0 are the NTCP on the i-th patient and on the patients who receive zero dose, respectively, and epi is equal to 1 if the complication has occurred and 0 otherwise. In the fitting process, the set of parameters of the model which maximizes the natural logarithm of L is determined. The confidence intervals for the obtained parameters were calculated using the bootstrap method [21], which consists of generating a large number of bootstrap samples from the original data. A bootstrap dataset is a series of data (DVH, number of fractions, occurrence of complication), each from a patient chosen randomly from the cohort of patients. Each bootstrap sample had the size of the patient cohort, fortyefive patients. A distribution of
4
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
Fig. 1. Radiographical changes after RT of the lung considered as the end-point in NTCP models: diffuse consolidation (a); patchy consolidation and ground-glass opacities (b).
values for each model parameter was obtained by repeating the fitting process on each bootstrap sample. The 95% Confidence Intervals (CIs) for each parameter were calculated from its distribution. In this analysis the CIs were calculated from 1000 bootstrap samples. Data fitting and statistical analysis was performed using MATLAB software (The Mathworks, Inc, Natick, MA). Models evaluation The models were ranked using the Akaike's second order information criterion (AIC) [22]. AIC is used for relative evaluation of models, and was previously used for selection of NTCP models [23]. The model with minimum AIC has the best trade-off between goodness of fit to clinical data and number of parameters. The models predictive power for the risk of acute severe RRLI in our patients was studied using different tests. The capability of the model to predict complication rate was evaluated using the area under receiver operating characteristic curve (AUC) curve [24]. 95% CIs for the AUC were calculated using the bootstrap method with 1000 bootstrap samples. The test was considered successful if the model reached AUCs, with their 95% CI, superior to the reference area of 0.5. The validity of the model was also assessed through the permutation test (Ptest), following the method by Xu et al. [25]. In this test, the endpoints of the patients are randomly exchanged between patients. Following this method, 1000 permutated datasets were generated. The model to be tested was fitted to each permutated dataset and the AUC was calculated from each obtained model on its permutated endpoint dataset. AUC of the original data should be outside the 95% confidence bound of the distribution of AUC from the permutated data, so that the test measures the probability of obtaining the observed prediction performance by chance. Models which passed both criteria, AUC CIs >0.5 and AUC > upper 95% band of the Ptest, were considered as usable for a reliable prediction of RRLI. The threshold doses for risk of RRLI were determined by finding the x-intercept to the tangent line of the curve at dose of TD50 [7] with the exception of the CV model which, because of a different shape of doseeresponse, does not allow extraction of threshold
dose with this method. The threshold in this case was assumed as the quantity of damaged lung volume, md, corresponding to an NTCP of 5%. In order to provide dose levels to classify patients at high risk of RRLI, we searched the NTCP values which correspond to the optimal points of their receiver operating characteristic (ROC) curve [24] and the associated dose levels. The ManneWhitney U-test test was used to assess the correlation between the dosimetric variables EUD and md, used in the models, and the occurrence of RRLI. The correlation of EUD with severity of RRLI was assessed using the Spearman's correlation coefficient. Results Patient results Patient and dosimetric data are shown in Table 1. Twenty-four (53.3%) patients experienced acute severe RRLI, 13 of which developed a diffuse consolidation pattern and 11 a patchy consolidation and GGO pattern. The results of statistical tests used to assess the correlation between patient variables and severe RRLI are shown in Table 2. None of the patient-related binary variables: age, sex, volume of the ipsilateral tumor, tumor laterality, lung lobe, administration of chemotherapy, corticosteroids, smoking status, were significantly related to occurrence of early severe RRLI. PTV size and average dose per fraction to PTV were also unrelated to risk of RRLI. The following variables related to dose distribution were significantly related to occurrence of severe RRLI: mean lung dose, V5, and V20. The number of fractions of the treatment was significantly related to occurrence of the side effect. Results of NTCP modelling In Fig. 2 aed, the modeled NTCPs with 68% confidence bands are plotted versus dose or damaged volume and compared with the observed incidence of the endpoint (% of patients with severe early RRLI) calculated from patients grouped in 4 bins of dose. The
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
analysis, all the NTCP models had the AUCs with their 95% CI superior to the reference area of 0.5. All the models passed the permutation test, as they had and AUC > upper 95% band of the Ptest. The EUDs calculated used for each model derived were correlated with the occurrence and severity of RRLI. The best model according to the Akaike's information criterion was the LEUD model.
Table 2 Univariate analysis of clinical and dosimetrical factors for predicting the development of early severe RRLI. For the binary variables (sex, tumor laterality, lower vs upper lung lobe, smoking status, administration of chemotherapy and corticosteroids) the Chi square test was used to compare incidence of the side effect on patients with and without presence of the factor considered. For variables that were not dichotomous, the ManneWhitney U test was used. Variable
p value
sex age Ipsilateral lung volume Tumor laterality (left/right lung) Tumor lobe (upper/lower) Chemotherapy Corticosteroids Smoke status Volume of PTV Number of fractions V20 V5 MLD Mean dose to PTV per fraction
0.530 0.554 0.882 0.841 0.143 0.071 0.062 1.000 0.052 0.046 0.035 0.013 0.010 0.199
Discussion Statistical analysis
parameter values obtained from fitting NTCP models to the patient data with 95% confidence intervals are reported in Table 3. Threshold EUDs for risk of severe RRLI were 5.9, 3.2, 7.8 Gy with the LEUD, LogEUD and RS models and md of 0.013 with the CV model. TD50s were 20.3, 18.3, 21.0 Gy with the LEUD, LogEUD and RS models and 55.9 with the CV model. Optimal values of EUD for predicting high risk of severe RRLI obtained were EUDs of 16.1, 15.2 and 18.1 Gy with the LEUD, LogEUD and RS models, and md of 0.071 with the CV model. Results of tests used to evaluate the goodness of fit and predictive power of models are shown in Table 4. From the ROC curve
The statistical analysis showed that severe RRLI in our patient cohort depends purely on variables of the RT treatment such as dose, irradiated volume, and fractionation schedule, as no correlation was found between RRLI and clinical and patient variables. In previous studies on radiographical lung CT changes, the statistical results on the influence between non dosimetrical variables have been mixed. Patient sex has never been found to be predictive for RRLI [2,14]. Age was generally not associated with increased density in CT [2,3], with the exception of one study [14]. Smoke status has been found to be associated with increase on CT density [3] but this result has not been confirmed by other studies [2,16]. Tumor location (upper vs lower and centra vs peripheral) was found to be associated with RRLI in one study [3] but the correlation was not confirmed [6]. In a dose response analysis the incidences of RRLI in the lower and upper tumor were fitted with two distinct dose response curves, but the differences between their parameters were not statistically significant [7]. No correlation has been found between RRLI and chemotherapy [3]. On the other hand, there is a substantial agreement among investigators on dose and irradiated volume as key factors for the occurrence and severity of RRLI. Dosimetric surrogates V7, V10, MLD and EUD were found to be correlated with lung CT changes in single
75
75
NTCP %
(b) 100
NTCP %
(a) 100
50
25
0 0
10
20 30 EUD (Gy)
40
0 0
50
(d)
75
10
20 30 EUD (Gy)
40
50
100
75
NTCP %
NTCP %
50
25
(c) 100
50
25
0 0
5
50
25
10
20 30 EUD (Gy)
40
50
0 0
0.1
0.2 μ
0.3
0.4
d
Fig. 2. NTCP using the LEUD (a) the LogEUD model (b) versus EUD calculated using Eq. (1), NTCP with the RS model versus EUD calculated using Eq. (A.6) (c), and NTCP with the critical volume population averaged model (d) versus the damaged volume, md (Eq. (A.7)). NTCP curves (solid lines) are shown with 68% confidence intervals (dashed lines). Dots indicate incidences of complication on patients with 68% CI bars.
6
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
Table 3 Parameters of NTCP models of obtained from fitting RRLI to clinical data. Model
a/b (value, 95% CIs)
TD50 (value, 95% CIs)
n/s/mer (value, 95% CIs)
m/k/g (value, 95% CIs)
Threshold values
Optimal values for prediction of RRLI
LEUD LogEUD RS CV
2.9 3.1 2.7 3.2
20.3 18.3 21 55.9
0.78 (0.45e2.19) 0.84 (0.43e1.91) 0.42 (0e1.13) mer ¼ 0.09 (0.01e0.28), smer ¼ 0.08 (0.01e0.23)
0.56 3.91 0.84 0.67
EUD ¼ 5.9 Gy EUD ¼ 3.2 Gy EUD ¼ 7.8 Gy md ¼ 0.013
EUD ¼ 16.1 Gy EUD ¼ 15.2 Gy EUD ¼ 18.1 Gy md ¼ 0.071
(0e3.9) Gy (0e4.7) Gy (0e38.1) Gy (0e15.5) Gy
(7.9e31.7) Gy (8.2e32.5) Gy (5.6e33.9) Gy (17.9e78.7) Gy
(0.29e1.38) (1.92e9.06) (0.29e2.07) (0.25e7.99)
Table 4 Indices and statistical tests used to evaluate the models of RRLI. The AIC test evaluates the goodness-of fit and the remaining tests are used to evaluate predictive capability of models for risk and severity of RRLI. Model
Aikike information criterion
Correlation of EUD and risk for RRLI
LEUD LogEUD RS CV
60.6 62.3 62.4 63.9
p p p p
¼ ¼ ¼ ¼
0.0119 0.0127 0.0135 0.0104
Correlation of EUD and severity of RRLI p p p p
fraction stereotactic RT [6]. The association of V20 with increase in CT density has been observed in several studies [3,5,26]. A correlation was found between CT number increases and larger PTV [2] which was not found in the present work. A possible explanation for the role of dose and irradiated volume in the occurrence of RRLI could be that a larger volume irradiated with high dose could produce more mediators and trigger a stronger tissue response [27]. Dose-response model can be corrected for other variables than dose and volume when they are related to the occurrence of the side effect, as the presence of risk factors such as comorbidities [11]. However, in our statistical analysis, no significant association was found between any clinical variable and the risk of RRLI. Consequently, we chose to model the dose response only from variables of the RT treatment: doses and volumes as represented on the DVH and number of fractions. Dose response models In the present study we chose to consider only the ipsilateral lung dose distribution assuming that the density changes in the lung assessed through CT examinations should not be affected by the behavior of the contralateral lung, as assumed in the previous work by Rancati et al. [10]. The radio-sensitivity of ipsilateral lung was described with threshold EUDs of 3.2e7.8 Gy and TD50 between 18.3 and 21 Gy with the LEUD, LogEUD and RS models. The CV model yielded a different value of TD50 ¼ 55.9 Gy as in this model TD50 describes the radio-sensitivity of a single functional subunit of the lung. The parameter md in the CV model indicates the volume of the lung that is damaged during irradiation. In the doseeresponse modeled with the CV model, the risk for RRLI has a threshold value of md ¼ 0.013 and reaches 50% at md ¼ 0.097. Considering that the reference volume in our model was set at the entire lung volume and this was 1811 cc average in our patients, these results indicate that 23.5 cc and 176 cc of damaged lung correspond to 5% and 50% risk of RRLI, respectively. The low values of the parameter for the effect of fractionation, a/b, obtained with all models show that the development of the most severe patterns of acute RRLI depends strongly on fractionation of the treatment. These results agree well with the value of 3 Gy that was previously used for modeling CT voxel density changes [7] and for CT pneumonitis [10].
¼ ¼ ¼ ¼
0.0100 0.0158 0.0273 0.0108
AUC (value, 95% CIs)
AUC from Ptest (value, 95% CIs)
0.72 (0.58e0.92) 0.71(0.59e0.88) 0.71(0.57e0.91) 0.72(0.60e0.92)
0.56 0.55 0.59 0.60
(0.41e0.67) (0.38e0.64) (0.43e0.65) (0.50e0.66)
The parameter TD50 represents the dose at which there is 50% of risk of complication and also the steepness of doseeresponse is at its maximum. Lower TD50 values of 16.1e18.2 Gy were found in a previous work, were dose responses for mild RRLI assessed from CT were derived using NTCP models LEUD, LogEUD and RS [10]. These different results can be explained with the less severe end-point compared with the present work, where severe patterns of early RRLI were considered as endpoints. Our threshold doses of 3.2e7.8 Gy obtained with different models are consistent with previously published clinical results. In the work of Aoki et al., it was observed that minimum dose for lung tissue change in their patients was 24 Gy median in four fractions [26]. Authors reported a median V20 of 5%. Assuming that the average lung volume was the same as in our study, and that the isodose of 24 Gy covered the same volume of the lung volume, this corresponds to EUD of 4.0 Gy calculated using the parameters of the LEUD model (Table 1). In the EUD based models, serial organs have small values for the volume effect parameter n and parallel organs have n close to 1. The values of n determined in the present work of 0.65 and 0.78 in the LEUD and LogEUD models respectively, also suggest that lung has a strong volume effect. The s parameter in the RS model is the ratio between serial and all subunits contained in the organ. Thus the value of s ¼ 0.42 in the RS model also suggests a predominance of parallel subunits. The results obtained in the work of Rancati et al. also indicated a strong volume effect for CT pneumonitis with n ¼ 1 in the LEUD model [10]. Our result of a slightly lower value of n ¼ 0.78 in the LEUD model can be explained with the use of a scoring system for RRLI which includes CT changes, which are related to local dose [1,2,5,7]. The values derived here are based on patterns 1 and 2 of lung injury only, therefore it cannot be exclude that mild patterns 3 and 4 have different parameters. CT density changes are influenced by the irradiated volume, possibly because a larger irradiated volume produces more mediators that trigger tissue response [27]. This effect manifests as a shift of dose-relationship for CT density increase to lower doses when higher volumes are irradiated [5]. These previous findings on radiation induced CT changes suggest that this side-effect depend on dose and irradiated volume as found in our analysis. The models obtained describe well the observed incidences of the most severe patterns of acute RRLI (Fig. 2aed) and have a significantly high capability of discriminating individual risk of RRLI
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
in our patients, as shown by the AUCs with their 95% confidence intervals and the results of the Ptest (Table 3). The model which provides best description of clinical data, according to the AIC, is the LEUD model. All models provided comparable values of AUC with overlapping 95% confidence intervals. Finally, EUDs calculated using model parameters were significantly related with risk of RRLI and its severity. In conclusion, the dose levels derived here can be used to predict the risk of RRLI. Patients with EUDs and md below the threshold levels can be considered to be reasonably safe from risk of s RRLI. The dose levels associated to optimal points in the ROC curve represent the dose levels with best capability of predicting the risk of side effects and correspond to relatively high NTCPs of (35.8, 38.9, 38.5, and 36.9% with the four NTCP models). Therefore they can be used for classification of patients at high risk of RRLI and be of support to the radiologist in the differential diagnosis of tumor recurrence from RRLI. Study limitations The study has a limited sample size. The lack of correlation between RRLI and patient and clinical factors could be likely due to the limited power with this sample size. The limited number of patients does not allow model validation, as full validation of a radiobiological model requires evaluation of predictive capability on a large patient dataset, possibly multicentric [25]. Dose received from imaging procedures was not considered in calculations. However, since dose for image guidance procedures in kVCT and Tomotherapy MVCT is in the order of a few cGy per fraction [28], while the prescribed dose per fraction in our study from 2.5 to 10 Gy, this approximation should only minimally affect the results. It should be also emphasized that the calculated values of EUD are derived from data of patterns of grade 1 and 2 of RRLI. These models should not be applied to milder grades (patterns 2 and 3) as they could have different parameter values. Lung injury is an effect which evolves in time. The CT increases assessed at follow-up in the intervals 0e3 months and 3e6 months have different doseeresponse curves [7]. Therefore a model to predict RRLI should include the effect of time over the development of RRLI, however, the limited number of patients at different followup time intervals in the present study prevents development of a time-dependent NTCP model. Future work will aim at collecting prospectively the data needed for validation of the model, and to study the effect of time after RT on doseeresponse. Conclusions Prospective patient data were analyzed in order to derive NTCP models for RRLI based on dose-volume data of the ipsilateral lung. The models describe well the incidence of side effect in our patient cohort and allowed derivation of dose levels to assess the risk for severe RRLI. The parameters extracted indicate that development of early severe patterns of RRLI depend on dose, irradiated volume and fractionation of the treatment. Cross validation is needed to prove prediction performance of the model outside the dataset from which it was derived. Appendix A A.1 The LymaneEUD model (LEUD)
1 NTCP ¼ F T ¼ pffiffiffiffiffiffi 2p
ZT
e
t2 2
dt
(A.1)
∞
and
T¼
EUD TD50 mTD50
(A.2)
where TD50 is the dose resulting in a 50% complication probability. Parameter m describes the slope of the NTCP curve at TD50. The function Ф in Eq. (A.1) is called probit function. A.2 The LogiteEUD model (LogEUD) The Logit formula describes the doseeresponse relationship for normal tissues through TD50 and k, the slope of the response curve at TD50 [9,11]:
1
NTCP ¼
TD50 EUD
1þ
(A.3)
k
A.3 The relative seriality model (RS) This model is based on the poisson model of cell survival and it accounts for the architecture of the organ through the parameter of relative seriality, s. According to the Poisson doseeresponse relationship, the probability of injury in a volume i receiving equivalent dose at 2 Gy/fraction, EQD2,i (see Eq. (1)) is defined as [19]:
PðEQD2;i Þ ¼ 2
exp
eg
EQD
1 TD 2;i 50
(A.4)
NTCP is calculated as
( NTCP ¼
1
Yh
s ivi 1 P EQD2;i
)1 s
(A.5)
i
Where vi is the fraction of the organ at risk, or the relative subvolume in the differential ipsilateral lung DVH, receiving dose EQD2,i, and g is the slope of the response curve at TD50. In this model, NTCP is not calculated from an equivalent uniform dose. However, a method to calculate EUD from NTCP with the RS model has been introduced [29] making it is possible to plot curves of NTCP and observed incidences of side-effects versus EUD as done in Fig. 2c:
eg lnð ðNTCPÞÞ EUD ¼ TD50 , eg lnðln 2Þ
(A.6)
A.4 The population averaged critical volume model (CV) In this NTCP model an organ is compromised when the fraction of functional subunits damaged md surpasses a critical fraction called critical volume [8]. In a population of patients, it is assumed that the critical volume is logelog normally distributed around a mean value mcr and inter-patient variability smcr. The mean damaged volume in the patient population is calculated as:
md ¼
X i
In the Lyman-EUD (LEUD) model, the probability of side-effect is calculated from the equivalent uniform dose, EUD (see Eq. (2)) using equations:
7
pffiffiffiffiffiffi EQD2;i vi F 2pg ln TD50
(A.7)
where TD50 and g are the parameters describing the probability of damaging a single functional subunit and may differ in value from
8
M. Avanzo et al. / Physica Medica 31 (2015) 1e8
the TD50 and m/g/k of the other models describing the probability of RRLI on patients. The population averaged NTCP is calculated using the same probit function used in Eq. (A.1), with the variable:
T¼
lnð lnðmd Þ þ lnðlnðmer ÞÞ s
(A.8)
where
sz
smer mer ,lnðmer Þ
(A.9)
References [1] Marks LB, Yu X, Vujaskovic Z, Small W Jr , Folz R, Anscher MS. Radiationinduced lung injury. Semin Radiat Oncol 2003;13:333e45. [2] Phernambucq EC, Palma DA, Vincent A, Smit EF, Senan S. Time and doserelated changes in radiological lung density after concurrent chemoradiotherapy for lung cancer. Lung Cancer 2011;74:451e6. [3] Jenkins P, Welsh A. Computed tomography appearance of early radiation injury to the lung: correlation with clinical and dosimetric factors. Int J Radiat Oncol Biol Phys 2011;81:97e103. [4] Linda A, Trovo M, Bradley JD. Radiation injury of the lung after stereotactic body radiation therapy (SBRT) for lung cancer: a timeline and pattern of CT changes. Eur J Radiol 2011;79:147e54. [5] Palma DA, van Sornsen de Koste J, Verbakel WF, Vincent A, Senan S. Lung density changes after stereotactic radiotherapy: a quantitative analysis in 50 patients. Int J Radiat Oncol Biol Phys 2011;81:974e8. [6] Kyas I, Hof H, Debus J, Schlegel W, Karger CP. Prediction of radiation-induced changes in the lung after stereotactic body radiation therapy of non-small-cell lung cancer. Int J Radiat Oncol Biol Phys 2007;67:768e74. [7] Lee S, Stroian G, Kopek N, AlBahhar M, Seuntjens J, El Naqa I. Analytical modelling of regional radiotherapy dose response of lung. Phys Med Biol 2012;57:3309e21. [8] Stavrev P, Niemierko A, Stavreva N, Goitein M. The application of biological models to clinical data. Phys Med 2001;17:1e12. [9] Rancati T, Fiorino C, Gagliardi G, Cattaneo GM, Sanguineti G, Borca VC, et al. Fitting late rectal bleeding data using different NTCP models: results from an Italian multi-centric study (AIROPROS0101). Radiother Oncol 2004;73:21e32. [10] Rancati T, Wennberg B, Lind P, Svane G, Gagliardi G. Early clinical and radiological pulmonary complications following breast cancer radiation therapy: NTCP fit with four different models. Radiother Oncol 2007;82:308e16. [11] Tucker SL, Liu HH, Liao Z, Wei X, Wang S, Jin H, et al. Analysis of radiation pneumonitis risk using a generalized Lyman model. Int J Radiat Oncol Biol Phys. 2008;72:568e74.
[12] Avanzo M, Stancanello J, Trovo M, Jena R, Roncadin M, Trovo MG, et al. Complication probability model for subcutaneous fibrosis based on published data of partial and whole breast irradiation. Phys Med 2012;28:296e306. [13] Allen Li X, Alber M, Deasy JO, Jackson A, Ken Jee KW, Marks LB, et al. The use and QA of biologically related models for treatment planning: short report of the TG-166 of the therapy physics committee of the AAPM. Med Phys 2012;39:1386e409. [14] Ma J, Zhang J, Zhou S, Hubbs JL, Foltz RJ, Hollis DR, et al. Regional lung density changes after radiation therapy for tumors in and around thorax. Int J Radiat Oncol Biol Phys 2010;76:116e22. [15] Mackie TR, Balog J, Ruchala K, Shepard D, Aldridge S, Fitchard E, et al. Tomotherapy. Semin.Radiat.Oncol 1999;9:108e17. [16] Trovo M, Linda A, El Naqa I, Javidan-Nejad C, Bradley J. Early and late lung radiographic injury following stereotactic body radiation therapy (SBRT). Lung Cancer 2010;69:77e85. [17] Vogelius IS, Westerly DC, Cannon GM, Bentzen SM. Hypofractionation does not increase radiation pneumonitis risk with modern conformal radiation delivery techniques. Acta Oncol 2010;49:1052e7. [18] Gay HA, Niemierko A. A free program for calculating EUD-based NTCP and TCP in external beam radiotherapy. Phys Med 2007;23:115e25. [19] Kallman P, Agren A, Brahme A. Tumour and normal tissue responses to fractionated non-uniform dose delivery. Int J Radiat Biol 1992;62:249e62. [20] Lagarias J, Reeds J, Wright M, Wright P. Convergence Properties of the NeldereMead simplex method in low Dimensions. SIAM J.Optim 1998;9:112e47. [21] Tucker SL, Liu HH, Wang S, Wei X, Liao Z, Komaki R, et al. Dose-volume modeling of the risk of postoperative pulmonary complications among esophageal cancer patients treated with concurrent chemoradiotherapy followed by surgery. Int J Radiat Oncol Biol Phys 2006;66:754e61. [22] Akaike H. A new look at the statistical model identification. Automatic Control IEEE Trans 1974;19:716e23. [23] Bakhshandeh M, Hashemi B, Mahdavi SR, Nikoofar A, Vasheghani M, Kazemnejad A. Normal tissue complication probability modeling of radiationinduced hypothyroidism after head-and-neck radiation therapy. Int J Radiat Oncol Biol Phys 2013;85:514e21. [24] Zweig MH, Campbell G. Receiver-operating characteristic (ROC) plots: a fundamental evaluation tool in clinical medicine. Clin Chem 1993;39:561e77. [25] Xu CJ, van der Schaaf A, Van't Veld AA, Langendijk JA, Schilstra C. Statistical validation of normal tissue complication probability models. Int J Radiat Oncol Biol Phys 2012;84:e123e9. [26] Aoki T, Nagata Y, Negoro Y, Takayama K, Mizowaki T, Kokubo M, et al. Evaluation of lung injury after three-dimensional conformal stereotactic radiation therapy for solitary lung tumors: CT appearance. Radiology 2004;230:101e8. [27] Chen Y, Williams J, Ding I, Hernady E, Liu W, Smudzin T, et al. Radiation pneumonitis and early circulatory cytokine markers. Semin Radiat Oncol 2002;12:26e33. [28] Bujold A, Craig T, Jaffray D, Dawson LA. Image-guided radiotherapy: has it influenced patient outcomes? Semin Radiat Oncol 2012;22:50e61. [29] Mavroidis P, Lind BK, Brahme A. Biologically effective uniform dose (D) for specification, report and comparison of dose response relations and treatment plans. Phys Med Biol 2001;46:2607e30.