Analytical modeling of depth-dose degradation in heterogeneous lung tissue for intensity-modulated proton therapy planning.

BACKGROUND AND
PURPOSE: Proton therapy may be promising for treating non-small-cell lung cancer due to lower doses to the lung and heart, as compared to photon therapy. A reported challenge is degradation, i.e., a smoothing of the depth-dose distribution due to heterogeneous lung tissue. For pencil beams, this causes a distal falloff widening and a peak-to-plateau ratio decrease, not considered in clinical treatment planning systems.
MATERIALS AND METHODS: We present a degradation model implemented into an analytical dose calculation, fully integrated into a treatment planning workflow. Degradation effects were investigated on target dose, distal dose falloffs, and mean lung dose for ten patient cases with varying anatomical characteristics.
RESULTS: For patients with pronounced range straggling (in our study large tumors, or lesions close to the mediastinum), degradation effects were restricted to a maximum decrease in target coverage (D 95 of the planning target volume) of 1.4%. The median broadening of the distal 80-20% dose falloffs was 0.5 mm at the maximum. For small target volumes deep inside lung tissue, however, the target underdose increased considerably by up to 26%. The mean lung dose was not negatively affected by degradation in any of the investigated cases.
CONCLUSION: For most cases, dose degradation due to heterogeneous lung tissue did not yield critical organ at risk overdosing or overall target underdosing. However, for small and deep-seated tumors which can only be reached by penetrating lung tissue, we have seen substantial local underdose, which deserves further investigation, also considering other prevalent sources of uncertainty.

Introduction

Lung cancer, which is the main cause of cancer-related mortality [1], is commonly treated with concomitant radio-chemotherapy. However, the application of radio-chemotherapy is limited by radiation-induced side effects in organs at risk close to the tumor. Proton therapy can improve overall survival by reducing the dose to the heart and healthy parts of the lung, as compared to photon therapy, due to the inverse depth-dose profile of the Bragg curve [2], [3]. This holds true especially for early-stage non-small-cell lung cancer, as well as large primary (diameter greater than 5 cm), centrally located, and recurrent tumors [2], [4].

Nevertheless, at present, proton therapy is not widely applied for treatment of lung cancer, which may be attributed to three main challenges: First, internal motion, mainly due to breathing, deteriorates the dose distribution [5], [6], [7]. Second, range uncertainties due to setup errors or defective conversion from x-ray computed tomography (CT) numbers to stopping power ratio can lead to underdosing of the tumor and overdosing of organs at risk [7], [8]. Third, degradation of the depth dose due to the heterogeneous fine structure in lung tissue yields fluctuation in the water-equivalent thickness over the length of the alveoli, which widens the advantageously sharp distal falloff. For clinical treatment planning systems, this degradation is currently not considered because the alveolar diameter is typically 100–500 µm, which is smaller than the voxel size of clinical CT images (≥1 mm3).

First approaches to account for degradation have been verified with experiments and Monte Carlo simulations of monoenergetic proton beams [9], [10], [11], [12]. For a few phantom and patient cases, it was shown that the degradation can lead to underdosing of the tumor and overdosing of organs at risk [11], [13], [14], [15].

Here, we describe the implementation of a model for dose degradation into an analytical dose calculation engine that is fully integrated into a treatment planning workflow. This enabled a first comprehensive study of the degradation effects in ten clinical cases.

Materials and methods

Model implementation

We used the open-source treatment planning toolkit matRad, which was validated against the clinical syngo RT planning software (Siemens Healthineers, Erlangen, Germany) [16]. MatRad’s source code is openly available at http://e0404.github.io/matRad/. MatRad facilitated a pencil-beam algorithm for dose calculation with a Gaussian parameterization of the depth dose. Degradation was implemented in matRad by convolving the pristine Bragg curve, which underlay the dose calculation, with a Gaussian filter. The Gaussian standard deviation σdeg depended on the depth of penetration into lung tissue tlung, which can be computed with matRad’s raytracing algorithm, and the modulation power Pmod

A detailed description of the degradation algorithm is presented in Supplementary Material A.

