我正在尝试从 lme 对象中获得对观察结果的预测。这应该很简单。然而,由于我在不同的试验中得到不同类型的错误,在我看来我错过了一些东西。我的模型如下:
model <- lme(log(child_mortality) ~ as.factor(cluster)*time +
my.new.time.one.transition.low.and.middle + ttd +
maternal_educ+ log(IHME_id_gdppc) + hiv_prev-1,
merged0,na.action=na.omit,method="ML",weights=varPower(form=~time),
random= ~ time| country.x,
correlation=corAR1(form = ~ time),
control=lmeControl(msMaxIter = 200, msVerbose = TRUE))
它运行良好,非常适合数据并且结果很有意义。现在为了得到预测,我尝试了以下方法:
test.pred <- data.frame(time=c(10,10,10,10),country.x=c("Poland","Brazil",
"Argentina","France"),
my.new.time.one.transition.low.and.middle=c(1,1,1,0),
ttd=c(0,0,0,0),maternal_educ=c(10,10,10,10),
IHME_id_gdppc=c(log(5000),log(8000),log(8000),log(15000)),
hiv_prev=c(.005,.005,.005,.005),
cluster=c("One Transition, Middle Income","One Transition,
Middle Income","One Transition, Middle Income","Democracy,
High Income"))
>
> predict(model,test.pred,level=0)
Error in X %*% fixef(object) : non-conformable arguments
如果我排除法国,并且只包括 cluster="OneTransition, Middle Income" 的国家,那么我会得到一个不同的错误
# create a toy data set
test.pred0 <-
expand.grid(time=20:29,country.x=c("Poland","Brazil","Argentina"))
z0 <-as.data.frame(cbind(my.new.time.one.transition.low.and.middle =
c(0,0,0,0,0,0,1,2,3,4), ttd=c(0,0,0,0,0,0,1,0,0,0),
maternal_educ=seq(from=10.0, to=12.0, length.out=10),
IHME_id_gdppc=log(seq(from=5000, to=8000, length.out=10)),
hiv_prev=rep(.005,10),
cluster=rep("One Transition, Middle Income",10)))
z <- rbind(z0,z0,z0)
test.pred <- cbind(test.pred0,z)
# check
head(test.pred)
> time country.x my.new.time.one.transition.low.and.middle ttd
> maternal_educ IHME_id_gdppc hiv_prev
> 1 20 Poland 0 0
> 10 8.51719319141624 0.005
> 2 21 Poland 0 0
> 10.2222222222222 8.58173171255381 0.005
> 3 22 Poland 0 0
> 10.4444444444444 8.64235633437024 0.005
> 4 23 Poland 0 0
> 10.6666666666667 8.69951474821019 0.005
> 5 24 Poland 0 0
> 10.8888888888889 8.75358196948047 0.005
> 6 25 Poland 0 0
> 11.1111111111111 8.80487526386802 0.005
> cluster
> 1 One Transition, Middle Income
> 2 One Transition, Middle Income
> 3 One Transition, Middle Income
> 4 One Transition, Middle Income
> 5 One Transition, Middle Income
> 6 One Transition, Middle Income
# run the predictions
predict(model,test.pred,level=0)
> Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
> contrasts can be applied only to factors with 2 or more levels
在此示例中,问题始终是由于 cluster="One Transition, Middle Income" 造成的。
我不明白为什么这是一个问题。如果我想让 predict() 工作,我必须包含模型中的所有变量,对吗?显然,模型调用中的输入数据不会包括所有情况下设置为相同值的因子。然而,如果我只想获得数据子集或新观察的预测,我可能只对某些因素始终设置为相同的情况感兴趣。是否有意义?在这种情况下如何获得预测?