2

伙计们，

我正在尝试使用 {nlme} 包中的 lme() 在 R 中拟合一个模型，该模型具有许多固定效果和附加随机效果（必须是 {nlme}，因为我还想稍后包含一个 AR(1) 相关矩阵） . 但是，该模型非常慢（由于大量的固定效应）。

现在，我们知道运行带有假人的模型在数值上但不是在计算上等同于在去均值数据上运行模型（在估计范围内）。因此，以下两个模型会产生相同的结果，但第一个模型运行得更快：

# Load library
library(nlme)   # For lme()
library(lfe)    # For demeanlist()

# Create reproducible dummy dataset
set.seed(1)
dat <- data.frame(y=rnorm(1000), x=rnorm(1000), cluster1=sample(LETTERS, 1000,  replace=TRUE), cluster2=sample(LETTERS, 1000, replace=TRUE))

# Create de-meaned dataset
dat_demean <- demeanlist(dat[,c("y", "x")], list(dat$cluster1))
dat_demean$cluster2 <- dat$cluster2 # Copy cluster2 column over to de-meaned data

# Model 1: Fixed effect on cluster1
orig <- lm( y ~ x + cluster1, data=dat)
deme <- lm( y ~ x - 1, data=dat_demean)

# This works as expected. Coefficients are exactly identical (SE are wrong and need to be fixed).
all.equal( coef(orig)["x"], coef(deme)["x"] )
R> [1] TRUE

当我们想在第二个分组变量上额外包含随机效应时，可以做类似的事情吗？以下是我的尝试，它不起作用。

# Model 2: Fixed effect on cluster1 and additional random effect for cluster2
orig <- lme( y ~ x + cluster1, random=~1|cluster2, data=dat)
deme <- lme( y ~ x -1, random=~1|cluster2, data=dat_demean)

# This does not work - coefficients are different
fixef(orig)["x"]
R>        x 
R> -0.005885881 
fixef(deme)["x"]
R>        x 
R> -0.006169549 
all.equal(fixef(orig)["x"], fixef(deme)["x"])
R> [1] "Mean relative difference: 0.04819478"

我错过了什么或我做错了什么？或者，将随机效应模型拟合到贬义数据上不是那么简单吗？

干杯

4

1 回答 1

0

我现在知道这是一个老问题，但我刚刚抓住了它。我相信您的估计在第二个示例中存在一定程度的差异，因为您没有对集群之间的方差进行建模。一般来说，如果簇间方差不太大，这种方法将提供类似的估计，除了线性 fe 方法之外，它还增加了一些关于随机效应的参数假设。

下面的示例parameters::demean使用方便（它有一个有用的命名约定）。

library(nlme)   # For lme()
library(lfe)    # For demeanlist()

# Create reproducible dummy dataset
set.seed(1)
dat <- data.frame(y=rnorm(1000), x=rnorm(1000), cluster1=sample(LETTERS, 1000,  replace=TRUE), cluster2=sample(LETTERS, 1000, replace=TRUE))

dat <- cbind(parameters::demean(dat, select = c("x"), group = "cluster1"), dat[,c("x","y","cluster1","cluster2")])

names(dat)
# [1] "x_between" "x_within"  "x"         "y"         "cluster1"  "cluster2" 
orig <- lm( y ~ x + cluster1, data=dat)
deme <- lm( y ~ x_within - 1, data=dat)

round(coef(deme)["x_within"]  - coef(orig)["x"], 3)
# x_within 
# 0 

rm(deme)
rm(orig)

## second example 
orig <- lme( y ~ x + cluster1, random=~1|cluster2, data=dat)
(origx <- fixef(orig)["x"])
# x 
# 0.002395727 
deme <- lme( y ~ x_within + x_between -1 , random=~1|cluster2, data=dat)
(origdeme <- fixef(deme)["x_within"])
#  x_within 
# 0.002395727 

round(origx  - origdeme, 3)
# x 
# 0
于 2022-01-04T22:27:39.873 回答