# Comparison of 2 glm logit models (anova.drc, F, Chisq, or Odds Ratio)

I would like to compare two logit models. One is a mixture of chemicals 'x' and 'y', the other is a model of just 'x'. I need help interpreting the 3 types of anova presented in the code/graphic below

However..

*A similar, but more complex situation: I would eventually like to compare independent tests of chemicals x and y, with a mixture of x and y. objective: look for possible synergistic effects (whether mixture is more potent, than either isolated chemical)

Would averaging the independent tests (x+y)/2 , making a model from averaged values, and then comparing this to a model made from the results of a true mixture of xy be folly?

Additional /supporting information can be found at bottom of post

#I prefer to use drc to create the glm
library(drc)

mod2 <- drm(probability ~ (dose), weights = total, data = mydata2, type ="binomial", fct=LL.2(names=c("Slope:b","ED50:e")))
mod10 <- drm(probability ~ (dose), weights = total, data = mydata10, type ="binomial", fct=LL.2(names=c("Slope:b","ED50:e")))

anova(mod2,mod10) #uses anova built into drc package
anova(mod2, mod10, test="F")
anova(mod2, mod10, test="Chisq")

#calculating Odd's Ratios (not included in output graphic below)
exp(coef(mod2))
exp(coef(mod10))


output:

Data:

mydata10 <-structure(list(dose = c(0, 0, 0, 3, 3, 3, 10, 10, 10, 27.5,27.5,     27.5, 50, 50, 50, 82.5, 82.5, 82.5), total = c(25L, 25L, 25L,     25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L,25L, 25L), affected = c(1, 0, 2.7, 3, 5.2, 5.3, 6.6, 7.3, 10.4,    14.3, 10.1, 14.7, 15.9, 12, 19.8, 20.4, 15.1, 19.4), probability = c(0.04,
0, 0.108, 0.12, 0.208, 0.212, 0.264, 0.292, 0.416, 0.572, 0.404,
0.588, 0.636, 0.48, 0.792, 0.816, 0.604, 0.776)), .Names = c("dose",
"total", "affected", "probability"), row.names = c(NA, -18L), class = "data.frame")

mydata2 <- structure(list(dose = c(0L, 0L, 0L, 3L, 3L, 3L, 10L, 10L, 10L,
30L, 30L, 30L, 55L, 55L, 55L, 85L, 85L, 85L), total = c(25L,
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L,
25L, 25L, 25L, 25L), affected = c(0, 0, 5.4, 4.2, 8.6, 9.6, 7.2,
13, 14, 17.2, 17.2, 16.8, 19.2, 15, 20.2, 21.8, 22.6, 16.2),
probability = c(0, 0, 0.216, 0.168, 0.344, 0.384, 0.288,
0.52, 0.56, 0.688, 0.688, 0.672, 0.768, 0.6, 0.808, 0.872,
0.904, 0.648)), .Names = c("dose", "total", "affected", "probability"
), row.names = c(NA, -18L), class = "data.frame")


Supporting info:

when comparing the models with anova.drc the log likelihoods are negative values, with a p of 0. When calculated with the argument anova(mod2, mod10, test="F"), the RSS (?) are positive values and the p is .1322. Which test should I use? Or would it be better to just compare the odd's ratios or some other method, e.g. (exp(coef(mod2)) - exp(coef(mod10)) ?

Plos1 article from drc author: Supporting Info Describes the use drc's anova with some examples

From the anova.drc Cran page it seems anova.drc should only be used to test nested models, but I am not certain of this - in the results below you can see it produces the same result as anova(mod2, mod10, test="Chisq").

From anova.drc Cran page"

Specifying only a single object gives a test for lack-of-fit, comparing the non-linear regression model to a more general one-way or two-way ANOVA model. If two objects are specified a test for reduction from the larger to the smaller model is given. (This only makes statistical sense if the models are nested, that is: one model is a submodel of the other model.)

stackoverflow recommends against use of r^2

stats.stackexchange thread recommended anova,test="Chisq" to compare 2 glm models

Apologies for strange formatting, Unable to post more than 2 links due to no reputation, apologies - I normally just use stackoverflow.

## 1 Answer

"I need help interpreting the 3 types of anova presented in the code/graphic below"

The test that you are trying to do is generally used for model selection on the same data set (i.e either mod2 or mod10). A significant p-value means you should not simplify your model from a 4 -parameter logistic model (LL.4) to a 2-parameter logistic model (LL.2). Because you've set the data as "Binomial", I'd imagine the default Chi-squared test would be appropriate.

" I would eventually like to compare independent tests of chemicals x and y, >with a mixture of x and y. objective: look for possible synergistic effects >(whether mixture is more potent, than either isolated chemical)"

This can be done in various ways. If you had a fully-crossed design, you could fit a GLM and look for an interaction. However, in your example you don't. A simple way to see if there is a difference is to compare ED50 values and see if the 95% CI overlap.

mod2.ed50=ED(mod2, 50, interval = 'delta')
mod10.ed50=ED(mod10, 50, interval = 'delta')


mod10.ed50 is sig. greater than mod2.ed50. Also look into Concentration Addition models.

A final note. You refer to you models as GLMs but I am not sure this is the case as I get a different model fit using the regular glm function. It might be that adding 'type="binomial" to the drm model just allows weighting of the data. I stand to be corrected.