 Technical advance
 Open Access
 Open Peer Review
 Published:
Comparing lethal dose ratios using probit regression with arbitrary slopes
BMC Pharmacology and Toxicologyvolume 19, Article number: 61 (2018)
Abstract
Background
Evaluating the toxicity or effectiveness of two or more toxicants in a specific population often requires specialized statistical software to calculate and compare median lethal doses (LD_{50}s). Tests for equality of LD_{50}s using probit regression with parallel slopes have been implemented in many software packages, while tests for cases of arbitrary slopes are not generally available.
Methods
In this study, we established probitlog(dose) regression models and solved them by the maximum likelihood method using Microsoft Excel. The z and χ^{2}tests were used to assess significance and goodness of fit to the probit regression models, respectively. We calculated the lethal doses (LDs) of the toxicants at different significance levels and their 95% confidence limits (CLs) based on an accurate estimation of log(LD) variances. We further calculated lethal dose ratios and their 95% CLs for two examples without assuming parallel slopes following the method described by Robertson, et al., 2017.
Results
We selected representative toxicology datasets from the literature as case studies. For datasets without natural responses in the control group, the slopes, intercepts, χ^{2} statistics and LDs calculated using our method were identical to those calculated using PoloPlus and SPSS software, and the 95% CLs of the lethal dose ratios between toxicants were close to those calculated using PoloPlus. For datasets that included natural responses in the control group, our results were also close to those calculated using PoloPlus and SPSS.
Conclusion
This procedure yielded accurate estimates of lethal doses and 95% CLs at different significance levels as well as the lethal dose ratios and 95% CLs between two examples. The procedure could be used to assess differences in the toxicities of two examples without the assumption of parallelism between probitlog(dose) regression lines.
Background
In toxicological, entomological and environmental studies, doses of toxicants that kill a defined proportion of organisms, e.g., the median lethal dose (LD_{50}) which kills 50% of the population, are typically used as indicators of acute toxicity. Comparing the activities of different toxicants in a specific population or determining the relative susceptibilities of different populations to a single toxicant are common research goals. The relative potency, which assumes that the regression lines of the two toxicants being compared are parallel, provides a convenient comparison of the toxicities of two toxicants [1].
However, in practice, many regression lines are not parallel, particularly those derived from bioassays of toxicants with different modes of action, or from sameaction toxicants administered to populations with different resistance levels. The 95% confidence limits (CLs) of a lethal dose ratio can be calculated by estimating the slopes and intercepts of two probit regression lines and constructing their variance and covariance matrices. The 95% CLs of this ratio indicate whether the lethal doses of the two toxicants are statistically different from one another [2]. PoloPlus software, developed by Robertson et al. [3], separately analyzes the data for each substance using probit or logit models based on the joint probability of all observations and calculates lethal dose ratios and their CLs at different significance levels. IBM SPSS provides solution to calculate the lethal doses with 95% CLs based on probit or logit models, and also the relative median potency (RMP) assuming that the two regression lines are parallel [4].
In this study, we calculated lethal doses and 95% CLs of toxicants at different significance levels, as well as the lethal dose ratio and its 95% CLs for two toxicants, from probitlog(dose) regression models constructed using the maximum likelihood method in Microsoft Excel. The effectiveness of this method was compared with that of PoloPlus and IBM SPSS.
Methods
Construction of probitlog(dose) regression models for a single toxicant or population
For a population treated with serial doses (i) of a toxicant, in which n subjects were treated and r subjects exhibited a characteristic response to each dose, the empirical proportion (p^{*}) of responders was given by
where i = 1 to k and k indicated the number of toxicant doses.
If the characteristic response occurred in the control group (natural response) with proportion C, the proportions of responders were corrected using the Abbott equation for each treatment dose [5]:
The corrected proportion (p_{i}) was then converted to a probit value (y_{i}) [1]:
which was calculated as y_{i} = NORM.S.INV(p_{i}) in Excel.
A provisional regression line between y_{i} and the logarithm of the dose (x_{i}) was established:
In the regression equation, i = 1 to m, where m is the number of toxicant doses at which the corrected proportion was not equal to 1 or 0. The intercept (α_{0}) and slope (β_{0}) could be calculated by the leastsquares procedure and were retrieved using the INTERCEPT(y_{i}, x_{i}) and SLOPE(y_{i}, x_{i}) functions, respectively, in Excel.
We then calculated the expected probits (Y) for all dose sets, included those where the corrected proportion was 1 or 0:
In Eq. (5), i = 1 to k.
We next calculated the expected response proportion (P_{i}) for each dose set [1].
where Φ(Y_{i}) returned the cumulative probability of the standard normal distribution corresponding to (Y_{i}), obtainable using the NORM.S.DIST (Y_{i}) function in Excel, and C was the natural response proportion, if one existed, in Eq. (2).
The working probit (y_{i}) was calculated from the following eq. [1]:
where
An optimized set of expected probits was then derived from the linear regression equation of working probits weighted on x_{i}, with each y_{i} being assigned a weight, n_{i}w_{i}, where w_{i} was the weighting coefficient. This was calculated as previously described [1]
where C was the natural response proportion in Eq. (2).
The slope β of the working probitlog_{10}(dose) regression equation was
The intercept α of the working probit regression equation was
where ȳ was the average of y and \( \overline{x} \) was the average of x:
The standard error of β was [1]
and the standard error of α was [6]
The χ^{2} statistic of the probit regression equation was [1]
The significance level p of the χ^{2} statistic was calculated as the righttailed probability of the chisquared distribution (CHISQ.DIST.RT) with k – 2 degrees of freedom (d.f.).
A significant χ^{2} statistic (p < 0.05) might indicate either that the population did not respond independently or that the fitted probitlog(dose) regression line did not adequately describe the doseresponse relationship in the test samples.
To get an optimal fit of the probitlog_{10}(dose) regression, we substituted α and β for α_{0} and β_{0} and repeated the calculations of Eq. (5) to Eq. (15) until a stable χ^{2} appeared, indicating convergence. This procedure was a maximum likelihood (ML) method [1].
The significance of the slope was assessed using the z test [7],
If the absolute zvalue was less than 1.96, the regression slope was not significant and the data were excluded from further analysis. Similarly, we might test the significance of the intercept (α).
The heterogeneity factor h of the regression equation was calculated to adjust for large χ^{2}. h was defined as [1]
If h < 1, the model provided a good fit to the data. Otherwise, standardized residuals were plotted to identify outliers or other possible causes of poorness of fit [8]. Each residual defined the difference between the observed r_{i} and the expected response number (n_{i}P_{i}) for each dose. The residuals were standardized by dividing them by their standard errors, \( \sqrt{n_i{P}_i\left(1{P}_i\right)} \). For models providing a good fit, the standardized residuals fell mostly between −2 and 2 [8]. Standardized residuals distributed randomly showed no systematic patterns or tendencies toward positive or negative sign.
Calculation of the lethal doses of toxicants or populations and their 95% CLs
In this step, we first calculated the logarithms of the doses (θ_{π}) at which levels of interest (π) gave the expected response proportion:
where y_{π} was the π^{th} percentile of the probit distribution curve calculated in Excel using NORM.S.INV(π) for the probit distribution. For example, if π = 10, 50, 90 and 99, y_{π} was calculated as − 1.282, 0, 1.282 and 2.326.
The π^{th} lethal dose was then calculated as
The standard error of θ_{π}, σ(θ_{π}), was given by [1]
The 95% CL of the LD_{π} was then given as
t_{0.05, k − 2} returned the twotailed inverse of the Student’s tdistribution at α = 0.05 with d.f. = k  2 [T.INV.2 T(0.05, k2)].
The g value could be calculated to adjust if the confidence limits were valid. The g value was given as [9]:
If p (χ^{2}) was less than 0.15, t = 1.96 and h^{*} = 1; otherwise, h^{*} = h and t = t_{0.05, k − 2} [4]. If g exceeded 1, the CLs for the LD_{π} did not have practical importance [1].
The above steps were repeated to determine all parameters for the second toxicant for the same population, or the same toxicant in the second population.
Comparison of lethal dose ratios of two toxicants or populations
If there were l toxicants or populations in the experiment, then we compared the LD_{π} values of the first (as a reference) to those of others. We first calculated the difference between the log(doses) yielding the expected response proportions (π^{th} percentile) for toxicants or populations 1 and j (j = 2 to l), θ_{π1j} = θ_{π1}  θ_{πj}. Its standard error was given by [2]
The ratio of the two lethal doses was then given as
and the 95% CLs were
If the 95% CLs of this ratio excluded 1.0, the lethal doses of the two toxicants or populations were significantly different; otherwise, there was no evidence to reject the null hypothesis of equal LDs [2].
Test for parallelism of the two regression equations
Although the above procedures did not assume equal slopes of the two regression lines, the specific LD_{π} level used depended on the parallelism of the regression lines. To examine parallelism of the two regression lines, we used the ztest [10]:
If the absolute zvalue exceeded 1.96, the two regressions were nonparallel; otherwise, they were parallel.
Case studies
The above procedures might be executed on an Excel (version 2010 or higher) spreadsheet (provided as an Additional file 1). To compare the results of the ML procedure in Excel with those of PoloPlus and SPSS, we extracted bioassay data from the literature: (1) chrysanthemum aphids dosed with Rotenone, Deguelin, and a mixture of these two substances [11], (2) three populations, Fairfax, Pixley and Schaefer, of the pest bug “Wicked Witch of the West” dosed with deltamethrin [12], and (3) two populations, BugRes and BugLab, of Godfather larvae dosed with pyrethroid [2] (Table 1).
Results
Slopes, intercepts and significance testing of probitlog(dose) models fitted to the example data
When we implemented the ML procedure to solve the probitlog(dose) equations for the three sample data in Excel, for the datasets in which there was no natural response (e.g., Rotenone, Deguelin, Mixture, Fairfax and Schaefer), the slope (β) and intercept (α) estimates of the converged probitlog(dose) regression were identical to those calculated using PoloPlus and SPSS (with two methods, SPSS^{1} and SPSS^{2}, to include the natural response proportion, C, by inputting the value of C and calculating it from the data, respectively) (Table 2). The standard errors of both β and α, calculated by Eq. (13) and Eq. (14), were close but not identical to those calculated using PoloPlus and SPSS (Table 2). When the data sets included natural responses (e.g., Pixley, BugRes and BugLab), β and α, as well as their standard errors, were close to those produced by Poloplus and SPSS. The results of our method and PoloPlus were closer to those calculated using SPSS^{1} method than those calculated using SPSS^{2} method (Table 2, Bold items).
The probitlog(dose) regression model assumes a linear relationship between the logarithm of serial doses and the probit of the response proportions. When ztests (this study and SPSS) or the tratios (PoloPlus) were used to evaluate the significance of the regressions, all z values and tratios for both β and α estimates calculated using all four methods exceeded 1.96 (Table 2), indicating that all regression parameters were significant. If the zvalue for the slope was less than 1.96, the regression model would be insignificant and the dataset should be excluded from further analysis.
Goodnessoffits of the probitlog(dose) regressions
While ztests evaluated whether a linear relationship existed between the probits and the log(dose), χ^{2} tests are usually used to test the goodnessoffit of the log(dose)probit regression model. For datasets that did not include natural responses, the χ^{2} and h values calculated in this study were identical to those calculated using PoloPlus and SPSS (Table 3). When the datasets included natural responses, the χ^{2} and h values were close to those produced by Poloplus and SPSS. Again, the results of our method and PoloPlus were closer to those calculated using SPSS^{1} method than those calculated using SPSS^{2} method (Table 3, Bold items).
For some datasets, χ^{2} was not significant but h was greater than 1 (Table 3). When standardized residuals were plotted against log(doses), one or more outliers were observed (outside the bounds of − 2 to 2) in the Schaefer and BugLab data. For the BugLab data especially, the standardized residuals were not distributed randomly and showed a tendency toward positive sign (Fig. 1), indicating that this data should be fitted using other models [13].
LD_{10}, LD_{50}, LD_{90} and LD_{99} estimates with 95% CLs
We further compared the LD_{π}s and their 95% CLs calculated using these four methods. For datasets that did not include natural responses, the LD_{π}s calculated in this study were identical to those calculated using PoloPlus and SPSS, and the 95% CLs of LD_{π}s calculated using our method were close but not identical to those produced by PoloPlus and SPSS (Table 4). For datasets that included natural responses, the LD_{π}s and their 95% CLs were close but not identical to those calculated using Poloplus and SPSS. The results of our method and Poloplus were closer to those calculated using SPSS^{1} method than those calculated using SPSS^{2} method (Table 4, Bold items).
Comparison of lethal dose ratios between two samples
For datasets that did not include natural responses, the LD_{π} ratios calculated using our method were identical to those calculated using PoloPlus and their 95% CLs were also close. For datasets that included natural responses, LD_{π} ratios and their 95% CLs calculated using our method were similar to those calculated using PoloPlus (Table 5, Bold items). The LD_{50} ratios and their 95% CLs calculated using our method were closer to those calculated using PoloPlus than to the relative median potency (RMP) calculated using SPSS (Table 5).
When judged by whether the 95% CLs of lethal ratios included 1.0, all methods reached the same conclusions for toxicity differences between two samples (Table 5).
Comparisons of two regression slopes
Parallelism between paired regression equations was examined using ztests. The conclusions of our method for the five regression pairs were identical to those arrived at by PoloPlus and SPSS, which used χ^{2} tests (Table 6).
Discussion
Many methods have been developed to calculate the lethal or effective doses of toxicants and their confidence limits. Probit analysis, developed by Bliss [14] and improved by Finney [11], is one such commonlyused method. To calculate the parameters of the probitlog(dose) regression, Finney suggested fitting the regression line by eye as precisely as possible and obtaining parameters, such as slopes and intercepts, of the provisional regression line at the first stage. Thereafter, one calculates the working probits Y, and repeats this process with the new set of Y values; when the iterations converge, this gives a precise estimate of the linear regression parameters [1]. In this study, we calculated slopes and intercepts for the provisional regression line by the leastsquares procedure, and calculated working probits and performed the iteration procedure (ML) using the popular software program, Microsoft Excel. We obtained similar results to those obtained using PoloPlus and SPSS.
Several software packages, such as PoloPlus and SPSS, might be used to calculate the lethal doses and 95% CLs at different significance levels, and even test the equality of the lethal doses. Such professional statistical programs are difficult to handle for common toxicologists and environmental ecologists, and are easily abused. Excel in the Microsoft Office Package is the most popular statistical program around the globe. As to the Excel spreadsheet developed in this study, the users are easily to trace the procedure which is used to solve the regression equations, and calculate the CLs of a lethal dose and also the lethal dose ratios. They may further redevelop it easily according to their request.
χ^{2} values were used as indicators of the goodnessoffit of the probitlog(dose) regressions as the iteration proceeded. The equations
or
could also be applied [1]. When there were no natural responses in the datasets, these two equations, along with Eq. (15), gave the same results when the iterations converged, and these results were identical to those produced by PoloPlus and SPSS. When the datasets included natural responses, Eq. (27) always gave the smallest χ^{2} value, Eq. (28) always gave the largest value, while Eq. (15) gave an intermediate value which was closer to the output of PoloPlus and SPSS (data not shown). During iteration for some datasets, the χ^{2} values produced from all these three equations might increase [1]. Most regression models converged after several iterations, and we reported the results after 20 iterations, as this was the default maximum used by SPSS.
Strictly speaking, the 95% CLs of LD_{π} were the values of x for which the boundaries of the fiducial band attained the relevant value of y_{π}. The exact CLs of θ_{π} could be calculated by constructing the variance matrices of the slope (var(β)) and intercept (var(α)) and their covariance (cov(α,β)) matrices as follow [1, 9]:
It has been theorized that, in practice, the method for determining 95% CLs of LD_{π} most often performed sufficiently good based on a trustworthy value for the variance of θ_{π} as Eq. (20) [1, 15]. It was suggested that 95% CLs of LD_{π} could be calculated using the formula \( {10}^{\theta_{\pi}\pm 1.96\sigma \left({\theta}_{\pi}\right)} \) [15]. The results of this equation were close to those calculated using Eq. (29) when the dose number (k) was large (e.g., close to 10), while the CLs were much narrow than those calculated exactly using Eq. (29) when k was small. By contrast, the results given by Eq. (21) were nearer to those calculated exactly at different levels of k. The 95% CLs of LD_{π} calculated using PoloPlus were often identical to those calculated using SPSS when there was no natural response, with some exceptions (e.g., the Mixture and Fairfax data; Table 4, italic brackets, although the g values were not large for both of these cases).
While it is common to find estimates of LDs obtained from probit analyses in the toxicology literature, it is less common to find a hypothesis test procedure to determine whether estimated differences between LDs are statistically significant [16]. Relative potency has been frequently used [1, 4], but this method assumes the regression lines being compared are parallel. When the regression lines were parallel, the LDs and their 95% CLs for two toxicants calculated from the two datasets simultaneously were similar to those calculated from the datasets separately. However, when the regression lines were not parallel, the LDs and their 95% CLs calculated from the two datasets simultaneously were quite different from those calculated from the datasets separately.
In cases where the data are suggestive of a trend toward significant differences between LD_{50}s, the use of nonoverlapping CLs for LD_{50} values has frequently been proposed as a criterion for assessing significance, while use of this criterion is thought to be conservative [17, 18]. An alternative method involves calculating the variances of θ_{π} using the deltamethod:
calculating the ratio of the LDs as in Eq. (24), then calculating the 95% CLs of the ratio as in Eq. (25) [2]. If the 95% CLs of the ratio include 1.0, the LDs of the two samples are not significantly different. We followed this procedure in this study, but we calculated the standard error of θ_{π} as in Eq.(20) by the maximum likelihood procedure. We obtained 95% CLs of the LD ratio similar to those obtained using PoloPlus.
Biologically, the slope of a probit or logit regression line represents the change in the proportion of responders per unit change in dose. Toxicological evidence suggested that the slope of a dose–response regression line reflected host enzyme activity [19]. Thus, nonparallel lines might indicate different modes of action of the two toxicants. Parallelism between regression pairs was essential for determining the level at which to compare the effects of two toxicants. Generally, there were three main categories of parallelism: (i) the two regression lines were statistically parallel (e.g., Fairfax vs Pixley; Fig. 2a); (ii) the two regression lines were not statistically parallel but did not cross within the dominant region (20–80%) of the response proportions (e.g., Rotenone vs Deguelin; Fig. 2b); and (iii) the two regression lines crossed around the median lethal dose (e.g., BugRes vs BugLab; Fig. 2c). In the first case, reporting the LD_{50}s of the two toxicants and their ratios was sufficient. In the second case, one should report both LD_{50}s and LD_{90}s (and/or LD_{10}s) and their ratios. In the third case, reporting the ratios of LD_{10}s, LD_{50}s, LD_{90}s is meaningless, but the significance of difference between the two slopes should be valid.
Conclusions
We successfully developed a method to calculate the lethal doses of a toxicant at different significance levels, and compare lethal dose ratios using probitlog(dose) regression by the ML procedure implemented in Microsoft Excel. Lethal doses calculated using this method at different significance levels, as well as lethal dose ratios with their 95% CLs, were identical or close to those calculated using PoloPlus and SPSS. When judged by whether the 95% CLs of the lethal ratios included 1.0, all methods reached the same conclusions regarding toxicity differences between two samples.
Abbreviations
 95% CLs:

