Open Access
How to translate text using browser tools
21 October 2015 A Simple Model to Quantify Radiolytic Production following Electron Emission from Heavy-Atom Nanoparticles Irradiated in Liquid Suspensions
Nathan Wardlow, Chris Polin, Balder Villagomez-Bernabe, Fred Currell
Author Affiliations +
Abstract

We present a simple model for a component of the radiolytic production of any chemical species due to electron emission from irradiated nanoparticles (NPs) in a liquid environment, provided the expression for the G value for product formation is known and is reasonably well characterized by a linear dependence on beam energy. This model takes nanoparticle size, composition, density and a number of other readily available parameters (such as X-ray and electron attenuation data) as inputs and therefore allows for the ready determination of this contribution. Several approximations are used, thus this model provides an upper limit to the yield of chemical species due to electron emission, rather than a distinct value, and this upper limit is compared with experimental results. After the general model is developed we provide details of its application to the generation of HO through irradiation of gold nanoparticles (AuNPs), a potentially important process in nanoparticle-based enhancement of radiotherapy. This model has been constructed with the intention of making it accessible to other researchers who wish to estimate chemical yields through this process, and is shown to be applicable to NPs of single elements and mixtures. The model can be applied without the need to develop additional skills (such as using a Monte Carlo toolkit), providing a fast and straightforward method of estimating chemical yields. A simple framework for determining the HO yield for different NP sizes at constant NP concentration and initial photon energy is also presented.

INTRODUCTION

It has been suggested that the chemical processes (particularly the creation of radiolysis products) in aqueous environments containing nanoparticles, are important to several areas of interest, including energy production, nuclear waste processing, radiation chemistry, chemical synthesis, radiotherapy, catalysis, sensing, nanotoxicity and nanomedicine (1).

In radiotherapy, the use of high-energy photons is one conventional approach in cancer treatment. The aim is for these photons to cause the excitation of electrons, which will interact, directly or indirectly, with the DNA in cancer cells and cause damage. The HO is responsible for between 50 and 70% of this damage (2).

Nanoparticles (NPs) that include heavy elements (gold, hafnium, etc.) in their composition have been shown to drastically increase radiation damage (39). The collective effect of several electrons emitted after a single inner shell ionization of one of the heavy elements contained in such nanoparticles has been modeled to explain the large increase in biological effect observed in vitro (10, 11). However, there are reasons to believe other mechanisms might be at work (1, 12, 13).

Very high yields of HO observed in radiolytic experiments (1, 12, 13) might contribute to the additional therapeutic efficacy (or “dose enhancement”) observed in vitro (14, 15). Using monoenergetic synchrotron radiation and the arguments about energetics, it was shown (13) that there must be an additional mechanism able to produce large yields of HO, possibly through radiolytic products created in the surrounding media diffusing to the nanoparticle where HO is produced.

To better evaluate the contributions to this type of chemical production, and also to aid other researchers in analyzing further measurements in the future, it is helpful to have a readily accessible model of the contribution from the electron emission process. The aim of this work was to produce a simple expression for chemical product yield due to the electrons produced by X-irradiated nanoparticles. As mentioned above, very high yields of HO have been observed experimentally, well in excess of that provided by the understood electron emission process. In this work a parameterization was developed, which can be used to quickly determine an approximate maximum value of chemical production by the electron emission process; any chemical production in excess of this must therefore be produced by an alternative process. In determining this maximum value, the result of this work will place an upper limit (or “upper bound”) on the electron emission process's chemical production. Thus, at any point in this treatment where there is uncertainty in the calculation method or where approximations are being made, the treatment errs on the side of overestimation, giving a larger yield than more accurate methods would produce. As a result, if there is still additional chemical product observed experimentally, then there can be no doubt that the excess originated from an alternative production pathway.

Given that X rays are a common medical imaging diagnostic tool, as the field of nanomedicine matures it will become increasingly important to determine any chemical effects of nanoagents exposed to ionizing radiation (16). Furthermore, these types of interaction have been linked with a large range of other applications (1).

To facilitate the greatest impact of the model described here, we keep the method general, both in terms of nanoparticle composition and chemical product.

Chemical Production Pathways

To clarify this idea, Fig. 1 shows two pathways for the production of an arbitrary chemical product by radiation exposure of a nanoparticle in some media. Pathway A is the electron emission process, where the radiation interacts with the nanoparticle to produce electrons, which goes on to interact with the surrounding media to produce the desired chemical product. Pathway B is radiolysis of the media, where the radiation directly interacts with the media to produce the desired product (or some intermediate product which itself reacts with the media to produce the desired product).

FIG. 1.

Production pathway for general NP/media radiolysis. The two reaction networks represent the G values for the production of an arbitrary product, P, in the medium, M, by electrons and photons, respectively.

i0033-7587-184-5-518-f01.tif

To outline a specific case addressed by the model presented, the different pathways able to produce HO in a NP and water solution are shown in Fig. 2.

FIG. 2.

The three HO production pathways for NP/water radiolysis.

i0033-7587-184-5-518-f02.tif

Until now, pathway A has been considered the main contributor for HO production in the presence of nanoparticles; again, this is the electron emission process, where the electrons produced by nanoparticle irradiation goes on to interact with the surrounding water to produce HO. Pathway B corresponds to water radiolysis and is part of conventional radiotherapy. The X rays interact with the water, resulting in HO via a complex network of reactions (not shown) (17). Pathway C represents some other interaction with the nanoparticle, where a different product of water radiolysis goes on to interact with the NP, eventually producing HO. A prime candidate is H2O2, which has been noted to have its decomposition catalyzed at nanoparticle surfaces (18). See sections 2.2 and 3.3 of Sicard-Roselli et al. (13) for a discussion of other potential candidates for this pathway for the particular case of gold nanoparticles (AuNPs).

The existence of pathway C was detected by subtracting an upper bound for the yield of pathway A and a literature value for the yield of pathway B from the measured yield (13). Thus, a better determination of the yield of pathway A, replacing the upper bound used by Sicard-Roselli et al. (13), will facilitate a better estimate of the yield due to pathway C.

The average number of product molecules from a given radiation energy input to a reaction is referred to as a G value. The aim is to determine a (relatively) simple expression for the G value of HO production in a system of known NP dimensions and X-ray energy.

THE MODEL

