Monte Carlo and analytic simulations in nanoparticle-enhanced radiation therapy

Abstract

Analytical and Monte Carlo simulations have been used to predict dose enhancement factors in nanoparticle-enhanced X-ray radiation therapy. Both simulations predict an increase in dose enhancement in the presence of nanoparticles, but the two methods predict different levels of enhancement over the studied energy, nanoparticle materials, and concentration regime for several reasons. The Monte Carlo simulation calculates energy deposited by electrons and photons, while the analytical one only calculates energy deposited by source photons and photoelectrons; the Monte Carlo simulation accounts for electron–hole recombination, while the analytical one does not; and the Monte Carlo simulation randomly samples photon or electron path and accounts for particle interactions, while the analytical simulation assumes a linear trajectory. This study demonstrates that the Monte Carlo simulation will be a better choice to evaluate dose enhancement with nanoparticles in radiation therapy.

Keywords:

• nanoparticle
• dose enhancement
• Monte Carlo simulation
• analytical simulation
• radiation therapy
• tumor cell
• X-ray

Introduction

Radiation therapy is one of the most commonly used methods for cancer treatments. A challenge for X-ray radiation therapy is that the therapeutic doses used can damage neighboring normal cells. Other than designing better beam techniques, a variety of radiosensitizers (including chemicals, ions, and nanoparticles) have been used to enhance the efficacy of radiation therapy. Because nanoparticles can be functionalized to target ligands expressed at tumor cell surfaces, nanoparticles of heavy atoms have been proposed as promising radiosensitizers for radiation therapy. Owing to high X-ray attenuation coefficient and biocompatibility, gold nanoparticles have been intensively studied as promising candidates to enhance radiation doses.^1–3 Secondary electrons and photons released when a photon is incident on a metallic surface are the cause of sensitization. Several interactions can take place between a photon and a metallic surface based on the energy of the photon. These include Compton scattering and photoelectric effects accompanied by photoelectrons and Auger electrons. Compton scattering and photoelectric effects dominate when metallic nanoparticles are used as imaging contrast agents.^4,5 As photoelectrons and Auger electrons are released due to the photoelectric effect, they travel through the cell and hydrolyze water molecules creating free radicals that interact with DNA causing DNA damage, and eventually cell death.^6,7

Analytical and Monte Carlo methods have been used to calculate the production of secondary electrons in terms of the dose enhancement factor (DEF).^8,9 The DEF is defined as the ratio of the average radiation dose absorbed by the tumor when it is loaded with nanoparticles to the dose absorbed without the nanoparticles. Previous analytic calculations were carried out for a single nanoparticle and scaled to account for the total number of nanoparticles attached to a cell surface. The probability of photoelectric interactions, energy balances, and particle transport equations are analytically solved to determine the amount of energy deposited within the cell. The simulation assumes that the nanoparticles are evenly distributed on the cell surface, each nanoparticle has the same probability of interacting with a photon and producing a photoelectron, there are no interactions between secondary electrons and neighboring nanoparticles, and for every electron that leaves the cell, another one enters creating a lateral equilibrium.¹⁰ Analytical methods yield solutions that are exactly reproducible between different runs and allow flexibility in changing parameters for any given application. In addition, computations take much less time, milliseconds, compared to Monte Carlo methods, which can take up to hours per run based on geometry complexity. However, analytical techniques are mostly limited to one-dimensional (1-D) space and errors are not easy to estimate.

The second method used to calculate DEF is a Monte Carlo simulation. The Monte Carlo method treats the electron or photon trajectory as a sequence of random steps. Each step is presumed to end in an interaction, chosen through another random sampling. The interactions can cause a change in direction, energy loss, or creation of secondary particles.^11,12 All types of physical interactions between matter and radiation can be modeled using Monte Carlo methods for 1-D, 2-D, and 3-D applications. However, high accuracy and resolution require long calculation times. In all cases, Monte Carlo methods suffer from high statistical uncertainties particularly near material boundaries and interfaces. For deep penetration problems or for regions where too few particles contribute to the solution, variance reduction techniques can be used to improve accuracy, but even then the statistical error may be unacceptable.¹³ Several software packages are available for radiation transport in matter such as GEANT4, MCNP5, and EGSnrc. MCNP5 is a radiation transport code distributed by the Radiation Safety Information Computational Center at Oak Ridge National Laboratory. Written in Fortran 90, this transport code uses a condensed history model.^14,15 EGSnrc is a free open-source software package distributed by the National Research Council of Canada. Written in C++, it has an exact boundary crossing algorithm that is used in conjunction with the PRESTA-II transport code.¹⁴ GEANT4 is another open-source package written in C++ originally intended for high-energy physics.¹⁶ Because of the broad range of applications of GEANT4, the user must supply a physics list that outlines the details of the model.^14,17 EGSnrc was decided upon due to it being free open-source software, with its ease of use for user input files, customizable material files, geometry visualization, and documentation.

An issue of these activities is that the results obtained in each group using different methods are largely different. Since each group models a different parameter, it is difficult to compare the soundness of each simulation method. Given the importance of simulation in guiding experiments, it is imperative to find a reliable simulation method with minimum requirements on empirical data. This study intends to compare an analytical simulation based on a 2010 paper¹⁰ and a Monte Carlo simulation based on a 2005 paper.¹⁸ This study aims to provide, for the first time, a head-to-head comparison between analytical solutions vs Monte Carlo simulations for radiation dose enhancement. It will look into how DEF differs with source energy, nanoparticle materials, and concentration for both analytical and the Monte Carlo simulation.

Methodology

A cell is modeled as a slab of tissue with a single nanoparticle attached to its surface (Figure 1). An analytical and a Monte Carlo simulation were compared for various source energies, nanoparticle materials, and concentrations. Endothelial cell lines of tumor vasculature are more likely to be exposed to nanoparticles. For both simulations, an endothelial cell was modeled as a 10×10×2 µm slab composed of 10% hydrogen, 11% carbon, 2% nitrogen, and 76% oxygen according to the International Commission on Radiation Units and Measurements. X-ray energies from 10 to 900 keV were chosen to model the low-energy photon beam. Nanoparticle concentrations of 5, 10, 15, and 20 mg/g of tissue, evenly distributed through the cell, were chosen based on their usage in literature. Both simulations were performed for gold, bismuth, silver, and copper.

Figure 1 Analytical simulation scheme, where a single nanoparticle is at the center of a cell modeled as a slab of tissue.

Analytical simulation

The analytical simulation modeled a single nanoparticle with radius of r, at the center of a cell. Dose enhancement was calculated using the equations presented in the 2010 Ngwa et al paper and the 2012 Hossain and Su paper.^9,10 For these equations, a base dose, D_w =2 Gy, chosen based on its clinical relevance, and a monoenergetic beam were used to calculate the photon fluence (Φ):

$\mathrm{\Phi }=\frac{D_{\mathrm{w}}}{E{\left(\frac{{\mathrm{\mu }}_{\text{en}}}{\mathrm{\rho }}\right)}_E}$(1)

and the number of photons interacting with a single nanoparticle, N_ph = Φ*π*r², where (µ_en/ρ)_E is the mass absorption of coefficient of water at a given photon energy E, and r is the nanoparticle radius. On the basis of the probability of a photoelectric interaction, P ≈ (µ_pe/ρ)_E*ρ_NP*d_NP, where d_NP is the average distance a photoelectron travels through the nanoparticle, ρ_NP is the nanoparticle density, and (µ_pe/ρ)_E is the photoelectric absorption coefficient of the nanoparticle at energy E. The total number of photoelectrons produced per nanoparticle is calculated using the following equation:

$N_{\text{pe}}=N_{\text{PH}}*P$(2)

Based on the energy of the photoelectron, E_pe = E − E_edge, the maximum distance traveled through the cell is:

$Rt=0.0431{(E+0.367)}^{1.77}-0.007$(3)

The energy deposited within the cell, per photoelectron is calculated using the following equation:

$E_{\text{ec}}=\frac{H_{\text{abc}}-C_{\text{xby}}}{2*H_{\text{abc}}}{\displaystyle \underset{r}{\overset{Rt}{\int }}\left[3.316{(Rt-x+0.007)}^{-0.435}+0.0055{(Rt-x)}^{0.33}\right]}dx,$(4)

