我正在尝试为配备 geeglm 函数（geepack 包）的 GEE 模型中的分类预测器获取模型预测的均值和 CI。该模型没有问题，但是当我试图估计模型预测的组均值时，我遇到了困难。这在其他软件包（即 SAS 和 SPSS）中相当容易做到，但我的意思是尝试在 R 中做到这一点（我也有点失望地发现没有直接的方法来获得分类预测器的整体测试除了分别拟合简化模型和完整模型然后比较它们；但另一方面，我能够找到一种使用 MESS 包计算 QIC 的简单方法）。无论如何，我四处搜索看看是否可能是用另一个包完成的（有一个名为 lsmeans 的包，但它似乎与 geepack 或 gee 包不兼容）。我很惊讶我无法找到像估计 R 中 GEE 模型的模型预测均值那样常见的东西，所以我想知道是否有人知道这个解决方案。我尝试联系 geepack 维护者，但没有得到答复。
为什么使用 GEE 而不是混合模型（可以使用 lmerTest 包计算混合模型的边际均值）？好吧，它的假设更少，并且对小样本更稳健。
doBy包中的 LSmeans 函数可能会有所帮助。
这是对小插图中示例的简单修改。
library(doBy)
library(geepack)
warp.gee <- geeglm(breaks ~ tension, id=wool, family=gaussian, data=warpbreaks)
LSmeans(warp.gee,effect="tension")
该emmeans软件包现在为 GEE 模型提供了估计的边际均值：
library(geepack)
library(emmeans)
warp.gee <- geeglm(breaks ~ tension, id=wool, family=gaussian, data=warpbreaks)
emmeans(warp.gee, ~tension)
tension lsmean SE df asymp.LCL asymp.UCL
L 36.38889 5.774705 Inf 25.07067 47.70710
M 26.38889 1.689200 Inf 23.07812 29.69966
H 21.66667 2.042753 Inf 17.66294 25.67039
plot(emmeans(warp.gee, ~tension), horizontal=FALSE, ylab="Estimated mean")
spind 包提供了一个类似于其他用途的预测功能predict：
https://www.rdocumentation.org/packages/spind/versions/2.1.3/topics/predict.GEE
# load packages
library(geepack)
library(spind)
n <- nrow(warpbreaks) # number of cases
trainIndex <- sample(1:n, n*.60) # random subset of 60% of cases (training)
testIndex <- setdiff(1:n, trainIndex) # rest of cases not in training (testing)
# model fit
warp.gee.fit <- geeglm(breaks ~ tension, id=wool, family=gaussian, data=warpbreaks[trainIndex,])
# predict test cases
predict(warp.gee.fit, warpbreaks[testIndex,])
```