95% confidence limits
 LD_{50} :

Median lethal dose
 ML:

Maximum likelihood
 RMP:

Relative median potency
References
 1.
Finney DJ. Probit Analysis. 3rd ed. Cambridge: Cambridge University Press; 1971.
 2.
Robertson JL, Jones MM, Olguin E, Alberts B. Bioassays with arthropods. 3rd ed. Boca Raton, FL: CRC Press, Taylor & Francis Group; 2017.
 3.
LeOra Software. PoloPlus, POLO for windows. Petaluma, CA: LeOra Software, 107 B St, 94952; 2007.
 4.
SPSS Inc., IBM Corp. IBM SPSS statistics 20.0. Chicago, IL: SPSS Inc.; 2011.
 5.
Abbott WS. A procedure of computing the effectiveness of a toxicant. J Econ Entomol. 1925;18:265–7.
 6.
Berkson J. Minimum χ^{2} and maximum likelihood solution in terms of a linear transform, with particular reference to bioassay. J Am Stat Assoc. 1949;44(246):273–8.
 7.
Walpole RE, Myers RH, Myers SL, Ye KY. Probability & Statistics for engineers & scientists, 9th ed. Boston: Prentic Hall; 2012. p. 135.
 8.
Preisler HK. Assessing insecticide bioassay data with extrabinomial variation. J Econ Entomol. 1988;81(3):759–65.
 9.