where H_abc = 2*π*R(0.5*t)^2, C_xby = 2*π*Rt*(Rt −0.5t), with t being the thickness of the cell. This equation assumes evenly distributed nanoparticles and equilibrium between secondary electrons leaving and entering the cell. This is then scaled by the total number of photoelectrons released, N_Pet = N_PE*N_NP, where N_NP is the total number of nanoparticles in the cell, to get the total energy deposited in the cell, E_ect = E_ec*N_pet. The DEF is then

$\text{DEF=}\frac{D_{\mathrm{w}}+D_{\text{NP}}}{D_{\mathrm{w}}},$(5)

where D_NP = E_ect/m_cell is the dose deposited from the photoelectrons.

Monte Carlo simulation

EGSnrc provided by the National Research Council of Canada was used as the basis for the Monte Carlo simulation. EGSnrc is an electron transport code that comes equipped with electron transport in matter, tracking for individual electrons and photons. Using Tutor7pp as a base, the code was modified to allow for photoelectric interactions, atomic relaxations resulting in Auger electrons, and Compton interactions. Tutor7pp is a standard code with the EGSnrc distribution that uses analog scoring to calculate energy deposition in each region along with energy reflected and transmitted through the geometry.¹⁴ The code was modified to tag and track photoelectrons, Auger electrons, Compton photons, or Compton electrons. If the tracked particle deposited its energy within the cell, its type and energy was outputted to a file and the energy was scored in a bin based on its type. After the simulation was complete, the code exported the energy deposited by Compton photons and electrons, Auger electrons, photoelectrons as well as the total energy deposited in the cell.

Four-component International Commission on Radiation Units and Measurements tissue was used as the base material for the cell, and the composition and density were altered for the different concentrations. Pegs4 was used to generate the material pegs file needed for the simulation. The cell was irradiated with photons using a parallel beam with a Gaussian distribution. The beam cross-section was 10×10 µm², had a full width at half maximum of 10, and varying mean energies. The global electron cutoff energy was set to 0.7 MeV, and the global photon cutoff energy was set to 0.01 MeV. DEF was calculated for each energy by dividing energy deposited in the cell with nanoparticles by energy deposited in the cell without nanoparticles for the same number of photon instances. The number of photon instances that provided a base dose of 2 Gy was used for each photon energy.

Results and discussion

Low-energy beams, ≤200 kV, with a penetration depth of ≤5 mm have limited applications and are mainly used for surface treatment of the skin. Medium energies 200 kV–1 MV (orthovoltage and supervoltage X-rays) are used for shallow treatments, while greater high-energy beams from 1 to 25 MV are used for treating deeply buried (>2 cm in depth) tumors.¹⁹

Figure 2A shows the results from analytical simulation, where DEFs as a function of source energy for copper, silver, gold, and bismuth nanoparticles are calculated. X-ray beams used for medical treatment are polyenergetic where the high-energy spectra have a low-energy component. At lower energies, X-ray matter interactions are dominated by the photoelectric effect, whereas the high-energy component interacts with matter by Compton scattering. The photoelectric effect is most effective until the photon energy from the external beam reaches a certain maximum that depends on the specific material under consideration. The interaction probability between the source photons and the nanoparticles, as well as the energy transfer coefficient in the medium, increases with decreasing source energy and increasing atomic number of the nanoparticles. The photoelectric cross-section varies as Z⁴ or Z⁵, where Z is the atomic number of the nanoparticle and is enhanced by an increased absorption by electron shells (K, L, M, etc) at low energies. The photoelectric effect is followed by the emission of an X-ray photon or an Auger electron. The Auger effect is mostly seen in Z<15 materials and is almost absent in materials with Z>60.¹⁹

Figure 2 (A) Analytical results for DEF as a function of energy for 20 mg of copper, gold, bismuth, or silver per gram of tissue. (B) Monte Carlo results for DEF as a function of energy for 20 mg of copper, gold, bismuth, or silver per gram of tissue.

Abbreviation: DEF, dose enhancement factor.

In materials with higher atomic numbers, such as gold and bismuth, there are two peaks corresponding to the K and M edge energies. Silver only has one peak, which corresponds to the L edge energy, while copper does not have these peaks. The peaks occur when the photon energy exceeds the edge energy and begins ionizing a different electron shell. Peaks near the edge energies support the theory that at low energies, photoelectrons and Auger electrons have a shorter trajectory in tissue, and therefore are more likely to deposit all of their energy in the cell. For high energies, the analytical simulation predicts a DEF of 1, ie, no enhancement, for all materials.

