How to calculate p-values in a glmer with multiple explanatory variables using likelihood ratio tests in R I would like to ask a quite general question about how to use likelihood ratio tests to calculate p-values in a glmer. The goal is to report these in a paper to show whether each of the explanatory variables has an effect on the response.
Let's say I have a final model like this:
m_final <- glmer(x ~ a + b + c, data = dat, family = family)
Now, I'd like to report if a, b and c affect the response. My understanding until now was that I simply create three new reduced models, each with ONE of the variables dropped. Then, I'd compare each of the reduced models to the full model in order to see whether the dropped variable has a significant influence. Like this:
m_final <- glmer(x ~ a + b + c + (1|d), data = dat, family = family)

m1 <- glmer(x ~ a + b, ...) # c dropped
anova(m_final, m1) # gives p-value, chisq, df for c

m2 <- glmer(x ~ a + c, ...) # b dropped
anova(m_final, m2) # gives p-value, chisq, df for b

m3 <- glmer(x ~ b + c, ...) # a dropped
anova(m_final, m3) # gives p-value, chisq, df for a

My first question is: Is this approach correct?
I talked to a colleague this morning, and he told me that for LRTs, you cannot only drop 1 variable at the time and compare these reduced models to the original model, but rather you have to drop them all after another, until you have no variables left. This logic confused me completely. So his approach would be:
m_final <- glmer(x ~ a + b + c + (1|d), data = dat, family = family)

m1 <- glmer(x ~ a + b, ...) # c dropped
anova(m_final, m1) # gives p-value, chisq, df for c

m2: <- glmer(x ~ a, ...) # c and b dropped, comparison now to model in which c was already dropped
anova(m2,m1) # gives p-value, chisq, df for b

m3 <- glmer(x ~ 1, ...) # c, b and a dropped, comparison now to model in which c and b were already dropped
anova(m3,m2) # gives p-value, chisq, df for a

If the second approach is correct, the sequence of dropping the terms would make a difference, right?
This whole discussion originally started because I'm working with some more complex models, where multiple interactions are present and for which I have to do LRTs. I'd also very much appreciate your advice in how to correctly compute likelihood ratio tests in a case like the following model:
m_final2 <- glmer(x ~ a*b + a*c + b*c + d + (1|r), data = dat, family = family)

The difficuly here are the interactions. d can easily be dropped, but not a,b and c. Nevertheless, I need to report p-values for all the variables alone and for the interactions.
Thank you very much for your help in advance. This might be a bit a silly question, but it really confuses me and all I can find is how to use LRTs for model selection purposes or for very simple examples.
All the best,
Asu
 A: The first approach would be very similar to a stepwise regression (see wikipedia) for variable selection, which is something that is widely done but generally not "good scientific practice". This is because the selection procedure does not lead to valid estimates and p-values, which you are actually interested in. See also this question.
The second approach proposed by your colleague is more sensible. You should decide (before seeing the data) which models you want to try beforehand. These should be "nested" within each other, that is, each subsequent model includes some additional variables compared to the previous one. Afterwards, you can run an ANOVA test in R that gives you the "best" model: anova(m3, m2, m1), which gives you two p-values for the model complexity steps m3-->m2 and m2-->m1 (where m3 is the baseline model with the intercept only). If both tests are significant, you take m1, if only m3-->m2 is significant, you take m2. Be careful, if m3-->m2 is not significant, but m2-->m1 is, you still have to go with m3.
I think you can report the p-values for the coefficents for the final chosen model.
