首页 > 解决方案 > 使用 `ordinal::clmm` 模型对新数据进行预测

问题描述

我有一些重复测量,有序响应数据:

dat <- data.frame(
  id = factor(sample(letters[1:5], 50, replace = T)),
  response = factor(sample(1:7, 50, replace = T), ordered = T),
  x1 = runif(n = 50, min = 1, max = 10),
  x2 = runif(n = 50, min = 100, max = 1000)
)

我已经建立了以下模型:

library(ordinal)

model <- clmm(response ~ x1 + x2 + (1|id), data = dat)

我有一些新数据:

new_dat <- data.frame(
  id = factor(sample(letters[1:5], 5, replace = T)),
  x1 = runif(n = 5, min = 1, max = 10),
  x2 = runif(n = 5, min = 100, max = 1000)
)

我希望能够使用该模型来预测每个级别dat$response发生的概率new_dat,同时还要考虑id

不幸的是predict()不适用于clmm对象。predict()确实适用于clmm2对象,但它忽略了包含的任何随机效果。

我想要实现的是类似于使用此代码在下面的图 3中所做的事情:

library(ordinal)

fm2 <- clmm2(rating ~ temp + contact, random=judge, data=wine, Hess=TRUE, nAGQ=10)

pred <- function(eta, theta, cat = 1:(length(theta)+1), inv.link = plogis){
Theta <- c(-1e3, theta, 1e3)
sapply(cat, function(j)
inv.link(Theta[j+1] - eta) - inv.link(Theta[j] - eta))
}

mat <- expand.grid(judge = qnorm(0.95) * c(-1, 0, 1) * fm2$stDev,
contact = c(0, fm2$beta[2]),
temp = c(0, fm2$beta[1]))

pred.mat <- pred(eta=rowSums(mat), theta=fm2$Theta)

lab <- paste("contact=", rep(levels(wine$contact), 2), ", ", "temp=", rep(levels(wine$temp), each=2), sep="")

par(mfrow=c(2, 2))

for(k in c(1, 4, 7, 10)) {
plot(1:5, pred.mat[k,], lty=2, type = "l", ylim=c(0,1),
xlab="Bitterness rating scale", axes=FALSE,
ylab="Probability", main=lab[ceiling(k/3)], las=1)
axis(1); axis(2)
lines(1:5, pred.mat[k+1, ], lty=1)
lines(1:5, pred.mat[k+2, ], lty=3)
legend("topright",
c("avg. judge", "5th %-tile judge", "95th %-tile judge"),
lty=1:3, bty="n")
}

除了,我的模型包含多个连续协变量(与二元因子相反)。

如何使用模型数据来预测每个级别dat$response发生的概率new_dat,同时还要考虑id

非常感谢。

标签: rpredictmixed-modelsordinal

解决方案


推荐阅读