Determining dose enhancement factors of high-Z nanoparticles from simulations where lateral secondary particle disequilibrium exists

Nanoparticles (NPs) containing high atomic number (high-Z) materials have been shown to enhance the radiobiological effectiveness of ionizing radiation. This effect is often attributed to an enhancement of the absorbed dose in the vicinity of the NPs, based on Monte Carlo simulations that show a significant local enhancement of the energy deposition on the microscopic scale. The results of such simulations may be significantly biased and lead to a severe overestimation of the dose enhancement if the condition of secondary particle equilibrium is not met in the simulation setup. This current work shows an approach to estimate a ‘realistic’ dose enhancement from the results of such biased simulations which is based on published photon interaction data and provides a way for correcting biased results.

In our paper (Rabus et al 2019) we presented a procedure to determine the dose enhancement factors (DEFs) for gold nanoparticles (GNPs) for an extended photon radiation field from simulations with narrow beam geometries. In the context of different work, it has become evident that there has been an error in the ROOT software routine used for processing the simulation data for the purpose of our article. More precisely, in the command line where the deposited energy obtained from the simulation was converted to the energy density in the spherical shells used for scoring, the volume of the spherical shells was implemented as (4/3) * 3.1416 * ((rmax * rmax * rmax) − (rmin * rmin * rmin)) Using integers in the factor 4/3 in the first bracket leads to an integer output, i.e. 1 instead of 1.3333… In consequence, the values shown as simulation results in figure 5 of our paper were enhanced by a factor of 4/3 and so were the values for the dose enhancement factors for an extended photon field in figures 6 and 7. The corrected figures are presented below.
Furthermore, we noted the following typos in the article: • A wrong sign in the denominator on the right-hand side of equation (4). The correct expression is In producing the data for figures 7 and 9, the correct equation as given above was used, however.
• In the paragraph on page 7 below equation (12), it should read '… is on the order of 1 /3 (µr) 2 . … z in equation (11) …' • In the Conclusion, the second recommendation should read '… as given in equation (8). …' Correcting the identified data processing error leads to a reduction of the predicted DEF for GNPs irradiated with an extended field of diagnostic x-rays by a factor of ¾. In consequence, the conclusions of our publication (regarding the magnitude of the DEFs for an extended beam and the limited radial range where significant DEFs occur) remain. Figure 5. Dose per fluence as a function of the radial distance from a gold nanoparticle (NP) of (a) 50 nm and (b) 100 nm diameter surrounded by water that is centrally irradiated by a 100 kVp photon beam with diameter of NP diameter + 10 nm for the case that the NP volume is filled with water (dark-grey squares) or gold (light-grey diamonds). The solid line is the value obtained for the absorbed dose in water estimated by taking into account only the primary photons and assuming secondary electron equilibrium (SEE) conditions. The light-grey circles indicate the estimated absorbed dose per fluence in the presence of the gold NP.

