首页 > 解决方案 > 用于惩罚样条的 R 中 lmer 中的指定协方差矩阵

问题描述

我一直在尝试使用惩罚样条和线性混合模型之间的连接来拟合 R 中的惩罚样条回归模型。虽然我对 R 函数非常熟悉,但lme使用它的竞争对手会lmer遇到一些困难。这是一个我想概括为的玩具示例lmer

require(nlme)

grid <- seq(0, 1, len = 100)
y <- rep(0, length(grid))
for(i in 1:length(grid)){
  y[i] <- sin(3*pi*grid[i])  + rnorm(1, 0, 1)
}
X <- cbind(rep(1, length(grid)), grid, grid^2, grid^3)
knots <- seq(0.01, 0.99, len = 100)
Z <- outer(grid, knots, "-")
Z <- Z*(Z>0)
Z <- Z^3
data.mixed <- data.frame(X, Z)
data.mixed$all <- rep(1 ,nrow(data.mixed))
fit.mixed <- lme(y~ X-1, random = list( all  = pdIdent(~Z-1) ), data = data.mixed )
curve(sin(3*pi*x), 0, 1, lwd = 2, ylim = c(-1.5, 1.5))
beta.hat <- fit.mixed$coefficients$fixed
u.hat <- unlist(fit.mixed$coefficients$random)
f.hat <- X%*%beta.hat + Z%*%u.hat
lines(grid, f.hat, lwd = 2, col = "red")

棘手的部分是 pdIdent 函数,它指定虽然我们有几个随机效应,但它们的协方差结构是单位矩阵的倍数;lme然后估计这个标量倍数。

如果有人熟悉lmer,我将不胜感激有关如何使用该功能完成此操作的一些提示。

先感谢您

标签: rstatisticslme4nlme

解决方案


推荐阅读