I'm reproducing an example from Generalized, Linear, and Mixed Models. My MWE is below:
Dilution <- c(1/128, 1/64, 1/32, 1/16, 1/8, 1/4, 1/2, 1, 2, 4)
NoofPlates <- rep(x=5, times=10)
NoPositive <- c(0, 0, 2, 2, 3, 4, 5, 5, 5, 5)
Data <- data.frame(Dilution, NoofPlates, NoPositive)
fm1 <- glm(formula=NoPositive/NoofPlates~log(Dilution), family=binomial("logit"), data=Data)
summary(object=fm1)
Output
Call:
glm(formula = NoPositive/NoofPlates ~ log(Dilution), family = binomial("logit"),
data = Data)
Deviance Residuals:
Min 1Q Median 3Q Max
-0.38326 -0.20019 0.00871 0.15607 0.48505
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 4.174 2.800 1.491 0.136
log(Dilution) 1.623 1.022 1.587 0.112
(Dispersion parameter for binomial family taken to be 1)
Null deviance: 8.24241 on 9 degrees of freedom
Residual deviance: 0.64658 on 8 degrees of freedom
AIC: 6.8563
Number of Fisher Scoring iterations: 6
Code
anova(object=fm1, test="Chisq")
Output
Analysis of Deviance Table
Model: binomial, link: logit
Response: NoPositive/NoofPlates
Terms added sequentially (first to last)
Df Deviance Resid. Df Resid. Dev Pr(>Chi)
NULL 9 8.2424
log(Dilution) 1 7.5958 8 0.6466 0.00585 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Code
library(aod)
wald.test(b=coef(object=fm1), Sigma=vcov(object=fm1), Terms=2)
Output
Wald test:
----------
Chi-squared test:
X2 = 2.5, df = 1, P(> X2) = 0.11
Estimated coefficients are perfectly matching with the results given in the book but SE's are far apart. Based on LRT test the slope is significant but based on Wald and Z-test slope coefficient is insignificant. I wonder if I miss something basic. Thanks in advance for your help.
Best Answer
The main problem is that if you're going to use the ratio as your response variable, you should be using the
weights
argument. You must have ignored a warning about "non-integer #successes in a binomial glm" ...The LRT and Wald test results are still quite different ($p$-values of $4 \times 10^{-4}$ vs. $7 \times 10^{-10}$), but for practical purposes we can go ahead say they're both strongly significant ... (In this case (with a single parameter),
aod::wald.test()
gives exactly the same p-value assummary()
.)The Wald vs profile confidence intervals are also moderately different, but whether CIs [shown below] of (0.7,2.5) (Wald) and (0.9, 2.75) (LRT) are practically different depends on the particular situation.
Wald:
Profile: