Environmental Health Perspectives 106, Supplement 6, December 1998
A Nonlinear Isobologram Model with
Box-Cox Transformation to Both Sides
for Chemical Mixtures
Ding G. Chen1 and Joel G. Pounds2
1Pacific Biological Station, Department of Fisheries and Oceans, Nanaimo, British Columbia, Canada;
2Institute of Chemical Toxicology, Wayne State University, Detroit, Michigan
Abstract
The linear logistical isobologram is a commonly used and powerful graphical and statistical tool for analyzing the combined effects of simple chemical mixtures. In this paper a nonlinear isobologram model is proposed to analyze the joint action of chemical mixtures for quantitative dose-response relationships. This nonlinear isobologram model incorporates two additional new parameters,
ymin and
ymax, to facilitate analysis of response data that are not constrained between 0 and 1, where parameters
ymin and
ymax represent the minimal and the maximal observed toxic response. This nonlinear isobologram model for binary mixtures can be expressed as:
In addition, a Box-Cox transformation to both sides is introduced to improve the goodness of fit and to provide a more robust model for achieving homogeneity and normality of the residuals. Finally, a confidence band is proposed for selected isobols, e.g., the median effective dose, to facilitate graphical and statistical analysis of the isobologram. The versatility of this approach is demonstrated using published data describing the toxicity of the binary mixtures of citrinin and ochratoxin as well as a new experimental data from our laboratory for mixtures of mercury and cadmium. -- Environ Health Perspect 106(Suppl 6):1367-1371 (1998).
http://ehpnet1.niehs.nih.gov/docs/1998/Suppl-6/1367-1371chen/abstract.html
Key words: chemical mixtures, dose-response modeling, maximum likelihood estimation, risk assessment, toxicology
This paper is based on a presentation at the Conference on Current Issues on Chemical Mixtures held 11-13 August 1997 in Fort Collins, Colorado. Manuscript received at EHP 17 February 1998; accepted 24 June 1998.
This study was supported in part by the Agency for Toxic Substances and Disease Registry cooperative agreement ATU541481.
Address correspondence to D.G. Chen, Pacific Biological Station, Department of Fisheries and Oceans, 3190 Hammond Bay Road, Nanaimo, British Columbia V9R 5K6, Canada. Telephone: (250) 756-7341. Fax: (250) 756-7053. E-mail: chend@dfo-mpo.gc.ca
Abbreviations used: BCTBS, Box-Cox transformation to both sides; CI, an interaction index for mixtures, proposed by Berenbaum; ED50, median effective dose; ED100p, dose or concentration to produce the effect in 100p% of the population; LDH, lactate dehydrogenase; MK2, monkey kidney cells; TBS, transformation to both sides.
Epidemiologists and toxicologists frequently conduct a type of toxicity experiment to study the joint action of chemical mixtures in which the biologic response to two or more chemicals is recorded. The issues of experimental design and appropriate analysis of the resultant data are very important and also controversial. The most recent references can be found in Ramamoorthy et al. (
1) and Kaiser (
2).
In the case of binary mixtures, the resultant observations consist of data in three-dimensional space (d1, d2, y), where d1 and d2 are the dose levels for the two chemicals and y the response. In this paper, special attention will be paid to the quantitative dose-response relationship where y is a continuous response, whereas a brief discussion will be given in "Discussion" for the qualitative dose-response relationship. The systematic part of the dose-response model would be
y=f(d1,d2,ß) [1]
based on the theory and knowledge from different applications of modeling with parameters of vector ß to be estimated from the dose-response relationship, f. The main objective of the statistical modeling is to estimate the parameters and make statistical inference, especially for the parameters describing the joint action of the mixtures. Therefore, an accurate model to describe the mechanism of mixtures is essential to statistically and toxicologically identify the interaction.
Carter et al. (3) applied the logistical model below for the mixture joint action and discussed the isobologram:
[2]
where p is the mortality rate and ß12 is used to describe the interaction. This model can be easily understood from statistical, biologic, and toxicologic points of view. For this model, however, the response is required to be from 0 to 1. Even if the response is the percentage (e.g., y ranging from 0 to 100), then the simple transformation of p=y/100 should be used. However, many experimental observations are not constrained in the range of 0 to 100, so Carter's model would present some difficulties.
Barton et al. (4) proposed a series of nonlinear models to characterize the chemical joint action. Specifically, the simple similar model is approximately equivalent to the dose/concentration additivity; the independent action model is approximately equivalent to the response additivity and the nonadditivity model is for the interaction. In these models two additional parameters, ymin and ymax, were introduced to overcome the difficulties of the response range not being uniformly within 0 to 100. Although these models could be easily understood by professional statisticians, the complexity of the computational code was a burden for many biologists and toxicologists who implemented these models.
Kodell and Pounds (5) considered an envelope for the mixtures' joint action and examined the dose additivity and response additivity for the dose and log dose scales. The results were extended by Razzaghi and Kodell (6) using a Box-Cox transformation (7). We suggest that this model is practical for characterizing toxicologic interactions.
In this paper, we expand the linear logistical isobologram of Carter et al. (3) to include the features of the Barton et al. (4) model, incorporating the two parameters ymin and ymax to overcome the difficulties in dealing with different data types in the analysis. Based on this nonlinear isobologram model, a confidence band is proposed for the observed isobol and the median effective dose (ED50) line for a statistically sound analysis of the isobologram. We also introduce a transformation to both sides (TBS) technique to achieve the homogeneity and normality of the residuals, where the Box-Cox transformation family as a special case is concentrated. Two real data sets demonstrate the application of these newly developed procedures.
Nonlinear Isobologram Model
We will concentrate on the binary mixture. For a binary chemical mixture, the observations from the experiment include a triplet such as (d1, d2, y), where d1 and d2 are the dose levels from two chemicals and y is the response from the mixtures. The nonlinear isobologram model for binary mixtures in the dose response [model 1] is proposed as
[3]
where
ymin is the parameter associated with the minimum response,
ymax is the parameter associated with the maximum response,
ß0 is the parameter associated with the effect of the control dose levels,
ß1 is the parameter associated with the effect of the first chemical,
ß2 is the parameter associated with the effect of the second chemical, and
ß12 is the parameter associated with the interaction of the two chemicals.
Simple mathematical manipulation can result:
[4]
This manipulation shows that Carter's logistic model [2] is the special case when parameters ymin and ymax are fixed to be 0 and 100, respectively.
Model [3] can be used to characterize the interaction of two drugs or chemicals. Following the discussion of an interaction index for mixtures proposed by Berenbaum (CI) (8), the interaction of chemicals can be claimed as synergism, antagonism, and additivity if the index, defined as:
[5]
is the <1, >1, =1, respectively. Here in the definition [5], d1 and d2 are the doses of chemical number 1 and number 2 in combination that gives a 100p% response (ED100p%). ED100p(1) and ED100p(2) are the respective single agent ED100ps.
In models [3] or [4], the ED100ps for each chemical can be easily obtained as
where p=(y-ymin)/(ymax-ymin). Therefore, with the combination of the interaction index definition [5] and model [4], and with little algebraic manipulation, the interaction index can be expressed as
It can be followed that the two chemicals interact synergistically, antagonistically, and additively if the parameter ß12 is greater than, less than, or not different than zero, respectively. This means that the CI can be also characterized by the interaction parameter ß12.
The nonlinear least square regression or maximum likelihood estimation can be used to fit the model to the observed data from the toxicologic experiment and to test the goodness of fit under the assumptions of homogeneity and normality of residuals. If the model fits the data appropriately, it is then possible to test the hypothesis of additivity (ß12=0) and also to get the confidence interval to test for synergism, antagonism, or additivity.
Analysis of Isobolograms
The analysis of isobolograms is a commonly used approach to identify the joint action of chemicals. This approach is conducted by comparing the observed isobols to the line of additivity, which is frequently constructed from the ED50 line associated with the mixtures. Specifically, if the observed isobol is bent under this line, synergism can be declared, and if the observed isobol is bent above the line, antagonism can be declared. When the observed isobol does not deviate from the line, the joint action can be claimed as additivity. In this sense, a statistically sound procedure should be developed and followed when the observed isobol is compared to the ED50 line to conclude synergism, antagonism, and additivity.
A general approach can be developed by Fieller's theorem (9) for any isobols under any proportion of response p, based on the nonlinear isobologram model proposed in "Nonlinear Isobologram Model." In fact, after the model [3] is fitted and tested to be appropriate, the model can be rearranged as in [4]. For any response y, the proportion of response p can be calculated as p=y-
min/max-min.
Therefore, the observed isobols d2Ip corresponding to ED100p for p can be constructed as
[7]
and its associated ED100p line d2Lp as
[8]
For the purpose of identifying the joint action for synergism, antagonism, and additivity, the observed isobol d2Ip for proportion 100p% can be selected to compare with the confidence band for the ED100p line d2Lp constructed by Fieller's theorem. If the observed isobol d2Ip is within the confidence band of d2Lp, the additivity is claimed. If the observed isobol d2Ip is not completely within the confidence band for d2Lp and is bent below the line, the synergism is claimed; if the observed isobol d2Ip is not completely within the confidence band of d2Lp and is bent above the line, the antagonism is claimed. Figure 1 illustrates the application of this procedure for data from "Binary Mixture of Cadmium and Mercury."
|
Figure 1. Isobologram for the cytotoxic interaction of cadmium and mercury to primary cultures of rat hepatocytes. The ED50 line (solid line) is the predicted ED50 with no interaction (ß12=0) and is bounded by the 95% confidence bands (---). ED50 isobol (dashed line) is the observed ED50 isobol.
|
In the application of any dose-response model [1], care should be taken for the assumptions of normality and homogeneity of the residuals; otherwise, the estimates and the statistical inference could be invalid. In the case of nonnormality and heteroscedasticity, an alternative model or appropriate transformation should be used to reanalyze the data. TBS is promising, particularly when the data come from a large range with a large quotient between the maximum and the minimum response, which is often the case for dose-response data structure. The TBS model takes the form:
[9]
where h (·, 
) is the transformation family with parameter associated with the transformation, and 
is assumed to have mean 0 and variance 1.
The Box-Cox (7) power transformation family is well known, with the form:
[10]
Barton et al. (4) used a log transformation to both sides for the nonlinear models to satisfy the homogeneity of variance assumption. It can be noted that the Box-Cox family [10] is continuous in and with 
0 being the logarithm transformation. Therefore, the log transformation in Barton et al. (4) is just a special case of the Box-Cox transformation family. It can be easily seen that if **1, the TBS model [10] is back to the untransformed model [3]. We can then consider the untransformed nonlinear model [3] and Barton's log transformation to both sides as the two extremes and by allowing in [10] to vary in [0,1], a whole family of transformation is obtained.
This continuous family depends on only a single parameter, , which is estimated from the data, as well as the vector of parameter ß. The criterion for determining the appropriate values of the transformation parameter is to minimize the sum of squares for the error for the nonlinear model based on the transformation, or a maximum likelihood estimation can be used to estimate the parameter along with the ßs.
Barton's Example
Barton et al. (4) proposed a series of nonlinear models for the effect of chemical mixtures. A study was conducted for the cytotoxic effect of citrinin and ochratoxin on renal cortical slices. Within the three models of simple similar action, independent joint action, and nonadditive model, the independent joint action model was indeed the best fit and therefore the most appropriate model. Further analysis suggested that the two chemicals acted slightly synergistically. Razzaghi and Kodell (6) reexamined these data using the procedures from Kodell and Pounds (5) with Box-Cox transformation to integrate the regression to dose and log dose. It was found that both the concentrate additivity and the response additivity provided almost the same fit; therefore these two chemicals are mathematically indistinguishable from the dose-response function. The experimental data can be found from Barton et al. (4) and were reproduced in Razzaghi and Kodell (6).
If the model [4] in Carter et al. (3) is used for the data, the parameter estimates and the associated statistical inference can be made under the linear regression setting. The values for ymax and ymin, 5 and 0.95, respectively, are selected based on the observed maximum and minimum values of the response (4.8 and 0.97, respectively). The R2=0.74 with p-value=0.0001 from the analysis of variance indicates that the model fit is highly statistically significant. The parameter estimates are ߈0=0.5997, ߈1=-0.0184, ߈2=-0.0056, and ߈12=-0.00001, respectively, with the p-value of t-test for ß12=0 to be 0.287, indicating that parameter ß12 is not statistically significantly different from zero. Based on these results and without checking the model assumptions, we would conclude that these two chemicals act additively. The test of normality with p-value of 0.0374 and a systematical pattern in the residual plot, however, indicate that basic assumptions of the model [4] are violated and therefore the parameter estimates and the statistical inference about ß12 may not be valid.
We then fit the nonlinear isobologram model [3] to the data. The fit is highly statistically significant, as the R2=0.9685 and p<0.0001. The Pearson's 
2 goodness-of-fit test for the standard residual for standard normal distribution yields p-value of 0.9391, confirming the normality and homogeneity assumptions. For this data, the Box-Cox transformation to both sides (BCTBS) is unnecessary because the assumptions of normality and homogeneity are satisfied.
In this situation, the parameter estimates are min=1.4612, max=5.3361, 
0=0.8783, 1=-0.0687, 2=-0.0080,
12=0.0001, and an asymptotic 95% confidence interval of (-0.00015, 0.00035) for ß12.
This analysis demonstrates that the chemical mixture of citrinin and ochratoxin may act toxicologically in slight synergism based on the positive parameter estimate, but in fact they are additive statistically and toxicologically. This analysis coincides with the results from Barton et al. (4) and Razzaghi and Kodell (6).
Binary Mixture of Cadmium
and Mercury
Cadmium and Hg are toxic metals to which human populations are often exposed concurrently in the workplace or from hazardous waste sites in the environment. To evaluate this exposure, we conducted an experiment in our laboratory. A rhesus monkey renal tubular epithelial cell line, LLC-monkey kidney (MK2), was selected as the principal cell system because the kidney is the major organ of trace element metabolism and storage, the results may be readily applied to in vivo animal and human studies, and the results can be readily extended to cultures of human tissues and other species facilitating the extrapolation of toxicant joint action principles from animal models to humans. The MK2 cells express many of the phenotypic characteristics of a tubular epithelial cell. Triplicate cultures of MK2 cells were exposed to Cd and Hg in a 1:2 mixture. After 24-hr treatment, the release of the soluble enzyme lactate dehydrogenase (LDH) is measured kinetically and expressed as the fraction of the total LDH activity released (response) for each culture. Specifically, the response is defined as follows:
where Lr is LDHreleased or the enzyme activity in a 200-µl aliquot, Lt is LDHtotal or the total activity following addition of 200 µl Triton X100, and Lb is LDHbackground or the endogenous activity of the serum-containing medium, and df is the dilution factor. The LDH release is a sensitive reliable measure of irreversible cell injury.
The data are given in Table 1 and the data analysis is summarized in Table 2.
Carter's model [4] is applied to this data from two aspects by using ymax=100 and ymin=0 or 1. The model with ymax=100 and ymin=1 and the model with ymax=100 and ymin=0 are noted in Table 2. The differences in parameter estimates and conclusions by changing just one parameter, ymin, from 0 to 1 are remarkable. Statistically this model is not robust or is too sensitive to model parameters, especially ymax and ymin. This indicates that when Carter's model is used, caution should be taken.
The nonlinear model [3] in Table 2 obviously has the advantage to estimate the parameters ymax and ymin itself. These parameters can be estimated from the data along with the ß parameters. When the nonlinear model is used the model is highly significant, but the test for assumptions failed, indicating that an alternative model or a transformation should be used.
We then apply the BCTBS to these data. The estimate for parameter is not statistically significantly different from zero, which is equivalent to a log TBS. In this case, the model is statistically significant and the test of assumptions cannot be rejected. The p-value for the test of ß12=0 is 0.059, indicating that the conclusion of additivity or slight synergism for chemicals Cd and Hg can be made. Data analysis is summarized in Table 2. The isobologram analysis supports the conclusion (Figure 1).
In Table 2, P(goodness of fit) is the p-value of Pearson's 2 goodness of fit test. To test that the standard residual of model fit is standard normal distribution, the parameter estimates are associated with its standard error (in parentheses).
The identification of synergistic or antagonistic interactions, e.g., environmental estrogens, may have profound impact on legislation, regulation, and distribution of scientific resources involving public and environmental health (10,11). Therefore, the accurate and appropriate characterization of toxicant interactions is a necessary and critical part of the hazard identification of simple mixtures (12).
Fundamentally, the characterization of chemical interactions consists of the rejection or acceptance of a single model or the selection of a particular model from among competing models. These decisions are based on statistical and/or graphical information involving regression equations, each with underlying implicit and explicit mathematical, statistical, and data structure assumptions (13). As such, the choice of model, criteria for goodness of fit, data transformation, and method of parameter estimation are important to support a mathematically and biologically sound decision. The statistical isobolographic approach is a useful and powerful tool to identify interactions. The application of isobolographic and other approaches to the analysis of interactions, however, often do not adequately consider the fitting of the data to the model (14,15).
In this paper, a nonlinear isobologram model was developed to extend the linear logistical isobologram of Carter et al. (3) to include the model features of Barton et al. (4), incorporating two parameters ymin and ymax to make the isobolographic analysis more versatile. This model [3] is superior to Carter's logistic model because the parameters ymin and ymax can be estimated from the data together with other model parameters. This overcomes the sensitivity of different choice of these parameters and therefore is robust on different specifications of model parameters. The BCTBS technique is recommended, especially when the model assumptions are violated. The newly developed model and its comparison to the model of Carter et al. (3) was implemented in SPLUS software (16). The new model, will soon be a module of the software named MAS (Mixture Analysis System), and will be used specifically for chemical mixtures research.
In the fitting of models [3] and [9], a scheme to find good starting values for the parameters is critical. This selection can be obtained from the following procedure. The initial values for parameters ymin and ymax can be chosen from the observed minimum and maximum values from the response. By using these two values for model [3], the initial values for ßi i=0, 1, 2, 12 can be obtained from linear regression [4]. The parameter * in the model BCTBS [9] can take any value from 0 to 1.
The confidence band for the ED50 line d2Lp can be derived from Fieller's theorem (9), which can be used visually to evaluate the effect of chemical mixtures. Specifically, this confidence band is plotted with the observed 50% isobol to see whether the observed 50% isobol is within the confidence band to make conclusions about synergism, antagonism, and additivity. Figure 1 illustrates the application of this confidence band structure for the observed ED50 isobol and the corresponding ED50 line.
We presented the model [3] for the interaction of only two chemicals. In fact, this model can be easily modified for the analysis of any number of drugs or chemical mixtures. Because the interaction of more chemicals is more toxicologically complicated to classify, the discussion would be complex. For example, in the mixtures of three drugs or chemicals, the nonlinear model [3] would be
[11]
where g(d1,d2,d3,ß)=ß0+ß1d1+ß2d2...[11]
In this situation, CI is defined as:
[12]
and can be also expressed as
where
and
It can be seen that whether the CI is greater than, less than, or equal to zero is dependent on the expression m2, which is a function of the parameters and also the magnitude of the three dose levels. In the case that all parameter estimates are statistically greater than, less than, or equal to zero, the conclusion of synergism, antagonism, or additivity of these three drugs or chemicals can be made.
It has been suggested that a good contribution could be made in this area by using a generalized linear model approach or the iteratively reweighted least square algorithm to fit the dose-response relationship. We have examined these approaches for the quantal dose-response relationship where the responses are qualitative, such as in mortality rates for rats (17).
In all, we are encouraged by the performance of this nonlinear isobologram model and believe that it will provide a plausible method for risk assessment, toxicology, epidemiology, and environmental research.
References and Notes
1. Ramamoorthy K, Wang F, Safe S, Norris JD, McDonnell DP, Gaido KW, Bocchinfuso WP, Korach KS. Potency of combined estrogenic pesticides. Science 275:405-406 (1997).
2. Kaiser J. Synergy paper questioned at toxicology meeting. Science 275:1879 (1997).
3. Carter WH, Gennings C, Staniswalis JG, Campbell ED, White KL. A statistical approach to the construction and analysis of isobologram. J Am Coll Toxicol 7:963-973 (1988).
4. Barton CN, Braunberg RC, Friedman L. Nonlinear statistical models for the joint action of toxins. Biometrics 49:95-105 (1993).
5. Kodell RL, Pounds JG. Assessing the toxicity of mixture of chemicals. In: Statistics in Toxicology (Krewski D, Franklin C, eds). New York:Gordon Breach, 1991;559-591.
6. Razzaghi M, Kodell RL. Box-Cox transformation in the analysis of combined effects of mixtures of chemicals. Envirometrics 3(3):319-334 (1992).
7. Box GEP, Cox DR. An analysis of transformations. J R Stat Soc Ser B 26:211-252 (1964).
8. Berenbaum MC. Criteria for analyzing interactions between biologically active agents. Adv Cancer Res 35:269-333 (1981).
9. Finney DJ. Probit Analysis. Cambridge, MA:Cambridge University Press, 1971.
10. National Academy of Sciences. Complex Mixtures. Methods for In Vitro Toxicity Testing. Washington:National Academy Press, 1988.
11. Sexton K, Beck BD, Bingham E, Brain JD, DeMarini DM, Hertzberg RE, O'Flaherty EJ, Pounds JG. Chemical mixtures from a public health perspective: the importance of research for informed decision making. Toxicology 105:429-491 (1995).
12. Mauderly JL. Toxicological approaches to complex mixtures. Environ Health Perspect 101(Suppl 4):155-165 (1993).
13. Jackson JE. Comparison of a class of regression equations. Am J Physiol 246:R271-R276 (1984).
14. Berenbaum MC. Synergy, additivism, and antagonism in immunosuppression. A critical review. Clin Exp Immunol 28:1-18 (1977).
15. Berenbaum MC. What is synergy? Pharmacol Rev 41:93-141 (1989).
16. SPLUS User's Manuals, Version 3.3 for Windows. Seattle, WA:StatSci Division Mathsoft, Inc. 1995.
17. Chen DG, Pounds, JG, Mumtaz MM. Unpublished data.
Last Updated: December 7, 1998