For illustrative purposes the model was developed with specific reference to production of HO radicals in a nanoparticle/water colloid (Appendix A). However, the treatment is completely general and can be applied to other media and products, provided the expression for the energy dependence of the G value for product formation by electrons is known and is relatively straightforward to express (see Appendix B). To determine the enhancement due to nanoparticles, first consider the amount of HO production due to pathway A (see Fig. 1) for a single irradiated NP. Initially, a model is developed for electron production inside a single NP of arbitrary dimensions surrounded by water to provide an expression for the number of HO produced per photon absorbed. This result is then used along with standard X-ray absorption data to develop means to determine the enhancement as a function of NP concentration.

There are several stages involved in the production of HO via pathway A. First, the photon is absorbed in the NP, producing a photoelectron. If this photoelectron is from an inner shell, the inner shell vacancy is stabilized through a mixture of radiative and Auger transitions with some number (possibly zero) of Auger electrons being emitted. The total energy of all the electrons emitted plus the energy of the photoelectron is less than the initial photon energy. Furthermore, each electron loses energy in transporting from the ionization site to the nanoparticle's edge. The remaining energy is used to drive the radiolytic reaction chain shown in Fig. 1.

As a first (somewhat simple) approximation that leads to an overestimation of the HO yield (as discussed in the Introduction), all of the emitted electrons are treated as a single high-energy photoelectron. This electron can be produced anywhere inside the NP volume, depending on where the photon was absorbed. The electron then has to travel some distance through the NP until it exits into the surrounding water, where it can interact to produce HO. On its outward journey, the electron will lose some of its energy through collisions with the electrons and nuclei inside the NP. After exiting the NP, some proportion of the electron's remaining energy is used to make HO, as dictated by the electron G value for HO production. If in fact multiple electrons are emitted (a combination of a photoelectron and one or more Auger electrons) then, on average, this ensemble will lose more energy while moving through the NP (as each electron is individually attenuated), leaving less energy to drive the production of HO. Thus, considering this ensemble of electrons as a single electron with an energy equal to the total energy of all the electrons in the ensemble provides an overestimate of the energy available to form HO. Since the HO production always increases with electron energy, this overestimate in turn leads to an overestimation in HO produced. The total combined magnitude of the overestimations used in this treatment can be up to the order of 40%, as will be shown in the Results section.

The purpose of this work was to piece these stages together in a useful way that is relevant to a range of nanoparticle materials and to develop a framework which can be used for similar analysis for alternative products, not just HO.

Electron Production and Attenuation

Consider the production of an electron by X-ray-induced photo processes somewhere inside a NP. To be able to produce HO by interacting with the surrounding water, this electron must first transport through the nanoparticle to some boundary. During this process it can lose energy. The loss of energy from electrons traveling through matter is reported by Cole (19) (see Appendix C):

i0033-7587-184-5-518-e01.gif

This expression can be used to determine the remaining electron energy after traveling along some path taking it to a distance rz (μm) from its starting point. Here E0 and E are the initial and final electron energies (keV) and ρe is the electron density of the NP material in electrons/m3. See Appendix C for the derivation of this expression. As with most methods of determining electron range or stopping power, this expression becomes inaccurate for low electron energies. A discussion of the accuracy of this method can be found in Appendix C. In regions where it loses accuracy, this method underestimates the electron energy loss and thus overestimates HO production, in keeping with the goal of providing an upper bound.

This electron is traced along all possible exit paths and an “exit energy” spectrum is found; this is effectively the solution of a 5-dimensional integral problem (3D for electron production point, another 2D for emission direction). This approach would lead to an unwieldy expression so a simplification in terms of the available symmetry is sought.

The easiest way to simplify this is to consider one point on the NP surface as an electron exit point, then step into the NP in equidistant shells from the surface point. Electrons “born” at any point in a shell will have the same exit energy if they move to the emission point (see 2D illustration in Fig. 3).The symmetry of this system enables the quick consideration of the full angular distribution of electrons. Every possible electron path is shown in Fig. 3: any initial electron position and direction has an analogous starting point and path towards the exit point by rotational symmetry. Due to the innate symmetry in the system, a consideration of all photoelectron position depth shells from one static emission point will give the same results as any other point on the surface, to a very good approximation. The approximation made here is because some of the electrons originating at the extremes of the equi-depth meniscus shown in Fig. 3 that have been counted as leaving at the specified emission point will instead leave earlier, due to their random walk between their origin and the emission point straddling the NP edge. However, this effect is very small except at very low energies. Taking a straight-line approximation here will underestimate the electron energy loss, thus overestimating the HO production, in keeping with determining an upper limit to HO yield. The overall effect of the random walks is included in the range expression, so these approximations are not outrageous. Using the above approach, the emitted photoelectron energy spectrum can be calculated for a NP of known size and material, with known X-ray energy.

FIG. 3.

Two-dimensional illustration of equi-distant initial photoelectron position shell.

i0033-7587-184-5-518-f03.tif

The volume of the shell is directly correlated to the probability of a photoelectron being emitted with the given attenuated energy for that shell; if the X ray is going to be absorbed by the NP then each location within the nanoparticle will have an equal probability of being the point where the photon is absorbed (since photon attenuation by the nanoparticle is negligible).

The above describes the treatment for a single photoelectron. In reality, X rays produce electrons inside nanoparticles, with a range of initial energies (see Fig. 4). Each of these electrons can be considered in turn in the model to give a more accurate idea of HO production, or a simpler system where the X ray gives all its energy to a single photoelectron can be considered to quickly give an upper bound to the yield with significantly less effort. It should be noted that this upper bound will be considerably less than the one used by Sicard-Roselli et al. (13) since the treatment presented here considers both the energy lost due to transport through the nanoparticle and the G value for the production of HO in water by electrons.

FIG. 4.

Electron production spectrum from gold irradiated by 40 keV X rays, generated in Geant4 using the Livermore Physics list for 10,000,000 histories. The large 10 keV+ peaks are photoelectrons, the lower energies are Auger electrons.

i0033-7587-184-5-518-f04.tif

HO Production

The amount of HO produced by electrons traveling through water is known (20, 21). The expressions used in the model for the production of HO were obtained from fitting the G value at nanoseconds after production from Fig. 6 of Yamaguchi et al. (20) and can be found in Appendix A. There are four expressions in total. One of these expressions is shown below in Eq. (2).

For electron energies greater than 0.5 keV but less than 25 keV,

