Complex mixed radiation fields exist in interplanetary space, and little is known about their late effects on space travelers. In silico synergy analysis default predictions are useful when planning relevant mixed-ion-beam experiments and interpreting their results. These predictions are based on individual dose-effect relationships (IDER) for each component of the mixed-ion beam, assuming no synergy or antagonism. For example, a default hypothesis of simple effect additivity has often been used throughout the study of biology. However, for more than a century pharmacologists interested in mixtures of therapeutic drugs have analyzed conceptual, mathematical and practical questions similar to those that arise when analyzing mixed radiation fields, and have shown that simple effect additivity often gives unreasonable predictions when the IDER are curvilinear. Various alternatives to simple effect additivity proposed in radiobiology, pharmacometrics, toxicology and other fields are also known to have important limitations. In this work, we analyze upcoming murine Harderian gland (HG) tumor prevalence mixed-beam experiments, using customized open-source software and published IDER from past single-ion experiments. The upcoming experiments will use acute irradiation and the mixed beam will include components of high atomic number and energy (HZE). We introduce a new alternative to simple effect additivity, “incremental effect additivity”, which is more suitable for the HG analysis and perhaps for other end points. We use incremental effect additivity to calculate default predictions for mixture dose-effect relationships, including 95% confidence intervals. We have drawn three main conclusions from this work. 1. It is important to supplement mixed-beam experiments with single-ion experiments, with matching end point(s), shielding and dose timing. 2. For HG tumorigenesis due to a mixed beam, simple effect additivity and incremental effect additivity sometimes give default predictions that are numerically close. However, if nontargeted effects are important and the mixed beam includes a number of different HZE components, simple effect additivity becomes unusable and another method is needed such as incremental effect additivity. 3. Eventually, synergy analysis default predictions of the effects of mixed radiation fields will be replaced by more mechanistic, biophysically-based predictions. However, optimizing synergy analyses is an important first step. If mixed-beam experiments indicate little synergy or antagonism, plans by NASA for further experiments and possible missions beyond low earth orbit will be substantially simplified.
INTRODUCTION
Background and Goals
Mixed radiation fields are important in radiobiology, especially where cancer risks from space travel are concerned [reviewed in ref. (1)]. A complex galactic cosmic ray (GCR) mixed radiation field is present outside of low earth orbit. This GCR field includes potentially significant contributions from high-linear energy transfer (LET) ions of high charge and energy (HZE), having high relative biological effectiveness (RBE) (2). Generally, the doses and dose-rates involved are low (3–5). However, chronic exposure over several years is worrisome, and effective astronaut shielding from many of the GCR components is not feasible [reviewed in refs. (5–7)].
Epidemiological data from high-LET exposures is as yet very sparse, especially at the relevant low doses and low dose rates. Indirect information on HZE risk, obtained from radiation studies of rodents, has therefore been considered [e.g., (8–13)]. One set of experiments, on murine Harderian gland (HG) tumorigenesis, has been influential in astronaut risk analysis [reviewed in ref. (6, 14, 15)]. These murine HG experiments started in the 1970s (16) and are still continuing (17).
Many of the murine HG experiments were performed with single-ion beams, that is, a single-ion isotope at a given initial energy enters the beam. The beam then inevitably encounters some material and thus is more complex at the HG (1, 5, 7). In some of these experiments the amount of material the beam encountered upstream of the target was increased by adding shielding material intended to modify the beam (18). However, all of our calculations will only use data from, and make predictions about, experiments where no such extra shielding is introduced.
Our group will perform corresponding mixed-beam murine HG experiments during the next several years starting in fall 2016, using recently developed capabilities (19) at the National Aeronautics and Space Administration (NASA) Space Radiation Laboratory (NSRL) at Brookhaven National Laboratory (BNL; Long Island, NY). These experiments will involve mixtures such as 60% H (at initial energy of 250 MeV/u) and 40% Si (at initial energy of 260 MeV/u). Here and throughout this article such percentages will refer to dose in Gy, and not to equivalent dose in Sv. Dose is used rather than fluence to facilitate the use of standard synergy analysis methods.
Synergy Analysis
In this article, we address how modern synergy analysis [reviewed in ref. (20)] can help plan such mixed-beam radiobiology experiments and interpret their results. Synergy analysis deals with mixtures of “agents” such as therapeutic drugs or single-ion radiation beams. It makes predictions for mixed agent dose-effect relationships based on individual dose-effect relationships (IDER) for each agent in the mixture, assuming no synergy or antagonism.
For radiation studies, work on predicting and experimentally estimating mixture effects has been actively pursued for almost 50 years (21–23). Many radiobiology mixture analyses (24–28) have simulated effects of chronic radiation exposure using two dose fractions separated by a time delay deemed comparable to, or larger than, relevant repair or relaxation times, while including theory or results also relevant to effectively simultaneous mixture exposures, without a time delay. Synergies between hadron therapy and other therapeutic agents have also been considered [reviewed in ref. (29)].
The literature on effects of agent mixtures, e.g., mixtures of therapeutic drugs in pharmacology (where drug–drug interactions can occur), of damaging chemicals in toxicology, of multiple stressors in evolutionary ecology, etc., is more extensive [reviewed in refs. (20, 30)]. This literature provided prediction methods that are potentially applicable to mixed radiation fields, including not only simple effect additivity but also various other methods now considered preferable to simple effect additivity in specific situations.
Current synergy analyses thus involve various mutually contradictory interpretations of “synergy”, corresponding to the different prediction methods. Interpretation inconsistencies are reviewed in refs. (20, 31, 32). To avoid terminological confusion, it is therefore always advisable to specify which of the possible meanings are being assigned to “synergy” and to its opposite, “antagonism”; an international agreement on how to implement this clarification is reviewed in ref. (33).
Default Hypotheses
This article emphasizes mixtures where each individual agent is a single-ion beam. Much of our current knowledge about murine HG tumorigenesis due to single-ion beams is summarized in their IDER. As in pharmacometrics and other fields (30, 34), this situation suggests the following question.
Suppose we know each IDER for a mixture and have essentially no further information. What is the most reasonable default hypothesis about the effect of the mixture?
Four points to note about this question are as follows.
1. The IDER should, as far as possible, be for experimental conditions appropriately matched to the conditions of the mixture experiment. For example, if, as in “representative beam” approaches (1, 7, 35) ∼20 gm/cm2 of shielding is used in single-beam experiments to help approximate the GCR field at the target (respectively is not used), the resulting single-beam IDER are appropriate for synergy analyses of a mixture experiment that uses such shielding (respectively does not use such shielding, as is the case in our planned experiments). Similarly, matching with regard to dose timing (e.g., acute irradiation, dose fractionation or chronic low dose rate) is also important. Recent modeling (6, 17) that provided updated IDER for single-ion murine HG tumorigenesis was based on acute exposure data. To use these IDER our planned mixed-beam experiments will also use acute exposures. In addition to enabling synergy analysis by matching dose timing, using acute exposures will also prevent severe difficulties that are encountered in low-dose-rate animal experiments at NSRL (7). However, a concomitant disadvantage is that acute exposures differ considerably from the chronic GCR dose rates encountered by astronauts above low-Earth orbit.
2. The phrase “default hypothesis” refers to predictions that can be rejected by statistical analyses of observed mixture effects. Thus, a default hypothesis must not only predict mixture effect size but also enable uncertainty analyses, such as estimation of 95% confidence intervals (CIs), from information on IDER uncertainties. Then the CI can be used to help determine if observed mixture effects above (or respectively below) the default prediction indicate statistically significant synergy (or respectively antagonism) relative to that default hypothesis.
3. Any general method using IDER to get default hypotheses on mixture effects must be able to handle heterogeneous IDER shapes, e.g., linear no-threshold (LNT), strictly concave (both shown in Fig. 1) and sigmoid (“S-shaped”, with a point of inflection) in one mixture.
4. Most of the default hypotheses that have been suggested in the literature emphasize using mathematical and computational manipulations of IDER rather than biologically or biophysically-based information in constructing default hypotheses [reviewed in (20, 36)]. As will be discussed, intentional avoidance of biological/biophysical arguments has important advantages, as well as clear disadvantages.
FIG. 1.
Convexity and concavity. Consider a dose-effect relationship E(d) with E(0) = 0, with positive slope, and with second derivative d2E/dd2. If d2E/dd2 ≥ 0 over the entire dose range [0, dmax] of interest, the relationship is convex. If the inequality is strict, that is, d2E/dd2 > 0, the slope is increasing, as shown, and the relationship can be referred to either as “strictly convex” or just as “convex”. If the inequalities are reversed then substituting “concave” for “convex” gives the appropriate terminology. Thus a strictly concave curve has decreasing slope. If d2E/dd2 = 0 at all doses, E(d) is linear no threshold (LNT).
In analyzing the above question, an N-beam irradiation with dose dj of single-ion beam j (j = 1, …, N) will be considered, with the putatively known IDER denoted by Ej(dj). Unless explicitly stated to the contrary, we will assume each IDER has the following properties in the dose range from 0 to the largest dose (dmax) considered: Ej(0) = 0 (i.e., background effect has been subtracted out); Ej(dj) is a twice continuously differentiable function; and the slope dEj/ddj is positive for dj > 0. Whether the second derivative is positive, zero or negative is important in describing IDER curve shapes (Fig. 1).
Examples of Dose-Effect Relationships
In radiobiology, linear-quadratic (LQ) IDER are frequently used. Consider irradiation for a time T or less with a dose rate r(t) and total dose 
. Then the LQ IDER is
1. α and β are non-negative constants and at least one of them is nonzero.
2. G is the generalized Lea-Catcheside functional [reviewed in (37)]. The equation for G is provided and discussed in the Supplementary Material, section S3.2 ( http://dx.doi.org/10.1667/RR14411.1.S1). G incorporates the effects of dose protraction such as fractionation, constant chronic dose rate during the time from 0 to T, variable chronic dose rate or any combination. Its functional form is based on assuming linear sub-lesion repair. It depends only on dose timing, not total dose. For a single acute dose, G = 1; otherwise 0 < G < 1. In the Supplementary Material proof is given that for fixed total dose and sufficiently low dose rate, G approaches zero and the quadratic term in Eq. (1) therefore drops out.
3. If β = 0, the IDER is LNT. Assuming instead α = 0 gives the pure quadratic case:
4. An LQ curve is convex, as defined and exemplified in Fig. 1. If βG > 0 the curve is strictly convex.
Many of the IDER used for analyzing nontargeted effects are instead concave (17, 38–40). The IDER used most commonly in pharmacometrics are Hill functions (30, 34, 41), reviewed in the Supplementary Material ( http://dx.doi.org/10.1667/RR14411.1.S1). Some Hill functions are concave. Others are sigmoid, having positive second derivative for small doses, negative second derivative for large doses and one point of inflection.
LNT IDER and Simple Effect Additivity
The simplest example of a default hypothesis on N-component mixture effects is simple effect additivity. For example, suppose all the IDER are well represented by LNT functions. Then Ej(dj) = αjdj, with αj > 0, and the simple effect additivity default prediction for the mixture is mixture effect,
For error analyses the simple effect additivity default hypothesis includes the assumption of statistical independence for the individual contributions in the sums (8). If the uncertainties of each IDER are known, then 95% CI for S(d1,…,dN) can be computed. Mixture CI are typically calculated by Monte Carlo simulation. In this LNT case, standard statistical formulas can also be used.
We hereafter write the simple effect additivity prediction as S, even if many different kinds of curvilinear IDER are involved, rather than only IDER that are LNT:
In Eq. (4) and from now on, it is assumed that the total mixture dose d determines each of the individual component doses, usually as a dose-independent fraction of the total mixture dose. As shorthand, we shall use S(d) to refer to the simple effect additivity default hypothesis as well as to the mixture dose-effect relationship, Eq. (4), that the hypothesis predicts.
Limitations of Simple Effect Additivity S(d)
In practice, S(d) is routinely used in planning experiments to make prospective power estimates. A default prediction of S(d) in the case that all IDER are LNT, with synergy or antagonism then judged by observed deviations from S(d), is generally accepted in radiobiology (8) and throughout biology. However, somewhat surprisingly, S(d) has often been found quite inappropriate as a default pharmacometric hypothesis if some of the component IDER are not LNT (20, 42). S(d) often treats the changes in slope that curvilinearity implies in an unrealistic way (36). One main argument has been that for curvilinear IDER, using S(d) gives flatly incorrect predictions for so-called “sham mixtures” of an agent with itself (20, 36), as exemplified in Fig. 2. Such problems with simple effect additivity S(d) have been discussed for more than a century (43).
FIG. 2.
An example where simple effect additivity is an underestimate. Consider a hypothetical case where a single-ion beam has purely quadratic IDER, E = d2 (black line). Regard the beam as a 50-50 mixture of two single-ion beams, both of which happen to have the same dose-response curve as the original beam. Then for total mixture dose d, each of the two beams contributes dose d/2 and thus has effect E/4. Using simple effect additivity S(d) therefore gives effect E/2 rather than the correct answer E. In the special case that all component IDER are LQ, one can correct simple effect additivity by using Eq. (5) instead. In general, some other method, such as incremental effect additivity, is needed to deal with mixtures having heterogeneously shaped IDER. For a mixture of components that are similar but not exactly identical, corresponding discrepancies arise. For example, suppose every IDER for a mixture is convex, as shown here. Then S(d) is likely to be an unrealistic underestimate (like the red curve). If every IDER is concave (Fig. 1) then S(d) is likely to be an unrealistic overestimate.
Therefore, the default hypotheses now favored in radiobiology and other fields often differ from S(d). One example is the biophysically motivated dual radiation action hypothesis D(d) on mixtures suggested by Zaider and others (24, 44–48). The hypothesis applies to mixtures only if each component has an LQ IDER, Eq. (1). For one acute dose (so that G = 1), the prediction D(d) differs from S(d) by using the square of a sum instead of a sum of squares for the LQ dose-quadratic term:
In the sham mixture example (Fig. 2), Eq. (5A) gives the correct result for the mixture effect, twice as large as the incorrect prediction given by simple effect additivity S(d).
In the Supplementary Material, sections S4 and S6 ( http://dx.doi.org/10.1667/RR14411.1.S1) some other common approaches to synergy are discussed, comparing them to each other and to Eq. (5A); the linear isobole approach commonly used in pharmacology (36) is emphasized. For reasons that are discussed in the Supplementary Material, section S4, none of the known alternatives to S(d) proved suitable for the mixed-ion beam HG tumorigenesis calculations in the current article. For example, D(d) in Eq. (5) could not be used, because the IDER we will use (6, 17) are not LQ; as will be discussed in the next section, some of these IDER include terms that describe nontargeted effects and cause the IDER to differ strongly from LQ dose responses at very low doses.
We will therefore use a new default hypothesis, incremental effect additivity I(d), defined in the Methods section. This approach borrows from the simple effect additivity and linear isobole approaches, while circumventing flaws in the linear isobole approach (described in the Supplementary Material, sections S4.2.1 and S4.3.3; http://dx.doi.org/10.1667/RR14411.1.S1) as well as the flaws in S(d) discussed above. I(d) is more suitable for our purposes, and perhaps in general, than other synergy analysis approaches.
Preview
While most murine HG tumors are not carcinomas, their dose-response curves have often been considered informative about carcinogenesis risks and about quality factors [(16), reviewed in refs. (17, 49)]. Our calculations compare simple effect additivity predictions S(d) with incremental effect additivity predictions I(d) for planned experiments on murine HG tumorigenesis after mixed-beam exposure at the NSRL. Our analysis is entirely in silico. It uses modeling (6, 17) on single-ion beam murine HG tumorigenesis based on data summarized in (15, 17). Most of the data is for acute irradiation and the IDER proposed by Cucinotta et al. and Chang et al. (6, 17) do not consider dose-rate effects. Therefore, in the rest of the current article we shall always assume a single acute dose unless explicitly stated otherwise. The calculations serve as a specific illustration of how current synergy analysis methods can model mixed-beam effects in biological systems. Frequently used abbreviations are summarized in Table 1.
TABLE 1
Frequently Used Terms, Acronyms and Symbols
MATHEMATICAL/COMPUTATIONAL METHODS
Software
We used the free, open-source computer language R (50), initially designed for statistical calculations but now rapidly gaining acceptance among modelers (51). Our customized source codes are available free of charge at https://github.com/rainersachs/NASA/ and at https://github.com/nopphons/NASA/.
Harderian Gland Tumorigenesis
We shall adopt the IDER and parameters used by Cucinotta and colleagues (6, 17, 39) in an analysis of murine HG tumorigenesis. They consider the data set as a whole, rather than analyzing each ion type separately, include analyses of parameter uncertainties and use information criteria to select optimal models. They approximate track-structure effects as being, for Z > 2, independent of Z when linear energy transfer (specifically, LET∞) L is fixed. A study utilizing track structure models that have been used for many other end points instead of merely using L may be available soon (40). However, the approximation that only L matters is adequate for our purposes here, i.e., reviewing current approaches to synergy analysis of mixtures, introducing I(d), illustrating its application to mixed-beam murine HG tumorigenesis studies and comparing it to S(d).
The HG analysis utilized by Cucinotta et al. (6, 17, 39) considers two alternative types of IDER. Formally, both types have maxima at comparatively high doses and then decrease at still larger doses. It is not clear from the published HG data whether there is such a decrease rather than just a leveling off. The decrease will not play any role in our analysis of planned mixture experiments, where no HZE doses >40 cGy are proposed.
The first, targeted effect (TE) IDER type, PTE(d), is
Here, the seven parameters (α0, α1, α2, β, λ0, λ1 and λ2) and the background incidence P0 were estimated from the published literature on single-ion experiments with ions of various LETs. The parameter values are given in table 6 of ref. (17), part of which is repeated with permission in modified form here as Table 2.
TABLE 2
Parameters Used
Alternative IDER, incorporating nontargeted effects [NTE; reviewed in (52)] were considered by others (6, 17, 39) for the same data and found to have better information criteria scores. These NTE IDER were based on a scenario where the high energy ions create, in addition to targeted effects, nontargeted effects that rise very sharply at very low doses, saturate at a dose in the mGy region (smaller than any dose used in the past or here in proposed murine HG experiments) and then gradually decrease as the dose increases. The NTE IDER is:
Here, κ(L) = κ1Lexp(–κ2L); Θ is a factor that is 0 for dose = 0, rises in the mGy region and is effectively 1 for all doses where data on heavy-ion murine HG tumorigenesis is available. The term in Eq. (7) involving κ(L) represents the nontargeted part of the total effect and for Z > 2 dominates at doses <1 cGy. The other term, which dominates at large doses, represents the targeted part of the total effect. It has the same form as in Eq. (6), but the fitted parameters are somewhat different (Table 2).
For Θ, in (17) a step function is used at d0 (Fig. 3), where d0 is set to 0.1 cGy with the understanding that any d0 in the mGy region, greater than 0 and less than 1 cGy, would give the same results for analyzing single-ion beam experiments that use doses >1 cGy (or in the controls, dose = 0). However, for our calculations on mixed-beam experiments, the infinite slope at d0 causes problems. We need a function Θ with finite slope. The criteria we used for choosing Θ were the following: 1. Θ should be very close to 1 for all nonzero doses used in the relevant published single-beam experiments; 2. It should involve only one adjustable parameter (replacing the adjustable parameter d0); 3. It should depend on dose and LET only through the average number H of direct hits to an HG nucleus (proportional to particle fluence); and 4. It should have the concave shape (Figs. 1 and 3) usually associated with nontargeted effects (53, 54) for end points where doses low enough to obtain data in the region where nontargeted effects are increasing from 0 to their saturating value are used.
FIG. 3.
NTE curves near d = 0. As an example, the NTE IDER, Eq. (7), for LET = 200 keV/μm is shown in the mGy region, with various assumptions about how the nontargeted part of the effect builds up, as governed by the factor Θ. The dashed line shows the assumption in reference (17), that Θ is a step function at d0 = 1 mGy. The red lines show the assumption of Eq. (8) with respective values, left to right, of φ = 105, 104, 3 × 103 and 103 (units Gy-μm/keV). Importantly, it is seen that whenever φ ≥ 3 × 103, the IDER slope in the figure at doses >4 mGy (which is approximately α Gy–1) is small compared to the average slope between 0 and 4 mGy.
Taking H = 150d/L when d is in Gy and L in keV/μm (17), a function that meets all these criteria is:
Here, as shown in Fig. 3, the inequality on the parameter φ corresponds approximately to the inequality d0 ≤ 0.1 cGy for a step function version of Θ. Extensive numerical explorations (Supplementary Material, section S4.3.2; http://dx.doi.org/10.1667/RR14411.1.S1) showed that our final results are insensitive to the value of φ, provided the inequality in Eq. (8) holds. The Supplementary Material, section S4.3.2 also discusses an intuitive interpretation of φ.
Figure 4 shows some of the resulting IDER: for protons (Z = 1); and for atomic nuclei with Z > 2 having the 3 LET values shown. Table 3 shows examples that have the LETs used in Fig. 4 at approximately the energies shown in the table. There are other Z > 2 ions that have the same LET for appropriate energies. Due to the aforementioned approximation in refs. (6, 17), only the LET will matter in our calculations for Z > 2, not the ion species.
FIG. 4.
Sample IDER HG tumorigenesis curves. High-energy protons and Z > 2 ions with the three LETs shown are used as examples of IDER in the dose and effect ranges of main interest in this article. Panel A represents Eq. (6) and panel B represents Eqs. (7) and (8). As shown, all curves in panels A and B have an approximately constant slope, apart from the drastic changes in the slope at mGy doses in panel B, shown in more detail in Fig. 3.
TABLE 3
Examples of Ions that have a Given L = LET∞
Default Predictions for Mixed Beams
Consider acute irradiation with a mixed beam of N ≥ 2 different radiation qualities. The dose proportions rj that the different qualities contribute to total dose 
obey the equations:
In our subsequent calculations, rj will always be independent of dose. This is a typical pattern for acute irradiation at the NSRL. Using Eqs. (6), (7) and (9) in the simple effect additivity prediction, (Eq. 4) gives:
S(d) will be compared to the incremental effect additivity prediction I(d), described below.
Incremental Effect Additivity I(d)
I(d) modifies S(d) predictions and linear isobole predictions. It uses small increments (i.e., derivatives) and “compositional inverses”. Compositional inverses are used in radiobiology when discussing RBE. They are needed when using effect, rather than dose, as the independent variable. They play a prominent role in computing isobole default hypotheses (Supplementary Material, sections S4.1 and S4.2.1; http://dx.doi.org/10.1667/RR14411.1.S1). The compositional inverse of a monotonically increasing function reverses the action of the function. For example, for x > 0,
, so the positive square root function is the compositional inverse of the squaring function; note that the compositional inverse of x2 is not x–2. As another example, exp[ln(x)] = x and ln[exp(y)] = y, so the functions exp and ln are compositional inverses of each other.
Suppose we have a mixture of N components and each IDER Ej(dj) has a compositional inverse function, denoted by Dj(Ej). Then we will define I(d) as a solution of the following first order, typically nonlinear, separable ordinary differential equation.
Here, d is the total mixture dose. The brackets containing the subscript reflect the following steps used for calculation: 1. Find the slope of the jth IDER curve as a function of individual dose dj. 2. Evaluate dj using the compositional inverse Dj with the argument of Dj being the effect I already present due to the influence of all the components acting jointly. Integrating the differential Eq. (11) using the initial value I(d = 0) = 0 defines the incremental effect additivity dose-response relationship I(d) for the mixture.
Using dj = Dj(I) in Eq. (11) instead of the seemingly more natural dj = Dj(Ej) is the key assumption made. It can be shown (Supplementary Material, section S4.3.1; http://dx.doi.org/10.1667/RR14411.1.S1) that using dj = Dj(Ej) would merely lead back to simple effect additivity S(d).
Equation (11) can be interpreted as follows, as the total mixture dose increases slightly, every mixture component contributes some incremental effect. The size of the incremental effect is determined in an appropriate way, by the state of the biological target, specifically by total effect already contributed by all the components collectively (and not by the dose the individual component has already contributed). In this way different components appropriately track changes of slope both in their own IDER and in the other IDER. A more detailed derivation of Eq. (11) is given in the Supplementary Material, section S4.3.1 ( http://dx.doi.org/10.1667/RR14411.1.S1).
To clarify the mathematical manipulations involved in Eq. (11), one can use a hypothetical illustrative example with purely quadratic IDER, Ej(dj) = βjdj2. This is one of the cases that is quite exceptional, where all the mathematical manipulations required to set up and solve Eq. (11) can be done with equations, rather than only by numerical simulations. The Supplementary Material, section S4.4.1 ( http://dx.doi.org/10.1667/RR14411.1.S1) shows proof that in this case incremental effect additivity I(d) gives:
Here, all the individual doses dj are given by Eq. (9) as linear functions of total dose d, so both I(d) and S(d) are themselves pure quadratic functions of d, but with different coefficients. Comparing Eq. (4) with G = 1 to Eq. (12) shows that the biophysically-based default hypothesis D(d) of dual radiation action and incremental effect additivity I(d) give the same predictions in this special case, even though I(d) is calculated by Eq. (11), which: 1. Can be applied even when the IDER are not LQ; and 2. Does not use biophysical arguments.
In integrating Eq. (11) it can sometimes happen that I becomes so large that it approaches, reaches and then exceeds the maximum Ej for a particular component. Then as I approaches maximum(Ej) from below, the component in question makes a smaller and smaller contribution to dI/dd, since the derivative of Ej at its maximum is zero. For values of I greater than maximum (Ej), the contribution of the jth component is taken here to be zero, as it was at maximum (Ej). This extra assumption makes incremental effect additivity applicable over dose and effect ranges sufficiently large for our calculations, as detailed in the Supplementary Material, sections S3.1.2 and S4.3.3 ( http://dx.doi.org/10.1667/RR14411.1.S1).
Calculations of Predicted HG Tumorigenesis for Mixed Beams
I(d) and S(d) were calculated by numerical simulations for six different ion mixtures, including four that are scheduled to be used in future experiments. The mixtures are described in Table 4. Subsequent figures and tables use the same labeling “A–F” (or A*–F* when discussing NTE IDER, in which case φ = 3 × 103 was used). Row F is a hypothetical mixture designed to illustrate certain differences between S(d) and I(d). Many other mixtures were analyzed. We were unable to find a systematic way to enlarge Table 4, even though we tried various methods. There is a “combinatorial complexity” problem involved: even after one has decided on ions and their energies, there is, for mixtures with many components, a bewildering variety of different possible percentage patterns (column 4, Table 4). Implications of this combinatorial complexity will be considered in the Discussion section.
TABLE 4
Mixtures
Because the examples in Table 4, or even many more examples, cannot adequately summarize the very large number of relevant mixtures, Supplementary Material, Section S6 ( http://dx.doi.org/10.1667/RR14411.1.S1) gives results for an illustrative and much simpler case, where a comparatively complete mathematical overview is attainable. Confining attention to mixtures of two radiation qualities both having LQ IDER, Supplementary Material, section S6 compares four different default mixture hypotheses in detail, using theorems, conjectures supported by strong numerical evidence and examples.
Confidence Intervals (95%)
Confidence intervals for mixture dose-effect relationships were calculated under the default hypotheses S(d) and I(d) using the following steps. Only uncertainties in the IDER, not any extra uncertainties in the mixture measurements, were taken into account. The IDER parameters were modeled as independent gamma distributions with the mean values in Table 2 and with a variance whose square root equals the errors listed there. The use of gamma distributions prevented the possibility of negative values. These assumptions allowed Monte Carlo sampling of the IDER parameter space. For each set of parameters, S(d) and I(d) were computed and at least 1,000 Monte Carlo runs gave numerical distributions from which the 2.5 percentile and 97.5 percentile could be estimated to get 95% CI. Supplementary Material, section S5 ( http://dx.doi.org/10.1667/RR14411.1.S1) considers the extent to which the assumption of parameter independence could influence the results of these calculations, giving one example that suggests that the CI calculated assuming parameter independence may be overestimates.
RESULTS
The results of calculating default predictions S(d) and I(d) for the 6 mixtures in Table 4 are shown in Fig. 5, which indicates various qualitatively different patterns for different mixtures and for TE vs. NTE IDER. Panels A and F* of Fig. 5 are extremes. In Fig. 5A, predictions are almost LNT, the two predictions almost coincide, and both predictions are nested between the IDER of the two components. On the other hand, in Fig. 5F* both predictions are markedly concave at very low doses and the two predictions differ very substantially, with S(d) predicting large effects (>20%) for a total mixture dose of only a few mGy and S(d) being higher than any of the component IDER. Implications of the patterns will be considered in the Discussion section. As detailed in the Supplementary Material, section S4.3.2 ( http://dx.doi.org/10.1667/RR14411.1.S1), all the qualitative differences can be explained intuitively in terms of two factors: different IDER shapes, such as the convex H curve and concave Fe curve in Fig. 5A; and the fact that I(d) uses effect, not dose, as the basic independent variable.
FIG. 5.
Comparing predictions of two default hypotheses. Panels A–F: Default predictions S(d) and I(d) are compared for the six mixtures shown in the corresponding rows of Table 4 when TE IDER are used; Panels A*–F*: Default predictions S(d) and I(d) are compared for the same six mixtures when NTE1 IDER are used. In panels A and B, the green IDER are labeled with the ion species and these identical curves are present in many of the other panels. For fuller descriptions of the eight ion and seven HZE mixtures, see Table 4. It is seen that: 1. The red curve for I(d) always lies nested between the curves for the individual mixture components, but in panels D*–F* the dotted S(d) curve is not everywhere nested; 2. In panels A–F and A*–C*, I(d) and S(d) are rather similar; and 3. There are pronounced differences between I(d) and S(d) in panels D*–F*.
Table 5 gives numerical information on the differences between S(d) and I(d) shown visually in Fig. 5. We defined a maximum relative excess ρ as follows.
Here dR is the dose at which the absolute value |R(d)| is a maximum for d in the closed interval [0.0001,dmax]; dmax is again the maximum dose considered (Table 4 and Fig. 5).
TABLE 5
Effect Magnitudes and Relative Differences
For example, in Table 5, column A* (for 60/40 H + Fe in the dose range 0 < d ≤ 100 cGy) ρ(d) had a maximum of 2.8% at 2.5 cGy and a minimum of –2.2% near 50 cGy. Since the absolute value of the maximum is larger than the absolute value of the minimum the value of ρ at 2.5 cGy is recorded in Table 5.
The maximum relative excess ρ must lie between –1 and +∞ by its definition and the fact that in the dose range of interest each IDER is positive with positive slope. Then ρ → ∞ means S(dR) >> I(dR). As shown in Fig. 5 and in Table 5 (bold text), the only drastic differences between S(d) and I(d) occur when NTE IDER are assumed and there are two or more HZE ions in a mixture, in which case ρ can be quite large at a few cGy, and be above 50% even at a few mGy.
To see how the differences between predictions compare to the uncertainties in each prediction due to uncertainties in the IDER, we used Monte Carlo sampling. Figure 6 shows 95% CI for the panels shown in Fig. 5. In most, but not all panels, the difference between the two predictions is considerably smaller than the 95% CI for either.
FIG. 6.
Confidence intervals for mixture effect predictions. Monte Carlo simulation 95% CI estimates for predictions on the mixtures in Table 4 are shown. Note the different dose scale in the last column, which leads to a visual underestimate of the steepness of the curves compared to the steepness in the other columns. In panel F* the blue shaded area has been truncated at 100% for plotting purposes. Except in panels D* and F*, the vertical difference between the two central value curves is small compared to the estimated 95% CI for either prediction.
DISCUSSION
Synergy analysis helps plan mixture experiments and interpret their results using IDER-based computations. It requires choosing a default hypothesis that defines what particular definition of synergy is being used (20). Simple effect additivity S(d), given by Eq. (4), is the most obvious default hypothesis, but is often not the best (20, 34, 36, 55, 56). In this work, we compared S(d) with incremental effect additivity I(d) as applied to mixed-beam experiments on murine HG tumorigenesis.
In this discussion, we review our results on I(d) vs. S(d), then outline the advantages and disadvantages of predicting and interpreting mixed agent effects using formal synergy analysis based on IDER rather than using more biologically-based approaches. Next, we address some differences between acute-dose mixed-beam murine HG tumorigenesis experiments without extra shielding compared to astronauts' exposure during prolonged travel above low-Earth orbit. Our conclusions are summarized in the final section of this article.
Review of Results
I(d) and S(d) give similar predictions unless nontargeted effects are important and a mixture contains a considerable admixture of two or more HZE ions.
Nontargeted Effect IDER
Panels D*–F* in Figs. 5, 6 and Table 5 give examples where, assuming NTE1 IDER, S(d) is much larger than I(d) at low doses [Supplementary Material, section S4.3.2 ( http://dx.doi.org/10.1667/RR14411.1.S1)], gives another example where the difference is even larger, with the red 95% CI ribbon for I(d) lying entirely below the blue 95% CI ribbon for S(d). Since NTE1 models were preferred in Chang et al. (17), and since realistic mixtures include contributions from many different HZE ions, cases where S(d) and I(d) differ drastically are important.
The reason for the large differences is that simple effect additivity S(d), as expected from its known drawbacks (20, 34, 36, 55, 56) and the pronounced changes of slope at low doses shown in Fig. 4B, gives unrealistic predictions. Specifically, S(d) does not take saturation of nontargeted effects (54) into account properly. For a mixture of N different ions S(d) predicts that the nontargeted part of the total effect rises rapidly to approximately
as the dose increases from 0 to approximately N/10 cGy. The sum can be very large if there are many ions in the mixture. Correspondingly, one can, by generalizing the arguments used in Fig. 2, construct sham mixtures of N identical components where the S(d) prediction, Eq. (14), directly contradicts the single-component HG tumorigenesis data in the literature.
In contrast, I(d) predicts that the nontargeted part of the mixture effect saturates at about the largest individual κ among all the ions in the mixture. This was determined by numerically analyzing examples such as the NTE1 IDER shown graphically in Fig. 5D*–F*. By treating total effect, rather than total dose, as the basic independent variable, I(d) automatically incorporates nontargeted effect saturation at a plausible level without using any of the biophysical arguments usually invoked to explain saturation in detailed bystander models. The fact that in Fig. 5D*–F* the S(d) and I(d) curves are approximately parallel from total doses of a few cGy to 50 cGy or more, indicates that the S(d) overestimates are due primarily to marked IDER curvature at very low doses. Intuitive explanations are given in the Supplementary Material, section S4.3.2 ( http://dx.doi.org/10.1667/RR14411.1.S1) for why S(d) is an overestimate whenever IDER are concave. The near parallelism means that the function Δ(d) = S(d) – I(d) is slowly varying in the experimentally relevant range d > 1 cGy.
Targeted Effect IDER
The total dose of relevant mixtures in space will presumably always include a contribution of substantially more than 50% (in Gy) from low-LET radiation, especially from high-energy protons (5, 19, 57). Results for mixtures incorporating this restriction and assuming TE IDER, which are virtually LNT or LQ in the dose and effect ranges of main interest rather than incorporating an additional nontargeted effect term, are shown in Fig. 5A–D. The I(d) and S(d) curves are visually quite close in Fig. 5A–D. The numerical estimate of the maximum relative excess (S(d) – I(d))/I(d) is <8%, which is small compared to the relative errors in the IDER parameters (Table 2) and well within the predicted 95% CI for mixture effects shown in Fig. 6. For TE IDER the fact that the predictions of I(d) and S(d) are so close increases the credibility of both predictions and also tends to validate the common procedure of using S(d) for experiment planning purposes. However in (16) TE IDER were not prefered.
Importance of High-Quality IDER
Information about synergy or antagonism is important. In practice, mixture experiments will be planned and interpreted using, explicitly or implicitly, synergy analyses (1). Our results, especially those on CI (Fig. 6), indicate that without high-quality IDER, synergy analysis is disabled.
Advantages and Disadvantages of Synergy Analysis
General Considerations. Ultimately, estimating cancer risk for astronauts due to radiation exposure will require biophysically-based knowledge of cancer etiology (Fig. 7, blue arrows). IDER-based synergy analysis default predictions of mixed-beam dose-effect relationships (Fig. 7, black arrow) are only temporary expedients (34, 36). Such expedients are much less reliable than predictions based on biophysical understanding. However, the expedients are typically orders of magnitude faster, cheaper and simpler (1, 20, 36).
FIG. 7.
A long hard road or a temporary shortcut. Eventually, but probably not soon, default predictions about GCR mixed field carcinogenesis based solely on IDER (far left arrow) will be replaced by biologically-based predictions that incorporate whatever synergy or antagonism actually occurs (blue arrows). For now, optimizing the far simpler shortcut is important.
Any IDER-based approach produces only a default hypothesis. If a mixture has substantial agent–agent interactions not encapsulated in the component IDER, synergy or antagonism in the sense of that hypothesis is involved. IDER-based analyses are needed to define synergy, but cannot predict it (22). If there is major synergy or antagonism, biophysical insights and/or multiple (expensive) mixture experiments are needed to clarify the situation (1, 7, 20, 36).
If, in fact, there is no major synergy or antagonism, then mixture results can be predicted from observed mixture component IDER. The severe combinatorial complexity problem mentioned earlier, that in general a mixture result cannot be extrapolated even to a mixture with the same components but in different proportions, is largely circumvented, a drastic simplification for mixtures containing many different components.
Using I(d). For upcoming mixed-beam murine HG tumorigenesis, predictions and interpretations using I(d), despite being just a temporary expedient (Fig. 7), have some substantial potential advantages. For example, even a single-ion beam will create a complex radiation field at the HG due to interactions with material in the beam, with any extra shielding introduced intentionally, with shielding due to mouse containers and with mouse self-shielding (1, 5, 7). The single-ion beam can nonetheless be considered as a single agent acting on the HG just as in pharmacometrics, a drug whose mode of action is complex and not fully understood can be considered as a single agent. The single-ion beam IDER can be used in standard synergy analysis (see the shortcut in Fig. 7) for a mixed beam matched with regard to shielding conditions. Then the synergy analysis automatically incorporates effects of secondary radiations due to the shielding, thereby including influences of shielding complications in the default mixture predictions. The long route in Fig. 7, involving biophysical understanding of HG tumorigenesis due to low-dose, low-dose-rate exposure to a complex “representative” radiation field (1, 5, 7), will not be available in the foreseeable future.
At times, like linear isobole analyses (36), I(d) can mimic biologically-based arguments. Examples include automatically incorporating NTE saturation, as discussed above, and the fact that in Eq. (12) I(d) gives an equation previously derived from mechanistic considerations of pairwise lesion interactions. That is gratifying in these special cases, but it is unfortunately not a general feature. For example, given two IDER that are LQ, not just pure quadratic, the biophysically based prediction Eq. (5) and the default prediction using I(d) in general differ, as shown in the Supplementary Material, section S6 ( http://dx.doi.org/10.1667/RR14411.1.S1).
Murine HG Tumorigenesis Experiments versus Astronaut Carcinogenesis. Current estimates of astronaut carcinogenesis risk use the predominantly low-LET data in the life span study of Japanese atomic bomb survivors, and use animal experiments that compare high LET with low LET; a wide variety of auxiliary data are also considered [reviewed in (6)]. Mixed field experiments designed to improve the risk estimates will almost inevitably use synergy analysis, explicitly or implicitly, both in experimental planning and in interpreting their results (1).
Our predictions on upcoming experiments (see panels B–E and B*–E* in both Figs. 5 and 6; Table 5) illustrate modern synergy analysis and its application to mixed-beam murine HG tumorigenesis. The calculations concern highly simplified situations where single-ion experiments have, over the years, led to comparatively well characterized IDER (TE and NTE1). A requisite for synergy analyses is that the exposure protocols in the upcoming experiments match those of the single-ion experiments as closely as possible (1).
Actual astronaut exposures differ substantially in a number of important ways. For example, the proposed experiments use one acute exposure, whereas, except for solar particle events, astronaut exposure involves chronic irradiation at much smaller dose rates. Specifically, 10 cGy in two years (1, 3, 5, 7) compared to 50 cGy in 10 min involves a factor of ∼2 × 10–6. Several additional discrepancies, such as marked differences in the composition of the radiation field at the target and in biological end point, are detailed in the Supplementary Material, section S2 ( http://dx.doi.org/10.1667/RR14411.1.S1).
The synergy analyses provided for simplified situations should help in understanding synergy, antagonism or the absence of both in experiments that are less idealized. The single-ion experiments that led to the IDER we used have, despite the discrepancies mentioned above, been influential in notions about astronaut risk (39). Corresponding mixture experiments potentially play an analogous role.
Conclusions
Overall, our main conclusions are the following.
First, and probably foremost, individual dose-effect relationships are important supplements to mixed-beam experiments, and should be obtained under matching experimental conditions. Accurate IDER are critical for determining the presence or absence of synergy or antagonism.
Synergy analysis using simple effect additivity S(d) is not optimal for experiments on murine HG tumorigenesis. If nontargeted effects are negligible and mixtures of high and low LET are being analyzed, S(d) predictions are numerically close to predictions using the preferred incremental effect additivity estimate I(d). However, S(d) gives unrealistic results in some other cases. In contrast, I(d) gave consistent and reasonable results in all cases considered. We therefore advocate the use of I(d).
In general, all synergy analysis approaches, including I(d) and S(d), have important disadvantages and also important advantages compared to approaches based on biophysical understanding (Fig. 7).
If accurate in silico IDER-based default predictions on mixed-beam effects can be found, NASA's planning and interpretation of mixed-beam experiments, and perhaps even NASA's planning for missions beyond low-Earth orbit, can be substantially simplified using such predictions.
SUPPLEMENTARY INFORMATION
Supplementary Material is available at http://dx.doi.org/10.1667/RR14411.1.S1.
ACKNOWLEDGMENTS
We thank F.A. Cucinotta and P.Y. Chang for highly informative discussions. This research was supported by NASA [grant no. NNJ16HP22I (RKS and EAB)] under U.S. Department of Energy contract no. DE-AC02-05CH11231. Additional support was provided by the Undergraduate Research Apprenticeship Program (URAP) at UC Berkeley (NS, AC and NH).