Figure 2B shows the results from the Monte Carlo simulation, where DEF varies as a function of energy for copper, silver, gold, and bismuth. For the energies of interest, the Monte Carlo simulation predicts a Gaussian distribution for the DEF, with a maximum at around 90 keV for gold and bismuth and 50 keV for silver and copper. The Monte Carlo simulation also predicts that at high energies the DEF levels out around 1.2 for copper and upward of 1.6 for gold. This is because the Z-dependent probability of the photoelectric effect gets smaller at the higher MeV ranges, where the Compton interaction starts to dominate. In Compton scattering, which is independent of Z, incoming X-rays eject an electron from the atom and scatters a photon. The scattered photon, depending on its energy, can further induce the photoelectric or Compton effect as well as generate secondary electrons. Secondary electrons generated by keV X-rays from the nanoparticle surface have lower energy, a shorter trajectory, and higher energy transfer rate to the surroundings, which in turn results in a higher DEF value.

Figure 3A–D shows a direct comparison between the analytical and Monte Carlo simulation for various materials at a concentration of 20 mg/g. For each material, the analytical and Monte Carlo simulations predict vastly different DEF values. One major difference is at high energies where the analytical simulation does not predict an enhancement (DEF =1), the Monte Carlo simulation does. One reason for this is in the limitations of the analytical simulation. The analytical simulation does not account for photoelectron interactions with neighboring nanoparticles. Since the Monte Carlo takes a random sampling for each trajectory step, it accounts for interactions with neighboring nanoparticles. These interactions could change the trajectory of the photoelectrons and Auger electrons, reduce the electrons energy, or cause another photoelectron or Auger electron to be ejected. In all the cases, the Monte Carlo simulation would continue to track its trajectory until all the energy is deposited within the cell. The analytical simulation on the other hand assumes a linear trajectory of the electron with a set probability of depositing the energy within the cell. The discrepancy at low energies may be due to the analytical simulation not accounting for recombination.

Figure 3 DEF results from the analytical simulation and Monte Carlo simulation for (A) gold, (B) silver, (C) bismuth, and (D) copper at a concentration of 20 mg/g.

Abbreviation: DEF, dose enhancement factor.

Finally, DEF was calculated as a function of energy for various concentrations of gold nanoparticles. Figure 4A shows the results from the analytical simulation, while Figure 4B shows the results from the Monte Carlo simulation. Both simulations predict an increase in DEF as the concentration increases, although the increase in the Monte Carlo simulation is larger than the increase in the analytical simulation for larger energies. At 10 keV, the Monte Carlo simulation does not predict concentration dependence. This may be due to several factors including recombination, which is accounted for in the Monte Carlo simulation.

Figure 4 The relationship between DEF and concentration for the (A) analytical simulation and (B) the Monte Carlo simulation.

Abbreviation: DEF, dose enhancement factor.

Conclusion

Nanoparticles of high atomic number elements could be used to enhance X-ray radiation therapy. Although both analytical and the Monte Carlo simulations predict an increase in dose enhancement in the presence of nanoparticles, the two methods do not give comparable results for few reasons.

1. The Monte Carlo simulation takes into account energy deposited by photoelectrons, Auger electrons, Compton photons and Compton electrons, and source photons, while the analytical one calculates energy deposited by source photons and photoelectrons.

2. The Monte Carlo simulation accounts for electron–hole recombination while the analytical simulation does not.

3. Monte Carlo simulation randomly samples the photon or electron path and accounts for particle interactions, while the analytical simulation assumes a linear trajectory.

This work demonstrates that for varying nanoparticle material, concentration, and X-ray beam conditions, the Monte Carlo simulations present a more comprehensive approach to calculate DEFs. Future works should involve improving the analytic simulation by adding in the Auger electron spectrum, taking into account the relaxation effects, and simulating discrete nanoparticles within the cell. Validation of the two simulations should be conducted through experiments carried out in phantoms or mice.

Disclosure

The authors have no conflicts of interest to disclose.

Acknowledgments

This project is supported by a Director’s New Innovator Award from National Institute of Health to Ming Su (1DP2EB016572).