如何解决有没有办法用 plm 或其他包从 panelAR 复制功能？
我正在尝试复制和扩展 2011 年的一项研究，该研究恰好用作演示之一 “panelAR”包中的示例。我不知道究竟是如何或为什么，但演示代码从原始研究的一个部分产生了完全相同的回归结果。其中一位作者在他们的网站上发布了他们的复制代码，但它在 Stata 中，所以我无法在“panelAR”演示代码中继续了解它如何完成与 Stata 代码相同的事情。
这里是原始 STATA code、article 和 data 的链接。
我已经能够成功地使用“panelAR”代码对我的新数据运行回归，但遗憾的是“panelAR”对象与“stargazer”不兼容，这是我的包用于制作我的格式化表格。
综上所述，有没有办法使用不同的面板数据包或包的组合来复制以下代码？我试过使用“plm”、“pcse”和“nmle”，但没有成功。
data(LupPon)
tibble(LupPon)
# A tibble: 858 x 14
country id year redist ratio9050 ratio5010 ratio9010 skew turnout fempar propind pvoc union unempl
<chr> <int> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 Australia 1 1960 NA NA NA NA NA NA NA NA NA NA NA
2 Australia 1 1961 NA NA NA NA NA NA NA NA NA NA NA
3 Australia 1 1962 NA NA NA NA NA NA NA NA NA NA NA
4 Australia 1 1963 NA NA NA NA NA NA NA NA NA NA NA
5 Australia 1 1964 NA NA NA NA NA NA NA NA NA NA NA
6 Australia 1 1965 NA NA NA NA NA NA NA NA NA NA NA
7 Australia 1 1966 NA NA NA NA NA NA NA NA NA NA NA
8 Australia 1 1967 NA NA NA NA NA NA NA NA NA NA NA
9 Australia 1 1968 NA NA NA NA NA NA NA NA NA NA NA
10 Australia 1 1969 NA NA NA NA NA NA NA NA NA NA NA
# … with 848 more rows
LupPon <- LupPon[!is.na(LupPon$redist),]
LupPon$redist.lag <- unlist(by(LupPon,LupPon$id,function(x){c(NA,x[,"redist"]
[1:(length(x[,"redist"])-1)])}))
LupPon$time <- unlist(by(LupPon,function(x) seq(1:nrow(x))))
out1 <- panelAR(redist ~ redist.lag + ratio9050 + ratio5010 + turnout + fempar + propind +
pvoc + union + unempl,data=LupPon,panelVar='id',timeVar='time',autoCorr='ar1',panelCorrMethod='pcse',rho.na.rm=TRUE,panel.weight='t-1',bound.rho=TRUE)
summary(out1
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 68 Avg obs. per panel 4.5333
Number of panels: 15 Max obs. per panel 9
Number of times: 9 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -3.26666 11.15944 -0.293 0.770776
redist.lag 0.50658 0.12652 4.004 0.000179 ***
ratio9050 3.81044 3.35976 1.134 0.261402
ratio5010 -4.76833 2.06327 -2.311 0.024405 *
turnout 0.09781 0.03644 2.684 0.009454 **
fempar 0.09134 0.05464 1.672 0.099973 .
propind 0.07253 2.54464 0.029 0.977360
pvoc 0.01860 0.03668 0.507 0.613909
union 0.08862 0.03736 2.372 0.021029 *
unempl 0.12415 0.13443 0.923 0.359580
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.8886
Wald statistic: -708.2314,Pr(>Chisq(9)): 1
** Redistribution models (Table 2)
preserve
keep if redist!=.
sort id year
by id: egen order = seq()
tsset id order
xtpcse redist l1.redist dvpratio9050 dvpratio5010 dvturnout dvfempar dvstddisp_gall dvpvoc dvunion dvunempl,pairwise cor(ar1)
predict pred if e(sample),xb
gen resid = redist-pred
egen stresid=std(resid)
gen outlier = 0 if e(sample)
replace outlier = 1 if abs(stresid)>1.5
编辑 3：以下是 R 和 STATA 中接下来 7 个回归模型的所有代码块
#### Regressions in R
# Removing outliers...
mod1.resid <- out1$residuals
index <- which(abs((mod1.resid-mean(mod1.resid))/sd(mod1.resid)) <= 1.5)
LupPon.nooutlier <- out1$model[index,]> out2 <- panelAR(redist ~ redist.lag + ratio9050 + ratio5010 + turnout + fempar + propind + pvoc + union + unempl,data=LupPon.nooutlier,bound.rho=TRUE)
The following units have non-consecutive observations. Use runs.analysis() on output for additional details: 12,15,16,17,4,6.
Panel-specific correlations bounded to [-1,1]
summary(out2)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 58 Avg obs. per panel 3.8667
Number of panels: 15 Max obs. per panel 8
Number of times: 9 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.57080 7.27261 0.078 0.9378
redist.lag 0.49404 0.07800 6.333 7.74e-08 ***
ratio9050 6.04188 2.81801 2.144 0.0371 *
ratio5010 -6.58628 1.32426 -4.974 8.82e-06 ***
turnout 0.06427 0.02554 2.516 0.0153 *
fempar 0.07852 0.03606 2.178 0.0344 *
propind -2.46670 2.05462 -1.201 0.2358
pvoc 0.01582 0.02327 0.680 0.4999
union 0.12558 0.01634 7.686 6.59e-10 ***
unempl 0.04132 0.10911 0.379 0.7066
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.931
Wald statistic: 2323.5873,Pr(>Chisq(9)): 0
out3 <- panelAR(redist ~ ratio9050 + ratio5010 + as.factor(id),bound.rho=TRUE)
Panel-specific correlations bounded to [-1,1]
summary(out3)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 77 Avg obs. per panel 5.1333
Number of panels: 15 Max obs. per panel 10
Number of times: 10 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 13.97114 9.87777 1.414 0.162412
ratio9050 14.05312 4.99515 2.813 0.006619 **
ratio5010 -8.13602 4.90541 -1.659 0.102420
as.factor(id)3 13.18767 1.50276 8.776 2.36e-12 ***
as.factor(id)4 -0.69241 2.81607 -0.246 0.806616
as.factor(id)5 11.97750 2.20502 5.432 1.07e-06 ***
as.factor(id)6 10.30933 1.94688 5.295 1.78e-06 ***
as.factor(id)7 -2.09143 1.43608 -1.456 0.150511
as.factor(id)8 -1.66623 1.03527 -1.609 0.112766
as.factor(id)9 -0.07301 2.12339 -0.034 0.972686
as.factor(id)12 6.05386 1.73534 3.489 0.000916 ***
as.factor(id)14 8.45693 1.95346 4.329 5.77e-05 ***
as.factor(id)15 13.59385 2.24826 6.046 1.03e-07 ***
as.factor(id)16 -12.92293 1.34996 -9.573 1.09e-13 ***
as.factor(id)17 -2.62601 1.37326 -1.912 0.060623 .
as.factor(id)18 -9.95612 2.26996 -4.386 4.74e-05 ***
as.factor(id)20 -13.69930 2.20810 -6.204 5.59e-08 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.8806
Wald statistic: 189246.9165,Pr(>Chisq(16)): 0
# Removing outliers...
mod3.resid <- out3$residuals
index <- which(abs((mod3.resid-mean(mod3.resid))/sd(mod3.resid)) <= 1.5)
LupPon.nooutlier <- out3$model[index,]> out4 <- panelAR(redist ~ ratio9050 + ratio5010 + as.factor(id),5,1]
summary(out4)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 68 Avg obs. per panel 4.5333
Number of panels: 15 Max obs. per panel 8
Number of times: 10 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 20.6537 6.0314 3.424 0.00122 **
ratio9050 9.8536 3.2130 3.067 0.00346 **
ratio5010 -7.7280 4.0575 -1.905 0.06248 .
as.factor(id)12 6.8209 1.5622 4.366 6.19e-05 ***
as.factor(id)14 7.1422 1.6633 4.294 7.87e-05 ***
as.factor(id)15 11.7269 1.3660 8.585 1.79e-11 ***
as.factor(id)16 -13.1042 1.3083 -10.016 1.22e-13 ***
as.factor(id)17 -2.3581 1.0988 -2.146 0.03664 *
as.factor(id)18 -8.6729 1.5719 -5.518 1.16e-06 ***
as.factor(id)20 -12.3829 1.4979 -8.267 5.57e-11 ***
as.factor(id)3 12.6117 1.1372 11.091 3.37e-15 ***
as.factor(id)4 -1.8655 2.2007 -0.848 0.40057
as.factor(id)5 12.7513 0.8727 14.612 < 2e-16 ***
as.factor(id)6 8.6724 0.8584 10.102 9.10e-14 ***
as.factor(id)7 -1.1486 1.0426 -1.102 0.27575
as.factor(id)8 -1.7659 1.0488 -1.684 0.09833 .
as.factor(id)9 0.6549 1.6795 0.390 0.69822
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.9517
Wald statistic: 5785762.7106,Pr(>Chisq(16)): 0
out5 <- panelAR(redist ~ redist.lag + ratio9010 + skew + turnout + fempar + propind + pvoc + union + unempl,1]
summary(out5)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 68 Avg obs. per panel 4.5333
Number of panels: 15 Max obs. per panel 9
Number of times: 9 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -14.73371 9.19697 -1.602 0.114585
redist.lag 0.49211 0.12412 3.965 0.000204 ***
ratio9010 -0.01548 1.13592 -0.014 0.989172
skew 10.17135 3.67271 2.769 0.007529 **
turnout 0.10182 0.03629 2.806 0.006819 **
fempar 0.08536 0.05333 1.601 0.114901
propind -0.06816 2.45060 -0.028 0.977905
pvoc 0.01991 0.03702 0.538 0.592875
union 0.09013 0.03607 2.499 0.015316 *
unempl 0.11177 0.13563 0.824 0.413280
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.8918
Wald statistic: -151.0228,Pr(>Chisq(9)): 1
# Removing outliers...
mod5.resid <- out5$residuals
index <- which(abs((mod5.resid-mean(mod5.resid))/sd(mod5.resid)) <= 1.5)
LupPon.nooutlier <- out5$model[index,]> out6 <- panelAR(redist ~ redist.lag + ratio9010 + skew + turnout + fempar + propind + pvoc + union + unempl,1]
summary(out6)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 58 Avg obs. per panel 3.8667
Number of panels: 15 Max obs. per panel 8
Number of times: 9 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -12.43089 6.18074 -2.011 0.0499 *
redist.lag 0.48096 0.07362 6.533 3.83e-08 ***
ratio9010 -0.16200 0.94572 -0.171 0.8647
skew 12.98571 2.58573 5.022 7.48e-06 ***
turnout 0.06363 0.02581 2.466 0.0173 *
fempar 0.07440 0.03485 2.135 0.0379 *
propind -2.37649 1.93445 -1.229 0.2252
pvoc 0.01183 0.02326 0.509 0.6134
union 0.12312 0.01525 8.073 1.71e-10 ***
unempl 0.05119 0.10653 0.480 0.6331
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.9346
Wald statistic: 2147.9426,Pr(>Chisq(9)): 0
out7 <- panelAR(redist ~ ratio9010 + skew + as.factor(id),1]
summary(out7)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 77 Avg obs. per panel 5.1333
Number of panels: 15 Max obs. per panel 10
Number of times: 10 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -4.6646 8.6392 -0.540 0.591242
ratio9010 1.3439 1.5360 0.875 0.385089
skew 24.4739 7.5166 3.256 0.001860 **
as.factor(id)3 12.3092 1.3360 9.214 4.32e-13 ***
as.factor(id)4 -0.0509 2.7927 -0.018 0.985518
as.factor(id)5 11.0080 2.1338 5.159 2.95e-06 ***
as.factor(id)6 9.0069 1.9432 4.635 1.97e-05 ***
as.factor(id)7 -2.6626 1.2938 -2.058 0.043947 *
as.factor(id)8 -1.6262 0.9011 -1.805 0.076137 .
as.factor(id)9 0.6049 2.1973 0.275 0.784038
as.factor(id)12 5.9046 1.6921 3.490 0.000913 ***
as.factor(id)14 7.9706 1.7490 4.557 2.60e-05 ***
as.factor(id)15 11.9357 2.3695 5.037 4.62e-06 ***
as.factor(id)16 -12.8997 1.5345 -8.406 9.96e-12 ***
as.factor(id)17 -2.1192 1.3775 -1.538 0.129196
as.factor(id)18 -9.3785 2.2897 -4.096 0.000128 ***
as.factor(id)20 -13.1480 2.2069 -5.958 1.45e-07 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.8874
Wald statistic: 63452.5759,Pr(>Chisq(16)): 0
# Removing outliers...
mod7.resid <- out7$residuals
index <- which(abs((mod7.resid-mean(mod7.resid))/sd(mod7.resid)) <= 1.5)
LupPon.nooutlier <- out7$model[index,]> out8 <- panelAR(redist ~ ratio9010 + skew + as.factor(id),1]
summary(out8)
Panel Regression with AR(1) Prais-Winsten correction and panel-corrected standard errors
Unbalanced Panel Design:
Total obs.: 67 Avg obs. per panel 4.4667
Number of panels: 15 Max obs. per panel 8
Number of times: 10 Min obs. per panel 1
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.8009 5.5166 0.689 0.49402
ratio9010 -1.5372 0.9529 -1.613 0.11298
skew 24.4161 5.8523 4.172 0.00012 ***
as.factor(id)12 6.1617 1.4439 4.267 8.80e-05 ***
as.factor(id)14 5.0717 1.4991 3.383 0.00140 **
as.factor(id)15 8.3799 1.3962 6.002 2.17e-07 ***
as.factor(id)16 -14.9084 1.2974 -11.491 1.22e-15 ***
as.factor(id)17 -0.7629 1.0720 -0.712 0.47999
as.factor(id)18 -5.5874 1.5338 -3.643 0.00064 ***
as.factor(id)20 -9.3915 1.5103 -6.218 1.00e-07 ***
as.factor(id)3 10.5130 1.0326 10.181 8.77e-14 ***
as.factor(id)4 1.4597 2.0448 0.714 0.47862
as.factor(id)5 10.6512 0.8865 12.015 2.35e-16 ***
as.factor(id)6 6.2242 1.0314 6.035 1.93e-07 ***
as.factor(id)7 -1.5296 0.7421 -2.061 0.04451 *
as.factor(id)8 -1.7578 0.7908 -2.223 0.03079 *
as.factor(id)9 3.7324 1.7773 2.100 0.04079 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
R-squared: 0.9683
Wald statistic: 1080793.3476,Pr(>Chisq(16)): 0
# Regressions in STATA
xtpcse redist l1.redist dvpratio9050 dvpratio5010 dvturnout dvfempar dvstddisp_gall dvpvoc dvunion dvunempl if outlier!=1,pairwise cor(ar1) hetonly
drop pred resid stresid outlier
xi: xtpcse redist dvpratio9050 dvpratio5010 i.id,xb
gen resid = redist -pred
egen stresid=std(resid)
gen outlier = 0 if e(sample)
replace outlier = 1 if abs(stresid)>1.5
xi: xtpcse redist dvpratio9050 dvpratio5010 i.id if outlier!=1,pairwise cor(ar1)
drop pred resid stresid outlier
xtpcse redist l1.redist dvratio9010 dvskew dvturnout dvfempar dvstddisp_gall dvpvoc dvunion dvunempl,xb
gen resid = redist -pred
egen stresid=std(resid)
gen outlier = 0 if e(sample)
replace outlier = 1 if abs(stresid)>1.5
xtpcse redist l1.redist dvratio9010 dvskew dvturnout dvfempar dvstddisp_gall dvpvoc dvunion dvunempl if outlier!=1,pairwise cor(ar1)
drop pred resid stresid outlier
xi: xtpcse redist dvratio9010 dvskew i.id,xb
gen resid = redist -pred
egen stresid=std(resid)
gen outlier = 0 if e(sample)
replace outlier = 1 if abs(stresid)>1.5
xi: xtpcse redist dvratio9010 dvskew i.id if outlier!=1,pairwise cor(ar1)
drop pred resid stresid outlier
restore
解决方法
如果您主要关心的是导出格式化的表格，那么用 texreg 替换 stargazer 比使用其他包复制结果要容易得多。根据我的经验，texreg 包支持 panelAR 对象。例如： screenreg(list(panelAR_model1,panelAR_model2)) 将在控制台中显示一个格式化的表，因此您可以检查输出，而 texreg(list(panelAR_model1,panelAR_model2)) 将生成一个 Latex 表。请参阅 texreg 文档以进行自定义。
版权声明：本文内容由互联网用户自发贡献，该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请发送邮件至 dio@foxmail.com 举报，一经查实，本站将立刻删除。