i0033-7587-184-5-518-e02.gif
where i0033-7587-184-5-518-ilm01.gif is the number of HO molecules produced per 100 eV of electron energy (E).

Note again at this point that this description is also applicable to other chemical products, provided their production can be described by a relatively straightforward expression for the G value as a function of energy; one simply replaces the expression used in Eq. (2) above.

The photoelectron energy spectrum of electrons leaving the nanoparticle is known from the treatment presented in the “Electron Production and Attenuation” section.

Combining these factors with the shell volume/probabilities described above allows for the calculation of the number of HO output for a range of initial electron energies. This is done at each point inside the NP, leading to an “average HO output” per electron of initial energy E (also known as a G value: often expressed in nM/J or molecules/100 eV). See Appendix D for further explanation of these steps.

NP Material and Size Dependence

To facilitate use in the general case, the HO yields described above were calculated for several different NP sizes and materials. The important material characteristic is the electron thickness (i.e., the product of the electron density and the distance traveled through the NP. See Appendix C.), as X rays are primarily attenuated by interaction with electrons.Electron density can be calculated using Eq. (A6) (Appendix C). This equation can be used to calculate the electron density of any linear combination of elements, so even complex compound-nanoparticles can be used with this model.

Application: Approximations and Varying Degrees of Accuracy

There are several ways to implement the parameterization once it has been determined, ranging from considering a single photoelectron to using the full electron production spectrum. These approximations produce a set of successively lower upper bounds on the HO production, but also require successively more effort to implement.

The simplest case is to assume an incident photon transfers all of its energy to a single electron, which is emitted without losing any energy i.e., assuming the electron was only just bound to the gold. In reality, many electrons are produced by a single photon, which causes more energy loss since multiple electrons are attenuated, so the single photoelectron approximation will give an upper bound to the number of HO produced.

The first approximation is the unattenuated single photoelectron approximation; this is for the single, unattenuated photoelectron described above, which has the same energy as the incoming photon. This gives the absolute maximum HO production, with no energy being lost due to attenuation through the NP. As no electron transport takes place, this simply uses the G value relations from Eq. (2).

The second approximation is the attenuated single photoelectron approximation, which assumes that a single photoelectron with the same energy as the incoming photon is produced somewhere within the NP, to which the parameterization [Eq. (3)] is applied to account for its loss of energy in passing through the particle. This is an improvement on the unattenuated single electron, but still neglects some of the energy loss due to energy division between many electrons and their subsequent transport.

The third approximation is the edge-crossing approximation, which uses a simplified version of the electron spectrum (not the full spectrum as shown in Fig. 4), applying the K-, L- and M-edge data (22). This is a general procedure that can be applied to any element using publicly available data and is described in Appendix E. In brief, for each photon energy, instead of using the full electron emission spectrum (e.g., Fig. 4), a simpler version of the spectrum is constructed based on publicly available data. To achieve this, the X-ray attenuation is treated as a set of continuous functions (straight lines on a log-linear scale) with sharp discontinuities at the NP material's X-ray absorption edges (see Fig. 5). Below each edge, the ionizations are treated as being from the previous shell, with a single photoelectron (initial photon energy-binding energy of the shell) and a single “super-Auger” electron with energy equal to the binding energy of the respective shell. The probability of each of these photoelectron/super-Auger pairs is weighted according to the competition between whichever shells have binding energy less than the initial X-ray energy. For example, for initial energies greater than 81 keV the M, L and K shells will all contribute to the simplified spectrum, resulting in three photoelectron/super-Auger pairs, i.e., the spectrum will contain six electrons of varying energies and ionization probabilities. Each of these electrons then has the parameterization [Eq. (3)] applied to it, providing a more accurate result than the previous two cases.

FIG. 5.

Illustration of the edge-crossing approximation using X-ray attenuation data. Here the simplified spectrum is calculated for a 40 keV photon.

i0033-7587-184-5-518-f05.tif

FIG. 6.

Emitted electron energy spectrum for 32.5 nm AuNP with 40 keV initial electron energy. Area under curve i.e., expected electron exit energy = 39.9978 keV.

i0033-7587-184-5-518-f06.tif

The final approximation uses the full electron spectrum: this applies the parameterization in Eq. (3) to the full electron emission spectrum (e.g., Fig. 4), and should be the most accurate. In this implementation, the emission spectrum was determined by Geant4 simulations. In these simulations, a phantom made of gold at very low density (about 10−4 g/cm3, equivalent to 10 nM of 32.5 nm diameter AuNPs) of size 1 × 1 × 10 mm was irradiated by 10 million mono-energetic 40 keV gamma rays using the Geant4 toolkit with the Livermore physics list. The de-excitation process was activated to take into account the Auger electrons produced in the phantom. The simulation calculates the energy spectrum of all electrons created in the phantom [both photoelectrons and Auger electrons (Fig. 4)]. All secondary particles were killed after their creation, with only their initial energies being retained to produce the final spectrum, producing a “clean”, just-created energy spectrum of electrons before transport through the nanoparticle has been considered. This enables the fresh spectrum to be “injected” into whichever size of nanoparticle is desired by using the parameterization [Eq. (3)] with each electron energy from the spectrum. Note that fresh spectra would need to be computed for each new NP material under consideration, however, if the earlier approximations are reasonably close to this method, they will be suitable as standalone methods of quickly calculating an upper bound for HO yield. The simulations were performed using CentOS 6.5 64-bit machines with twin 8-core 2.5 GHz Xeon CPUs and 64 Gb of RAM, and took between 5 and 25 min to run, depending on photon energy. Computing the HO yield from an input electron spectrum with the parameterization takes a few seconds, whereas running the electron transport through the NP in Geant4 takes significantly longer than the spectra generation for each individual nanoparticle size. Completing a study of 100 sizes would take many days using Geant4, while the parameterization gives an upper bound (with reasonable accuracy) in a few minutes.

Recall that these terms are arranged in order of accuracy, and each improvement in accuracy should provide a successively lower value for the upper bound.

It is noteworthy that all of these approximations ignore the interaction of the outgoing electron(s) with the plasmonic structure as they leave the nanoparticle. However, on average this will give as much energy to the electron as it removes and hence should have a negligible effect on the chemical yield of any products in the liquid phase.

RESULTS

Application to Nanoparticles

Some results of the model for the particular case of HO radical generation by gold nanoparticles are shown below. The energy spectrum of electrons exiting the NP is shown in Fig. 6, for initial electron energy of 40 keV in a 32.5 nm diameter AuNP. For the majority of energies and material densities there is very little attenuation, as the nanoparticles are so small.