Introduction
Since the pioneering work of Hainfeld et al (2004), who demonstrated that better tumour control could be achieved by injecting mice with gold nanoparticles (NPs) prior to ionizing radiation exposure, so-called radioenhancing effects have been observed in vitro and in vivo for a variety of NPs containing high atomic number (high-Z) elements, as recently reviewed in the literature (Cui et al 2017, Kuncic and Lacombe 2018. A widely accepted assumption in studying photon irradiation of high-Z NP, is that the number of secondary electrons emitted by NPs is higher than in soft tissue and that a higher energy deposition within and around NPs is expected, due to the enhanced photo-absorption. These effects are usually quantified by the socalled dose enhancement factor (DEF), which is essentially the ratio of absorbed dose to water with and without the presence of the NP. In the latter case (i.e. NP absence), the volume otherwise occupied by the NP is filled with water. Among others, studies by Regulla et al (1998Regulla et al ( , 2000 showed that the presence of a 150 µm thick gold layer produces a physical enhancement of the absorbed dose to water that can be of several orders of magnitude at a distance of some micrometres from the metal/tissue interface. Moreover, Regulla et al (2000) could demonstrate a corresponding enhancement of biological damage to cells in contact with the metallic surface. In order to better understand the observed dose enhancement and radiobiological efficiency at the gold/tissue interface, the authors performed Monte Carlo simulations with the PARTRAC code to compare with their experimental results.
Unlike for these experiments, a direct measurement of the dose enhancement in the presence of metallic NP is challenging due to their inhomogeneous distribution in cells and in tissue, so that the scientific community is still relying on theoretical estimations of the local energy deposition. To this purpose, many Monte Carlo studies on spherical NPs were carried out in the last years. The majority of these investigate the radial dependence of the dose enhancement obtained around a single NP (for example, see Chow et al (2012aChow et al ( , 2012b, Lin et al (2014), Yahya Abadi et al (2014) and Incerti et al (2016)), while a smaller amount of work was developed up to now to estimate dose enhancement or secondary-electron production in the vicinity of small NP clusters (see Gadoue Determining dose enhancement factors of high-Z nanoparticles from simulations where lateral secondary particle disequilibrium exists 2 H Rabus et al et al (2018) and references therein) or in tissue-equivalent materials that are loaded with a given NP concentration (Jones et al 2010, Douglass et al 2013, Zygmanski et al 2013, Martinov and Thomson 2017, Zabihzadeh et al 2018. There have also been studies investigating the impact of NPs on nanodosimetric quantities or radiation damage yields (Dressel et al 2017(Dressel et al , 2019. The common challenge with such Monte Carlo simulations is the small target size of NPs, which typically have diameters of up to a few tens of nanometers. Consequently, techniques for variance reduction often need to be applied in order to obtain results with sufficiently low statistical uncertainties. The seemingly straightforward approach of using pencil beams or beam geometries confined to the size of the NP is a common pitfall that potentially biases results. This problem was already addressed by Zygmanski et al (2013) and Zygmanski and Sajo (2016), who well described the complexity of numerical studies of radiation transport involving multiple scales. In their work, Zygmanski and Sajo point out the importance of an accurate scaling of a typical 'macroscopic beam' used in radiotherapy to a 'microscopic' dimension that is suitable for particle-track simulation.
The problem of using such confined beam geometries is that there is a lack of lateral secondary particle equilibrium, which affects the scoring of absorbed dose (Zygmanski and Sajo 2016). In addition, for the case of electron radiation interacting with NPs, Rogers (2013) claimed that the enhanced dose deposition reported by Chow et al (2012b) was irrelevant since the effects of an extended primary beam were not considered by their use of a pencil beam. Rogers also concluded that in realistic scenarios there would be 'no advantage using gold NPs in an electron beam, unlike the case for photon beams'. When considering photon beams with small diameters (i.e. of the order of the NP diameter), the outcome of such a simulation will not only lead to an overestimation of the DEF, but also of the spatial range of such dose enhancement (Zygmanski et al 2013, Zygmanski andSajo 2016). This is based on the fact that when the beam diameter is significantly smaller than the range of the secondary electrons, the following three effects are unaccounted for: (1) In the case of irradiation with an extended photon field, electrons generated by photon interactions outside the confined simulated beam can also deposit energy in the shell segments used for scoring. Hence, there will be a severe underestimation of the real absorbed dose to water in both cases, i.e. with and without NPs present. If this underestimation is corrected for, this will significantly decrease the DEF.
(2) Due to the confined nature of the beam used in the simulation, the dose contribution from Compton or Rayleigh scattered photons interacting with the NP is underestimated. In a more realistic (i.e. extended) radiation beam, most of these photons would be produced outside the 'confined beam'. Including this contribution, would increase the absorbed dose in both the absence and presence of NPs. Nevertheless, a net increase in the DEF might be expected due to the NP. (3) The presence of NPs may lead to a radiation sink effect (Brivio et al 2017), i.e. an absorption of electrons produced outside the volume covered by the narrow photon beam in the simulation. The resulting reduction of energy deposition in the NP's 'shadow' (i.e. directly beyond it) may be expected to reduce the DEF. On the other hand, high-energy electrons produced outside the confined beam and interacting with the gold nanoparticle may also produce additional secondary electrons by impact ionization on the AuNP which would lead to an enhancement of the DEF. Although, considering the uncertainties involved in the corrections of the first two effects, this third correction may be irrelevant.
In the current work, an approach is presented to obtain corrections for the first two aforementioned effects by using literature data on photon interaction cross sections. The approach is applied exemplarily to the most frequently studied case of gold NPs but would also work for other high-Z materials. Furthermore, a few 'quality assurance' checks are recommended to verify uncompromised results when simulating radiation effects of NPs.

Simulation setup
The simulation data used in this work were produced in the frame of an intercomparison exercise of radiation transport and track structure Monte Carlo codes, which was recently conducted within the European Radiation Dosimetry Group (EURADOS) (Rühm et al 2018) by Task Group 7 'Internal Microdosimetry' of Working Group 7 'Internal Dosimetry' (Li et al 2019). The simulation geometry was such that a spherical gold NP surrounded by liquid water was irradiated by a parallel photon beam of circular cross section with a diameter 10 nm larger than that of the NP. The NP was at a distance of 100 µm from the circular source that emitted photons of a given spectral distribution and with a constant fluence (number of photons per area) over the source area. The exercise comprised simulations for two different NP diameters (50 nm and 100 nm) and three different photon energy spectra produced by an x-ray tube with a tungsten anode that was operated with acceleration voltages of 50 kV, 60 kV, and 100 kV, respectively. The radiation was filtered by 3.9 mm Al and 0.8 mm Be. The spectra were calculated using the SpekCalc code (Poludniowski et al 2009) and were provided in 500 eV energy bins for energies between 10% and 100% of the respective maximum possible x-ray energy.
For each of the six combinations of photon spectrum and NP size, simulations were performed with and without the NP. In the latter case (i.e. without NP), the volume otherwise occupied by the NP was filled with water. In both simulations, the energy deposition was scored in spherical shells around the NP and then normalized to the volume of the respective shell. For radial distances (from the NP center) between r np (radius of the NP) and r np + 1 µm, the radial bin size (or thickness of the spherical shell) was 10 nm. For radial distances between r np + 1 µm and r np + 50 µm, a radial bin size of 1 µm was used. The results from both simulations were normalized to the number of source photons and the bin-wise ratio of the results with and without the NP was reported as the DEF. The spectral flux of electrons ejected from the gold NP per primary photon was also reported.
In the intercomparison exercise, the narrow-beam geometry was intentionally chosen to allow for a better test of variations between codes. The data presented in this paper were among those that gave the largest DEFs, which were obtained from simulations performed with Geant4 version 10.02p01 (Agostinelli et al 2003) using the low-energy extensions of Geant4-DNA (Incerti et al 2010, Bernal et al 2015. The simulations inside the gold NP used the G4EmStandardPhysics mode. The Livermore models were used for photon interactions within the whole geometry as well as for electron transport inside the gold NP. A low-energy limit of 10 eV was used for electron transport, even though the cross-sections calculated with these models at such low energies do not give a good description of electron interactions in gold. Nevertheless, there are currently no other models in Geant4 that are well adapted for describing electron transport in gold for kinetic energies below 1 keV. (In fact, testing and benchmarking the validity of an artificial extension of the low energy limit of the Livermore models for Gold was the objective of the respective participant in the exercise.) The simulation of electron transport in liquid water surrounding the gold NP was performed using option7 of the Geant4-DNA cross sections. With this option, the whole energy range of electrons generated in the photon interactions for this particular problem was covered (50 kVp, 60 kVp and 100 kVp x ray spectra), where the cross sections for inelastic and elastic electron scattering of options 4 (E < 10 keV) and 2 (E > 10 keV) are combined. For an explanation of the physical models used in options 2 and 4, the reader is referred to Incerti et al (2018). For each combination of NP size and photon energy spectrum, the simulations were performed with 10 8 primary photons emitted from the source.
It should be noted that the simulation geometry used in these simulations leads to a bias in the determination of the DEF for several reasons. The two most important are the following: (1) The simulation uses an incident plane-parallel photon beam that is 10 nm wider than the diameter of the target NP. Secondary electrons that can contribute to the energy deposition in the scoring volumes may have longer range than 10 nm. For instance, 560 eV electrons have an approximate range of 10 nm in water. In contrast, as an example, Compton-scattered electrons from 30 keV photons (near the mean energy of 100 kVp spectrum) have a range up to 300 nm. Therefore, the selection of only a 10 nm additional layer of water beyond the NP perimeter excludes the contribution of a significant number of Compton electrons. (2) Spherical shell scoring volumes were used up to a radial distance of 50 µm. Therefore, beyond a radius of r np + 10 nm, all scoring shells included volumes of water unirradiated by the primary beam. This likely leads to a significant error in the computation of the DEF.

Approach to correct for the deficiencies of the simulation
The approach used in this paper is aimed at obtaining more realistic estimates for the radial dependence of the DEF around a high-Z NP from simulations that suffer from the deficiencies mentioned above by applying (a posteriori) corrections on the simulation results. It is based on the following assumptions: 1. For a setup ensuring secondary electron equilibrium (SEE) in the volume around the NP, simulations in the absence of the NP would yield a constant value of absorbed dose and this value should agree (within uncertainties) with the value derived from known photon interaction data. This is basically an assumption regarding the accuracy of the simulation code. It allows one to derive an upper limit on the average dose enhancement (see section 2.4).
2. With the NP present, the absolute deficiency in absorbed dose in each water volume around the NP (due to the lack of SEE) is essentially the same as for the case when the NP is absent (i.e. when the NP volume is filled with water). This assumption implies that the NP's effect on the contribution of the secondary electron field in its vicinity from electrons that originate from photon interactions outside the confined photon beam does not significantly alter the dose distribution around the NP. Using the simulation results together with the input photon spectrum, one can derive a correction for the DEF (see section 2.5) 3. The ratio of the contribution to the absorbed dose around the NP that results from interactions of scattered photons with the NP and the absorbed dose around the NP that is obtained when considering only primary photons is approximately given by the ratio of absorbed dose to gold obtained for the two spectral photon fluence distributions (scattered photons and primary photons) under SEE. This assumption forms the basis of the considerations presented in sections 2.6-2.8. It is justified by the fact that the photon energy absorption coefficient deviates only slightly from the photon energy transfer coefficient (Hubbell and Seltzer 2004). Hence, the absorbed dose to gold under SEE can be used as an estimate of the kinetic energy of electrons released by photon interactions in the gold NP. This assumption also implies neglecting the potentially different fraction of electron energy absorbed inside the NP owing to the different photon and, thus, secondary electron energy spectra. This difference would lead to an overestimation of the dose to water emerging from the NP.
4. The radiation sink effect (Brivio et al 2017), i.e. electrons losing more energy when crossing the NP than when the NP volume is filled with water (so that the energy deposited of the 'transmitted' electrons in the 'shadow' behind the NP is reduced), is negligibly small compared to the other corrections.
The approach is based on the relation between absorbed dose and known photon interaction cross sections (see section 2.3) under the condition of SEE. SEE means that the fluence of spectral secondary electrons leaving a volume is equal to the spectral fluence of secondary electrons entering the volume. Under SEE conditions, the absorbed dose in the volume is equal to the collision kerma, i.e. the kinetic energy per unit mass that is transferred to secondary electrons and deposited by non-radiative interactions of these electrons (ICRU 2011). It should be noted that the dose enhancement due to NPs is, in fact, due to the absence of SEE for the volume of the NP. The lack of SEE that is corrected with the approach presented here is, however, referring to the volume (filled mostly with water) that is passed by the photon beam in the simulation setup described in section 2.1.

Photon interaction data and absorbed dose
A photon beam propagating in matter is attenuated due to interactions with the atoms and molecules in the material. The rate of photon interactions per path length is given by the (total) attenuation coefficient µ, which is the product of the total scattering cross section and the number density of scatterers. µ has the dimension of a reciprocal length and depends on the photon energy. As the total scattering cross section is the sum of the cross sections for the different interaction processes, namely elastic (Rayleigh) scattering, incoherent (Compton) scattering, photoabsorption, and pair production, each of them can also be associated with an attenuation coefficient. In Compton scattering, a fraction of the photon energy is transferred to an electron. In photoabsorption and pair production, the photon disappears and its energy minus the binding energy is transferred to an electron or its full energy is transferred to an electron-positron pair, respectively. The fraction of the photon energy that is on average imparted to the material when a photon interaction occurs multiplied by the rate of photon interactions per path lengths is the energy absorption coefficient µ en . The energy is imparted to matter by the ionizing interactions of the charged secondary particles (electrons and positrons) released in photon interactions, and it differs from the energy transferred by the amount of energy lost in radiative scattering (bremsstrahlung) of the charged secondary particles.
Photon interaction data are often reported after normalization to the mass density ρ. The rationale behind this is that for mixtures, the mass attenuation and energy absorption coefficients can be obtained as sum of the respective quantities for the components weighted by the mass fraction of the respective components in the mixture.

Relation between energy absorption coefficient and absorbed dose
If the condition of secondary charged particle equilibrium is fulfilled, the absorbed dose D can be obtained from the energy absorption coefficient µ en by the following relation: where E is the photon energy, ρ is the mass density, and dΦ dE (E) is the spectral photon fluence (particles per area per energy interval).
For a small volume V of matter traversed by a photon beam, the average number of photon interactions in that volume, n, which is called event frequency in the field of microdosimetry (ICRU 1983), is given by where m is the mass contained in the volume. It should be noted that this equation refers to photon interactions and is valid for independent of prevalence of charged particle equilibrium.

Photon interaction data of gold
The mass energy absorption coefficients µ en /ρ in the energy range between 1 keV and 20 MeV were obtained for water and gold from the XAAMDI (X-Ray Attenuation and Absorption for Materials of Dosimetric Interest) data base (Hubbell and Seltzer 2004). The data are shown in figure 1, where the y axis shows the product of photon energy E and µ en /ρ as this combination is relevant for dosimetry (see equation (1)). Figure 1 also shows the curve for a mixture of gold and water with a mass fraction for gold of 2 × 10 −4 as this fraction is considered a typical value for radiological applications of gold NPs (Liu et al 2010, Kim et al 2012. The relative enhancement of E × µ en /ρ for such a mixture as compared to pure water is shown in the inset. Under SEE conditions, this ratio is equal to the ratio of absorbed dose in a mixture of gold and water to the absorbed dose in water. This ratio constitutes an upper limit for the DEF in water medium (see section 2.4).
Data for the total mass attenuation coefficient µ/ρ of gold and water as well as the mass attenuation coefficients for coherent (Rayleigh) and incoherent (Compton µ inc /ρ) scattering and photoelectric absorption in water were obtained from the XCOM data base (Berger et al 2010) over an energy range between 1 keV and 20 MeV. In addition to the standard energy grid, which is the same as for listings in the XAAMDI data base, XCOM can provide results for additional energy values. Hence, additional data were obtained in 1 keV steps in the range from 1 keV to 100 keV.
The minimum in E × µ en /ρ that occurs for water around 60 keV (see figure 1) coincides with the maximum in the cross section for Compton scattering (Berger et al 2010). In fact, incoherent (Compton) scattering is the largest contribution to the photon interaction cross sections for water in the energy range between 30 keV and 20 MeV and accounts for more than 80 % (and up to almost 100%) in the energy range from 50 keV to 800 keV (see figure 2). For gold, on the other hand, the maximum Compton scattering cross section occurs at around 85 keV (Berger et al 2010), but in this energy range (up to 500 keV) the photon interaction cross section is dominated by photoelectric absorption (figure 2).
In photoelectric absorption, the photon energy is almost completely transferred to the photoelectron (except for the binding energy which is in the order of 10 eV). In Compton scattering, the energy transferred to electrons is given by where E is the photon energy, ϑ is the scattering angle of the photon and m e c 2 = 511 keV is the rest energy of the electron. Hence, for photon energies significantly smaller than 511 keV, only a small fraction of the photon energy is transferred to electrons and, in water, this energy is mostly deposited within several micrometres (Berger et al 2005) while the major fraction of the energy is carried away by the scattered photon (having a mean free path in water in the order of centimetres). Hence, the one to two orders of magnitude difference in E × µ en /ρ between gold and water is in the energy range from 50 keV to 250 keV.

Upper limit for the average dose enhancement around a NP
As was mentioned in section 2.3, the normalization of the mass energy absorption coefficient with the mass density offers the advantage that the mass energy absorption coefficient for a mixture of materials can be obtained as a weighted sum of µ en /ρ of the components using the mass fractions of the components as weighting factors. In this section, an estimate for the radial dependence of the DEF is derived by considering the dose enhancement of a homogenous mixture of gold and water within a sphere where the mass fraction of gold is the same as when a spherical gold NP is in the center of this sphere and surrounded by water. If the NP radius is r np , the radius of the sphere is r sp , and the mass densities of gold and water are ρ Au and ρ w , respectively, then the mass fraction x Au is given by Under SEE conditions, the energy E deposited in a volume filled with and surrounded by a mixture of gold and water with mass fraction x Au of gold is the sum of the energies m w K col,w and m Au K col,Au deposited by electrons that originate from a photon interaction with water molecules or gold atoms, respectively.
Here, m w and m Au are the masses of water and gold in the considered volume, while K col,w and K col,Au are the collision kermas (ICRU 2011) defined by where (µ en /ρ) w and (µ en /ρ) Au are the mass energy absorption coefficients of water and gold, respectively. Under SEE, the collision kermas K col,w and K col,Au are equal to absorbed dose to water in water, D w , and absorbed dose to gold in gold, D Au , respectively. Therefore, the absorbed dose in the considered volume of a mixture of gold and water is given by It should be noted that equation (7) is based on the additivity of absorbed energy and does not imply additivity of absorbed dose.
From equation (7), an average DEF DEF mix for a mixture of gold and water with respect to water only is obtained as If we now consider a spherical volume with radius r sp that is filled with a mixture of gold and water and that is surrounded by water, SEE between this volume and its surrounding may be violated due of the higher interaction cross sections of the gold atoms, if the sphere dimensions are smaller than the range of the secondary electrons. Therefore, the average DEF for the spherical volume is smaller than the value obtained by equation (8). If gold atoms in the considered sphere are not uniformly mixed with water, but rather concentrated in a NP, the average dose in the water surrounding the NP is smaller than the absorbed dose in the case of the mixture of gold and water, because the energy deposited in gold is deposited within the NP, i.e. outside this water-filled volume. Therefore, the average dose enhancement given by equation (8) is larger than the average DEF in the sphere of radius r sp around the NP center which is given by and constitutes an upper limit for the average DEF in the NP case.
2.5. First order estimate for absorbed dose to water in the vicinity of the nanoparticle Using equation (1) to estimate the value of the absorbed dose to water in the absence of the NP, the spectral fluence of photons in the respective volume has to be known. In a first order approximation, this spectral fluence, dΦ dE , can be obtained by considering the attenuation of the primary photon beam, i.e. the fluence of photons that have not undergone collisions after travelling the distance z from the photon source.
Here, µ is the total mass attenuation coefficient, Φ is the total fluence (photons per area, summed over all photon energies) at the photon source (z = 0), and dΦ rel dE is the relative spectral distribution (photons per energy interval) that satifies the normalization condition In the 5 keV to 100 keV photon energy range, µ has values below 5 × 10 −3 µm −1 . If e −µz is averaged over a spherical shell of radius r, the relative deviation from the value at the sphere's center is on the order of ⅓ (µr)2. Hence, for the spherical shells around the NP volume used for scoring the energy deposition in the simulations, z in equation (12) can be set equal to the distance of the NP volume center from the photon source plane.
Inserting equation (11) in (1) and taking into account that the photon energy spectrum is given as a histogram gives: Here Φ rel,j is the relative frequency of photons in the energy bin of width ∆E j centered at E j .

Estimating the contributions of first generation Compton scattering
The dominance of Compton scattering for water in the photon energy range above 30 keV (see section 2.3.2) implies that in an extended photon field, scattered photons originating from primary photons whose trajectory passes the NP at large lateral distances (in the order of centimetres) may also reach the volume around the NP location and contribute to the absorbed dose. The spectral fluence dΦc1 dEc of photons that have undergone only one Compton scattering process before reaching the NP can be estimated in the following way: where r and cosϑ are spherical polar coordinates with respect to the NP center and azimuthal symmetry has been considered. In these two equation, E p and E c are the photon energies before and after the Compton scattering process, E p,m is the maximum photon energy in the primary photon spectra, r sp is the radius of the spherical NP, r max is the maximum radial distance of the points in space where a Compton interaction event occurs, z 0 is the distance of the NP from the plane of the primary photon source, n w is the number density of molecules in water and d 2 σ c /dE c dΩ is the double differential cross section for Compton scattering. The 'transfer function' T E c |E p given by equation (15) is the total probability that a photon of energy E c that originates from a Compton scattering event of a photon of energy E p reaches the NP. To obtain this probability, one has to integrate over all radial distances and polar and azimuth angles defining points where this Compton scattering interaction occurs. The terms under the integrals are the probability that a primary photon reaches the interaction point without a preceding interaction, e −µ(Ep)(z0−r cos θ) , the probability (per unit length) for a Compton scattering interaction in which the scattered photon is emitted towards the NP, n w d 2 σc dEcdΩ ∆Ω, and the probability that the Compton photon reaches the NP without interaction on its way, e −µ(Ec)r .
In first order approximation, the dependence on scattering angle of the double differential cross section for Compton scattering can be assumed to be the same as for Compton scattering on free electrons, so that d 2 σc dEcdΩ is given by where µ c is the attenuation coefficient for Compton scattering, i.e. the probability per unit length that a Compton interaction occurs, d 2 σ KN /dE c dΩ is the double differential cross section described by the Klein-Nishina formula and σ KN is the respective total cross section. Furthermore, m e c 2 and r e are the rest energy of the electron and the classical electron radius, respectively. δ is Dirac's delta function, which establishes the one-to-one correspondence between cos θ and E c for a given E p , where the Compton scattering angle θ c is given by: Using equations (17) and (18) and the value ∆Ω = r 2 sp π/r 2 .of the solid angle covered by the NP (seen from the point of the Compton interaction) transforms equation (15) into where the single differential Klein-Nishina cross section is given by If we consider an infinitely extended primary photon field, the upper limit of the remaining integral in equation (19) is infinity for cos θ c 0 and z 0 / cos θ c otherwise, such that Here T 2 = e −µ(Ec)z0/ cos θc) for cos θ c > 0 and 0 otherwise.
rel,i ∆Ei of photons that underwent a single Compton scattering event before reaching the NP with an energy in the energy bin of width ∆E i centered at E i , was obtained from where the sum extends over all energy bins of width ∆E j centered at E j with E j E i .

Estimating the contributions of first generation Rayleigh scattering
Rayleigh scattering of photons in water contributes by as much as 12% to the total attenuation coefficient. As photons do not lose energy in Rayleigh scattering, in an extended photon field, also Rayleigh scattered photons originating from primary photons whose trajectory passes the NP at large lateral distances may reach the volume around the NP location and contribute to the absorbed dose.
The spectral fluence dΦr1 dEp of first generation Rayleigh scattered photons interacting with the NP can be estimated in a similar approach to that for Compton scattering. As the photon energy is conserved, the overall expression simplifies to: where r and cosϑ are spherical polar coordinates with respect to the NP center and azimuthal symmetry has been taken into acount. Again, r sp is the radius of the spherical NP, z 0 is the distance of the NP from the plane of the primary photon source, n w is the number density of molecules in water, dσ R /dΩ is the differential cross section for Rayleigh scattering and ∆Ω = r 2 sp π/r 2 is the solid angle covered by the NP as seen from the point of the scattering interaction. The terms under the integrals are the probability that a primary photon reaches the interaction point without a preceding interaction, e −µ(Ep)(z0−r cos θ) , the probability (per unit length) for a Rayleigh scattering interaction in which the scattered photon is emitted towards the NP, n w dσR dΩ ∆Ω, and the probability of the scattered photon reaching the NP without interaction on its way, e −µ(Ep)r .
Analogous to equation (16), the differential cross section dσR dΩ for Rayleigh scattering can be approximated by where µ R is the attenuation coefficient for Rayleigh scattering, dσ T /dΩ is the differential Thomson cross section σ T is the respective total cross section and r e is the classical electron radius. Using equations (24) and (25) and the value of the solid angle ∆Ω transforms equation (23) into The upper limit of the radial integral in equation (26) is infinity for cos θ c 0, and z 0 / cos θ otherwise, such that where S 2 = e −µ(Ep)z0/ cos θ) for cos θ > 0, and 0 otherwise.
rel,i ∆Ei of photons that underwent one single elastic scattering event before reaching the NP with an energy in the energy bin of width ∆E i centered at E i , was obtained from

Estimating the contribution of higher generation Compton and Rayleigh scattered photons
Calculating the contributions of photons undergoing multiple scattering interactions is not as straightforward as for the first generation, and thus would require a Monte Carlo approach. To get a rough estimate of the contributions of photons arriving at the position of the NP after having undergone two or more scattering events, one can consider the production rates of higher-order generations of photons. The probability for production of a Compton photon of energy E c by incoherent scattering of a photon of primary energy E p is given by The probability for Rayleigh scattering is then given by Consequently, the spectral fluence dΦn dEc of photons of energy E c that have undergone n scattering interactions is given by the recursive relation with the starting condition In equation (31), E p,m is the maximum energy of the primary photon spectrum. Equation (31) follows from the consideration that a photon of energy E c in the nth generation of scattered photons can either result from a Rayleigh scattering event of a photon of the same energy E c (occuring with a probability given by equation (30)) or from a Compton scattering event (occuring with a probability given by equation (29)) of a photon from the preceding generation that has an energy E p E c .
For all considered primary photon energy spectra used in the starting condition given in equations (32) and (31) was repeatedly applied to obtain the photon energy spectra of the first 20 generations of scattered photons. Figure 3(a) shows the development of the photon fluence spectra with the contribution from different generations for the case of the 100 kVp spectrum. It can be seen that after about ten scattering processes the cumulative spectral fluence distribution converges. For the other two photon spectra, 50 kVp and 60 kVp, the effect of the higher-order generations of scattered photons is smaller, as can be seen in figure 3(b) which shows the ratio R 20,1 of the sum of photon spectral fluences of the first 20 generations of scattered photons and of the first generation of scattered photons.
For the 100 kVp spectrum, the additional fluence from the higher-order generations results in a total fluence that is as much as a factor of five higher than the fluence of the first generation. For the 50 kVp and 60 kVp spectra, on the other hand, the increase in fluence is less than a factor of three owing to the combined cross sections for Compton and Rayleigh scattering being significantly smaller in these energy ranges.
In order to estimate the fluence contribution dΦ all /dE from all generations of scattered photons reaching the NP, the sum of the fluences of the first generation of scattered photons (calculated based on equations (14) and (23)) were multiplied with the ratios shown in figure 3 Excel VBA programs developed for implementing equation (31) were also used to calculate the secondary photon spectra for a Co-60 photon spectrum calculated with BEAMnrc (Rogers et al 1995). In this case, the spectrum was binned in 7 keV steps, and pair production was taken into account by adding two 511 keV annihilation photons for each pair production event. The results are shown in figure 3(c).
This approach for estimating the contribution of higher order scattered photons using equation (34) involves the assumption that the ratio of the sum of spectral fluences from all generations of scattered photons to the spectral fluence of the first generation of scattered photons at the NP location is the same as would be obtained if all scattering events had an isotropic distribution of the secondary photons. To get an estimate of the uncertainty introduced by this assumption, the fluences of photons that reach the NP after either one single Rayleigh or one single Compton scattering event were also calculated under the condition that the distribution of these scattered photons is isotropic. For this, the double differential cross section given in equation (16) was replaced by and the single differential cross section given in equation by where µ c and µ r are the attenuation coefficient for Compton and Rayleigh scattering, respectively, dσ KN /dE c and σ KN are the single differential and total Klein-Nishina cross sections, and n w is the number density of molecules in liquid water.
Using these two equations in conjunction with equations (14), (15) and (23), the fluence of photons reaching the NP after single scattering event (Rayleigh or Compton) is obtained as: where the asterisk indicates the unrealistic assumption of isotropic scattering.

Energy distribution and estimated range of secondary electrons produced in photon interactions
In the iterative procedure described in the preceding section, an estimate of the energy spectrum of secondary electrons that are produced by photon interactions in water was obtained as a by-product. For Compton scattering, it was assumed that the difference of the photon energies before and after the scattering interaction, E p − E c , is fully transferred to the electron, i.e. the binding energy of the electron was neglected. In addition, secondary electrons are also produced by photoelectric absorption. As pair production is not possible in the energy range used in this study, the probability for photoelectric absorption is given by For energies higher than about 543 eV, the photoabsorption cross section of water is dominated by the partial cross section for photoabsorption on the oxygen 1s orbital, which accounts for about 95% of the total cross section (Cullen et al 1997). The 1s core hole of oxygen has 99.5% probability for non-radiative (Auger) decay which leads to further emission of an electron of about 500 eV energy (Krause 1979). To simplify the analysis, it was assumed that photoabsoption in water always occurs in the O 1s orbital and always lead to emission of a 500 eV Auger electron. The energy of the photoelectron was consequently taken as the photon energy minus the 1s binding energy of oxygen, which is about 543 eV.
In order to obtain an estimated range distributions of the electrons released in photon interactions, the binning of electron energies was done in a such a way that the energy bin boundaries corresponded to equidistant bins for the logarithm of the electron range in the continuous slowing-down approximation (CSDA) with 25 bins per decade. To determine the energy bin boundaries, data for the energy dependence of the electron CSDA range were obtained from the online ESTAR data base using the standard energy grid for energies between 10 keV and 100 keV (Berger et al 2005). The data were plotted in an EXCEL diagram of electron energy as a function of the CSDA where the power function was used as the trendline. The resulting dependency of energy on electron range R was obtained as 5.897 keV × (R/µm) 0.5686 with a maximum deviation of less than 1% over the whole range of data points retrieved from the ESTAR data base. This power function was then implemented to calculate the energy bin boundaries for CSDA ranges between 100 nm and 100 µm. The first energy bin was used as underflow bin so that all electrons with energies below its upper bin bounday were counted in this bin.
It should be noted that this approach requires an extrapolation of the established energy-range relation to energies below 10 keV. Furthermore, it should be noted that the CSDA approximation is not appropriate for electrons with energies below 1 keV (with a CSDA range on the order of 50 nm). However, the purpose of this approach was to estimate the fraction of electrons with ranges exceeding 1 µm (i.e. of energies exceeding about 6 keV) which is not affected by the inappropriateness of the CSDA at electron energies below 1 keV, because all these electrons are scored in the first bin (upper bin boundary corresponding to about 1.7 keV).

Range of secondary electrons
The frequency distribution of the secondary electrons produced by photon interactions in water is shown in figure 4 as a function of the electron range in the continuous slowing down approximation (CSDA). The data includes electrons from the first 20 generations of scattered photons as described in section 2.8. The graph in figure 4 is in the form of a so-called lethargy plot, where the x axis is logarithmic, and the y axis is the product of the frequency per bin width multiplied by the x value of the bin centre. This representation has the advantage that the areas of the columns are proportional to the relative contribution of the respective bin to the average value over all bins.
In the distribution of electron ranges, the two broad structures represent the energy spectrum of Compton electrons and the higher-energy distribution of photoelectrons. As expected, these two components of the range distribution shift towards larger ranges when the photon energy spectrum shifts towards higher energies. The additional intense peak in the first bin of the electron range is the distribution of Auger electrons which have energies around 500 eV.

Estimated absorbed dose to water and resulting DEF
As an example, figure 5 shows the Geant4 simulation results obtained for both 50 nm and 100 nm NP sizes with a 100 kVp x-ray spectrum. These results are compared with the absorbed dose per fluence obtained from equation (13), this is taking only the primary photons into account and assuming SEE. The latter assumption and the negligible attenuation of the primary photon beam over the small dimensions considered here, implies that the absorbed dose in the absence of the gold NP is a constant. As can be seen from the 50 nm NP results, this value of absorbed dose is two to four orders of magnitude higher than that obtained in the simulation with a photon beam collimated to the size of the NP. For the 100 nm NP, this difference amounts to between one and three orders of magnitude. The larger cross section of the primary photon beam enhances the part of the electron range distribution (see figure 4) that is properly considered in the scoring of imparted energy.
Generally, the absorbed dose to water for SEE in the absence of the NP is also significantly larger than the extra dose coming from the NP. The light-grey circles in figure 5 indicate the sum of the constant absorbed dose to water and the difference between the simulation results obtained with and without the NP. These values are a first-order estimate for the absorbed dose to the water shells surrounding the 50 nm NP. As can be seen for the 100 kVp x-ray spectrum, the presence of the gold NP leads to significant dose enhancement only within about 100 nm from its surface. Figure 6 shows a comparison of the DEFs obtained for the simulation setup following correction of the absorbed dose to water for SEE (calculated with equation (13)). The simulation results suggest DEFs higher than 300 for both NP sizes within the first 10 nm spherical shell around the NP surface. Furthermore, after an initial drop within the first 200 nm, the DEF saturates at values of about 15 and 35 for the 50 nm NP and 100 nm NP respectively, with a small decrease out to 1 µm from the NP surface. Significant DEF values are also observed at larger distances from the NP. For example, the DEF exceeds a value of two for distances up to 10 µm from the NP surface.
On the contrary, the DEF values are considerably lower if SEE is considered. In fact, with SEE taken into account DEFs of about 10 and 17 can be observed for the 50 nm NP and 100 nm NP, respectively.
In the case of SEE, the DEF for both the 50 nm and 100 nm NP drops rapidly within the first 100 nm from the NP surface to a value below 1.2 and 1.8, respectively (see inset of figure 6). Similar results were observed for other radiation qualities considered in this work, as shown in tables 1 and 2 for the 50 nm and 100 nm NPs, respectively. While the general trend of DEF with radiation quality is maintained under SEE conditions, the total DEF values close to the NP surface are smaller by as much as a factor of 30. Furthermore, the radial distribution can be seen  . Dose per fluence as a function of the radial distance from a gold nanoparticle (NP) of (a) 50 nm and (b) 100 nm diameter surrounded by water that is centrally irradiated by a 100 kVp photon beam with diameter of NP diameter +10 nm for the case that the NP volume is filled with water (dark-grey squares) or gold (light-grey diamonds) obtained with Geant4. The solid grey line is the value obtained for the absorbed dose estimated by taking into account only the primary photons and assuming SEE conditions. The light-grey circles indicate the estimated absorbed dose per fluence in the presence of the gold NP.
to drop more rapidly for the SEE case, even though DEFs depart from unity by a few percent at 1 µm from the 100 nm NP surface.
As illustrated in figure 7, also the average dose enhancement in a water shell around the NP obtained from the simulations is significantly higher than the upper limit established in section 2.4 for respective water shell thicknesses of about 40 nm and 80 nm for the 50 nm NP and 100 nm NP, respectively.

DEF modification due to Compton or Rayleigh scattered photons
Based on equations (22) and (28), the contribution of Compton and Rayleigh scattered photons to the photon fluence are shown in figure 8 for the three different radiation qualities. In order to obtain an estimate of the additional contribution to the fluence from photons reaching the NP after more than one scattering interaction, ratios of the sum of the fluences for all generations and the fluence of the first generation of scattered photons (shown in figure 3(b)) were considered. Thus, for each primary photon spectrum, the sum of the first generation photons reaching the NP after one Compton or Rayleigh scattering was multiplied by the difference between the respective ratio and unity. The resulting estimated contribution of higher-order generations of scattered photons to the spectral fluence are also shown in figure 8.
As can be seen in figure 8, the highest primary photon energies are heavily suppressed in the single Compton scattered spectra. This is due to the fact that these energies can only be obtained for small scattering angles, and hence only a fraction of the photon beam cross section contributes to the spectrum.
The spectral features appearing in the Compton scattered 100 kVp spectrum (figure 8(c)) can be explained as follows: The large increase in fluence with decreasing photon energy around 72.5 keV is due to the contribution of Compton-scattered photons from the intense 84.5 keV characteristic line, which are scattered at angles of 90° or greater. Therefore, the effective source region where these Compton photons are produced extends to infinity. The reduction in fluence with decreasing photon energy at around 63.5 keV is due to the fact that this energy corresponds to the minimum energy that a primary photon of 84.5 keV can retain after Compton scattering. Similar arguments can be used for the other features observed in the Compton spectrum as well as their relation with the other characteristic lines in the primary spectrum.
The sums of the spectral fluences of first-generation scattered photons as obtained by equations (22) and (28) as well as the spectral fluences for all generations of scattered photons estimated by equation (34) were used to calculate via equation (1) for gold and water the dose per fluence under SEE which is identical to the collision kerma, i.e. the kinetic energy transferred to electrons by photon interactions and later deposited in the material by non-radiative electron interactions. The ratio of the respective resulting integral value (from equation (1)) to Figure 6. DEFs for the 100 kVp photon energy spectrum calculated from the radial distribution of imparted energy (see figure 5). The diamonds represent the results from the simulation with circular photon beam sizes of 60 nm and 110 nm diameter for the 50 nm and 100 nm diameter spherical gold NPs, respectively. The circles show the results obtained for a larger photon field ensuring SEE, where the absorbed dose to water was estimated based on equation (13) (solid line in figure 5). The inset shows a close-up of the latter data in the 100 nm to 500 nm range. Figure 7. Average DEF in a spherical water shell around a gold nanoparticle (NP) with increasing outer radius of the spherical shell. The data apply to the 100 kVp photon energy spectrum and are calculated from the radial distribution of imparted energy (see figure 5). The diamonds represent the results from the simulation with Geant4 for circular photon beam sizes of 60 nm and 110 nm diameter for the 50 nm and 100 nm diameter spherical gold NPs, respectively. The circles show the results obtained for a larger photon field ensuring secondary electron equlibrium (SEE), where the absorbed dose to water was estimated using equation (13) (solid line in figure 5). The solid lines are the upper limits according to equation (8). Table 1. Comparison of the results for the DEF in the vicinity of a 50 nm gold nanoparticle (NP) for different distances from the NP surface. (The values in column 1 give the upper bin boundary, i.e. 10 nm corresponds to a spherical shell between 25 nm and 35 nm radius.) The values obtained from the Geant4 simulations with a photon beam confined to the NP diameter +10 nm for the three radiation qualities (50 kVp, 60 kVp and 100 kVp) are given in the first three columns. The last three columns show the results for the case that the dose is estimated from equation (1) Table 3. Ratios of the dose per fluence calculated with the combined spectrum of primary and scattered photons to the value calculated with the primary photon spectra only. The last column gives the resulting overall enhancement factor for the deviation of the DEF from unity.
Relative enhancement of dose/fluence with respect to only primary photon spectrum Resulting enhancement factor of (DEF-1) Water Gold 50 kVp 1st generation (equations (21) and (27) (21) and (27) (21) and (27) the value obtained for the the NP with the primary photon spectrum is the factor by which the deviation of DEF from unity changes when scattered photons are considered in addition to the primary photon spectral fluence. Table 3 summarizes the resulting relative increase in the dose-to-fluence ratios in water and gold due to the estimated spectral fluence of first and higher-order scattered photons that reach the volume occupied by the NP (filled with water in the absence of the NP). The last column in table 3 gives the ratio of the values in the two preceding columns which is the correction factor to be applied to the deviation of the DEF from unity assuming that the excess dose in the surrounding water due to photon interactions in the gold NP scales with the expected dose per fluence for gold under SEE conditions (i.e. the collision kerma).
For the three x-ray spectra considered in this work, the first generation of scattered photons causes an increase in the dose to water in the absence of the NP that amounts to 215%, 224%, and 241%, respectively, of the absorbed dose per fluence found for the primary spectrum only. When the estimated contribution from second and higher-order generation photons is also included, the absorbed dose to water in absence of the gold NP is enhanced by about 3.7, 4.3 and 6.6 times the value from the primary radiation for the respective 50 kVp, 60 kVp and 100 kVp primary spectra. If the angular dependence of the scattered photons is ignored and the photon spectrum of the first generation of scattered photons is calculated according to equation (37), then the relative enhancement compared to the primary spectrum is significantly smaller and indicates that the predicted dose per fluence may be understimated by about 30% when assuming isotropic emission of the scattered photons in the case of the 50 kVp spectrum. For the higher energetic spectra, this underestimation decreases to below 10% for the 100 kVp spectrum.
For the case of gold, the relative increase in the dose per fluence is found to be almost the same as for water with only a 8% larger increase for gold with the 100 kVp spectrum. For this spectrum, the correction factor obtained when considering only the first generation of (assumed) isotropically scattered photons amounts to 0.96 as opposed to 1.03 if the proper angular distribution is considerd. The 7% difference between the two values may be considered as an estimate for the uncertainty imposed by the approach used to obtain the spectral fluence for all generations of scattered photons (equation (34)).
It should be noted that the results reported in figure 8 are independent of the NP size because the fluence is normalized to the NP cross sectional area and the solid angle covered by the NP as seen from an interaction point is proportinal to this cross section. Within the approximation established by assumption 3 in section 2.2, this independence on NP dimensions is also true for the results reported in table 3. Table 3 also contains the corresponding information for the case of a Co-60 spectrum as found at one of the reference irradiation fields of PTB used for detector calibration (Chofor et al 2007). For this photon energy spectrum the source to NP distance was taken to be 1 cm instead of 100 µm such as to assure longitudinal SEE owing to the large range of the secondary electrons in this case. Here, the first generation of scattered photons Figure 9. Comparison of the upper limits for the average DEF in a spherical water shell around a gold nanoparticle (NP) with increasing outer radius of the spherical shell for a 100 kVp x-ray photon energy spectrum and for a Co-60 gamma spectrum of an irradiator facility. The curves were calculated using equation (8). that reach the NP lead to an enhancement of the deviation of the DEF from unity by about 15% when the proper angular distribution of scattered photons is considered. If alternatively the first generation scattered photons are assumed to have an isotropic distribution, this increase turns out to be 30%. If the difference is again taken as an estimate for the uncertainty of the approach presented in section 2.8, this uncertainty amounts to 15%. Considering all generations of scattered photons according to equation (34) leads to an enhancement of the deviation of the DEF from unity by a factor of about 6.5.
To explore the impact of this large enhancement factor for the deviation of the DEF from unity found for the scattered photons from the Co-60 spectrum the upper limits according to equation (8) for the 50 nm and 100 nm gold NPs are compared in figure 9 with the respective limits for the 100 kVp spectrum. For both NP sizes, the upper limits for the deviation of the average DEF from unity are by a factor of about 250 lower than for the 100 kVp spectrum, so that even with the 6.5 times increase due to scattered photons, the DEF for the Co-60 will still deviate by a factor of 40 less from unity than the DEF for the case of the 100 kVp x-ray spectrum.

Electron range and photon field size
As can be seen from figure 4, for the photon energy spectra considered in the simulations, a significant fraction of electrons released in photon interactions in water have ranges exceeding 100 nm. For the 50 kVp spectrum, this is the case for 73% of the released electrons, while for the 60 kVp and 100 kVp spectra the respective fractions are 75% and 83%, respectively. Furthermore, the fraction of electrons with ranges exceeding 10 µm amounts to 26%, 23% and 17% for the respective 50 kVp, 60 kVp and 100 kVp spectra. Hence, for obtaining lateral SEE, also for these comparatively low-energetic photon energies, a diameter of the photon beam in the order of 50 µm to 100 µm would be needed (see figure 4). This beam diameter would then ensure lateral SEE in the spherical shells around the NP position, which is in good agreement with the recommendations of Zygmanski and Sajo (2016), that are based on the maximum electron range in water.
For higher-energy primary photons, the photoelectron energies will increase linearly with photon energy, thus worsening the situation. Even more important is the increase in the range of electrons released in Compton scattering, as it follows from equation (3) that the upper limit of possible electron energies may increase to values close to those of photoelectrons. For the principal gamma photons of Co-60 radiation, with energies of 1.1732 MeV and 1.3325 MeV (Helmer and van der Leun 2000), the maximum possible Compton electron energies are 963 keV and 1119 keV, respectively, with an average energy of about half these values. The ensuing CSDA range of the maximum-energy Compton electrons is then between 400 µm and 500 µm (Berger et al 2005).
These considerations suggest that only few simulation studies performed to determine the DEF around NPs were able to obtain SEE conditions, and thus in most cases DEFs are generally overestimated.
The problem that simulations face is that to obtain results with small statistical uncertainty a sufficiently large number of photons need to interact with the NP. If the simulation is set up in such a way that the photons are emitted as a parallel beam starting uniformly from a source area of dimensions large enough to assure SEE, most of the primary photons will simply pass the NP without interacting. For a 100 nm-diameter NP and a circular photon beam of 100 µm diameter, the probability for a primary photon trajectory to cross the NP is only 10 −6 . The simulation results shown in figure 5 required 6-12 h of CPU time, where 10 8 primary photons originated from a circular area which was only 10 nm larger in diameter than the NP. In order to obtain the same simulation statistics in such a straight-forward simulation approach, one would have to calculate 10 14 primary photons which would require about 1000 years of CPU time (based on current computer technology) and is, hence, unfeasible. Fully exploiting the spherical symmetry of the problem and applying variance reduction techniques, would significantly reduce the CPU time required. However, correct application of variance reduction techniques requires a level of expertise as they should be applied with care.
Whilst discussion of these issues is beyond the scope of this paper, a potentially better way to obtain DEFs from radiation transport simulation involving NPs is to perform a multi-scale simulation, as suggested by Zygmanski et al (2013). In such a simulation, a broad photon beam is transported in water up to a depth that allows for charged particle equilibrium, where a phase-space file with the same angular and energy distribution is then used to transport a 'small' photon beam towards a single GNP. Alternatively, one could first simulate the spectral fluence of photons and secondary electrons as a function of depth and lateral displacement with respect to the starting point of the primary photons. As such simulations could be done for an ideal pencil beam, scoring the lateral displacement would be necessary to account for the finite size of the photon field. In the second step, these fluence spectra would be used to simulate the energy deposition on the region of interest in the presence and absence of the NP.
Such a two-step simulation approach would also rely, to some extent, on the assumption that the NP has only a minor influence on the spatial distribution of energy imparted by secondary electrons produced outside the region of interest. Therefore, the approach presented in this paper offers a conceptionally easier alternative to extract a more realistic estimate of the DEF from simulations that were performed with a confined photon beam size.

Contribution of scattered photons
The findings reported in table 3 suggest that, at least for the low-energy x-ray spectra considered in the simulations, taking into account the first generations of scattered photons that reach the NP will not significantly influence the value of the DEF. It should be noted that the fluence of these extra photons was derived in sections 2.6 and 2.7 for the case of an infinite lateral extension of the photon beam. For a photon beam of finite size, the fluence of photons that reach the NP region after one scattering interaction will be lower than the one-to twofold increase seen in figure 8. However, the principal finding of a negligible influence on the DEF is still expected as long as the relative increase in dose due to interactions of the scattered photons with the gold NP can be estimated based on the mass energy absorption coefficient of gold. Since surface effects of the NPs may lead to an increased number of low-energy electrons emitted from the NP (Verkhovtsev et al 2015), this requirement may no longer be fulfilled. However, as these surface effects are generally not included in the cross sections used in Monte Carlo simulations, they need to be considered as an additional contribution to the uncertainty of DEF values obtained from simulations.
The additional fluence contribution due to higher-order generations of scattered photons (also shown in figure 8) indicate that these photons would further increase the spectral fluence by a factor of three to four. This contribution is estimated by scaling the fluence of photons scattered once prior to reaching the NP with the ratio of the cumulative fluence of all generations of Compton photons and the fluence of the Compton photons produced by interactions of the primary photon beam. Owing to the large photon mean free paths, the different generations of photons may be produced in different spatial regions, so that the fluence distribution may significantly vary with distance from the source of primary photons. Although, even for the first generation of scattered photons, those reaching the NP are mostly backscattered photons. This is an obvious artefact of the simulation setup since in the energy range where Compton scattering dominates photon interaction in water, i.e. above 30 keV (see figure 2), the photon mean free path is greater than 2.5 cm (Berger et al 2010), which is much larger than the distance between the photon source and NP. Nevertheless, one may expect that the estimate of the fluence by multiply-scattered photons may be reasonably accurate within 10%. Considering this, together with the aforementioned uncertainty due to potential surface effects, the potential increase in the DEF due to scattered photons seen in table 3 is insignificant. Analysing the comparison exercise reported by Li et al (2019), which shows a large spread of different results obtained with different codes by different modellers, one could attribute an uncertainty in the DEF in the order of 60% for the deviation of the DEF from unity, particularly for the highest DEF values at distances from the NP surface in the nanometer range (Li et al 2019). This is to some extent related to the microscopic size of the scoring volumes (see section 4.3), and moreover reflects the state of the art of Monte Carlo simulation on the nanometer scale (Villagrasa et al 2019).
For the low-energy x-ray spectra considered in the EURADOS exercise, the results shown in figures 5 and 6, tables 1 and 2 clearly indicate that the DEFs derived by taking the ratio of results from simulations in a narrowbeam geometry significantly overestimate the real dose enhancement by orders of magnitude, particularly at distances from the NP exceeding 100 nm. Figure 7 demonstrates that such narrow-beam geometries tend to overestimate DEFs even beyond the theoretical upper limit (section 2.4) in this distance range. However, the values shown in table 3 indicate that the DEF in figure 6 is a good approximation of the expected DEF in the case of the full secondary photon spectra for the primary x-ray spectra.
To which extent the finding for the low-energy x-ray spectra can be transferred to higher photon energies, particularly to megavoltage photon beams from medical linear accelerators.
For MeV photons, the term in brackets in the formula for this upper limit of the DEF (equation (8)) is much smaller than that for the x-ray spectra considered here. In fact, while this factor (which multiplies the mass fraction of gold) amounts to about 150 and about 270 for the 50 kVp and 100 kVp spectra, respectively, it is only about 3.7 for an idealized Co-60 gamma spectrum and about 0.5 for the real Co-60 spectrum considered in figure 3(c). This implies that the upper limit for MeV photon radiation is about two orders of magnitude smaller than for the 100 kVp x-ray spectrum as is illustrated in figure 9 for the Co-60 spectrum where the upper limit of the average DEF from unity is by a factor of about 250 lower than for the 100 kVp spectrum. Even though the spectral fluence from the different generations of scattered photons converges slower than for the lower energy x-ray spectra considered in the simulations (see figure 3) and an enhancement factor for the deviation of the DEF from unity as large as 6.5 is found when the flueunce of higher order scattered phtotons is taken into account, the DEF for the Co-60 photon spectrum will still deviate by a factor of 40 less from unity than the DEF for the case of the 100 kVp x-ray spectrum. This implies that even though scattered photons result in a large enhancement of the extra energy deposition from electrons emerging from the NP, the overall DEF will still be significantly deviating from unity only in a small region of few 100 nm around the NP for the Co-60 spectrum, and it is quite plausible that this will also hold for other spectra of MeV photons. One reason for this is that the difference between the photon interaction coefficients of gold and water is much smaller in the MeV photon range than in the several 10 keV range (see figure 1). The other reason is that these higher-energy photon spectra produce higher energy secondary electrons that have ranges of the order of several 100 µm. Hence, their additional contribution to the absorbed dose is distributed over larger volumes, which pushes the estimated DEF closer to unity.
It is also worth noting that the upper limits shown in figure 7 do not require the NP to be located in the centre of the spherical sampling volume. As long as the conditions of SEE are fulfilled, they can also be applied for the case of a spatial distribution of NPs. Further discussion of this point is beyond the scope of this paper.

Microdosimetric aspects
As also illustrated in figure 7, the DEF will generally decrease if the radial bin size increases. This is a triviality, as the DEF is monotonically decreasing with increasing radial distance (figure 6) and averaging over a larger bin volume reduces the value in this bin.
In this context, it is important to note that the volumes of the first 10 nm spherical bins are approximately equal to that of spheres with 60 nm and 90 nm diameter for 50 nm-and 100 nm-diameter NPs, respectively. For bins extending to radial distances of 1 µm from the NP surfaces, the equivalent sphere diameters are 630 nm and 640 nm, respectively. This implies that the imparted energy per mass in these volumes, i.e. the specific energy z according to the terminology of microdosimetry (ICRU 1985), is a stochastic quantity, whose values have a statistical distribution reflecting the inchoate nature of radiation interaction. As is well known from the literature (Rossi andZaider 1996, Lindborg andWaker 2017), the expectation value of these distributions (the frequencymean specific energy z F ) for small target volumes is the absorbed dose, and the variance of these distributions depends both on the value of absorbed dose and the target volume dimensions.
An estimate of the standard deviation of the specific energy distribution, which would be interpreted in simulations as the numerical accuracy, can be obtained by considering the mean number of ionizations ν occurring in such a target volume, V t , for a value of absorbed dose, D, which can be obtained as Here W ≈ 25 eV (ICRU 1979) is the mean energy per ionization and ρ is the mass density of liquid water. For the first 10 nm thick spherical shells around the 50 nm and 100 nm NPs and an absorbed dose of 1 Gy, the mean number of ionizations are 0.27 and 0.91, respectively, if the NP volume is filled with water (i.e. NP absent) and a factor of about 10 to 14 (table 1) or 17 to 24 (table 2) higher for the NP present. For the 10 nm-thick spherical shell whose outer boundary is at 100 nm from the NP surface, the respective values would be about 4.3 and 6.3 without NP and by a factor of about 1.2 to 1.3 (table 1) and 1.9 to 2.3 (table 2) higher with NP. For the 10 nm-thick spherical shells whose outer boundary is at 1 µm from the NP surface, the values would be about 314 and 329 with and without NP.
Owing to the statistical correlations of the ionizations within a charged particle track, the statistics of the number of ionizations is a compound Poisson process. Nevertheless, the standard deviation for the case of a Poisson process, √ ν, may still be a good approximation for the standard deviation of the number of ionizations occurring in the different spherical shells and, hence, the relative standard deviation of the number of ionizations, s (ν) /ν , as well as of the specific energy, s (z) /z F , are both approximately given by Then, for the assumed average value of absorbed dose of 1 Gy, the relative standard deviation of the specific energy distribution would be higher than 100% for the first spherical shell without NP and of the order of several 10% with NP. This would also be the order of magnitude for the spherical bin ending at 100 nm from the NP surface with and without NP in place, while for the1 µm spherical bin the relative standard deviation is only a few percent.
Hence, if the numerical simulation is carried out such that a certain value of relative uncertainty u rel is achieved, one can estimate the value of absorbed dose for which the specific energy would have such a small relative standard deviation as For the first radial bin around the 50 nm and 100 nm NPs with a relative numerical uncertainty of 1%, the respective dose values are about 27 kGy and 8 kGy, respectively. For the 10 nm-thick spherical shell whose outer boundary is at 1 µm from the NP surface, the respective dose values are about 24 Gy and 23 Gy for the respective 50 nm and 100 nm NPs. It should be noted, however, that these values apply in the case of SEE. If the simulations are performed with the confined beam geometry, the absorbed doses corresponding to a relative numerical uncertainty of 1% may even be one to two orders of magnitude larger (see figure 5).
These considerations illustrate that a small numerical uncertainty can only be achieved when simulating an irradiation exposure for absorbed doses that are three or four orders of magnitude higher than what is generally used in therapeutic applications of ionizing radiations, or two to three orders higher than doses typically used in diagnostic applications. For realistic values of absorbed dose, the simulation results would have a standard deviation that may significantly exceed the mean values, such as the case in microdosimetry. Under such circumstances, the validity and usefulness of the concept of a DEF may be questionable.

Conclusions
In this paper we have demonstrated that radiation transport simulations performed with a narrow-beam geometry that does not assure SEE can result in a substantial overestimation of not only the spatial range around high-Z NPs where absorbed dose is significantly enhanced, but also the magnitude of the DEF. For the 50 kVp to 100 kVp x-ray spectra used in a recent EURADOS exercise, the maximum DEF in 10 nm spherical shells around spherical gold NPs of 50 nm and 100 nm diameter was between 10 and 20 rather than several hundred as suggested by published simulation results. Furthermore, significant dose enhancement only occurred within the first 100 nm to 200 nm around the gold NP rather than extending to several micrometers. This work also showed that scattered photons interacting with the gold NP do not significantly change these observations. It can be expected that these finding will also be the same for NPs of high-Z materials other than gold.
Owing to the intrinsic challenge of the simulations due to the small size of the NP and the required photon field sizes for assuring SEE, many DEFs reported in the literature may be overstimations arising from the lack of SEE due to the simulation setup. This work therefore presents an approach to remedy such deficiencies of the simulation setup in order to obtain more realistic estimates of the dose enhancement around high-Z NPs.
As a guideline for future assessment of DEFs, the following recommendations outline how to proceed with the analysis of simulation results in terms of quality assurance: 1. Verify that for the simulations without the NP, the absorbed dose values in the different microscopic scoring volumes agree with each other within statistical uncertainties (e.g. three times their standard deviations). 2. If the absorbed dose values obtained from the simulations reveal a systematic trend such as that in figure 5, the condition of SEE is most likely not fulfilled. The obtained DEF should therefore be checked against the theoretical upper limit as given in equation (6) for spherical NPs and spherical shells as scoring volumes. The transfer to other geometries of the NP and the scoring volumes is straightforward.One should keep in mind that close to the NP the upper limit is generally much larger than the true DEF. 3. To compensate for a lack of SEE, a more realistic estimate of absorbed dose to water in the absence of the NP can be obtained by equation (13). As a first order estimate, it may be sufficient to simply take the photon interaction data in the energy bin centers. This avoids the extra programming effort for the integration over the bin. For gold NP, the required interpolations for the photon interaction data can be easily implemented using the generalized spline functions given in the appendix. 4. Simulation results should preferentially not only be reported as absorbed dose per primary photon but as absorbed dose per fluence of primary particles, because this is a physically more meaningful quantity as any real source has a finite lateral extension in which case the source strength is characterized by the fluence of particles emitted and not by their total number. Furthermore, the value of fluence as well as the size of the photon field should be stated. In the case of ideal pencil beams and a spherical NP, one could, for example, normalize the number of particles to the geometrical cross section of the sphere. Furthermore, the approximate dose value corresponding to the simulated fluence should also be determined and stated.
CV performed the simulations using Geant4, and WL calculated the primary photon energy spectra for the EURADOS exercise. All co-authors carefully reviewed and edited the manuscript.

Appendix
For the procedures described in the sections 2.3 through 2.6, interpolating function were obtained for the photon energy dependence of the mass attenuation and mass energy absorption coefficients of gold and water. For the total mass attenuation coefficients and the mass energy absorption coefficients of water and gold as well as the attenuation coefficient for Compton scattering in water, µinc/ρ, the logarithm of the quantity of interest (in cm 2 /g) as a function of the logarithm of the photon energy (in keV) was fitted by generalized spline functions, i.e. by a set of polynomials defined in m subintervals of the considered abscissa range.
with x = log 10 (E/keV) and the conditions The interpolation for the quantity under consideration µ x is then given by For the mass energy absorption coefficients, the spline fit was performed on the data retrieved from the XMAADI data base in the energy range between E min = 1 keV and E max = 20 MeV. For the total mass attenuation coefficients of gold and water and the attenuation coefficients for coherent and incoherent scattering on water, data obtained from the XCOM data base were used and fitted in the base in the energy range between E min =1 keV and E max = 20 MeV. Cross checking showed that the data of the two data bases for µ/ρ agreed within 0.05% for energy values of the standard grid.
In general, continuity of the function and its first two derivates (slope and curvature) at the interval boundaries was required as boundary conditions as was zero curvature at the maximum considered energy. Using VBAprogrammed macros, the set of coefficients c jk minimizing the mean square deviation between the data points and the fit curve was calculated by solving the ensuing set of linear equations with boundary conditions. The total number of subintervals, their boundaries and the degrees n j of the polynomials in each interval was varied manually in an Excel 2016 spreadsheet such as to obtain a maximum deviation between the fit curve and each data point of below 0.02% (corresponding to a maximum deviation between fit curve and the data value of the quantity of interest of below 0.05%). This kind of deviation is much smaller than the estimated uncertainties of the photon interaction data in the data bases, which amount to a few percent (Andreo et al 2012).
In the case of the mass energy absorption coefficient of gold, the interval boundaries were fixed to the logarithms of the energies corresponding to the K, L and M absorption edges of gold, and the requirement of continuity of the function at the interval boundaries was abandoned to allow for the step-like rise of the energy absorption at these energies.
The parameters determining the different intervals and the coefficients obtained by this minimization procedure are listed in tables A1-A6. These parameters were applied in further data analysis by making use of the EXCEL built-in function POWERSERIES. As pair production is not possible in the photon energy range of the simulations described in section 2.1, the mass attenuation coefficient for photoeffect could be simply obtained as the difference between the total and the sum of the coefficients for coherent and incoherent scattering.
To evaluate the integrals over the energy bins in equations (13), (22) and (28), the generalized spline interpolations with the parameters given in this appendix were used to calculate the quantities µ/ρ, µ en /ρ, µ inc /ρ and µ R /ρ at the support points for the numerical integration which was based on the Newton-Cotes 3/8 rule.