如何在 lm_robust() 后获得具有聚集标准误差的边际效应?
How to get the marginal effects after lm_robust() with clustered standard errors?
我正在 运行 进行回归分析,其中包含按年份分类的标准误差。使用 Stata 很容易做到这一点,但我必须使用 R 来做到这一点,所以我 运行 它使用 estimatr
包中的 lm_robust()
函数。问题是我现在必须得到一些变量的边际效应,但我做不到,我猜这是因为集群标准错误。我按照 lm_robust()
手册上的内容进行操作,我看到他们只将 margins 包中的 margins 命令用于其他功能而没有聚集标准错误......有没有人知道我如何获得和绘制边际效应?
set.seed(42)
library(fabricatr)
library(randomizr)
dat <- fabricate(
N = 100, # sample size
x = runif(N, 0, 1), # pre-treatment covariate
y0 = rnorm(N, mean = x), # control potential outcome
y1 = y0 + 0.35, # treatment potential outcome
z = complete_ra(N), # complete random assignment to treatment
y = ifelse(z, y1, y0), # observed outcome
# We will also consider clustered data
clust = sample(rep(letters[1:20], each = 5)),
z_clust = cluster_ra(clust),
y_clust = ifelse(z_clust, y1, y0)
)
然后当我 运行 使用 lm_robust()
函数进行回归时:
library(estimatr)
lmout_cl <- lm_robust(
y_clust ~ z_clust + x,
data = dat,
clusters = clust
)
最后,我尝试获得利润...
library(margins)
mar_cl <- margins(lmout_cl)
但这会导致错误:
Error in attributes(.Data) <- c(attributes(.Data), attrib) :'names' attribute
[1] must be the same length as the vector [0]
问题是 estimatr::lm_robust()
生成一个 "lm_robust"
对象,目前 margins()
似乎不支持该对象。我们可以改用 miceadds::lm.cluster()
——并获得与 Stata 相同的聚类标准误差。
library(miceadds)
lmout_cl <- lm.cluster(y_clust ~ z_clust + x, data=dat, cluster=dat$clust)
这会产生一个包含两个元素的列表,其中正常 lm
对象存储在第一个元素中,具有聚类标准误差的方差-协方差矩阵存储在第二个元素中(参见 str(lmout_cl)
) :
> names(lmout_cl)
[1] "lm_res" "vcov"
margins()
现在可以指定为margins(model=model, vcov=vcov)
,所以我们说:
mar_cl <- with(lmout_cl, margins(lm_res, vcov=vcov))
屈服
> mar_cl
Average marginal effects
stats::lm(formula = formula, data = data)
z_clust x
0.6558 1.444
和
> summary(mar_cl)
factor AME SE z p lower upper
x 1.4445 0.3547 4.0728 0.0000 0.7494 2.1396
z_clust 0.6558 0.1950 3.3633 0.0008 0.2736 1.0379
具有集群标准错误。
与Stata的比较
R
foreign::write.dta(dat, "dat.dta") # export as Stata data to wd
Stata
. use dat, clear
(Written by R. )
. quietly regress y_clust z_clust x, vce(cluster clust)
. mfx
Marginal effects after regress
y = Fitted values (predict)
= .67420391
------------------------------------------------------------------------------
variable | dy/dx Std. Err. z P>|z| [ 95% C.I. ] X
---------+--------------------------------------------------------------------
z_clust*| .6557558 .19498 3.36 0.001 .273609 1.0379 .5
x | 1.444481 .35466 4.07 0.000 .749352 2.13961 .524479
------------------------------------------------------------------------------
(*) dy/dx is for discrete change of dummy variable from 0 to 1
.
我们可以清楚地看到——在这样做时,R 在集群标准误差和边际效应方面与 Stata 产生的结果相同。
对于此错误表示歉意,该错误阻止 margins()
在 estimatr
版本 0.10 和更早版本中使用具有非数字簇的 lm_robust()
对象。这是通过 estimatr::lm_robust()
和 margins::margins()
处理模型中的变量的内部方式创建的。
此错误已得到解决,因此您在 estimatr
内有两个解决方案。
让我先生成数据。
library(fabricatr)
library(randomizr)
dat <- fabricate(
N = 100,
x = runif(N),
clust = sample(rep(letters[1:20], each = 5)),
y_clust = rnorm(N),
z_clust = cluster_ra(clust),
)
获取最新版本estimatr
(v0.11.0)
https://declaredesign.org/r/estimatr上的dev版本修复了这个bug,下个月左右会上CRAN
install.packages("estimatr", dependencies = TRUE,
repos = c("http://r.declaredesign.org", "https://cloud.r-project.org"))
library(estimatr)
lmout_cl <- lm_robust(
y_clust ~ z_clust + x,
data = dat,
clusters = clust
)
library(margins)
mar_cl <- margins(lmout_cl)
使用 CRAN 版本 estimatr
(v0.10.0)
的数字簇
CRAN 上 estimatr
现有版本的解决方法是使用数字簇而不是字符簇
dat <- fabricate(
N = 100,
x = runif(N),
clust = sample(rep(1:20, each = 5)),
y_clust = rnorm(N),
z_clust = cluster_ra(clust),
)
install.packages("estimatr")
library(estimatr)
lmout_cl <- lm_robust(
y_clust ~ z_clust + x,
data = dat,
clusters = clust
)
mar_cl <- margins(lmout_cl)
我正在 运行 进行回归分析,其中包含按年份分类的标准误差。使用 Stata 很容易做到这一点,但我必须使用 R 来做到这一点,所以我 运行 它使用 estimatr
包中的 lm_robust()
函数。问题是我现在必须得到一些变量的边际效应,但我做不到,我猜这是因为集群标准错误。我按照 lm_robust()
手册上的内容进行操作,我看到他们只将 margins 包中的 margins 命令用于其他功能而没有聚集标准错误......有没有人知道我如何获得和绘制边际效应?
set.seed(42)
library(fabricatr)
library(randomizr)
dat <- fabricate(
N = 100, # sample size
x = runif(N, 0, 1), # pre-treatment covariate
y0 = rnorm(N, mean = x), # control potential outcome
y1 = y0 + 0.35, # treatment potential outcome
z = complete_ra(N), # complete random assignment to treatment
y = ifelse(z, y1, y0), # observed outcome
# We will also consider clustered data
clust = sample(rep(letters[1:20], each = 5)),
z_clust = cluster_ra(clust),
y_clust = ifelse(z_clust, y1, y0)
)
然后当我 运行 使用 lm_robust()
函数进行回归时:
library(estimatr)
lmout_cl <- lm_robust(
y_clust ~ z_clust + x,
data = dat,
clusters = clust
)
最后,我尝试获得利润...
library(margins)
mar_cl <- margins(lmout_cl)
但这会导致错误:
Error in attributes(.Data) <- c(attributes(.Data), attrib) :'names' attribute
[1] must be the same length as the vector [0]
问题是 estimatr::lm_robust()
生成一个 "lm_robust"
对象,目前 margins()
似乎不支持该对象。我们可以改用 miceadds::lm.cluster()
——并获得与 Stata 相同的聚类标准误差。
library(miceadds)
lmout_cl <- lm.cluster(y_clust ~ z_clust + x, data=dat, cluster=dat$clust)
这会产生一个包含两个元素的列表,其中正常 lm
对象存储在第一个元素中,具有聚类标准误差的方差-协方差矩阵存储在第二个元素中(参见 str(lmout_cl)
) :
> names(lmout_cl)
[1] "lm_res" "vcov"
margins()
现在可以指定为margins(model=model, vcov=vcov)
,所以我们说:
mar_cl <- with(lmout_cl, margins(lm_res, vcov=vcov))
屈服
> mar_cl
Average marginal effects
stats::lm(formula = formula, data = data)
z_clust x
0.6558 1.444
和
> summary(mar_cl)
factor AME SE z p lower upper
x 1.4445 0.3547 4.0728 0.0000 0.7494 2.1396
z_clust 0.6558 0.1950 3.3633 0.0008 0.2736 1.0379
具有集群标准错误。
与Stata的比较
R
foreign::write.dta(dat, "dat.dta") # export as Stata data to wd
Stata
. use dat, clear
(Written by R. )
. quietly regress y_clust z_clust x, vce(cluster clust)
. mfx
Marginal effects after regress
y = Fitted values (predict)
= .67420391
------------------------------------------------------------------------------
variable | dy/dx Std. Err. z P>|z| [ 95% C.I. ] X
---------+--------------------------------------------------------------------
z_clust*| .6557558 .19498 3.36 0.001 .273609 1.0379 .5
x | 1.444481 .35466 4.07 0.000 .749352 2.13961 .524479
------------------------------------------------------------------------------
(*) dy/dx is for discrete change of dummy variable from 0 to 1
.
我们可以清楚地看到——在这样做时,R 在集群标准误差和边际效应方面与 Stata 产生的结果相同。
对于此错误表示歉意,该错误阻止 margins()
在 estimatr
版本 0.10 和更早版本中使用具有非数字簇的 lm_robust()
对象。这是通过 estimatr::lm_robust()
和 margins::margins()
处理模型中的变量的内部方式创建的。
此错误已得到解决,因此您在 estimatr
内有两个解决方案。
让我先生成数据。
library(fabricatr)
library(randomizr)
dat <- fabricate(
N = 100,
x = runif(N),
clust = sample(rep(letters[1:20], each = 5)),
y_clust = rnorm(N),
z_clust = cluster_ra(clust),
)
获取最新版本estimatr
(v0.11.0)
https://declaredesign.org/r/estimatr上的dev版本修复了这个bug,下个月左右会上CRAN
install.packages("estimatr", dependencies = TRUE,
repos = c("http://r.declaredesign.org", "https://cloud.r-project.org"))
library(estimatr)
lmout_cl <- lm_robust(
y_clust ~ z_clust + x,
data = dat,
clusters = clust
)
library(margins)
mar_cl <- margins(lmout_cl)
使用 CRAN 版本 estimatr
(v0.10.0)
CRAN 上 estimatr
现有版本的解决方法是使用数字簇而不是字符簇
dat <- fabricate(
N = 100,
x = runif(N),
clust = sample(rep(1:20, each = 5)),
y_clust = rnorm(N),
z_clust = cluster_ra(clust),
)
install.packages("estimatr")
library(estimatr)
lmout_cl <- lm_robust(
y_clust ~ z_clust + x,
data = dat,
clusters = clust
)
mar_cl <- margins(lmout_cl)