The modulation power determines the magnitude of the degradation of the Bragg peak as a function of the alveolar diameter and is motivated in detail in Baumann et al. [11]. In agreement with previous publications [9], [11], [15], [17] and consistent with in-house measurements of a lung phantom and a swine lung sample [12], we used a representative value of Pmod,rep = 256 µm as a most probable scenario. Additionally, the maximum reported modulation power [9] of Pmod,max = 750 µm served as a worst case scenario.

Validation

The implementation of our degradation model was verified with in-house shoot-through experiments using proton beam kinetic energies of 100 MeV to 135 MeV and corresponding Monte Carlo simulations [12]. A single Bragg curve was measured twice in a water phantom (pristine curves) and once with a 30-mm slab of lung phantom (Gammex lung 455, Sun Nuclear Corporation, serial no. 45564732, relative electron density 0.30) inserted in front of the water phantom (degraded curve). For each energy, the integral depth-dose curves were normalized by the area under the pristine curves for energy conservation reasons. Both the pristine and the degraded Bragg curves calculated with matRad matched the measurements and Monte Carlo simulations (Fig. 1); the average of the dose difference (relative to the peak) between matRad, measurements, and simulation did not exceed 0.045%.

Fig. 1

Model validation. Integrated pristine depth-dose curves for pure water (right peak) and integrated degraded depth-dose curves (left peak) were measured, Monte Carlo (MC)–simulated, and calculated with matRad. The degraded curves were obtained behind a 30 mm-thick lung phantom (Gammex lung 455, Sun Nuclear Corporation, serial no. 45564732). The proton energy was 108.88 MeV. The pristine and degraded depth doses were normalized to the area under the curve for energy conservation.

Treatment planning

The impact of the degradation was investigated in ten patient cases that were previously treated with photon IMRT at Heidelberg University Hospital11 (target volume properties and beam configurations in Table B1 in Supplementary Material B). All patients consented to the anonymous use of their data for research purposes as part of their treatment agreement with Heidelberg University Hospital. Original treatment plans using protons were generated with syngo RT planning based on the PTV without robust optimization. The cases were grouped by tumor size and location:

Group L: target volumes in the middle of the left or right lung.

Group C: centrally located tumors, i.e., close to the mediastinum.

Group W: widespread tumors that covered a large part of the lung(s).

Beams for group L traveled completely through lung tissue, while beams for groups C and W traveled partly through lung tissue, i.e., parts of the tumor were reached by protons that did not undergo degradation. More specifically for group W, beams traveled through lung tissue only within or behind the planning target volume (PTV).

Original DICOM treatment plans and dose distributions from syngo RT planning were imported into matRad on a (2 mm)3 grid, and then recalculated with the heterogeneity correction turned off (“homogeneous lung”) and turned on (“heterogeneous lung”). In either case, pencil-beam grids and intensities from the original treatment plan generated with syngo RT planning were used. The “homogeneous lung” plans were used as reference to account for dose differences arising from heterogeneity only, rather than those arising from different planning systems.

Dosimetric effects

The recalculated plans were compared by means of transversal dose distributions, absolute dose differences, depth-dose curves, dose-volume histograms (DVHs), quality indicators (QIs) of the clinical target volumes (CTVs) and PTVs for tumor coverage (D95, minimum dose received by 95% of the PTV), and the median dose (D50). As DVHs and QIs did not consider single dose points but were affected by volume averaging, additional contours were added around the PTV, where the strongest degradation effects occurred: The “PTV ring 1” was a 10 mm thick hollow sphere outside of the PTV, while the “PTV ring 2” was a 10 mm thick hollow sphere outside of the “PTV ring 1”. For both ring contours, the changes in the mean dose were measured.

To investigate the mean lung dose, the healthy part of the lung was separated from the CTV by “Ipsilateral Lung excluding (w/o) CTV” and “Both Lungs w/o CTV”.

According to the International Commission on Radiation Units and Measurements, the target underdose was defined as the subvolume, Vunder, of the PTV receiving less than 95% of the prescribed dose, thus Vunder = VPTV – V95 [18]. VPTV is the volume of the delineated PTV and V95 is the volume within the PTV receiving at least 95% of the prescribed dose.

