通过多次回归提取所选变量的系数、标准误差和 p 值

Question

我有以下假设数据。

library(data.table)

city <- c(1,1,1,1,2,2,2,2,2,3,3,3,3,3,3,3)                                       
village <- c(1,2,3,4,1,2,3,4,5,1,2,3,4,5,6,7)                              
village_status <- c(1,0,1,0,1,1,1,0,0,1,1,1,1,0,0,0)
y <- c(1,1,0, 1, 0, 1, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1)
x1 <- c(2,5,4,4,3,4,5,3,2,1,3,5,1,5,4,5)
x2  <- c(21,33,46,8,19,30,20,2,34,38,19,89,35,67,60,37)
x3 <-  c(23,28,30,15,7,18,29,27,14,22,24,10,12,6,17,10)
datei <- data.table(city, village, village_status, y, x1, x2, x3)

我想做的是：

1. 我想通过根据城市洗牌村庄来随机化村庄状态 1000 次
2. 我想在每次回归中使用随机村庄状态进行 1000 次回归。
3. 我想提取随机村庄状态的系数、标准误差和p值并将它们存储为data.table

对于第一点和第二点，我已经做了以下工作：

1. 我对村庄状态进行了洗牌，并将随机村庄状态放入新列中（randomvil1，randomvil2，...，randomvil1000）

n <- 1000
datei[, paste0("randomvil", 1:n) := replicate(n, sample(village_status), simplify=F), city]

2. 我使用 lapply 进行了 1000 次回归

library(estimatr)

# Extract the randomized village columns
varvil0 <- c(paste0("randomvil", 1:n)) 
varvil <- datei[,..varvil0]


fit <- lapply(varvil, function(randomvil) lm_robust(y~ randomvil + randomvil:x1 + x2 + x3 + factor(city), data=datei, se_type= "stata", clusters= city))

现在，我想提取 randomvil 和 randomvil:x1 的系数、标准误差和 p 值并将它们存储为 data.table。我期望的是与此类似的东西。

|               | Estimate     | Std. Error |    Pr(>|t|)  |
| ------------- | ------------ | --------   | --------     |
| randomvil1    | -0.945474623 | 0.474268   |   0.1843847  |
| randomvil1:x1 | -0.004905517 | 0.012388   |   0.7303777  |
| randomvil2    | -0.005449813 | 0.022198   |   0.828959   |
| randomvil2:x1 | -0.341242368 | 0.167598   |   0.1786814  |

有人可以帮我吗？谢谢你。

Answer 1

要提取相关系数，您可以这样做

do.call(what = "rbind",
        args = lapply(
          X = fit, 
          FUN = \(x) coef(summary(x)) |>  
            {\(x) x[c("randomvil", "randomvil:x1"), c("Estimate", "Std. Error",  "Pr(>|t|)")]}() 
          )
        )

一种简单、简单且可能容易出错的重命名系数的方法可能是

c(rbind(paste0("randomvil", 1L:1000L), paste0("randomvil", 1L:1000L, ":x1")))

大家一起：

result <- do.call(what = "rbind",
                  args = lapply(
                    X = fit, 
                    FUN = \(x) coef(summary(x)) |>  
                      {\(x) x[c("randomvil", "randomvil:x1"), c("Estimate", "Std. Error",  "Pr(>|t|)")]}() 
                  )
                  )

row.names(result) <- c(rbind(paste0("randomvil", 1L:1000L), paste0("randomvil", 1L:1000L, ":x1")))

这给出了

                 Estimate Std. Error   Pr(>|t|)
randomvil1     0.76351541  0.2216514 0.07492725
randomvil1:x1 -0.09627251  0.1202819 0.50745202
randomvil2    -0.64983837  0.9565918 0.56700775
randomvil2:x1  0.15961059  0.3407142 0.68555231
randomvil3    -1.02139308  0.7079342 0.28585879
randomvil3:x1  0.33022175  0.1870911 0.21960329

如果您需要进一步帮助，请告诉我。