在不平衡面板上使用 plm 对杂波样本进行随机效应估计

Random Effects Esimation on cluter sample using plm on a unbalanced panel

我通常在广义最小二乘框架内工作,估计 Wooldridge 的 Introductory (2013) 调用的随机效应和固定效应模型对由个人和时间索引的纵向数据的影响维度.

我一直在使用包中plm()中的可行GLS估计来估计随机效应模型——一些统计数据文学术语混合模型plm() 函数采用 index 参数,我在其中指示个人和时间索引。但是,我现在面对一些数据,其中每个人在每个时间点都有多个度量,即多好的分组结构。

我发现可以使用 包中的 lmer() 来建立这样的模型,但是我对行话的差异以及可能性感到有点困惑框架,我想知道是否正确指定了模型。由于我不熟悉框架和这个术语,我担心我可能会忽略更多内容。

我可以使用 lmer() 复制我常用的 plm() 模型,但我不确定如何添加分组。我试图在下面说明我的问题。

我找了一些和我的数据有点像的数据来说明我的情况。首先是一些需要的包,

install.packages(c("mlmRev", "plm", "lme4", "stargazer"), dependencies = TRUE)

然后是数据

data(egsingle, package = "mlmRev")

egsingle 是一个不平衡的小组,由 1721 名学童组成,分为 60 所学校,跨越五个时间点。这些数据最初与 HLM 软件包(Bryk、Raudenbush 和 Congdon,1996)一起分发,但可以在 软件包中找到,详情请参阅 ? mlmRev::egsingle

一些轻数据管理

dta <- egsingle
dta$Female <- with(dta, ifelse(female == 'Female', 1, 0))

这是数据的片段

dta[118:127,c('schoolid','childid','math','year','size','Female')]
#>     schoolid   childid   math year size Female
#> 118     2040 289970511 -1.830 -1.5  502      1
#> 119     2040 289970511 -1.185 -0.5  502      1
#> 120     2040 289970511  0.852  0.5  502      1
#> 121     2040 289970511  0.573  1.5  502      1
#> 122     2040 289970511  1.736  2.5  502      1
#> 123     2040 292772811 -3.144 -1.5  502      0
#> 124     2040 292772811 -2.097 -0.5  502      0
#> 125     2040 292772811 -0.316  0.5  502      0
#> 126     2040 293550291 -2.097 -1.5  502      0
#> 127     2040 293550291 -1.314 -0.5  502      0

以下是我如何使用 plm()

设置 随机效应模型 而没有 schoolid
library(plm)
reg.re.plm <- plm(math~Female+size+year, dta, index = c("childid", "year"), model="random")
# summary(reg.re.plm)

我可以像这样重现这些结果lme4

require(lme4)
dta$year <- as.factor(dta$year)
reg.re.lmer <- lmer(math~Female+size+year+(1|childid), dta)
# summary(reg.re.lmer)

现在,从阅读 Bates (2010) 的第 2 章“lme4:混合效应建模 使用 R” 我相信这就是我指定模型的方式,包括集群级别,schoolid

reg.re.lmer.in.school <- lmer(math~Female+size+year+(1|childid)+(1|schoolid), dta)
# summary(reg.re.lmer.in.school)

但是,当我查看结果时,我不太相信我实际上已经正确指定了它(见下文)。

在我的实际数据中,重复的措施是在个人内部,但我认为我可以以这个数据为例。对于如何进行的任何建议,我将不胜感激。也许是对 notation/terminology 的工作示例的引用,与 Wooldridge (2013) 中使用的示例相差不远。而且,我如何逆向工作并编写 reg.re.lmer.in.school 模型的规范?

