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
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(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
pchisqmethod and the
lrtestmethod are equivalent for doing loglikelihood test? As I dunno how to use
lrtestfor 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
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  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
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 `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.