Monte Carlo modelling of tumour growth, hypoxia and radiotherapy in head and neck squamous cell carcinoma.

Abstract

Tumour hypoxia is the inadequate supply of oxygen in living tissue. Hypoxia is a major problem in the treatment cancer with ionising radiation because of the associated increase in radioresistance of hypoxic tumour cells. This effect can cause up to a three fold increase in the radiation dose required to kill the hypoxic cells compared to well oxygenated cells. Many locally advanced head and neck tumours exhibit hypoxia to some degree, and there is direct evidence that hypoxic tumour sub-volumes and their associated mean oxygenation levels have a direct influence on local tumour control after radiotherapy (Nordsmark 2005). Currently, head and neck cancer radiotherapy local control rates lie at approximately 80% for early stage disease, but reduce significantly (often below 50%) for locally advanced tumours. Efforts to improve these statistics through dose and fractionation modifications in randomised clinical trials have been made in recent decades using alternate fractionation schedules, but the average prognosis has not improved significantly. The effects of tumour reoxygenation during fractionated radiotherapy can assist in re-sensitising previously hypoxic tissue; however the complex dynamics and patient dependent characteristics of this phenomenon make the benefits difficult to quantify. Head and neck cancers, specifically head and neck squamous cell carcinoma (HNSCC), have also been shown to experience the phenomenon of accelerated repopulation during fractionated radiotherapy. Accelerated repopulation enhances cellular proliferation as a response to the trauma caused by treatment, and contributes to the low HNSCC local control rates after radiotherapy. The modelling work developed for this report was undertaken to better understand the mechanisms and quantitative effects of HNSCC cellular kinetics and tumour oxygenation during growth and radiotherapy. The goal of individualising treatment planning for this disease was the motivation for developing the model. A key aim was to produce an end product to be used as an efficient and user-friendly radiobiological tool for the input on tumour specific properties such as tumour oxygenation and reoxygenation onset time, to investigate their effects on cell kill during radiotherapy. To this end, a Monte Carlo model, named HYP-RT (for HYPoxic-RadioTherapy simulation), was developed. HYP-RT simulates the tumour cell division process according to epithelial proliferative hierarchy, starting from a single stem cell. Monte Carlo methods were used to simulate the probabilistic nature of the biological and radiobiological mechanisms and parameters incorporated into the model, e.g. the distribution of cell cycle times (normal or exponential) and oxygenation levels (normal or log-normal), and the randomised methods of cell kill and oxygenation increase during treatment. Probabilistic methods were also used to make decisions during cell division, as to the type of daughter cell products that would emerge after the division of a mother cell. After the growth of a 10⁸ cell tumour, an algorithm was developed to model the effects of fractionated radiotherapy. This algorithm was designed to simulate the oxygen dependent radiosensitivity of individual tumour cells, as well as the effects of gradual reoxygenation and accelerated repopulation (through loss of stem cell division asymmetry). Both reoxygenation and accelerated repopulation could be onset at varying times after the start of treatment. Experimental animal work using HNSCC (FaDu cell line) xenografts was undertaken during this research, and showed that reoxygenation occurred very late in an accelerated radiation schedule (40 Gy in 2 weeks), indicating the need to investigate a range of reoxygenation onset times in the model (0 to 3 weeks). Dynamic cell data was stored in a pre-allocated vector (the Cellarray) containing just over 10⁸ object elements, with each element representing one tumour cell. This enabled efficient random access to the data. Linked list methods were used to chronologically order cells in the Cellarray based on their times of division. Model efficiency was paramount during model development, to ensure convenient use of the model for the current work and potential future research. Using linked list methods, the goal of a one hour maximum computation time to grow and treat a tumour was successfully achieved. The model source code was written with the FORTRAN 95 programming language (complier v7.1.0, Lahey Computer Systems Inc.), within the Visual Studio (2003, Microsoft Corporation) framework. Two additional graphical user interface programs were developed using the JAVAprogramming language (Java SE Development Kit 6.17), to 1) read in and interpret tumour data files, and 2) allow for the input of key tumour parameters before a simulation (or batches of simulations) and iteration over multiple parameter sets. Cellular data and key algorithm parameters were written out to file at regular intervals (1000 hours by default), during tumour growth and before and after every dose fraction during treatment, for retrospective analysis. This data included the tumour pO₂ distribution, the instantaneous tumour growth rate and the number of cells of various types comprising the tumour. Simulation results showed that tumour growth rate was strongly dependent on the percentage of stem cells in the tumour (modelled to be approximately 1% during growth). Incorporating a “moderately” hypoxic oxygen distribution increased tumour doubling times significantly, from 37 days for oxic tumours up to 65 days for moderately hypoxic tumours. This was attributed to the effects of oxygen dependent cell cycle slowing, cellular quiescence and necrosis. Simulated conventional radiotherapy (5x2 Gy/wk) required on average an extra 16 Gy in total to achieve tumour control for moderately hypoxic compared to well oxygenated tumours. The effects of both accelerated repopulation and reoxygenation significantly altered the total doses required for tumour control, with accelerated repopulation effects dominating model outcomes. Accelerated repopulation and reoxygenation were found to be dependent on one another, making simulations of every combination of onset time for each effect necessary during model analysis. During accelerated repopulation, a dose per fraction of 2.5 to 3.0 Gy was required to control the extra cell growth in an otherwise 2 Gy per fraction schedule. This equated to an extra 5 Gy being requiring to maintain tumour control for every week that the onset accelerated repopulation was brought toward the start of treatment. The benefits of reoxygenation reduced as the time of onset was delayed, with +1 Gy required to maintain tumour control for every week that reoxygenation was delayed. Conventional fractionation simulation results had good agreement with standard Linear Quadratic theory, for the dose required to control well oxygenated tumours. However, comparison results were mixed for more complex cases involving hypoxic tumours with and without accelerated repopulation. When modelling altered fractionation schedules, simulation outcomes in terms of the total doses required for tumour control, agreed well with the prescriptions from published clinical trials. The most beneficial schedule, based on predicted total dose as well as biological effective doses (BED’s) calculations for normal tissues, was the 10x1.1 Gy/week schedule (Pinto et al. 1991). However, there were up to 30 Gy differences in total dose and BED results when simulating specific sets of tumour parameters for the same radiation schedule, highlighting the need for individualisation of treatment planning to improve the therapeutic ratio. Four newly designed altered schedules were also simulated with the HYP-RT model. Results showed that using a concomitant boost at the beginning, rather than at the end of treatment, or using a “less aggressive” continuous hyper-accelerated radiotherapy (CHART) schedule (compared to the UK CHART schedule) may have potential therapeutic benefits compared to existing clinical schedules. Altering the oxygen enhancement ratio (OER) curve based on dose per fraction for the altered fractionation schedules, changed model results significantly for hyperfractionated schedules (up to 20 Gy). This highlighted the critical nature of the OER curve in predictive radiobiological tumour models. In summary, the current research has involved the development, analysis and use of an efficient Monte Carlo tumour growth and radiotherapy model (HYP-RT). The model simulates a biologically plausible epithelial cell hierarchy, a large number of individual cells, tumour hypoxia, and the dynamics of reoxygenation and accelerated repopulation during radiotherapy. The user can input the desired oxygen distribution to describe the degree of tumour hypoxia as well as and the times of onset of treatment related effects, among many other cellular parameters. The model provides quantitative results regarding the total dose required to control a tumour, for a given fractionation schedule and tumour parameter set. It is hoped that computer models such as HYPRT will be used in the near future as a tool to aid in the individualisation of radiotherapy planning, based on specific tumour experimental/imaging information, to improve prognosis for patients with HNSCC.