Fieller EC. Some problems in interval estimation. J R Statist Soc. 1954;B16:175–85.
 10.
Clifford CC, Petkova E, Haritou A. Statistical models for comparing regression coefficients between models. Am J Sociol. 1995;100:1261–93.
 11.
Finney DJ. Probit analysis. Cambridge. England: Cambridge University Press; 1952.
 12.
Robertson JL, Russell RM, Preisler HK. Bioassays with arthropods. 2nd ed. Boca Raton, FL: CRC Press; 2007.
 13.
Robertson JL, Preisler HK. Bioassays with arthropods. Boca Raton, FL: CRC Press; 1992.
 14.
Bliss CI. The determination of the dosethe proportion responding curve from small numbers. Quart J Pharma Pharmacol. 1938;11:192–216.
 15.
Hayes WJ, Kruger CL. Haye’s principles and methods of toxicology. 6th ed. Boca Raton: CRC Press; 2014.
 16.
Jeske DR, Xu HK, Blessinger T, Jensen P, Trumble J. Testing for the equality of EC50 values in the presence of unequal slopes with application to toxicity of selenium types. J Agr Biol Environ Stat. 2009;14(4):469–83.
 17.
Schenker N, Gentleman JF. On judging the significance of differences by examining overlap between confidence intervals. Am Stat. 2001;55:182–6.
 18.
