我正在尝试在 R 中模拟 GLM 函数，包括随机不确定性。我将基于公式的方法与基于 R 的simulate()函数进行了比较，得到了不同的结果。不确定我（可能是我而不是 R）做错了什么。
我首先创建一个模拟队列：
set.seed(1)
library(MASS)
d <- mvrnorm(n=3000, mu=c(30,12,60), Sigma=matrix(data=c(45, 5, 40, 5, 15, 13, 40, 13, 300), nrow=3))
d[,1] <- d[,1]^2
拟合模型：
m <- glm(formula=d[,1]~d[,2] + d[,3], family=gaussian(link="sqrt"))
在线性范围内预测（使用公式）
p_linear <- m$coefficients[1] + m$coefficients[2]*d[,2] + m$coefficients[3]*d[,3]
比较两种方法的预测，包括反向链接函数。它似乎是相似的：
sum((predict(object=m, type="response")-p_linear^2)^2)==0
将随机部分添加到两个预测中：
sd_residuals <- sd(p_linear^2 - d[,1])
p <- p_linear^2 + rnorm(n=1000, mean=0, sd=sd_residuals)
p_simulate <- simulate(m)$sim_1
将两个预测（基于公式和simulate()基于公式）与图中的原始数据进行比较：
par(mfrow=c(1,2), mar=c(4,2,2,1))
plot(d[,1], p, col="blue", pch=20, cex=0.5, xlim=c(-500, 3000), ylim=c(-500, 3000))
abline(lm(p~d[,1]), col="blue")
points(d[,1], p_simulate, col="red", cex=0.5)
abline(lm(p_simulate~d[,1]), col="red")
abline(a=0, b=1)
...并比较两个预测：
plot(p, p_simulate, col="green", cex=0.5, xlim=c(-500, 3000), ylim=c(-500, 3000))
abline(lm(p_simulate~p), col="green")
abline(a=0, b=1)
似乎基于公式的预测反映了与原始数据相似的随机不确定性。基于 r 的simulate()函数方法似乎比原始数据表示的随机不确定性更小。基于公式的方法是否正确？如果不是，我应该如何适应它？（出于兴趣：是什么导致两种方法之间的差异？）
最终，我想在数据之外使用 GLM 模型作为预测的输入值。我也想在Excel中应用它。因此，我正在寻找一种基于公式的方法，而不是使用像simulate().
欢迎阅读任何有关阅读预测材料的提示，包括 R 中 GLM 中的随机性。我用谷歌搜索了很多，但是很难在 R 中找到基于公式的示例，而不是家庭/链接是高斯/身份。