Note that due to the y-axis being probability/photoionization, the integral of the product of this probability multiplied by the energy is equivalent to the expected value of the exit energy of the photoelectron. In the majority of cases this is very close to the initial electron energy due to there being very little attenuation.

The dependence of HO radical production on initial electron energy is shown in Fig. 7. For each energy shown in Fig. 7, an electron spectrum such as those shown in Fig. 6 was calculated. Each of these spectra was treated according to the process described in Appendix D to give the number of HO radicals per photoionization, shown in Fig. 7A. This is then converted into a G value by dividing by the initial photon energy (which is equal to the electron energy at this level of approximation), shown in Fig. 7B. The code used to implement this model is available in a previously published report (23).

FIG. 7.

Panel A: HO yield energy dependence up to 40 keV for a 32.5 nm AuNP. Panel B: HO G value energy dependence up to 40 keV for a 32.5 nm AuNP.

i0033-7587-184-5-518-f07.tif

Electron Transport Contribution

As shown in Fig. 6, there is very little energy loss for energetic photoelectrons as they travel through the nanoparticle. As previously mentioned, the model becomes inaccurate at very low-electron energies and very high material densities. The model currently overestimates electron penetration at these extremes, so is in keeping with the goal of establishing an upper limit for HO production.

Electron transport through nanoparticles has recently been studied elsewhere [e.g., Casta et al. (24)]. Unfortunately however, since the published results are in an extreme low-energy regimen, direct comparison is not possible. This further underscores the importance of future work to improve the model's accuracy at both extremely low energies and high material densities.

General Parameterization

Each part of Figs. 6 and 7 represent a single NP size and electron density. To obtain as general a solution as possible, similar spectra were calculated for many different NP sizes and electron densities, as described in Appendix D. These data points were then used to produce a polynomial fit. The purpose of this was to provide a single equation to allow quick estimation of HO yields without running through the detailed calculations above.

Due to the complexity of the data produced, a very large order polynomial would be required to fit all of the data. Alternatively, splitting the data in half yields two parameterizations, one for each set of energies on either side of 2.5 keV. These parameterizations both have 15 coefficients and extend to the 4th order in NP diameter, photon energy and electron density. It is noteworthy that the most important terms all have NP diameter and electron density raised to the same power, so it can be calculated as a polynomial in two variables: E and Dρe. This reduction of variables follows directly from the observation that the product Dρe is a measure of the “thickness” of electrons passed through by the photoelectron on its way out of the nanoparticle. The parameterizations are shown below.

For energies less then 2.5 keV,

i0033-7587-184-5-518-e03.gif
and for energies greater than 2.5 keV,
i0033-7587-184-5-518-e04.gif
where N is the number of molecules of HO, D is the NP diameter in nm, E is the initial photoelectron energy in keV and ρe is the NP electron density in electrons/m3. When Eq. (3) is referred to in the text, Eq. (3.1) should be used for photoelectron energy less than 2.5 keV and Eq. (3.2) when it is greater than 2.5 keV. The resulting fit is a very good approximation to the HO production given by the model. Table 1 shows a comparison between the model described above and its parameterization [as given in Eq. (3)] for various initial electron energies, nanoparticle sizes and nanoparticle densities. None of these condition sets were used to determine the parameters of Eq. (3). These results clearly show that in most cases the parameterization is correct to within 2%, extending to 10% in certain extreme conditions such as low energies. A custom calculator using the above parameterization is available (23).

TABLE 1

Comparison of the Number of HO Molecules Produced for Various Conditions

i0033-7587-184-5-518-t01.eps

For energies above ~50 keV, Eq. (3) is very slowly converging and significant fractional errors result. This can be seen by looking at the trend in the terms that are independent of the product Dρe, for each successive term the coefficient decreases by a factor of around 100, where the power of E increases by 1. This means each successive term is still around E/100 = 50/100 = 50% as large as the previous term, so the equation converges very slowly. For lower energies this factor is significantly reduced, so the fractional errors are much smaller. This is not a concern for the Dρe terms as their contribution is minimal.

Using the Parameterizations with the Various Electron Spectrum Approximations

Table 2 shows the HO production per photon for several photon energies in 32.5 nm AuNPs, using approximations of varying accuracy.

TABLE 2

HO Production per Photon for 32.5 nm AuNPs, Using Approximations of Varying Accuracy

i0033-7587-184-5-518-t02.eps

These show that the attenuated single photoelectron approximation overestimates HO production by around 30% in each instance compared to the full electron spectrum method; this approximation is then of the correct order of magnitude and so is acceptable for use to quickly establish the upper bound for pathway A.

Concentration Dependence

The above discussion has only considered a single X ray being absorbed by a single nanoparticle. In reality, there are many nanoparticles in the solution and there is only a small chance of the photon being absorbed by them. If the photon does not interact with the NPs at all, then pathway B (see Fig. 1) is the only pathway in effect. As such, the HO production enhancement comes from the small proportion of photons that are absorbed by the NPs. A “wrapper” needs to be put on the model described above to include the consideration of this.

Calculating the ratio of weights of NP to water enables the determination of the fraction of energy absorbed by the NPs. To determine the mass of NP material within one liter (l) of solution, see Appendix F.

From Appendix F,

i0033-7587-184-5-518-e05.gif
where cNP is the molar concentration of NPs in M/l, the NP diameter D is in nm and the density ρNP is in kg/m3.

As there are 1,000 g of water in one liter, the weight-to-weight ratio of NP to water in solution is given by

i0033-7587-184-5-518-e06.gif

Using this with the mass-energy attenuation coefficients allows the determination of the fraction absorbed by the NPs.

i0033-7587-184-5-518-e07.gif
where i0033-7587-184-5-518-ilm02.gif accounts for the approximate photo-ionization cross-section for the most recently open shell, described in the “Edge-Crossing Approximation” section.

The above f will be required for each given energy. Combining this with the calculated or measured G values for water alone and for NPs allows the determination of the relative yields due to each component.

Nanoparticles Composed of Mixtures of Elements

If the nanoparticle is a mixture of two or more elements, obtain i0033-7587-184-5-518-ilm03.gif for each component, multiply by their respective mass fractions then sum them:

i0033-7587-184-5-518-e09.gif
where ρ is the partial density. Use this summation in place of i0033-7587-184-5-518-ilm04.gif from before.