Two supplemental metrics specific for changes in the spread-out Bragg peak (SOBP) distal dose falloff were applied. The distal falloff, z80-20, measured the water-equivalent distance that the depth-dose curve needed to drop from 80% to 20% of the prescribed dose behind the high-dose plateau of the SOBP. It was calculated along the direction of each pencil beam which reached at least 95% of the prescribed dose.

Results

Target dose

Table 1 provides a comprehensive overview of the QI changes for the target volumes and ring structures. For centrally located and widespread tumors (groups C and W), relative QI changes were less than 0.8%. The largest changes occurred for tumors fully surrounded by lung tissue (group L), where the maximum QI change was 1.6% in PTV ring 1. In all cases, the D50 value of the PTV decreased.

Table 1

Relative changes in quality indicators of the clinical target volume (CTV), planning target volume (PTV), the ring contours, the underdosed volume, Vunder, (homogeneous lung), and its relative increase, ΔrelVunder, as defined in Section 2.4. Negative values indicate a decrease by introducing degradation, i.e., ΔD = (D (hetero) − D(homo))/D(homo)·100%.

Group	Patient ID	CTV ΔD95 (%)	PTV ΔD95 (%)	PTV ΔD50 (%)	PTV ring ΔDmean (%)		Vunder (%)	ΔrelVunder (%)
					ring 1	ring 2
L	S01	−1.1	−1.3	−1.0	−1.1	−0.7	10.8	20.4
L	S02	−0.7	−1.4	−0.7	−0.8	0.2	15.4	25.3
L	S06	0.3	0.1	−1.5	−1.6	−0.5	20.2	26.2
C	S07_1	0.5	0.3	−0.2	−0.4	−0.4	56.8	1.2
C	S07_2	0.2	0.4	−0.1	−0.5	−0.4	9.3	−7.1
C	S08	−0.1	0.4	−0.1	−0.1	0.0	36.8	1.1
C	S04	−0.1	−0.1	−0.1	−0.2	0.1	28.5	2.1
W	S03_2	−0.1	−0.7	−0.1	−0.3	0.3	6.2	15.1
W	S03_3	−0.2	−0.5	−0.1	−0.2	0.3	8.2	12.4
W	S05	−0.5	−0.8	0.0	−0.1	0.2	15.9	3.8

An overview of the underdosed volume Vunder within the PTV in the case of homogeneous lung and its relative increase, ΔrelVunder, are given in Table 1. The relative increase of the underdosed volume was smallest for group C (2% at the maximum); for one case, the underdose even decreased as cold spots of the inhomogeneous dose distribution were filled up during degradation. For group L, i.e., tumors fully surrounded by lung tissue and PTVs being comparably small (see Table B1 in Supplementary Material B), the underdosed volume increased by 20–26%.

The most pronounced dose changes occurred in patients S06 and S01. An exemplary transversal CT slice of patient S06 with the dose and dose difference overlay is depicted in Fig. 2. It illustrates the broader distal dose falloff along the beam directions with heterogeneity correction turned on. Local dose differences behind the target volume were up to 15% of the prescribed dose.

Fig. 2

Two-dimensional dose distribution of patient S06 for homogeneous lung tissue (a) and corresponding dose difference distribution by introducing heterogeneity correction (b). In (b), the blue areas indicate volumes that received lower doses with heterogeneity correction than without. The prescribed dose was 11 Gy (RBE) with the use of three coplanar treatment field orientations (gantry at 40°, 300°, and 340°; couch at 0°).

Distal falloffs

The distal falloff values, z80-20, as calculated along each pencil beam and their respective median value per treatment field were compared between individual treatment fields among the cases (median values in Table B2 in Supplementary Material B). The median falloff values were between 7 and 64 mm in groups C and W, and between 4 and 9 mm in group L for both the homogeneous and heterogeneous lung.

By applying heterogeneity correction, the median distal falloff did not increase by more than 0.4 mm except for one pencil beam, which increased by 3.8 mm.

For each pencil beam, the absolute change in the distal falloff was associated to the pristine falloff, i.e., to the z80-20 value for the homogeneous lung plan. Small z80-20 values, i.e., steep pristine falloffs, demonstrated a larger increase than larger z80-20 values. This relationship was observed for all proton beams that traveled through the lung tissue. For cases with very large pristine falloff values, e.g. due to opposing beams, the beam paths through the lung tissue did not cause a substantial additional widening of the already wide falloffs.

