Dose reconstruction technique using non-rigid registration to evaluate spatial correspondence between high-dose region and late radiation toxicity: a case of tracheobronchial stenosis after external beam radiotherapy combined with endotracheal brachytherapy for tracheal cancer

Purpose

For early-stage tracheal cancer, high-dose-rate endobronchial brachytherapy (HDR-EBT) can be a treatment choice as a boost for external beam radiotherapy (EBRT) or as a definitive therapy [1,2,3]. High-dose-rate endobronchial brachytherapy is also used successfully as a palliative treatment for symptomatic airway obstruction by malignant tumors [4,5]. However, in some cases with longer follow-up, HDR-EBT leads to severe complications, which are sometimes lethal, such as radiation-induced bronchitis, bronchial stenosis, and fetal hemoptysis [6,7]. Although these complications might reflect the high-dose region resulted from the small distance between the applicator and tracheobronchial mucosa [8], further investigation is needed on how high-dose volume topographically correlates with the location of late radiation injury. However, there exists technical challenge because the location of small high-dose volumes in fractionated radiotherapy can be affected by intra- and inter-fractional organ motion, and the anatomical features of an organ post-treatment will be different from those pre-treatment due to the sequelae of radiation treatment [9].
For the purpose of spatiotemporal dose assessment of HDR-EBT in combination with EBRT, one important technique is non-rigid image registration, which enables the estimation of dose summation based on the corresponding structures across the image data sets. Previously, we have demonstrated an image-processing framework using surface-based, non-rigid registration for assessing the spatially cumulative dose in consecutive brachytherapy treatments [10]. The core of our approach is based on the point set registration using a Gaussian mixture model (GMM-REG), which was originally introduced by Jian and Vemuri [11]. Because GMM-REG demonstrated the biomechanically natural registration with high accuracy and validity [12], we considered that it can be a plausible technique for transforming planned dose distribution to the post-treatment organ.
Here, we present a case of tracheobronchial radiation stenosis after the definitive therapy using HDR-EBT and EBRT for early-stage tracheal cancer. By using the proposed dose reconstruction technique, the spatial association between the dosimetric profile and the site of the radiation morbidity with a severe tracheobronchial stenosis was evaluated.

Material and methods

Patient data

A 78-year-old man developed a second primary tracheal cancer 5 years after a right lower lobectomy for early stage lung carcinoma. Bronchoscopy revealed a tumor localized at the left lateral wall of the middle part of the trachea and histology confirmed a squamous cell carcinoma. The patient was then treated with HDR-EBT (24 Gy in 6 Gy fractions) followed by EBRT (40 Gy in 2 Gy fractions) with a curative intent. Follow-up bronchoscopies revealed a patent airway with granulomatous mucositis for 1 year after the treatment, with no evidence of biopsy-proven malignancy. Twenty-four months after the treatment, at the age of 80, the patient was admitted to our department because of progressive dyspnea. Computed tomography (CT) detected a tracheobronchial stenosis with aspiration pneumonia of the right lung. Massive pneumonia infiltrated the right lung completely due to the stenosis of the right primary bronchus. Bronchoscopy revealed that the surface of the tracheal mucosa was irregular with inflammation, and bled easily on contact. It also revealed a small ulcer above the circumferential tracheal stenosis. Because the tissue examination of the tracheal mucosa was consistently negative for malignancy before the admission, we diagnosed late radiation-induced mucositis with severe circumferential stenosis.

Treatment protocol

Our HDR-EBT treatment protocol has been described in detail elsewhere [6]. In short, a 12 Fr Malecot-type catheter with two wings was placed from the middle part of the trachea to the right primary bronchi during bronchoscopy under local anesthesia. Based on the CT images with 2 mm slice thickness, treatment planning was done on Oncentra (Nucletron, an Elekta company, Elekta AB, Stockholm, Sweden). A dose of 6 Gy was prescribed at 1 cm from the catheter. High-dose-rate endobronchial brachytherapy was performed using the Microselectron 192Ir-HDR (Nucletron). Subsequently, EBRT was delivered using a three-dimensional conformal technique with linear accelerator (Clinac iX, Varian Medical System, Palo Alto, CA, USA). Treatment planning was based on the CT images of 3 mm slice thickness taken with a CT simulator (AquilionTM LB, Toshiba Medical Systems, Tokyo, Japan). The EBRT fields included the primary tumor with adequate margins of 2 cm.

Image-processing framework

