

  
Dear Geoff,



   Thank you so much for your kind reply. Actually, the problem I'm
finding is more like this:

      fm.1 <- lmer(Y ~ A*B + C)

      fm.2 <- lmer(Y ~ A*B)

      anova(fm.1, fm.2)

   

   But when your full model includes C*B. So, I have to remove C*B
from the model (as C is involved in it) to test the effect of C, or that's what
I understood.



   Anyway, let's see if someone else replies.



   Thank you very much for your help again.



--- El jue, 20/10/11, Geoff Brookshire <broog731@newschool.edu> escribió:

De: Geoff Brookshire <broog731@newschool.edu>
Asunto: Re: Calculating main effects
Para: "Iker Vaquero Alba" <karraspito@yahoo.es>
Fecha: jueves, 20 de octubre, 2011 14:45

Hi Iker,

Actually, that post of mine was meant to be more a request for
information than giving it. No one responded, so I'm not sure what
people's opinions on that are. I've done some reading in the meantime
that's moved me towards one of those options, though. I recommend
reading Baayen, Davidson, & Bates (2008) in the Journal of Memory and
Language to get a better picture of what these models are doing and
some recommendations for how to test significance values. That issue
of JML has a couple other articles about mixed-effects models, too.

My impression, to get to your question, is that the best way of
looking at significance for specific terms is to leave everything you
can in the model except the term of interest. So if you have the
factors A, B, and C, and want to know the influence of C, i think I'd
do:

fm.1 <- lmer(Y ~ A*B + C)
fm.2 <- lmer(Y ~ A*B)
anova(fm.1, fm.2)

That way you leave in everything to do with A and B. I choose this
option not because i know whether or not it gives the most accurate
p-values (i have no idea if this is the case, or if p-values can be
described as 'accurate'), but because it seems to be the closest
conceptually to what's going on in more basic tests. The pvals.fnc
option is a little different. It would just be to run something like
this:

fm.full <- lmer(Y ~ A*B*C)  # This is your full model, with every term
of potential interest
library(languageR)
pvals.fnc(fm.full)

This function doesn't give you wald chi-sq stats, but bayesian
confidence intervals (and their associated p-values) using MCMC
sampling and a non-informative prior. I'm avoiding that option because
I don't understand everything involved in Bayesian analysis, so I'm
afraid I'll misuse it. I know this function wasn't working for a
while, but it works for me now. If you want to take this route, update
your packages and it should be fine, at least for normally distributed
data (I think).

It's possible that some people would recommend something else.

Good luck!

geoff

On Mon, Oct 17, 2011 at 1:26 PM, Iker Vaquero Alba <karraspito@yahoo.es> wrote:
>
>    Dear Geoff:
>    I´m so sorry to bother you again with the same issue, but I have posted a
> new mail in the list, asking for help about calculating the main effects
> involved in higher order interactions of a model. I just wanted to know if I
> have understood the methods you suggest, adapted for my own data. I would
> really appreciate if you could take a pair of minutes to reply, just to make
> sure that I'm doing it right, as I am getting lower p-values compared to
> making the simplification straight away, and I don't know wether that is an
> acceptable result. Obviously, I appreciate very much the single fact that
> you are reading this, regardless of wether you reply or not.
>
>    I copy and paste the last post, if that's OK for you:
>
>    I have my minimum adequate model, which is:
>
>
>>s.hatchmodel1.103<-lmer(prophatch1~sex+inslarge+tlength+briventral+brithr+sex:brithr+tlength:brithr+briventral:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>    and I want to test the effect of sex, tlength, briventral and brithr when
> they're not taking part in any interactions. The first possible method you
> suggest is using the pvals.fnc from languageR on the full model. I don't
> know wether you mean the results of the Wald's test of the summary table,
> but if so, that's not a possibility since lmer does not give p-values for
> models with a Gaussian distribution, but only for GLMMs.
>
>    According to the second possibility, it would be simply like this (for
> example, to test brithr):
>
>    >brithrmodel1<-lmer(prophatch1~brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>    >brithrmodel2<-lmer(prophatch1~1+(1|site/pair)+(1|year),REML=FALSE)
>
>    >anova(brithrmodel1,brithrmodel2)
>
>
>    And according to the third suggestion, the one I'm really interested in
> understanding (as you say it gives the most accurate pvals): basically, for
> calculating the significance of brithr, I have to remove all the
> interactions in which brithr is taking part (in this case, all of them), and
> test the model with the ramaining terms PLUS brithr against the same model
> without brithr, like this:
>
>>brithrmodel1<-lmer(prophatch1~sex+inslarge+tlength+briventral+brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>>brithrmodel2<-lmer(prophatch1~sex+inslarge+tlength+briventral+(1|site/pair)+(1|year),REML=FALSE)
>
>>anova(brithrmodel1,brithrmodel2)
>
>
> For briventral, I only have to remove one of the 2-way interactions:
>
>
>>briventralmodel1<-lmer(prophatch1~sex+inslarge+tlength+briventral+brithr+sex:brithr+tlength:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>
>>briventralmodel2<-lmer(prophatch1~sex+inslarge+tlength+brithr+sex:brithr+tlength:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>>anova(briventralmodel1,briventralmodel2)
>
> Something similar for tlength:
>
>
>>tlengthmodel1<-lmer(prophatch1~sex+inslarge+tlength+briventral+brithr+sex:brithr+briventral:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>>tlengthmodel2<-lmer(prophatch1~sex+inslarge+briventral+brithr+sex:brithr+briventral:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>
>>anova(tlengthmodel1,tlengthmodel2)
>
> And the same for sex:
>
>
>>sexmodel1<-lmer(prophatch1~sex+inslarge+tlength+briventral+brithr+tlength:brithr+briventral:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>>sexmodel2<-lmer(prophatch1~inslarge+tlength+briventral+brithr+tlength:brithr+briventral:brithr+(1|site/pair)+(1|year),REML=FALSE)
>
>
>>anova(briventralmodel1,briventralmodel2)
>
>
> Is that right? I hope I have understood it properly.
>
> Thank you so much!!
>

	[[alternative HTML version deleted]]