# library(stargazer)
stargazer::stargazer(reg.re.plm, reg.re.lmer, reg.re.lmer.in.school, type="text")
#> =====================================================================
#>                                    Dependent variable:               
#>                     -------------------------------------------------
#>                                           math                       
#>                                panel                   linear        
#>                               linear                mixed-effects    
#>                                 (1)                (2)        (3)    
#> ---------------------------------------------------------------------
#> Female                        -0.025              -0.025     0.008   
#>                               (0.046)            (0.047)    (0.042)  
#>                                                                      
#> size                        -0.0004***          -0.0004***  -0.0003  
#>                              (0.0001)            (0.0001)   (0.0002) 
#>                                                                      
#> year-1.5                     0.878***            0.876***   0.866*** 
#>                               (0.059)            (0.059)    (0.059)  
#>                                                                      
#> year-0.5                     1.882***            1.880***   1.870*** 
#>                               (0.059)            (0.058)    (0.058)  
#>                                                                      
#> year0.5                      2.575***            2.574***   2.562*** 
#>                               (0.059)            (0.059)    (0.059)  
#>                                                                      
#> year1.5                      3.149***            3.147***   3.133*** 
#>                               (0.060)            (0.059)    (0.059)  
#>                                                                      
#> year2.5                      3.956***            3.954***   3.939*** 
#>                               (0.060)            (0.060)    (0.060)  
#>                                                                      
#> Constant                     -2.671***          -2.669***  -2.693*** 
#>                               (0.085)            (0.086)    (0.152)  
#>                                                                      
#> ---------------------------------------------------------------------
#> Observations                   7,230              7,230      7,230   
#> R2                             0.735                                 
#> Adjusted R2                    0.735                                 
#> Log Likelihood                                  -8,417.815 -8,284.357
#> Akaike Inf. Crit.                               16,855.630 16,590.720
#> Bayesian Inf. Crit.                             16,924.490 16,666.460
#> F Statistic         2,865.391*** (df = 7; 7222)                      
#> =====================================================================
#> Note:                                     *p<0.1; **p<0.05; ***p<0.01

经过研究 Robert Long's great answer on stats.stackexchange 我发现模型的正确规格是嵌套设计,即 (1| schoolid /childid)。然而,由于数据的编码方式(unniqe childidschoolid 中),交叉设计或规范,即我上面使用的 (1|childid)+(1|schoolid),产生相同的结果。

这是使用与上述相同数据的插图,

data(egsingle, package = "mlmRev")
dta <- egsingle
dta$Female <- with(dta, ifelse(female == 'Female', 1, 0))

require(lme4)
dta$year <- as.factor(dta$year)

重新运行交叉的 design-model, , reg.re.lmer.in.school, 以进行比较

reg.re.lmer.in.school <- lmer(math~Female+size+year+(1|childid)+(1|schoolid), dta)

这里是嵌套结构

reg.re.lmer.nested <- lmer(math~Female+size+year+(1| schoolid /childid), dta)

最后使用惊人的 包比较两个模型,

# install.packages(c("texreg"), dependencies = TRUE)
#  require(texreg)
texreg::screenreg(list(reg.re.lmer.in.school, reg.re.lmer.nested), digits = 3)
#> ===============================================================
#>                                    Model 1        Model 2      
#> ---------------------------------------------------------------
#> (Intercept)                           -2.693 ***     -2.693 ***
#>                                       (0.152)        (0.152)   
#> Female                                 0.008          0.008    
#>                                       (0.042)        (0.042)   
#> size                                  -0.000         -0.000    
#>                                       (0.000)        (0.000)   
#> year-1.5                               0.866 ***      0.866 ***
#>                                       (0.059)        (0.059)   
#> year-0.5                               1.870 ***      1.870 ***
#>                                       (0.058)        (0.058)   
#> year0.5                                2.562 ***      2.562 ***
#>                                       (0.059)        (0.059)   
#> year1.5                                3.133 ***      3.133 ***
#>                                       (0.059)        (0.059)   
#> year2.5                                3.939 ***      3.939 ***
#>                                       (0.060)        (0.060)   
#> ---------------------------------------------------------------
#> AIC                                16590.715      16590.715    
#> BIC                                16666.461      16666.461    
#> Log Likelihood                     -8284.357      -8284.357    
#> Num. obs.                           7230           7230        
#> Num. groups: childid                1721                       
#> Num. groups: schoolid                 60             60        
#> Var: childid (Intercept)               0.672                   
#> Var: schoolid (Intercept)              0.180          0.180    
#> Var: Residual                          0.334          0.334    
#> Num. groups: childid:schoolid                      1721        
#> Var: childid:schoolid (Intercept)                     0.672    
#> ===============================================================
#> *** p < 0.001, ** p < 0.01, * p < 0.05