See Appendix G for a working example for 32.5 nm gold nanoparticles.

NP Size Dependence

Combining the concentration wrapper above with the parameterizations in Eq. (3) allows the evaluation of the HO yield for a certain concentration of a certain nanoparticle at a certain initial photoelectron energy, as shown in Appendix G. As the effect of photoelectron energy has already been shown in Table 2, and NP concentration has a linear effect on HO yield [Eqs. (5) and (6)], it is now of interest to study the effect that NP size has on HO yields.

From Eqs. (5) and (6),

i0033-7587-184-5-518-e10.gif
which can be combined with the individual NP HO yields, N(E, Dρe) from Eq. (3) to give the total HO yield from all NPs:
i0033-7587-184-5-518-e11.gif

Figure 8 shows the HO yield from 1 nM concentration of several sizes of gold NPs, comparing the three main approximations detailed above and in Table 2. The three approximations display very similar trends. Figure 9 compares the single attenuated photoelectron approximation for gold and silver NPs on a log–log scale, clearly showing the relationship between NP size and HO yield. This is discussed further below.

FIG. 8.

HO yield of all three approximations for initial photoelectron energy of 20 keV incident on 1 nM of gold NPs.

i0033-7587-184-5-518-f08.tif

FIG. 9.

HO yield of the single unattenuated photoelectron approximations for initial photoelectron energy of 20 keV incident on 1 nM of gold and silver NPs.

i0033-7587-184-5-518-f09.tif

DISCUSSION

This theoretical model has been developed from a series of previously reported relationships, all of which were based on experimental results: electron transport (19), G value of HO (20) and the intersection of two spheres, which was first solved by the early Greek mathematicians. As such, the mathematical models used in this work include many physical parameters not expressly stated within the formulae given – this does not mean the physical parameters are forgotten; rather they are factored in to the constants contained in each formula.

We used these results from the literature since up until now there have been no reports of measurements that isolate each individual mechanism and observe it in absence of the others. The data/formulae produced have been treated as phenomenological and analyzed as such, to obtain descriptions of the relationships involved and thus determine the above fits to the observed trends.

All of the results discussed here give an upper bound for the HO produced via this mechanism – we have not included HO – HO recombination or other scavenging effects, which would further reduce the yield via this process.

A trend comparison of the results of the model for gold, the maximum possible energy from conservation of energy (when the photoelectron is produced with maximum initial energy on the GNP surface and immediately escapes) and some experimental observations are shown in Fig. 10, which illustrates how these results pull pathway A's upper bound down. See Sicard-Roselli et al. (13) for further discussion.

FIG. 10.

Comparison of model results with observed experimental results (32.5 nm AuNPs) showing how the improved estimate of the component due to pathway A in turn leads to an improved estimate of the yield due to pathway C. This figure is used with permission and derived from Fig. 6 of Sicard-Rosellin et al. (13), and has been modified to highlight the magnitude of the contributions of the various pathways.

i0033-7587-184-5-518-f10.tif

Here we show the use of the model in providing a bound for the emitted electron-driven HO production through a re-analysis of the data of Sicard-Roselli et al. (13), using our new model, which considers 20 keV X rays producing electrons in a solution of 32.5 nm AuNPs at 1 nM concentration. Pathway B as reported by Sicard-Roselli et al. (13) gives 386 HO per 20 keV X-ray photon. At this concentration, the photon absorption efficiency by NPs is 0.025.

To create one HO by simply breaking water bonds costs 5 eV. The absolute maximum number of HO produced would be by a single electron produced by ionization on the edge of the NP by a 20 keV photon. This gives 4,000 HO per photon absorbed by a NP, or 100 additional (486 total) HO per photon in the system [the upper bound (Fig. 10)]. Using the model we obtain 936 HO per photon absorbed by a NP (24 extra, 410 total, HO per photon in the system). Note that the new estimate of the upper bound for the emitted electron contribution is about a factor of 8 lower than the previous estimate.

Experimentally, 854 HO per photon absorbed in the system is observed in these conditions, so using our improved model, we can say that at least 444 HO per photon absorbed are produced by pathway C, i.e., the ratio of contributions of pathways C:A is about 18.5 whereas the less sophisticated (larger) upper bound provided by simple energy conservation arguments for the ratio of contributions is only 4.5. It is clear that while pathway A is providing only a small amount of the observed HO, it must be precisely modeled to estimate the relative contribution of other pathways accurately.

Comparing with Experiment: Silver Nanoparticles

It is clear from the data shown in Fig. 10 and the discussion that there is significant HO production enhancement by pathway C for gold NPs. Silver nanoparticles have been briefly studied to determine whether these gold nanoparticles have some special property that other materials do not have.

Experimental Method

The silver NPs used are the commercially available 30 nm Biopure citrate coated spherical silver nanoparticles (product no. AGCB30-25M; nanoComposix Europe, Prague, Czech Republic) (28). These were washed by three cycles of centrifugation to remove most of the citrate, following the method published elsewhere (13).

The samples were irradiated at Beamline B16 (Diamond Light Source, Didcot, UK) by an 11.8 keV X-ray beam at 1 Gy/s. Beam size was typically 5 × 5 mm. Sample delivery and collection were performed using the hanging drip method described previously (13). Fluorescence analysis was performed using a microplate reader, also described by Sicard-Roselli et al. (13).

Analysis

Figure 11 shows some preliminary results of HO yields for 32.5 nm silver nanoparticles irradiated with 11.8 keV X rays at the Diamond Light Source, and clearly demonstrates that the trend observed in gold is not present, and that the upper bound from the characterization [Eq. (3)] falls roughly within the error bars for the data; it appears that pathways A and B provide the observed HO yield in the presence of silver nanoparticles. Again, the HO production due to pathway A is quite small at 4 nM of NPs pathway A provides a 2.3% increase in yield over pathway B alone.

FIG. 11.

Comparison of model results with observed experimental results (32.5 nm AgNPs) showing that the estimate of the component due to pathway A provides all the additional HO. Silver NPs were irradiated with 11.8 keV X rays on Beamline B16 at Diamond Light Source.

i0033-7587-184-5-518-f11.tif

Gold Versus Silver

Since this model quantifies pathway A alone, it cannot be validated in full, although it is clear from experiment that pathway A and pathway B together account for the HO production observed in silver. It has been shown that for gold nanoparticles, a large proportion of the HO yield is not accounted for by pathway A and pathway B together. In cases such as this the model presented here is useful to put a lower bound on the production from pathway C.

These results suggest there is more to pathway C than a Z dependence. Further investigation of the silver nanoparticles is required, along with other materials. Future work will involve studies of hafnium oxide nanoparticles to determine whether gold or silver are anomalies and the development of a diffusion model to attempt to mechanistically define pathway C.

NP Size Dependence

Figures 8 and 9 show a clear cubic relationship between NP size and HO yield for a given NP material and initial photoelectron energy. Therefore, Eq. (8) needs to be evaluated only once to enable the determination of the HO yield from pathway A for many NP sizes to a very good approximation. See Appendix H for an example calculation. From Appendix H, for 1 nM of AuNPs with initial photoelectron energy of 20 keV,

i0033-7587-184-5-518-e12.gif
where D is the NP diameter in nm.

By simply calculating the multiplication constant, an equivalent expression can be easily obtained for any photoelectron energy, NP concentration and NP material, saving multiple evaluations of Eq. (8). Once this constant is known for a given scenario, the HO yield of NPs can be quickly estimated for any NP size by evaluating the equivalent of Eq. (9).

CONCLUSION

In summary, a method has been found to quickly calculate an estimate of the amount of HO produced by the electron-driven component of radiation-nanoparticle interaction, which is sufficiently accurate to allow a ready reckoning of this contribution. The size dependence on the HO yield from the electron-driven component was briefly investigated, which is shown to be cubic in nature, thus it enables quick calculation of the yields for other sizes using Eqs. (8) and (9), as shown in Appendix H. From comparison with experimental results for gold nanoparticles, it is clear that in some cases this electron-driven pathway is not responsible for the majority of the observed HO production enhancement, confirming the existence of a new pathway (13). In other studies (Fig. 11) for different nanoparticles, the entirety of the enhancement appears to be accounted for by the electron-driven component. Further work is needed to determine when and why the new pathway becomes important.

Appendices

APPENDIX A

Expressions for the G Value of HO

The following expressions were obtained by fitting to the (near-) instantaneous G value data shown in Fig. 6 of Uehara and Nikjoo (21). Due to the complex nature of the relationships, several expressions are required for various energy ranges.

For electron energies less than 0.5 keV,

i0033-7587-184-5-518-e13.gif

For electron energies greater than 0.5 keV but less than 25 keV,

i0033-7587-184-5-518-e14.gif

For electron energies greater than 25 keV but less than 300 keV,

i0033-7587-184-5-518-e15.gif

For electron energies greater than 300 keV,

i0033-7587-184-5-518-e16.gif
where i0033-7587-184-5-518-ilm05.gif is the number of HO molecules produced per 100 eV of electron energy (E).

APPENDIX B

Formulation for General Product

The model primarily determines the electron production within a nanoparticle of any material, eventually providing an output electron spectrum, where the electrons have been attenuated through the nanoparticle volume and escaped with a known energy.

The G value for HO used throughout this work has a straightforward relationship with electron energy and so can be readily applied to this output electron spectrum.

Similarly, other well known, easily characterized G values could be applied at this point, requiring a minimal change to the code used to generate the characterization data. More complex G values can also be used, but it will take more computer time to generate the characterization data under these circumstances.

APPENDIX C

Electron Density Dependence

From the previously published data of Cole (20) we have:

i0033-7587-184-5-518-e17.gif

Here, r is the attenuation distance in μm and ΔE is the attenuated energy in keV (collodion is a plastic used by Cole).

For different materials, r depends linearly on the electron density. The more electrons, the shorter the distance required to attenuate the initial electron by the given energy. Thus, we can determine r for any nanoparticle material by comparing the electron density of the material with that of collodion. Collodion's electron density is 1/1.106 of the electron density of water.

i0033-7587-184-5-518-e18.gif
where ρe is the electron density.

Generalizing for any material, we then have

i0033-7587-184-5-518-e19.gif
where Z denotes the atomic number dependence of the electron density.

Accuracy of Electron Transport Method

The accuracy of electron range calculations is inherently uncertain at low energies. Checking the calculations from Eq. (1) against NIST's ESTAR database (26) gives good agreement in the energy ranges for which ESTAR is accurate. For low-Z materials this method is accurate to within a few percentages down to 5 keV electron energies and overestimates the range (within an order of magnitude) at lower energies. For higher Z materials, this method agrees with NIST's data to within a few percentages above 20 keV. At lower energies it again overestimates the range (within an order of magnitude). For these lower, uncertain energies, the range is so small (e.g., 4.3 nm for a 1 keV electron in gold) that the fraction escaping from a nanoparticle of clinically relevant size is miniscule. These overestimates of the range will provide more electrons exiting the NP than would be seen in practice, in keeping with the “upper bound approximation” aim described above. As the electron density includes Z dependence, and the majority of electron energy loss is via collisions with electrons in the nanoparticle, the electron density contains all required information for a good approximation of electron penetration range.

Electron Density Calculation

The electron density of a material is calculated as follows:

i0033-7587-184-5-518-e20.gif
where ρ is the density in kg/m3, NA is Avogadro's constant (atoms/M), Z is the atomic number of the material (electrons/atom) and Ar is the relative atomic mass of the material in g/M. We multiply by 1,000 to convert from kg to g.

APPENDIX D

Steps in Determining Energy and Chemical Yield Spectra

From Appendix C,

i0033-7587-184-5-518-e21.gif
gives the distance rZ an electron is required to move through a material of electron density ρe(Z) to suffer an energy loss ΔE.

To obtain the emitted energy spectra, consider a single possible electron exit point, as shown in Fig. 3. Take incremental steps in the final emitted energy E from 0 up to the initial (beam) energy. The attenuation distance rZ is calculated for each of these, using ΔE = E0E. This gives a data set of emitted energies and the depth into the NP at which an electron must be “born” to reach the exit point and escape the NP.

Each of these depths rZ corresponds to a section of a shell within the NP, centered on the exit point, as shown in Fig. 3. Due to the innate symmetry in the system, and as the full set of shells covers the entire nanoparticle volume, the individual volumes of these shells are analogues of the probability that an electron born somewhere within the NP will have the emission energy related to that specific shell. Calculating the (normalized) volumes of each energy/depth shell will then result in a probability spectrum for final electron emission energy, first shown in Fig. 5. Note that the volumes of these shells are dependent on the NP diameter, so these spectra will vary with NP size, electron densities and initial photoelectron energy.

Combining this probability spectrum with the G value relationship for the chemical product [in this case i0033-7587-184-5-518-ilm06.gif for hydroxyl] gives the average yield of product for photoelectrons of initial energy E0 in the specified NP.

When this procedure is performed for many initial photoelectron energies it will result in a production spectrum (versus energy) for the NP in question (shown in Fig. 6).

Data Collection

To obtain all the data used to produce the parameterization [Eq. (3)], the spectra described above were calculated for many NP sizes and electron densities.

For a given NP size, electron density was varied in incremental steps from 1.5 × 1029 electrons/m3 to 8 × 1030 electrons/m3, a range spanning all the solid elements, and the HO production spectrum calculated as described above for each electron density.

This entire process was then carried out for NP sizes in 5 nm steps from 5 to 300 nm, resulting in a large data set in three variables, NP diameter, NP electron density and initial photoelectron energy, each of which increases in steps of equal size to ensure the discrete population is as representative as possible. The parameterizations in Eq. (3) were calculated from this data set.

APPENDIX E

Edge-Crossing Approximation

The mass attenuation data shows the cross-section that is the sum of all principle photon interactions (23). In this energy range, the dominant interaction is the photoelectric effect (27). The mass attenuation coefficient, when plotted on a log-linear (coefficient vs. energy) scale, shows a series of straight lines sloping downwards with discontinuities (absorption edges) found when a new orbital becomes available for photoionization. Below 11.9 keV photon energy, the ionization is treated as all being from the M shell with the result being a single photoelectron and a single Auger electron, the photoelectron being assigned an energy of 2.21 keV less than the photon energy (the binding energy of the M shell). The remaining energy (i.e., the binding energy) is then considered to manifest as a single Auger electron. The effect of transport of both of these electrons through the nanoparticle is then accounted for using Eq. (3). Since only 2 electrons are considered (any Auger cascade being represented by a single electron), the effect of energy loss during transport through the nanoparticle is overrepresented (but to a lesser degree than in previous approximations), resulting in an underestimate of the chemical yield.

Between the L and K shells there are considered to be two ionization processes in competition, one from the M and one from the L shell. The contribution from the M-shell ionization was determined by extrapolating the attenuation coefficient from below the L edge (this is easily achievable as the trend is a straight line on a log-linear graph). The difference between the extrapolated value and the tabulated value was then ascribed to L-shell ionization. Using these two values as weighting factors, the relative probability of M- and L-shell ionization was determined. In each case, two electrons were then considered, one being the photoelectron and one accounting for the Auger cascade as described for the sub-L-shell case above. Again, the transport though the nanoparticle for each electron is accounted for through the use of Eq. (3).

Above the L edge the procedure can be repeated considering three ionization processes in competition, from the M, L and K shells. Note that as the dependence of the attenuation coefficient between edges is a straight line with a downward slope on a log-linear graph, the small approximation introduced by the extrapolation process decreases as one goes to higher energy away from an edge, as does the neglect of the sub-structure in the L and M edges.

APPENDIX F

Nanoparticle Material Density in Solution

Determination of the number of NP atoms in one liter:

i0033-7587-184-5-518-e22.gif
where the NP diameter D is in nm;
i0033-7587-184-5-518-e23.gif
where ANP is the atomic mass number of the NP and the density ρ is in kg/m3.

To obtain the number of grams of NP material in one liter, the above is multiplied by i0033-7587-184-5-518-ilm07.gif, to give

i0033-7587-184-5-518-e24.gif
i0033-7587-184-5-518-e25.gif
i0033-7587-184-5-518-e26.gif

APPENDIX G

Example: 32.5 nm Gold Nanoparticles

The density of gold is 19,300 kg/m3.

For 20 keV photons, we have 835 HO/photon absorbed in the NP (full spectrum parameterization, Table 2). At this energy, we have 484 nM of HO/J in water (25).

At 20 keV, i0033-7587-184-5-518-ilm08.gif for gold is 65.22 cm2/g, i0033-7587-184-5-518-ilm09.gif for water is 0.5503 cm2/g (26).

From the Sicard-Roselli's experiment (13), 1 nM of gold NPs give 442.5 nM/J = 854 HO/photon. Then

i0033-7587-184-5-518-e27.gif

This shows that around 2.5% of the energy is absorbed in the gold NP material, with the rest being absorbed in the water. When combined with the G values, this gives a total HO yield of 407 HO per photon (i.e., 21 HO from gold on top of the 386 HO per photon in water alone). Compared with experimental results of 854 HO per photon (13), it is clear that this mechanism is not providing the observed HO enhancement. The new pathway suggested in (13) must provide the majority of the enhancement:

i0033-7587-184-5-518-e28.gif
i0033-7587-184-5-518-e29.gif

APPENDIX H

Example: 1 nM of Gold NPs with 20 keV Initial Photoelectron Energy

For 100 nm NPs, Eq. (8) evaluates to 674.055 HO from NPs. For 200 nm NPs, Eq. (8) evaluates to 5392.4401 HO from NPs. The equation of a line that becomes linear on a log–log plot is y = A · x b.

Taking the log of both sides: ln(y) = ln(A) + b · ln(x).

Evaluating for both sizes and subtracting gives

i0033-7587-184-5-518-e30.gif

b = 3 shows that the relationship between HO and NP diameter is cubic. Since y = A · x b, A = y/x b = 674.055/1003 = 6.74 × 10−4.

With these two constants, the HO produced by 1 nM of gold NPs of any size for this energy can be easily evaluated:

i0033-7587-184-5-518-e31.gif
where D is the NP diameter in nm.

ACKNOWLEDGMENTS

The authors acknowledge Diamond Light Source for time on Beamlines B16, I15 and the side laboratories no. 67 under proposals MT9104, EE8481 and MT8482. The authors also thank the teams from both beamlines and other support staff who made these measurements possible. Samples were also analyzed using equipment in the Research Complex at Harwell (Didcot, UK). Again, the authors thank the support staff who made this work possible. The authors also thank the Engineering and Physical Sciences Research Council for funding the corresponding author's research. The authors thank one of the referees of this manuscript who suggested a systematic study of NP size, as this led to the discussion and formulation of Eqs. (8) and (9).

REFERENCES

1.

NN Cheng Z Starkewolf RA Davidson A Sharmah C. Lee Chemical enhancement by nanoparticles under X-ray radiation. J Am Chem Soc 2012; 134:1950–3. Google Scholar

2.

CV. Sontag Free-radical-induced DNA damage and its repair. A chemical perspective. 1st ed. Berlin: Springer Verlag; 2010. Google Scholar

3.

KT Butterworth SJ McMahon FJ Currell KM. Prise Physical basis and biological mechanisms of gold nanoparticle radiosensitization. Nanoscale 2012; 4:4830–8. Google Scholar

4.

SJ McMahon MH Mendenhall S Jain FJ. Currell Radiotherapy in the presence of contrast agents: a general figure of merit and its application to gold nanoparticles. Phys Med Biol 2008; 53:5635–51. Google Scholar

5.

KT Butterworth JA Wyer M Brennan-Fournet CJ Latimer MB Shah FJ Currell et al . Variation of strand break yield for plasmid DNA irradiated with high-z metal nanoparticles. Radiat Res 2008; 170:381–7. Google Scholar

6.

Y Zheng L. Sanche Gold nanoparticles enhance DNA damage induced by anti-cancer drugs and radiation. Radiat Res 2009; 172:114–9. Google Scholar

7.

DB Chithrani S Jelveh F Jalali M van Prooijen C Allen RG Bristow et al . Gold nanoparticles as radiation sensitizers in cancer therapy. Radiat Res 2010; 173:719–28. Google Scholar

8.

WN Rahman N Bishara T Ackerly CF He P Jackson C Wong et al . Enhancement of radiation effects by gold nanoparticles for superficial radiation therapy. Nanomedicine 2009; 5:136–42. Google Scholar

9.

T Kong J Zeng X Wang X Yang J Yang S McQuarrie et al . Enhancement of radiation cytotoxicity in breast-cancer cells by localized attachment of gold nanoparticles. Small 2008; 4:1537–43. Google Scholar

10.

SJ McMahon WB Hyland MF Muir JA Coulter S Jain KT Butterworth et al . Biological consequences of nanoscale energy deposition near irradiated heavy atom nanoparticles. Nature Scientific Reports, 2013. ( http://bit.ly/1J2upQWGoogle Scholar

11.

SJ McMahon WB Hyland MF Muir JA Coulter S Jain KT Butterworth et al . Nanodosimetric effects of gold nanoparticles in megavoltage radiation therapy. Radiother Oncol 2011; 100:412–6. Google Scholar

12.

M Misawa J. Takahashi Generation of reactive oxygen species induced by gold nanoparticles under x-ray and UV irradiations. Nanomedicine 2011; 7:604–14. Google Scholar

13.

C Sicard-Roselli E Brun M Gilles G Baldacchino C Kelsey H McQuaid et al . A new mechanism for hydroxyl radical production in irradiated nanoparticle solutions. Small 2014; 10:3338–46. Google Scholar

14.

S Jain DG Hirst JM. O'Sullivan Gold nanoparticles as novel agents for cancer therapy. Br J Radiol 2012; 85:101–13. Google Scholar

15.

JF Hainfeld FA Dilmanian DN Slatkin HM. Smilowitz Radiotherapy enhancement with gold nanoparticles. J Pharm Pharmacol 2008; 60:977–85. Google Scholar

16.

ST Stern JB Hall LL Yu LJ Wood GF Paciotti L Tamarkin et al . Translational considerations for cancer nanomedicine. J Controlled Release 2010; 146:164–74. Google Scholar

17.

V. Arena Ionizing radiation and life. 1st ed. St. Louis: Mosby; 1971. Google Scholar

18.

W He YT Zhou WG Wamer X Hu X Wu et al . Intrinsic catalytic activity of Au nanoparticles with respect to hydrogen peroxide decomposition and superoxide scavenging. Biomaterials 2013; 34:765–73. Google Scholar

19.

A. Cole Absorption of 20-eV to 50,000 eV electron beams in air and plastic. Radiat Res 1969; 38:7–33. Google Scholar

20.

H Yamaguchi Y Uchihori N Yasuda M Takadaand H. Kitamura Estimation of yields of OH radicals in water irradiated by ionizing radiation. J Radiat Res 2005; 46:333–41. Google Scholar

21.

S Uehara H. Nikjoo Monte Carlo simulation of water radiolysis for low-energy charged particles. J Radiat Res. 2006; 47:69–81. Google Scholar

22.

JH Hubbell SM. Seltzer Tables of X-ray mass attenuation coefficients and mass energy absorption coefficients from 1 keV to 20 MeV for elements Z = 1 to 192 and 48 additional substances of dosimetric interest. NIST Report No. NISTIR 5632. National Institute of Standards and Technology Database 1989. ( http://www.nist.gov/pml/data/xraycoef/index.cfmGoogle Scholar

23.

Radiotherapy involving targeted nanoparticles. Centre for Plasma Physics, Queen's University Belfast. 2008. ( http://www.radphys.orgGoogle Scholar

24.

R Casta JP Champeaux P Cafarelli P Moretto-Capelle M. Sence Model for electron emission of high-Z radio-sensitizing nanoparticle irradiated by X-rays. J Nanopart Res 2014; 16:2480. Google Scholar

25.

JWT Spinks RJ. Woods An introduction to radiation chemistry. 2nd ed. New York: John Wiley & Sons; 1976. Google Scholar

26.

MJ Berger JS Coursey MA Zucker J. Chang Stopping-power and range tables for electrons, protons and helium ions. NIST, Physical Measurement Laboratory. ( http://1.usa.gov/1iiLhOBGoogle Scholar

27.

EB. Podgorsak Basic radiation physics. In: EB Podgorsak editor. Radiation oncology physics: a handbook for teachers and students. Vienna: IAEA; 2005. p. 1–7. Google Scholar

28.

nanoComposix Europe, 2015. ( http://bit.ly/1FBq6fgGoogle Scholar
©2015 by Radiation Research Society.
Nathan Wardlow, Chris Polin, Balder Villagomez-Bernabe, and Fred Currell "A Simple Model to Quantify Radiolytic Production following Electron Emission from Heavy-Atom Nanoparticles Irradiated in Liquid Suspensions," Radiation Research 184(5), 518-532, (21 October 2015). https://doi.org/10.1667/RR14059.1
Received: 24 February 2015; Accepted: 1 September 2015; Published: 21 October 2015
Back to Top