Payton ME, Greenstone HH, Schenker N. Overlapping confidence intervals or standard error intervals: what do they mean in terms of statistical significance? J Insect Sci. 2003;3:34–9.
 19.
Kuperman AS, Gill EW, Riker WF. The relationship between cholinesterase inhibition and druginduced facilitation of mammalian neuromuscular transmission. J Pharmacol Exp Ther. 1961;132:65.
Funding
We appreciate the support of the National Key Research and Development Program of China (2017YFD0201206), and the WIV “OneThreeFive” strategic programs (Y602111SA1).
Availability of data and materials
Additional file 1 is available via a link: https://figshare.com/s/f94393f752fcc15faea7.
Author information
Affiliations
Contributions
XS and CL conceived and designed the study. CL edited the Excel file, analyzed the data and prepared the manuscript. XS revised the manuscript. Both authors read and approved the final manuscript.
Corresponding author
Correspondence to Xiulian Sun.
Ethics declarations
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Competing interests
The authors declare that they have no competing interests.
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Additional file
Additional file 1:
Calculation of LDs and their ratios.xlsx (344 kb), which requires Microsoft Excel 2010 or higher. It is available via a link: https://figshare.com/s/f94393f752fcc15faea7. (XLSX 1493 kb)
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Received
Accepted
Published
DOI
Keywords
 Toxicity
 Probit regression
 Lethal dose ratio
 Maximum likelihood