Cholesky 变换 lmer 中的残差
Cholesky transformed residuals in lmer
我目前正在使用线性混合模型 (LMM) 并希望执行一些残差分析。由于线性混合模型的残差 "are correlated and do not necessarily have constant variance",它们不足以检查模型的假设(Fitzmaurice 等人,2011 年,第 266 页)。因此有人建议(Waternaux 等人,1989 年)使用 Cholesky 变换残差来研究模型假设,另请参见 Houseman 等人。 (2004) 和 Santos Nobre & da Motta Singer (2007)。
也就是说,假设我们有一个形式为 Yi = Xiβ + Zi 的 LMM bi + Ɛi,其中β是固定效应参数的向量,bi 是正态分布随机效应的向量,其中 Xi 和 Zi 是设计矩阵。
那么,残差ei = Yi - Xiβ是相关的,他们观察到的对应物使用估计的 β 系数而不是未知的真实值。我们将这些残差表示为 ri。通过找到 ri 的协方差矩阵的 Cholesky 分解来完成对残差的去相关。即令Cov(ei) = Σi,令Lii' 是基于 ri 的 Σi 估计的 Cholesky 分解。然后将 Cholesky 变换残差定义为 r*i = Li-1r我。这些是我试图从 lme4 (Bates et al. 2015) 包中获取的残差,但到目前为止没有成功。
但是,使用 Fitzmaurice 等人主页上的数据。 (2011),我已经设法使用 nlme 包(而不是 lme4 包)从书中的案例研究中重新创建残差。代码如下。
require(foreign) # To read dta-files
require(nlme) # The nlme package used to estimate the LMM.
require(mgcv) # Needed for the extract.lme.cov function
fat <- read.dta("https://content.sph.harvard.edu/fitzmaur/ala2e/fat.dta") # We read the data
# Add spline function for time, with a knot at age of menarch, denoted time = 0 in the data.
fat$timeknot <- fat$time * (fat$time > 0)
# Estimate the LMM with nlme
nlme.model <- lme(pbf ~ time + timeknot, data=fat, random = ~ 1 + time + timeknot | id)
# Calculate the residuals
raw.residuals <- residuals(nlme.model, level=0)
# Cholesky residuals
est.cov <- extract.lme.cov(nlme.model, fat) # We extract the blocked covariance function of the residuals
Li <- t(chol(est.cov)) # We find the Cholesky transformation of the residuals. (The transform is to get the lower trangular matrix.)
cholesky.residuals <- solve(Li) %*% raw.residuals # We then calculate the transformed residuals.
hist(cholesky.residuals) # We plot the residuals
因此,我的问题是,是否有使用 lme4 包提取残差的估计协方差矩阵的简单方法?或者直接得到 Cholesky 变换后的残差?
感谢任何帮助或指点。关于,
/菲尔
参考文献:
Bates, D.、Mächler, M.、Bolker, B. 和 Walker, S.(2015 年)。使用 lme4 拟合线性混合效应模型。 统计软件杂志, 67(1), 1-48.
Fitzmaurice, G.M., Laird, N.M., & Ware, J.H. (2011)。 应用纵向分析,第二版。 John Wiley & Sons.
Houseman, E.A., Ryan, L.M., & Coull, B.A. (2004)。 Cholesky 残差用于评估具有相关结果的线性模型中的正态误差。 美国统计协会杂志n, 99(466), 383-394。
Santos Nobre, J. 和 da Motta Singer, J. (2007)。线性混合模型的残差分析。生物统计学杂志:生物科学数学方法杂志, 49(6), 863-875。
Waternaux, C., Laird, N. M., & Ware, J. H. (1989)。纵向数据分析方法:血铅浓度和认知发展。 美国统计协会杂志, 84(405), 33-41。
经过大量挖掘,我发现这个问题的答案已在本网站的其他地方 post 编辑,而我在 [=15= 之前找不到 post ]我自己的。 (虽然我搜索过!)
可以在这里找到解决方案:
我目前正在使用线性混合模型 (LMM) 并希望执行一些残差分析。由于线性混合模型的残差 "are correlated and do not necessarily have constant variance",它们不足以检查模型的假设(Fitzmaurice 等人,2011 年,第 266 页)。因此有人建议(Waternaux 等人,1989 年)使用 Cholesky 变换残差来研究模型假设,另请参见 Houseman 等人。 (2004) 和 Santos Nobre & da Motta Singer (2007)。
也就是说,假设我们有一个形式为 Yi = Xiβ + Zi 的 LMM bi + Ɛi,其中β是固定效应参数的向量,bi 是正态分布随机效应的向量,其中 Xi 和 Zi 是设计矩阵。
那么,残差ei = Yi - Xiβ是相关的,他们观察到的对应物使用估计的 β 系数而不是未知的真实值。我们将这些残差表示为 ri。通过找到 ri 的协方差矩阵的 Cholesky 分解来完成对残差的去相关。即令Cov(ei) = Σi,令Lii' 是基于 ri 的 Σi 估计的 Cholesky 分解。然后将 Cholesky 变换残差定义为 r*i = Li-1r我。这些是我试图从 lme4 (Bates et al. 2015) 包中获取的残差,但到目前为止没有成功。
但是,使用 Fitzmaurice 等人主页上的数据。 (2011),我已经设法使用 nlme 包(而不是 lme4 包)从书中的案例研究中重新创建残差。代码如下。
require(foreign) # To read dta-files
require(nlme) # The nlme package used to estimate the LMM.
require(mgcv) # Needed for the extract.lme.cov function
fat <- read.dta("https://content.sph.harvard.edu/fitzmaur/ala2e/fat.dta") # We read the data
# Add spline function for time, with a knot at age of menarch, denoted time = 0 in the data.
fat$timeknot <- fat$time * (fat$time > 0)
# Estimate the LMM with nlme
nlme.model <- lme(pbf ~ time + timeknot, data=fat, random = ~ 1 + time + timeknot | id)
# Calculate the residuals
raw.residuals <- residuals(nlme.model, level=0)
# Cholesky residuals
est.cov <- extract.lme.cov(nlme.model, fat) # We extract the blocked covariance function of the residuals
Li <- t(chol(est.cov)) # We find the Cholesky transformation of the residuals. (The transform is to get the lower trangular matrix.)
cholesky.residuals <- solve(Li) %*% raw.residuals # We then calculate the transformed residuals.
hist(cholesky.residuals) # We plot the residuals
因此,我的问题是,是否有使用 lme4 包提取残差的估计协方差矩阵的简单方法?或者直接得到 Cholesky 变换后的残差?
感谢任何帮助或指点。关于,
/菲尔
参考文献:
Bates, D.、Mächler, M.、Bolker, B. 和 Walker, S.(2015 年)。使用 lme4 拟合线性混合效应模型。 统计软件杂志, 67(1), 1-48.
Fitzmaurice, G.M., Laird, N.M., & Ware, J.H. (2011)。 应用纵向分析,第二版。 John Wiley & Sons.
Houseman, E.A., Ryan, L.M., & Coull, B.A. (2004)。 Cholesky 残差用于评估具有相关结果的线性模型中的正态误差。 美国统计协会杂志n, 99(466), 383-394。
Santos Nobre, J. 和 da Motta Singer, J. (2007)。线性混合模型的残差分析。生物统计学杂志:生物科学数学方法杂志, 49(6), 863-875。
Waternaux, C., Laird, N. M., & Ware, J. H. (1989)。纵向数据分析方法:血铅浓度和认知发展。 美国统计协会杂志, 84(405), 33-41。
经过大量挖掘,我发现这个问题的答案已在本网站的其他地方 post 编辑,而我在 [=15= 之前找不到 post ]我自己的。 (虽然我搜索过!)
可以在这里找到解决方案: