Environ. Sci. Technol. 2008, 42, 1308–1314
A General Approach to Modeling Biphasic Relationships W I L L I A M N . B E C K O N , * ,† C A R Y P A R K I N S , ‡ ALEXEY MAXIMOVICH,§ AND ANGELA V. BECKON4 U.S. Fish and Wildlife Service, 2800 Cottage Way, Suite W-2605, Sacramento, California 95825, Inspironix, 2100 Watt Avenue, Suite 130, Sacramento, California 95825, Saint Petersburg State University, Universitetskaya nab. 7–9, Saint Petersburg 199034, Russia, and Pomona College, 170 East Sixth Street no. 848, Claremont, California 91711-7004
Received May 15, 2007. Revised manuscript received September 24, 2007. Accepted October 22, 2007.
Biphasic relationships can be found throughout the sciences, especially in the dose-response relationships of pharmacology, toxicology, agriculture, and nutrition. Accurate modeling of biphasic dose-response is an essential step in establishing effective guidelines for the protection of human and ecosystem health, yet currently-used biphasic mathematical models lack biological rationale and fit only limited sets of biphasic data. To model biphasic relationships more closely over wider ranges of exposures, we suggest a simple, general, biologically reasonable modeling approach leading to a family of mathematical models that combine log-logistic functions: at least one for the upslope and one for the downslope of the biphasic relationship. All parameters employed are meaningfully interpretable. These models can be used to test for the presence of biphasic effects, and they simplify to a standard log-logistic model in the special case where no biphasic effect can be detected. They offer the promise of improvement in assessment of the safety and efficacy of pharmaceuticals and nutrients as well as in determination of the toxicity of contaminants. Additionally, they may be useful in modeling nonmonotonic cause-effect relationships in other scientific disciplines.
Introduction Science has been handicapped by a lack of good mathematical models for describing a common kind of causeeffect relationship, one in which the dependent variable is maximal or minimal at an intermediate level of the independent variable. Such relationships are known as “biphasic” because they include both a rising and a falling phase. Biphasic relationships are depicted by hill-shaped or valleyshaped curves. They can be found throughout the sciences, particularly describing phenomena related to homeostasis, equilibrium, stability, and optimization. Examples emerge from disciplines as diverse as physiology (e.g., temperature affecting growth rate (1)), ecology (e.g., disturbance affecting species diversity (2)), psychology (e.g., intensity of punishment affecting rate of training (3)), and economics (e.g., tax rate affecting tax revenue 4, 5). * Corresponding author phone (916) 414-6597; fax (916) 414-6712 ; e-mail:
[email protected]. † U.S. Fish and Wildlife Service. ‡ Inspironix. § Saint Petersburg State University. 4 Pomona College. 1308
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 42, NO. 4, 2008
Biphasic relationships are especially common in the fields of nutrition, agriculture, pharmacology, and toxicology. Indeed, anything that is has a positive effect on an organism in small or moderate amounts has a negative effect in sufficiently large amounts (6), thus exhibiting a biphasic “dose-response” relationship. Well known examples affecting human health include exposure to sunlight (7, 8) and oxygen, and dietary intake of calories, sodium, iron, zinc, copper, iodine, fluorine, and selenium. Even substances that are generally viewed as toxic may, at very low doses, induce some sort of stimulatory response, at least temporarily (9–11). Our perspective and much of the language we use here derive from the disciplines of toxicology and pharmacology. Nevertheless, we recognize that the traditional dose-response terminology of toxicology and pharmacology is too restrictive because the general approach we suggest may be applicable to some cause-effect relationships in other disciplines as well. To facilitate these broader applications, here we use the term “affecter” to denote an independent variable in any relationship between independent and dependent variables. Basics of Dose-Response Modeling. In any population, individuals differ in sensitivity to any given affecter, such as a chemical that affects their physiology. The sensitivity of each individual may be quantified by a sensitivity threshold, the minimum exposure to the affecter sufficient to cause some given level of effect. If the probability of an individual having a sensitivity threshold at exposure or “dose” x is given by some function f(x), then the expected proportion y of individuals exhibiting at least the given level of effect is the cumulative function: y ) F(x) )
∫ f(x′) dx′ x
0
(1)
This dose-response relationship applies to positive measures of positive effects (e.g., the efficacy of pharmaceuticals or the growth induced by fertilizers), as well as to negative measures of negative effects (e.g., mortality caused by toxic environmental contaminants). If a positive measure, such as survival, is used for a negative effect, such as the toxicity of a contaminant, then the dose-response relationship is inverted: y ) 1 - F(x)
(2)
In toxicological and pharmacological applications, and perhaps in other spheres, log-normal distributions (normal distributions of responses to logarithm-transformed doses) are good candidates for f(x). Log-normal distributions characterize phenomena throughout living systems as well as the physical environment (12). Contaminant concentrations are well described by log-normal distributions (13), and an assumption of log-normal frequency distributions of sensitivity thresholds underlies the widely-used probit analysis of dose-response relationships (14). A log-normal probability distribution of sensitivity thresholds f(x) yields a cumulative log-normal dose-response relationship F(x) (15) that is cumbersome and not suited for nonlinear regression procedures. However, the cumulative log-normal function is closely approximated by the simple log-logistic model (16, 17), which can be expressed in terms of parameters that are readily interpretable and useful in the fields of toxicology and pharmacology: y ) Min +
Max - Min 1 + ( ⁄ x)β
(3)
where the independent variable x is the “dose,” that is, the amount, concentration, or intensity of the affecter; the 10.1021/es071148m CCC: $40.75
2008 American Chemical Society
Published on Web 01/15/2008
dependent variable y is the expected amount of the effect or proportion of the population experiencing the effect; parameter Min is the minimum effect or proportion affected; parameter Max is the maximum effect or proportion; parameter is the dose at the inflection point of the doseresponse curve, i.e., the ED50, the dose at which the effect is reduced by 50% or half the population is affected; β is a parameter related to the maximum slope of the curve, which occurs at dose . Often, Min is 0 and Max is 1, or the data may be normalized to meet these conditions; then the loglogistic equation can be simplified to a two-parameter version: y)
1 1 + ( ⁄ x)β
(4)
Forms of the log-logistic equation (eqs 3 and 4) are widely used in pharmacology and toxicology to model both upward sloping (β > 0) and downward sloping (β < 0) dose-response relationships (eqs 1 and 2 respectively). This broad applicability is consistent with log-normal distributions of sensitivity thresholds for both positive and negative effects. The discussion above is framed in terms of dichotomous outcomes (above or below an effect threshold). This approach is particularly suitable for quantal (all-or-none) effects, such as survival. However, the log-logistic equation that fits such a conception also has been used successfully in modeling count effects, such as fecundity, and continuous effects, such as growth (e.g., refs 16 and 17). The wide utility of the loglogistic equation is fortuitous for model development, but the reasons behind it remain to be elucidated. Existing Biphasic Models. The log-logistic equation yields a monotonically increasing or decreasing function. Without modification it cannot represent the combination of positive and negative effects that is characteristic of biphasic dose-response relationships. To accommodate biphasic relationships, Brain and Cousens (18) modified the standard log-logistic equation by the addition of a parameter, φ, representing the slope of a hypothesized linear doseresponse relationship at low doses: y)ω+
R - ω + φx 1 + ( ⁄ x)β
(5)
In this formulation, the meanings of the original parameters of eq 3 are altered or lost. The parameter R (replacing Max) is the effect level asymptotically approached as the dose approaches 0, and parameter ω (replacing Min) is the effect level asymptotically approached as the dose approaches infinity. The parameter cannot be interpreted as ED50. Therefore, to facilitate determination of EDx (the dose at which x percent of the population is affected), Van Ewijk and Hoekstra (19) and Schabenberger et al. (17) reexpressed the Brain-Cousens model in terms of ED50 and EDx, respectively. The Brain-Cousens equation, in the original or a reparameterized version, has become the most widely used model of biphasic dose-response relationships (20, 21). It adequately fits many data sets but poorly represents dose-response relationships across a broad range of doses from substantial deficiency at low doses to toxicity at high-doses. It is particularly ill suited for affecters that are essential to the test organism (R ) 0), such as selenium (for animals) and other essential nutrients. In both reparameterizations of the Brain-Cousens model (19, 17), the baseline from which the EDx is calculated is the zero dose asymptote R; therefore, these models cannot be used for essential affecters. Additional models to represent biphasic dose-response relationships have been proposed, including polynomials (22–24), piecewise functions (25),a modification of the Brain-Cousens model in which the φx term is replaced with φ × exp(-1/xa) (26), and others (27, 28). All of these biphasic models suffer
from lack of biological rationale, and all of them fail to closely fit data that span a broad range of doses. To fill this void, we suggest the following approach. Proposed General Approach. We suggest that the positive and negative effects of biphasic affecters have distinct sensitivity threshold probability distributions, f(x). This is a logical inference if positive and negative effects are caused by different mechanisms or by different consequences of a single mechanism. It appears that generally one or the other of these conditions is met. For example, the mechanism by which high concentrations of fertilizers “burn” plants is quite different from the way lower concentrations provide nutritional benefits. Similarly, the positive and negative effects of ultraviolet radiation on human health involve quite different mechanisms, enabling the synthesis of vitamin D at low levels (7), but causing skin damage and cancer at high levels (8). Both the therapeutic and toxic effects of warfarin on mammals are due to a single anticoagulant mechanism, inhibition of the synthesis of vitamin K, but small doses of warfarin can prevent pathological clotting within the vessels of the circulatory system, whereas higher doses result in hemorrhage from those vessels. These antithrombotic and hemorrhagic effects have distinct sensitivity threshold distributions (29). If we choose a positive measure to represent both the positive and negative effects of a biphasic affecter, and we designate the distributions of thresholds of sensitivity for positive and negative effects p(x) and n(x) respectively, then in accordance with eqs 1 and 2, the positive and negative sides of the biphasic dose-response relationship would be described separately by the equations y ) P(x) and y ) 1 N(x) (Figure 1A and B). Assuming that the positive and negative mechanisms are statistically independent, we treat this as a case of conditional probability. Therefore, we multiplicatively combine the two sides into a single biphasic relationship (Figure 1C) yielding y ) P(x)(1 - N(x))
(6)
in which, for any given dose x, y is the total proportion of the population at or above the given level of effect. This is the proportion of the population for which dose x is sufficient for the effect to exceed the given threshold level on the lowdose side of the relationship, but not so great that the effect is reduced back below the same given threshold on the highdose side of the relationship. Equation 6 represents the simplest case of biphasic doseresponse relationships, one due to a single positive and a single negative mechanism. A biphasic affecter may cause positive and negative effects by multiple mechanisms, each with a distinct distribution of sensitivity thresholds. The approach suggested here can be elaborated to accommodate these possibilities (see Supporting Information). Because the dose-response curve is the integral of the probability distribution of sensitivity thresholds (eq 1), it follows that the implicitly assumed distribution of sensitivity thresholds underlying a dose-response curve is the first derivative of the dose-response curve (with respect to dose). In the Brain-Cousens family of models, the portion of each model representing the low-dose effect is linear. The first derivative of a linear function is a constant. Therefore, these models are implicitly founded on the assumption that the distribution of sensitivity thresholds for the low-dose effect of a biphasic affecter is constant with respect to dose. Such a straight horizontal distribution of sensitivity thresholds seems generally unlikely. In contrast, the wide success of the log-logistic model in (separately) representing positive as well as negative effects suggests that log-normal distributions often well describe the distributions of sensitivity thresholds underlying the low-dose positive effects as well as the highdose negative effects of biphasic affecters. Therefore, in the VOL. 42, NO. 4, 2008 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
1309
FIGURE 2. Illustration of the meanings of the parameters of the general biphasic model for hill-shaped curves (eq 8). the steepness of the falling (negative) slope, and parameter Dn approximates the dose at the midpoint of the falling slope. Equation 7 represents the common circumstances in which y is a positive measure of effect, such as survival; y approaches 0 at sufficiently great exposures to the affecter; and the affecter is essential to the expression of the effect so y also approaches 0 at very low exposures. Generalizing eq 7 to accommodate nonzero low and high dose asymptotes and a maximum y other than 1 yields R-ω+ y)ω+
Max - R 1 + (Up ⁄ x)βUp
1 + (Dn ⁄ x)βDn (βUp >0, βDn 0, βDn < 0) (7)
in which the parameter βUp represents the steepness of the rising (positive) slope, parameter Up approximates the dose at the midpoint of the rising slope, parameter βDn represents 1310
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 42, NO. 4, 2008
where parameters R and ω are defined as in eq 5, and the parameter Max is the theoretical maximum that would be approached asymptotically by the rising component of the equation in the absence of the descending component (or vice versa); the other parameters are the same as in eq 7 (Figure 2). In cases where positive and negative sensitivity threshold distributions overlap appreciably, the maximum value of y that is reached by the fitted model, a value we designate Peak, falls somewhat short of Max (see Supporting Information). Max is a good estimate of Peak, and Up and Dn are good estimates of ED50s only where both slopes are steep relative to their separation from each other (as in Figure 1C). If one of the slopes approaches 0, Peak approaches Max/ 2, and if both slopes approach 0, Peak approaches Max/4. In such situations, Up and Dn lose their interpretive utility, but useful ED50s based on Peak may be determined graphically. Parameters are defined above (eqs 7 and 8) for positive measures of effects, such as survival, which have hill-shaped curves. However, with some redefinition of parameters, these equations apply just as well for negative measures of effects, such as mortality, which have trough-shaped curves. For example, the generalized model (eq 8) adjusted for troughshaped curves is: Min - R + y ) R+
ω - Min 1 + (Up ⁄ x)βUp
1 + (x ⁄ Dn)βDn (βUp >0, βDn
