REML vs ML stepAIC

I feel overwhelmed after attempting to dig into the literature on how to run my mixed model analysis following it up with using AIC to select the best model or models. I do not think my data is that complicated, but I am looking for confirmation that what I have done is correct, and then advise on how to proceed. I am unsure if I should be using lme or lmer and then with either of those, if I should be using REML or ML.

I have a value of selection and I want to know which covariates best influence that value and allow for predictions. Here's some made up example data and my code for my test that I am working with:

ID=as.character(rep(1:5,3))
season=c("s","w","w","s","s","s","s","w","w","w","s","w","s","w","w")
time=c("n","d","d","n","d","d","n","n","n","n","n","n","d","d","d")
repro=as.character(rep(1:3,5))
risk=runif(15, min=0, max=1.1)
comp1=rnorm(15, mean = 0, sd = 1)
mydata=data.frame(ID, season, time, repro, risk, comp1)
c1.mod1<-lmer(comp1~1+(1|ID),REML=T,data=mydata)
c1.mod2<-lmer(comp1~risk+(1|ID),REML=T,data=mydata)
c1.mod3<-lmer(comp1~season+(1|ID),REML=T,data=mydata)
c1.mod4<-lmer(comp1~repro+(1|ID),REML=T,data=mydata)
c1.mod5<-lmer(comp1~time+(1|ID),REML=T,data=mydata)
c1.mod6<-lmer(comp1~season+repro+time+(1|ID),REML=T,data=mydata)
c1.mod7<-lmer(comp1~risk+season+season*time+(1|ID),REML=T,data=mydata)


I have ~19 models that explore this data with various combinations and up to a 2 way interaction terms, but always with ID as a random effect and comp1 as my dependent variable.

• Q1. Which to use? lme or lmer? does it matter?

In both of these, I have the option to use ML or REML - and I get drastically different answers - using ML followed by AIC I end up with 6 models all with similar AIC values and the model combinations simply do not make sense, whereas REML results in 2 of the most likely models being the best. However, when running REML I cannot use anova any longer.

• Q2. is the main reason to use ML over REML because of use with ANOVA? This is not clear to me.

I am still not able to run stepAIC or I do not know of another way to narrow down those 19 models.

• Q3. is there a way to use stepAIC at this point?
• For Q2, ML is necessary because comparisons using REML are not valid when the fixed effects change. A possible useful related question is here: stats.stackexchange.com/a/16015/3601 – Aaron left Stack Overflow Oct 24 '12 at 15:57
• @Aaron I had looked at that question before, but was still confused. Use REML only "works" when the random effect changes? I obviously do not understand enough of ML vs REML. Thanks though, that helps with one of my questions! – Kerry Oct 25 '12 at 10:35
• Yes, that is correct. When comparing models, REML should only be used if the models have the same fixed effects. Answer expanded below. – Aaron left Stack Overflow Oct 25 '12 at 14:27

Q1. Which to use? lme or lmer? does it matter? Either is fine. They will give you same fits. lme will give you p-values, and lmer won't, but that's more than I want to get into here. The most famous reference is one of Doug Bates's posts to the R-help mailing list here.