Dose-volume histograms

Fig. 3 shows the DVHs for cases S08 and S06. S08 is representative for group C with no differences between recalculations with the homogeneous and heterogeneous lung. For case S06, the hot spots of almost 13 Gy (RBE) within the target volumes shrunk upon heterogeneity correction. Changes in D50 values were consistent with previous observations. Within the ring contours, the dose was redistributed as the high-dose volume slightly decreased, while the low-dose volume slightly increased. The dose within the lung contours did not change in any DVH.

Fig. 3

Dose-volume histograms of patients S08 (a) and S06 (b, c). Histograms for homogeneous lung tissues are indicated by solid lines and histograms for heterogeneous lung tissues by dotted lines. (c) Hot spots in target dose decreased upon heterogeneity correction.

Mean lung dose

The mean lung dose within the ipsilateral lung excluding the CTV decreased (up to −1.2%) for all cases. The degraded dose falloff pushed the dose partially from the lung into the tissues behind the lung. The mean dose to both lungs excluding the CTV remained constant or decreased in all but one case (increase by 0.1%). For this case, S03_3, the broader dose falloff of a beam coming from 90° led to a higher dose in the contralateral lung, which outweighed the mean-dose decrease in the ipsilateral lung.

Worst case assessment

The dose distribution in the worst case assessment (Pmod,max = 750 µm) showed similar but larger effects than for Pmod,rep = 256 µm. Hence, we focused on patient S06 where the degradation had the strongest effect on the dose distribution. Table 2 displays the corresponding QI changes. D50 decreased by 4.4% by introducing heterogeneity correction, while maximum dose differences remained comparable to those for Pmod,rep = 256 µm. The underdosed target volume increased by 119%. The dose to the lung did not change by applying degradation with the maximum modulation power.

Table 2

Relative changes of quality indicators (QIs) for patient case S06 by introducing degradation with maximum modulation power (Pmod,max = 750 μm) and representative modulation power (Pmod,rep = 256 μm). For definitions of volumes and QIs, see Section 2.4.

Volume	Single voxel	PTV	PTV	PTV	PTV ring 1	PTV ring 2
QI	Max. dose difference	ΔD50	Vunder heterogeneous	ΔrelVunder	ΔDmean	ΔDmean
Relative change for Pmod,max	−27%	−4.4%	44%	119%	−4.9%	−2.2%
Relative change for Pmod,rep	−23%	−1.5%	25%	26%	−1.6%	−0.5%

Discussion

This is the first comprehensive study investigating the clinical relevance of dose degradation effects in intensity-modulated proton therapy for lung cancer. Our study comprised ten patient cases with large to small tumors located centrally within the thorax or more toward its periphery. In agreement with previous work [10], [11], [15], the overall effects were moderate, especially for larger tumors where volume averaging effects within the target may compensate negative effects regarding local underdose that may be observed for smaller tumors.

For the ten cases, the dosimetric changes within organs at risk were negligible. As apparent from the PTV ring volume and lung dose analysis, there was no evidence supporting the previously debated overdose to organs at risk close to the tumor [9], [11], [12], [19]. The distal falloffs in our study were in the order of centimeters. Consequently, a median increase of Δz80-20 < 0.5 mm can be neglected, particularly in consideration of a dose grid resolution of (2 mm)3.

Phantom simulations, as in Fig. 1, with varying thicknesses of lung tissue demonstrated a noticeable impact of the actual spot modulation resulting from treatment plan optimization on the distal dose falloff. The large impact of optimization of the spot weights and positions on the falloff width is shown in Fig. 4. The position of the last spot of the spread-out Bragg peak relative to the target boundary, which was dependent on the geometrical distance, zgeo, of penetrated lung tissue in this experiment, caused variations in the z80-20 value of ± 1 mm. (The geometrical distance zgeo refers to the spatial distance that the protons travel through lung as opposed to the water-equivalent distance which accounts for the smaller density of lung tissue.) A high weight on the highest-energy Bragg peak led to a steep pristine distal dose falloff, as seen in Fig. 4b. In contrast, more similar weights on all single Bragg peaks led to a wider pristine dose falloff, as seen in Fig. 4c. Therefore, the median distal falloff increase of less than 0.5 mm when introducing heterogeneity correction is of minor relevance.