The analysis was performed using software developed in-house, which was written in Python using VTK/ITK library and published modules of GMM-REG [11]. A single observer (KK) contoured the inner surface of the tracheal bronchus on every CT image at the treatment planning (the 1st to 4th fraction of HDR-EBT and EBRT) and at post-treatment evaluation for the radiation stenosis. Delineation, which was bounded by the glottis superiorly and the level of second carina inferiorly, was performed using treatment planning software (Oncentra). Then, the contours were exported as binary masks from digital imaging and communication in medicine, and these binary masks were smoothened via Gaussian filter with a sigma value of 2.0 for reducing artifacts. A marching cubes filter was applied to produce surface mesh, which was represented by points that were connected with each other to form triangles. For refining the triangles which defined the surface, a quadric clustering filter with the division spacing of 2.0 mm was used as the decimate modifier to reduce the number of vertices of the mesh.
Next, non-rigid registration based on GMM-REG was guided by the point set distributed on the surface mesh. The theoretical background has been described in detail in the original article written by Jian and Vemuri [7]. Briefly, a Gaussian mixture is generally defined as

p(x) = ki = 1i(x|µi, i),

where ωi is the weight factor associated with the Gaussian component density (x|µi, i), µi is the mean vector and Σi is the covariance matrix. The scale factor σ represents the size of the spherical covariance matrix of the model, which effects the magnitude of transformation.
GMM-REG treated the problem of point set registration as that of aligning two Gaussian mixtures by minimizing the L2 distance between them. Thin-plate spline (TPS) function was selected for the non-rigid registration, which added a penalty term to the final cost function to control the strength of regularization. Therefore, the registration method between the model set M and the scene set S finds the warping function u of the TPS function, which minimizes the following cost function in the presence of the weight parameter λ for the bending energy:

E = ∫ {gmmσ(S) – gmmσ(u(M))}2dx + λBending (u)

The process was driven by a numerical optimization algorithm called the annealing step to minimize local minima. We defined the total number of steps for annealing as five. The transformation parameters of each step was as follows:

(σ1, σ2, σ3, σ4, σ5) = (0.6, 0.3, 0.1, 0.05, 0.02),

(λ1, λ2, λ3, λ4, λ5) = (1.0, 1.0 × 10–2,1.0 × 10–4, 1.0 × 10–6, 0.0).

Finally, the transformation function computed vectors that connected the points on one tracheobronchial surface to another tracheobronchial surface, which enabled the estimation of cumulative dose on the post-treatment tracheobronchial stenosis. For expressing the equivalent dose in a fractional dose of 2 Gy (EQD2) in HDR-EBT, / of a linear quadratic model was assumed as 3 Gy for tracheobronchial mucosa [13].

Validation for transformation

To measure the geometric accuracy of the transformation, surface distance error (SDE) and surface coverage error (SCE) were used [14]. Surface distance error was defined as the mean distance between the transformed points and the fixed surface. Surface coverage error represented the coverage of the fixed surface by the transformed points and was derived by calculating the standard deviation (SD) of the distances from each vertex of the fixed structure to the nearest transformed point. Transformation with small SDE and small SCE indicated a high geometric accuracy.

Homogeneity index

For accessing the homogeneity of the dose distribution, we calculated the homogeneity index (HI) as follows:

HI = Dmax/Dp

where Dmax is the maximum point dose to the inner surface of the tracheal bronchus and Dp is the prescribed dose.

Quantification of tracheobronchial stenosis

To quantify the severity of the stenosis, the luminal area of the tracheobronchial stenosis was evaluated as an intersection between the tracheobronchial surface and the plane, which was orthogonal to the principal axis of the organ. The principal axis was obtained using principal component analysis for the points on the structure surface. Each corresponding plane was computed based on the deformation fields for each tracheal bronchus at the time of treatment planning (the 1st to 4th fraction of HDR-EBT and EBRT), and the mean area of intersection was defined as the pre-treatment luminal area (Figure 1A). The local dose at each plane was calculated to be the average dose irradiated to the points constituting the outer line of intersection of the post-treatment tracheobronchial structure. The relationship between the local dose and the ratio of luminal areas was calculated from the level of the plane 100 mm above the first carina (inferiorly) to the plane just below the glottis (superiorly). Spearman’s rank correlation coefficient was calculated for the association between the irradiated dose and the ratio of luminal areas. Statistical significance was set at a two-sided p-value of less than 0.05. Data are presented as the mean ± standard deviation unless otherwise specified.

Results

Transformation results

An example of the deforming model point set to the scene point set on the basis of GMM-REG was shown in Figure 2. When the structure of post-treatment tracheobronchial stenosis was transformed to each tracheobronchial surface at the treatment planning (the 1st to 4th fraction of HDR-EBT and EBRT), SDE and SCE was 0.18 ± 0.01 mm and 0.39 ± 0.01 mm, respectively.

Mapping reconstructed dose distribution

