Ошибка в модели прогнозирования с LME

Question

Я пытаюсь делать прогнозы от моего LME-модели:

mod<-lme(height~direction+time+distnest+loop+twc+twc:direction,random=~1|bird_id/FT_no,data=dat,correlation=corAR1(0.5))  

К сожалению, у меня есть вложенные случайные и 5 фиксированных эффектов и на самом деле не знают, как обращаться с ними все.

Все фиксированные эффекты являются числовыми, за исключением «направления» и «времени», которые являются двухуровневыми факторами («земля» или «море», «день» или «ночь»).

Для предсказания модели я попытался следующие:

newdat<-data.frame(loop=seq(min(dat$loop),max(dat$loop),length=100), 
       direction=factor("land",levels(dat$direction)), 
       time_code=factor("1",levels(dat$time)), 
       distnest=mean(dat$distnest), 
       twc=mean(dat$twc)) 

newdat$pred<-predict(mod,newdata=newdat,level=0) 
plot(dat$loop,dat$height,pch=16,las=1,cex.lab=1.2) 
lines(newdat$loop,yhat,lwd=2) ### for plotting one of the fixed effects 
newdat$predse<-predict(mod,newdat,se.fit=TRUE)$se.fit 

Но при использовании se.fit = TRUE есть ошибка: «Ошибка в predict.lme (mod3, newdat, se.fit = TRUE): невозможно определить группы для желаемых уровней на 'newdata'

Не работает ли se.fit для lmes? Опущение уровня = 0 не сработало. Код неправильный?

Я также попытался код с glmm.wikidot:

newdat<-expand.grid(direction=c("land","sea"),time_code=c("1","2"),loop30=c(0.08,1),distne st=c(0.01,99.43),twc=c(-56.88744,57.93735)) 
newdat$pred<-predict(mod3,newdat,level=0) 
Designmat <- model.matrix(eval(eval(mod3$call$fixed)[-2]), newdat[-ncol(newdat)]) 
predvar <- diag(Designmat %*% mod3$varFix %*% t(Designmat)) 
newdat$SE <- sqrt(predvar) 
newdat$SE2 <- sqrt(predvar+mod3$sigma^2) 
library(ggplot2) 
pd <- position_dodge(width=0.4) 
g0 <- ggplot(newdat,aes(x=loop30,y=pred,colour=direction))+ 
geom_point(position=pd) 
g0 + geom_linerange(aes(ymin=pred-2*SE,ymax=pred+2*SE), position=pd) 
## prediction intervals 
g0 + geom_linerange(aes(ymin=pred-2*SE2,ymax=pred+2*SE2), position=pd) 

Но сюжет совершенно неправильно. Может ли кто-нибудь помочь мне с правильным кодом? Заранее большое спасибо.

С наилучшими пожеланиями, Анна

Answer 1

Для предсказания вам нужно все элементы в RHS вашей модели ... _as_named_:

direction + time + distnest + loop + twc 
...and... bird_id ...and... FT_no 

На данный момент вы неверно названы переменную time в time_code и не включили ни одну из ваших переменных «случайных эффектов». Вот соответствующий бит из справочной страницы nlme :: lme: «Все переменные, используемые в моделях с фиксированными и случайными эффектами, а также факторы группировки, должны присутствовать в кадре данных. Если они отсутствуют, возвращаемые значения возвращаются».