Fig. 4

(a) Distal falloff values, z80-20, of heterogeneous and homogeneous lung tissue, and their difference, Δz80-20. The phantom setup was a proton beam impinging on a 30 mm-thick water-equivalent chest wall, a lung-phantom slab (relative electron density 0.30 and variable geometrical distance, zgeo), and a water-equivalent cubic target with an edge length of 80 mm. The plans were optimized for homogeneous lung tissue and recalculated with heterogeneity correction on. The mean and standard deviation values of more than 1500 rays are shown. The connection lines are for visual guidance only. (b, c) Different optimizations of the distal end of a spread-out Bragg peak (SOPB). Colored solid lines represent the single Bragg peaks composing the integrated depth dose of the SOBP. A high weight on the highest-energy Bragg peak leads to a steep pristine dose falloff (b), more similar weights on all single Bragg peaks lead to a wider pristine dose falloff (c).

For some cases, the mean dose to the “PTV ring 2” increased by up to 0.3%. Nevertheless, there were no organs at risk other than healthy lung located directly behind the target volume in the treatment plans.

We observed a clinically meaningful underdose of the target volume only in patient geometries with small target volumes which can only be reached by penetrating long distances of lung tissue relative to the chest wall.

An increase in the D95 value, as observed in four cases (Table 1), is caused by filling up cold spots within the target when the dose distribution is smeared through degradation. A decrease in the D95 value can occur when an initially homogeneous dose is degraded and the 95%-isodose is pushed beyond the target boundary.

All simulations were acquired using a pencil-beam dose calculation algorithm, which is still the dominating standard for particle therapy treatment planning. As such, the degradation approach would be applicable not only for intensity-modulated proton therapy but also for passive scattering systems. The general limitations of pencil-beam algorithms within heterogeneous lung tissue, however, are well known and clinical treatment planning for lung cancer should be based on a Monte Carlo dose calculation engine [20], [21]. However, in the present study the specific algorithm enabled an isolated investigation of the degradation effects. In addition, we focused on relative effects, rather than the absolute dose values. Therefore, the pencil-beam model does not compromise the overall validity of the conclusions. It may rather be considered as a worst case estimate as pencil-beam dose calculation algorithms generally predict steeper dose gradients within the lung tissue than Monte Carlo dose calculations.

Robust treatment planning could become an important alternative to conventionally optimized plans for proton beam therapy. Robust plans have generally less in-field modulation to compensate for potential mismatches between the geometry during planning and treatment. They are thus more resilient for lung patient treatments [22], which potentially also extends towards degradation.

For a comprehensive clinical assessment, an isolated investigation of the degradation effects may only provide partial guidance. Also, other sources of uncertainty, particularly intrafractional motion [5], [7], [23], interfractional anatomical changes [24], and setup uncertainties [25] need to be considered, as more severe dosimetric changes are suggested in literature than observed in the context of our study. Many proton centers do not have access to high-quality 3D imaging for patient positioning during treatment which further increases the setup uncertainties for lung tumors [26].

In conclusion, this work presents a comprehensive study of the degradation effects on clinical lung cancer proton treatment plans. For larger and centrally located tumors, the degradation effects were moderate. For small tumors surrounded by lung tissue, however, problems regarding local target underdose might occur. Due to pronounced anatomical variations within lung cancer, our study does not allow for a general recommendation regarding the consideration of degradation effects for clinical practice. Instead, we advocate further research in consideration of other prevalent sources of uncertainty in lung cancer, such as breathing motion or patient setup – especially for small tumor volumes that can only be reached by penetrating substantial amounts of heterogeneous lung tissue. Overall, our study provides promising results as related work has shown that intensity-modulated proton therapy is particularly favorable for large, central, and recurrent tumors, where we see less of a problem for safe clinical application of intensity-modulated proton therapy due to degradation effects.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.