[R] mixed model random interaction term log likelihood ratio test

[Ben Bolker]

seatales <ssphadke <at> uh.edu> writes:



> 
> Hello,
> I am using the following model
> 
> model1=lmer(PairFrequency~MatingPair+(1|DrugPair)+(1|DrugPair:MatingPair),
> data=MateChoice, REML=F)
> 
> 1. After reading around through the R help, I have learned that the above
> code is the right way to analyze a mixed model with the MatingPair as the
> fixed effect, DrugPair as the random effect and the interaction between
> these two as the random effect as well. Please confirm if that seems
> correct.


  You should probably send this sort of question to the 
r-sig-mixed-models mailing list ...

  You probably want (MatingPair|DrugPair) rather than
(1|DrugPair:MatingPair).   
Whether REML=FALSE or REML=TRUE depends what you want
to do next.

> 
> 2. Assuming the above code is correct, I have model2 in which I remove the
> interaction term, model3 in which I remove the DrugPair term and model4 in
> which I only keep the fixed effect of MatingPair.
> 
> 3. I want to perform the log likelihood ratio test to compare these models
> and that's why I have REML=F. However the code anova(model1, model2, model3,
> model4) gives me a chisq estimate and a p-value, not the LRT values. How do
> I get LRT (L.Ratio) while using lmer?

  The chi-squared values are the differences in deviance (-2 log likelihood)
between the respective papers of models, which under the null hypothesis
of the LRT will be chi-squared distributed.  In other words, these
*are* the LRT test statistics.
> 
> 4. I am under the impression after reading a few posts that LRT is not
> usually obtained with lmer but it is given if I use lme (the old model).

   I don't know what you mean by this.
   The main difference between lmer and lme in the testing/inference
context is that lme is willing to guess at "denominator degrees of freedom"
to perform conditional F-tests.

> 
> 5. I could not find how to input the random interaction term while using
> lme? Is it the following way? Would someone please guide me to some already
> existing posts or help here?

  = ran
> 
> lme(PairFrequency~MatingPair, random=~(1|DrugPair)+(1|DrugPair:MatingPair),
> data=MateChoice, method="ML")...is this the right way? would lme give me
> loglikelihood ratio test values (L.ratio)?
> 

  See above.

> Thanks a lot. I hope someone can help. Most posts I have found deal with
> nesting but there is absolutely no nesting in my data. 
> 
> Sujal P.
> p.s: If it matters how data is arranged, then I have one vector called
> MatingPair which has 3 levels and another vector DrugPair which also has 3
> levels. The PairFrequency data is a count data and is normally distributed.
> The data are huge, hence I am not able to post it here.

  It is probably unwise to estimate DrugPair as a random effect if
it only has three levels.


> 
> Also, here is what I mean by getting chisq value rather than L.Ratio:

  See above.

> Data: MateChoice
> Models:
> model2: PairFrequency ~ MatingPair + (1 | DrugPair)
> model3: PairFrequency ~ MatingPair + (1 | DrugPair:MatingPair)
> model1: PairFrequency ~ MatingPair + (1 | DrugPair) + (1 |
> DrugPair:MatingPair)
> 
>        Df    AIC    BIC  logLik  Chisq Chi Df Pr(>Chisq)   
> model2  5 274.90 282.82 -132.45                            
> model3  5 282.44 290.36 -136.22 0.0000      0    1.00000   
> model1  6 276.90 286.40 -132.45 7.5443      1    0.00602 **
> ---
> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1
> ‘ ’ 1 
> 
> --
> View this message in context:
http://r.789695.n4.nabble.com/mixed-model-random-interaction-term-log-likelihood-ratio-test-tp3448718p3448718.html
> Sent from the R help mailing list archive at Nabble.com.
> 	[[alternative HTML version deleted]]
> 
>