在 R 模拟中估计线性回归和线性混合模型中的偏差

Question

我想运行模拟来估计线性模型和线性混合模型中的偏差。偏差是 E(beta)-beta，其中 beta 是我的 X 和 Y 之间的关联。

我从正态分布生成 X 变量，从多元正态分布生成 Y。

我了解如何通过模拟计算 E(beta)，即所有模拟的 beta 估计总和除以模拟总数，但我不确定如何估计真实的 beta。

meanY <- meanY + X*betaV

这就是我生成 meanY（betaV 是效应大小）的方式，然后用于生成多变量 Y 结果，如下所示。

Y[jj,] <- rnorm(nRep, mean=meanY[jj], sd=sqrt(varY))

我知道如何根据模拟计算 E(beta)，即所有模拟的 beta 估计值之和除以模拟总数，但我不确定如何估计真实的 beta。

据我有限的了解，真正的 beta 不是从数据中获得的，而是从我设置固定 beta 值的设置中获得的。

根据我生成数据的方式，我如何估计真实的贝塔值？

Answer 1

有两种模拟偏差的方法。我将使用线性模型举一个简单的例子。线性混合模型可能会使用类似的方法，但我不确定它是否适用于广义线性混合模型（我只是不确定）。

使用简单线性模型时，估计偏差的一种简单方法是 'choose' 从哪个模型估计偏差。比方说 Y = 3 + 4 * X + e。我选择了 beta <- c(3,4)，因此我只需要模拟我的数据。对于线性模型，模型假设为

Observations are independent

Observations are normally distributed

The mean can be described as by the linear predictor

使用这 3 个假设，模拟固定设计很简单。

set.seed(1)
xseq <- seq(-10,10)
xlen <- length(xseq)
nrep <- 100
#Simulate X given a flat prior (uniformly distributed. A normal distribution would likely work fine as well)
X <- sample(xseq, size = xlen * nrep, replace = TRUE)
beta <- c(3, 4) 
esd = 1
emu <- 0
e <- rnorm(xlen * nrep, emu, esd)
Y <- cbind(1, X) %*% beta + e
fit <- lm(Y ~ X)
bias <- coef(fit) -beta

>bias
 (Intercept)            X 
0.0121017239 0.0001369908

表示有小偏差。为了测试这种偏差是否显着，我们可以执行 wald-test 或 t-test（或将过程重复 1000 次，并检查结果的分布）。

#Simulate linear model many times
model_frame <- cbind(1,X) 
emany <- matrix(rnorm(xlen * nrep * 1000, emu, esd),ncol = 1000)
#add simulated noise. Sweep adds X %*% beta across all columns of emany
Ymany <- sweep(emany, 1, model_frame %*% beta, "+")
#fit many models simulationiously (lm is awesome!)
manyFits <- lm(Y~X)
#Plot density of fitted parameters
par(mfrow=c(1,2))
plot(density(coef(manyFits)[1,]), main = "Density of intercept")
plot(density(coef(manyFits)[2,]), main = "Density of beta")
#Calculate bias, here i use sweep to substract beta across all rows of my coefficients
biasOfMany <- rowMeans(sweep(coef(manyFits), 1, beta, "-"))

>biasOfMany
  (Intercept)             X 
 5.896473e-06 -1.710337e-04

在这里我们看到偏差减少了很多，并且改变了 betaX 的符号，这让我们有理由相信偏差是微不足道的。

改变设计将允许人们使用相同的方法研究交互偏差、异常值和其他东西。

对于线性混合模型，可以执行相同的方法，但是在这里你必须设计随机变量，这需要更多的工作，以及 lmer 据我所知，模型不适合 Y.

的所有列

但是 b（随机效应）可以被模拟，任何噪声参数也可以被模拟。但是请注意，由于 b 是包含单个模拟结果的单个向量（通常是多元正态分布），因此必须 re-run 每个 b 模拟的模型.基本上，这将增加 re-run 模型拟合过程的次数，以便获得对偏差的良好估计。

在 R 模拟中估计线性回归和线性混合模型中的偏差

Estimating bias in linear regression and linear mixed model in R simulation

beta

simulation

r

linear-regression

mixed-models