We reviewed the planned dose distribution on the organ surface at the time of treatment planning. As shown in Figure 3, the surface dose of EBRT was homogenously distributed (HI = 1.1); on the other hand, that of HDR-EBT varied unevenly from fraction to fraction (HI = 3.7 ± 0.2). Unexpected hot spots irradiated by more than 20 Gy per fraction resulted from the small distance between the catheter and the tracheobronchial inner surface (Figure 4), as the minimum distance between them was 2.4 ± 1.1 mm for the four fractions of HDR-EBT. The dose-surface histogram of the HDR-EBT was calculated (Figure 5) and the area of the tracheobronchial mucosa irradiated by more than 20 Gy was 94.7 ± 82.7 mm2 (0.6 ± 0.5%).
According to the non-rigid transformation function, we mapped cumulative dose during the course of the radiation therapy on the surface of the pre-treatment (the structure at the time of the treatment planning of EBRT) (Figure 6) and post-treatment tracheal bronchus (Figure 7). For each fraction of HDR-EBT (24 Gy in 6 Gy fractions) and the total EBRT dose (40 Gy in 2 Gy fractions), the fractional dose was normalized to EQD2 with / of 3. On the pre-treatment tracheal bronchus, the area with cumulative dose more than 100 Gy3 and 200 Gy3 was 23.3 × 102 mm2 (18.6%) and 5.7 × 102 cm2 (4.5%), respectively. Figure 7 shows that the subvolume of the tracheobronchial surface irradiated by more than 100 Gy3 developed stenosis, and the severity of the stenosis apparently corresponded with the cumulative local dose. When quantitating the structural deformation of radiation stenosis by the ratio of luminal areas, the subvolumes irradiated by more than 100 Gy3 showed severe stenosis with a more than 50% decrease in the luminal area (Figure 1B). The association between the local dose and the ratio of luminal area was statistically significant (p < 0.0001). In particular, small volumes irradiated by more than 200 Gy3 correlated with bronchoscopically-proven ulceration and the circumferential stenosis leading to a severe decrease in the luminal diameter (Figure 8).

Discussion

In our study, we have presented a novel dose reconstruction technique to investigate the topographic correspondence between high-dose volume and the site of late radiation injury represented by tracheobronchial stenosis. Although the post-treatment tracheobronchial surface showed a large anatomical deformation, the measured low-distance error implies that the proposed approach with GMM-REG successfully registered the structures. As a result, the estimated cumulative local dose on the tracheobronchial surface significantly corresponded with the severity of radiation stenosis. The present study is of importance because the surface-based, non-rigid registration can be a core of the technique to predict the location-by-location difference in the severity of late radiation injury.
While there is variation in dose-fractionation regimens and treatment intentions among previous studies, the crude incidence of symptomatic radiation-induced bronchitis is 4-12% in HDR-EBT with EBRT [7,15,16,17,18]. Speiser and Spratling proposed a grading system to classify radiation-induced bronchitis, which ranges from mild mucosal inflammatory response with swelling (Grade 1) to major fibrosis with circumferential stenosis (Grade 4) [7]. According to the grading scale, our presented case represented the most severe phase of Grade 4 radiation tracheobronchitis in small regions, which were cumulatively irradiated by more than 200 Gy3. Severe inflammatory response with mild fibrosis of the tracheobronchial wall, which corresponds to Grade 3, was apparently associated with mucosa irradiated by more than 100 Gy3.
It is necessary for radiation oncologists to avoid excessive irradiation to the tracheobronchial mucosa because once radiation-induced tracheobronchitis becomes progressive, it never heals spontaneously, and often develops into fatal hemoptysis. For the present case, we prescribed a dose of 6 Gy at 1 cm from the catheter without enough spatial consideration of the distance from the tracheobronchial surface; thus, this may have been a primary cause of incidental high-dose volume.
There have been few reports on the spatial correspondence between high-dose volume and the site of radiation injury. By using rectosigmoidoscopy and rectal morbidity assessment, Georg et al. reported the locational correlation of dosimetric parameters of gynecological brachytherapy with clinical changes [19]. Although they emphasized the clinical relevance of DVH parameters, such as D2cc, D1cc, and D0.1cc, it might be an alternative approach to focus on the surface dose for reporting when non-rigid registration is used, because non-rigid registration sometimes produces artifacts by condensing or stretching volumes [20]. The present study also implies that the surface dose can be used in the prediction of radiation morbidity in the walls of a hollow organ.
Although the present study proposed a non-rigid registration-based assessment for the dose reconstruction of a previously irradiated dose to the post-treatment organ with large deformation, the technique has so far been applied in only one case of tracheal cancer treated with HDR-EBT and EBRT. General validation for other treatment sites with other organs at risk, different inter-fraction variations, and more inhomogeneous and complex dose distributions is still a big challenge.

Conclusions

