### Likelihood ratio test in R

• Suppose I am going to do a univariate logistic regression on several independent variables, like this:

``````mod.a <- glm(x ~ a, data=z, family=binominal("logistic"))
mod.b <- glm(x ~ b, data=z, family=binominal("logistic"))
``````

I did a model comparison (likelihood ratio test) to see if the model is better than the null model by this command

``````1-pchisq(mod.a\$null.deviance-mod.a\$deviance, mod.a\$df.null-mod.a\$df.residual)
``````

Then I built another model with all variables in it

``````mod.c <- glm(x ~ a+b, data=z, family=binomial("logistic"))
``````

In order to see if the variable is statistically significant in the multivariate model, I used the `lrtest` command from `epicalc`

``````lrtest(mod.c,mod.a) ### see if variable b is statistically significant after adjustment of a
lrtest(mod.c,mod.b) ### see if variable a is statistically significant after adjustment of b
``````

I wonder if the `pchisq` method and the `lrtest` method are equivalent for doing loglikelihood test? As I dunno how to use `lrtest` for univate logistic model.

@Gavin thanks for reminding me, as comparing with stackoverflow, I need to spend more time to "digest" the answer before deciding whether the answer is appropriate or not, anyway, thanks again.

I would not recommend using waldtest from lmtest. Use the aod package for model testing. Its far more straightforward. https://cran.r-project.org/web/packages/aod/aod.pdf

`epicalc` was removed (source). An alternative could be `lmtest`.

10 years ago

Basically, yes, provided you use the correct difference in log-likelihood:

``````> library(epicalc)
> model0 <- glm(case ~ induced + spontaneous, family=binomial, data=infert)
> model1 <- glm(case ~ induced, family=binomial, data=infert)
> lrtest (model0, model1)
Likelihood ratio test for MLE method
Chi-squared 1 d.f. =  36.48675 , P value =  0
> model1\$deviance-model0\$deviance
[1] 36.48675
``````

and not the deviance for the null model which is the same in both cases. The number of df is the number of parameters that differ between the two nested models, here df=1. BTW, you can look at the source code for `lrtest()` by just typing

``````> lrtest
``````

at the R prompt.

thanks, and I just found that i can use glm(output ~ NULL, data=z, family=binomial("logistic")) for creating a NULL model, and so i can use the lrtest afterwards. FYI, thanks again

@lokheart `anova(model1, model0)` will work too.

@lokheart `glm(output ~ 1, data=z, family=binomial("logistic"))` would be a more natural null model, which says that `output` is explained by a constant term (the intercept)/ The intercept is implied in all your models, so you are testing for the effect of `a` after accounting for the intercept.

Or you can do it "manually": p-value of the LR test = 1-pchisq(deviance, dof)

Content dated before 6/26/2020 9:53 AM
• {{ error }}