我使用以下包:
library(aod)
library(MASS)
library(ggplot2)
我正在关注以下链接中的示例R代码:http://stats.idre.ucla.edu/r/dae/logit-regression/
这是我的GLMM的代码
str(data1)
flocation <- factor(data1$location)
fID <- factor(data1$ID)
GLMM1 <- glmmPQL(presence ~ water + location + turbidity + temperature +
sp.cond, random = ~ 1|fID, family = binomial, data = data1)
summary(GLMM1)
我根据不同的位置和水位做出预测,同时保持温度和涡轮恒定
newdata1 <- with(data1,
data.frame(water = water,
temperature = mean(temperature),
turbidity = mean(turbidity),
sp.cond = mean(sp.cond),
flocation = flocation))
newdata1$water.levelPred <- predict(GLMM1, type = "response")
newdata1
为了获得置信区间,我使用了下面的代码
newdata2 <- cbind(newdata1, predict(GLMM1, newdata = newdata1, type = "link",
se = TRUE))
newdata2 <- within(newdata2, {
PredictedProb <- plogis(fit)
LL <- plogis(fit - (1.96 * se.fit))
UL <- plogis(fit + (1.96 * se.fit))
})
运行置信区间代码后,我收到以下错误:
predict.lme中的错误(object,newdata,level = level,na.action = na.action):无法在'newdata'上评估所需级别的组plogis(fit)中的错误:找不到对象'fit'
为什么会这样?
因为我无法通过这一步,所以我在使用下面的代码编制具有预测概率的CI时遇到问题:
在ggplot2中绘图
ggplot(newdata2, aes(x = water, y = water.levelProb)) + geom_ribbon(aes(ymin = LL, ymax = UL, fill = flocation), alpha = 0.2) + geom_line(aes(colour = flocation),size = 1)+facet_wrap(~flocation)+xlab("Water Depth (m)")+ylab("Predicted Probability")+theme_bw()