[R-sig-ME] How to include a quadratic term in an lmer model

[Emmanuel Charpentier]

Le jeudi 04 juin 2009 à 17:28 +0200, Linda Mortensen a écrit :
> Hi list members,
>  
> I have recently e-mailed this list asking for some advice on how to
> use mixed-effects models on ordinal responses (see posts entitled "How
> to use mixed-effects models on multinomial data"). This query concerns
> the same data set, but since the topic is a different one, I post the
> query under a different header. 
> Following several list members advice, I'm using a linear
> mixed-effects model to analyse the data described in the earlier
> posts, so I'm still working within this model framework. But when
> trying to decide on the best-fitting model, I have run into a problem:
> In the data set, there is a curved relationship between one of the
> explanatory variables (i.e., the serial position of items in a list)
> and the response variable. A model that includes both a linear and a
> quadratic term for this variable would most likely describe this
> relationship better than a model that includes only a linear variable.
> But when I try to include the quadratic term in the model, using the
> formula "lmer(y ~ x + I(x^2)", I get the following error message:    
> "Error in `contrasts<-`(`*tmp*`, value = "contr.treatment") : 
>   contrasts can be applied only to factors with 2 or more levels"
> And the following warnings:
> "In addition: Warning messages:
> 1: In Ops.factor(cposition, 2) : ^ not meaningful for factors
> 2: In Ops.factor(cposition, 2) : ^ not meaningful for factors"
> 
> Judging from the error message, the problem seems to be that I have
> coded the serial position variable using the (default) contrast coding
> system.

So, your "x" is coded as a factor ? If so, x^2 isn't even defined.

To speak of a "curvature" *supposes* an a priori ordering of your "x"
variable. Two possibilities :

- code it as an *ordered* factor,
- code it as a continuous variable.

The second possibility entails more assumptions about x than the first
(additivity, fixed distances between levels, etc ...).

You might want to use poly(x,2) rather than x+I(x^2) : poly(x,order)
will use orthogonal polynomials, ensuring independence of the
coefficients of the linear and quadratic term estimators. However, the
coefficients have no easy interpretation, and getting back to the
original value might no be *that* easy. Therefore, you may use poly(x,2)
for testing purposes, and x+I(x^2) to get the numerical values :-). I'm
sure a better solution can be found by looking examples and R-help for
poly().

>          The serial position variable has five levels (positions 1
> through 5), so I don't understand why R is complaining about the
> number of levels. - There are 5 levels, so the serial position factor
> has "2 or more levels". Is it the coding system that is the problem?
> Can someone explain to me what I'm doing wrong here? Any suggestions
> are greatly appreciated.
>  
> Linda