The image-processing framework based on GMM-REG enabled to correlate the cumulative surface dose of HDR-EBT and EBRT with the location of radiation injury, by applying it to one case of tracheal cancer. This technique seems to be a promising tool that might be applied to other treatment sites in the future.

Acknowledgements

This work has been financially supported by Grants-in-Aid for Young Scientists (B) from the Ministry of Education, Culture, Sports, Science, and Technology of Japan (Grant number 15K19837), Health and Labor Sciences Research Grants for Commission (H26-applied-general-036), and the Practical Research for Innovative Cancer Control from Japan Agency for Medical Research and development, AMED.

Disclosure

Authors report no conflict of interest.

References

1. Furuta M, Tsukiyama I, Ohno T et al. Radiation therapy for roentogenographically occult lung cancer by external beam irradiation and endobronchial high dose rate brachytherapy. Lung Cancer 1999; 25: 183-189.
2. Marsiglia H, Baldeyrou P, Lartigau E et al. High-dose-rate brachytherapy as sole modality for early-stage endobronchial carcinoma. Int J Radiat Oncol Biol Phys 2000; 47: 665-672.
3. Kawamura H, Ebara T, Katoh H et al. Long-term results of curative intraluminal high dose rate brachytherapy for endobronchial carcinoma. Radiat Oncol 2012; 7: 112.
4. Barber P, Stout R. High dose rate endobronchial brachytherapy for the treatment of lung cancer: current status and indications. Thorax 1996; 51: 345-347.
5. Delclos ME, Komaki R, Morice RC et al. Endobronchial brachytherapy with high-dose-rate remote afterloading for recurrent endobronchial lesions. Radiology 1996; 201: 279-282.
6. Hara R, Itami J, Aruga T et al. Risk factors for massive hemoptysis after endobronchial brachytherapy in patients with tracheobronchial malignancies. Cancer 2001; 92: 2623-2627.
7. Speiser BL, Spratling L. Radiation bronchitis and stenosis secondary to high dose rate endobronchial irradiation. Int J Radiat Oncol Biol Phys 1993; 25: 589-597.
8. Carvalho Hde A, Gonçalves SL, Pedreira Jr. W et al. Irradiated volume and the risk of fatal hemoptysis in patients submitted to high dose-rate endobronchial brachytherapy. Lung Cancer 2007; 55: 319-327.
9. Tanderup K, Nesvacil N, Pötter R et al. Uncertainties in image guided adaptive cervix cancer brachytherapy: impact on planning and prescription. Radiother Oncol 2013; 107: 1-5.
10. Kobayashi K, Murakami N, Wakita A et al. Dosimetric variations due to interfraction organ deformation in cervical cancer brachytherapy. Radiother Oncol 2015; 117: 555-558.
11. Jian B, Vemuri BC. Robust Point Set Registration Using Gaussian Mixture Models. IEEE Trans Pattern Anal Mach Intell 2011; 33: 1633-1645.
12. Jian B, Vemuri BC. A Robust Algorithm for Point Set Registration Using Mixture of Gaussians. Proc IEEE Int Conf Comput Vis 2005; 2: 1246-1251.
13. Dale RG. The application of the linear-quadratic dose-effect equation to fractionated and protracted radiotherapy. Br J Radiol 1985; 58: 515-528.
14. Vásquez Osorio EM, Hoogeman MS, Bondar L et al. A novel flexible framework with automatic feature correspondence optimization for nonrigid registration in radiotherapy. Med Phys 2009; 36: 2848-2859.
15. Speiser BL, Spratling L. Remote afterloading brachytherapy for the local control of endobronchial carcinoma. Int J Radiat Oncol Biol Phys 1993; 25: 579-587.
16. Gauwitz M, Ellerbroek N, Komaki R et al. High dose endobronchial irradiation in recurrent bronchogenic carcinoma. Int J Radiat Oncol Biol Phys 1992; 23: 397-400.
17. Zajac AJ, Kohn ML, Heiser D et al. High-dose-rate intraluminal brachytherapy in the treatment of endobronchial malignancy. Work in progress. Radiology 1993; 187: 571-575.
18. Ozkok S, Karakoyun-Celik O, Goksel T et al. High dose rate endobronchial brachytherapy in the management of lung cancer: response and toxicity evaluation in 158 patients. Lung Cancer 2008; 62: 326-333.
19. Georg P, Kirisits C, Goldner G et al. Correlation of dose-volume parameters, endoscopic and clinical rectal side effects in cervix cancer patients treated with definitive radiotherapy including MRI-based brachytherapy. Radiother Oncol 2009; 91: 173-180.
20. Andersen ES, Noe KØ, Sørensen TS et al. Simple DVH parameter addition as compared to deformable registration for bladder dose accumulation in cervix cancer brachytherapy. Radiother Oncol 2013; 107: 52-57.