Stata：动态面板数据模型与xtabond2应用
👇 连享会 · 推文导航 | www.lianxh.cn
🍎 Stata：Stata基础 | Stata绘图 | Stata程序 | Stata新命令 📘 论文：数据处理 | 结果输出 | 论文写作 | 数据分享 💹 计量：回归分析 | 交乘项-调节 | IV-GMM | 时间序列 | 面板数据 | 空间计量 | Probit-Logit | 分位数回归 ⛳ 专题：SFA-DEA | 生存分析 | 爬虫 | 机器学习 | 文本分析 🔃 因果：DID | RDD | 因果推断 | 合成控制法 | PSM-Matching 🔨 工具：工具软件 | Markdown | Python-R-Stata 🎧 课程：公开课-直播 | 计量专题 | 关于连享会
作者：姚旭生 (天津大学)
邮箱：xsyao@tju.edu.cn
目录
1. Nickell 偏差
2. 动态面板估计方法
3. 实证案例
4. 参考资料
5. 相关推文
温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。或直接长按/扫描如下二维码，直达原文：
为了处理内生性问题，我们常常使用工具变量法。但是实际操作中，工具变量往往难以获取，我们可以通过运用自然实验法、处理效应模型、引入固定效应模型、广义矩估计等替代方法进行解决。如果你的数据是面板数据，那么使用动态面板模型（Dynamic Panel Data，DPD）将会是不错的选择。动态面板模型是解释变量中包含被解释变量滞后项的模型，可以理解为过去会在某种程度上影响未来的模型，通常使用广义矩估计方法进行参数估计。
1. Nickell 偏差
在正式引入动态面板估计方法之前，我们首先简要说明一下为什么之前介绍的静态面板估计方法失效了。在 DPD 模型的场景下，特别是“小 T,大 N”的场景下，之前单项固定效应模型出现了严重的问题。正如 Nickell (Econometrica, 1981)所述，这是因为在每个变量中减去平均值的做法导致了解释变量和误差项之间存在相关性。假设初始模型如下，
一阶差分之后表示为，
可见，差分后的解释变量和差分后的误差项（即一阶移动平均过程，或 MA（1））：前者包含 ，而后者包含 。因此，固定效应模型是不一致的，存在“动态面板偏差”，需要引入其它估计方法。
值得注意的是，我们这里所针对的模型主要是“小 T,大 N”的情形。对长动态面板，可以使用偏差矫正的 LSDV 方法进行估计。
2. 动态面板估计方法
通过观察差分后的方程，可以使用一个简单的工具变量（如，）进行估计。这种思想是由 Anderson 和 Hsiao (AH) 提出的，其 Stata 命令是xtivreg, fd。但是，上面提到的工具变量方法并没有利用样本中所有可用的信息，为此 Arellano 和 Bond (AB) 于 1991 年提出了使用更多工具变量的差分 GMM 估计方法。其考虑如下的方程，
其中， 是严格外生的解释变量， 是前定解释变量，包括 的滞后项和内生解释变量，它们都有可能和 和未观察到的个体效应相关。差分 GMM 的思路是对基本模型进行一阶差分以去除固定效应的影响，并用一组滞后的解释变量作为差分方程中相应变量的工具变量。其 Stata 命令是xtabond，xtdpd，xtabond2。
然而，AB 的估计方法被证实被解释变量的滞后项往往是一阶差分项的弱工具变量，特别在当变量接近随机游走的场景下。而这会产生向下的大的有限样本偏差，为此，Blundell 和 Bond (1998) 提出了系统广义矩估计即系统 GMM 估计方法,其 Stata 命令是xtdpdsys，xtabond2。他们在前人的基础上调整了滞后水平和滞后差分项。系统 GMM 可以提高估计的效率，并且相对来说具有更好的有限样本性质。此外，还可以估计出不随时间变化的变量的系数。更一般的动态面板数据模型如下，
使用 GMM 时需要进行两个假设检验。由于 DPD 估计是基于工具变量方法，因此在应用时进行 Sargan-Hansen 过度识别检验结果尤为重要。另外一个重要的检验是误差项的自相关检验，具体而言误差项的一阶差分应存在一阶自相关而不存在二阶自相关。
3. 实证案例
为了说明上述几个动态面板估计方法的效果，我们参照 Christopher F Baum，使用了 Arellano 和 Bond 的原始数据集，这一数据集可在 Stata 中通过 webuse 获取。这是一个英国 140 家公司 1976-1984 年的年度数据，该数据集是一个非平衡面板数据。在他们的论文中，被解释变量是 n ,是公司的雇佣人数，存在两期滞后。其他变量包括当前和滞后的工资水平 w ，当前、滞后一期和滞后两期的资本存量 k,当前、滞后一期和滞后两期的公司部门产出 ys，所有变量都取对数表示。为了捕捉商业周期的影响，一组时间哑元变量也被加入到模型中。
. webuse abdata.dta,clear
. set linesize 80
. set cformat %5.4f
. set sformat %4.2f
. set pformat %4.3f
. sum n nL1 nL2 w wL1 k kL1 kL2 ys ysL1 ysL2
Variable | Obs Mean Std. Dev. Min Max
---------+---------------------------------------------------
n | 1,031 1.056002 1.341506 -2.263364 4.687321
nL1 | 891 1.083518 1.338469 -2.095571 4.687321
nL2 | 751 1.107716 1.333478 -2.079442 4.687321
w | 1,031 3.142988 .2630081 2.081577 3.8118
wL1 | 891 3.132166 .2639638 2.081577 3.8118
---------+---------------------------------------------------
k | 1,031 -.4415775 1.514132 -4.431217 3.852441
kL1 | 891 -.4131872 1.501461 -4.431217 3.852441
kL2 | 751 -.392113 1.486371 -4.431217 3.852441
ys | 1,031 4.638015 .0939611 4.464758 4.85488
ysL1 | 891 4.651039 .0923352 4.464758 4.85488
---------+---------------------------------------------------
ysL2 | 751 4.67097 .085 4.490414 4.85488
如果我们在估计这个模型时忽略了它的动态面板性质，我们仅使用了基于面板聚类标准误的回归命令 regress：regress n nL1 nL2 w wL1 k kL1 kL2 ys ysL1 ysL2 yr*, cluster(id)
这种方法一个最明显的问题就在未考虑公司层面的异质性。为此，我们使用静态面板数据模型进行回归：xtreg n nL1 nL2 w wL1 k kL1 kL2 ys ysL1 ysL2 yr*, fe cluster(id)
固定效应估计将受到 Nickell 偏差的影响，考虑到可用的时间序列较短，这种偏差可能会很严重。
regress n nL1 nL2 w wL1 k kL1 kL2 ys ysL1 ysL2 yr*, cluster(id)
est store OLS
xtreg n nL1 nL2 w wL1 k kL1 kL2 ys ysL1 ysL2 yr*, fe cluster(id)
est store FE
local mm "OLS FE"
esttab `mm', mtitle(`mm') compress star(* 0.1 ** 0.05 *** 0.01) b(%10.3f) nogap
------------------------------------
(1) (2)
OLS FE
------------------------------------
nL1 1.045*** 0.733***
(20.17) (12.28)
nL2 -0.077 -0.139*
(-1.57) (-1.78)
w -0.524*** -0.560***
(-3.01) (-3.51)
wL1 0.477*** 0.315**
(2.78) (2.20)
k 0.343*** 0.388***
(7.06) (6.82)
kL1 -0.202*** -0.081
(-3.10) (-1.49)
kL2 -0.116*** -0.028
(-3.22) (-0.65)
ys 0.433** 0.469***
(2.42) (2.74)
ysL1 -0.768*** -0.629***
(-3.05) (-3.04)
ysL2 0.312** 0.058
(2.36) (0.44)
------------------------------------
N 751 751
------------------------------------
t statistics in parentheses
* p<0.1, ** p<0.05, *** p<0.01
在原来的 OLS 回归中，滞后因变量与误差正相关，使其系数向上偏置，所以真实的估计系数应该小于 1.045。由于被解释变量滞后一期的组内离差与误差项的组内离差存在负相关系，固定效应的估计结果会产生一个严重向下的偏误系数，从而低估了被解释变量滞后一期系数，所以真实值应该大于 0.733。被解释变量滞后一期的真实系数应该介于 0.733-1.045 之间，使用 OLS 和 FE 回归是一种有用的检验。
下面，我们使用 Anderson–Hsiao 的估计方法来进行参数估计：ivregress 2sls D.n (D.nL1 = nL2) D.(nL2 w wL1 k kL1 kL2 ys ysL1 ysL2 yr1979 yr1980 yr1981 yr1982 yr1983 )
估计结果如下：
-----------------------
(1)
A-H
-----------------------
D.nL1 2.308
(1.17)
D.nL2 -0.224
(-1.25)
D.w -0.810***
(-3.10)
D.wL1 1.422
(1.21)
D.k 0.253*
(1.75)
D.kL1 -0.552
(-0.90)
D.kL2 -0.213
(-0.89)
D.ys 0.991**
(2.14)
D.ysL1 -1.938
(-1.35)
D.ysL2 0.487
(0.96)
-----------------------
N 611
-----------------------
t statistics in parentheses
* p<0.1, ** p<0.05, *** p<0.01
虽然这些结果应该是一致的，但它们相当令人失望。被解释变量滞后一期系数在 OLS 和 FE 之外。接着，我们使用xtabond2命令来进行 GMM 估计。xtabond2命令语法格式：xtabond2 depvar varlist ,gmm(varlist,lag(# #) collapse) iv(varlist) nolevel twostep robust noconstant small or
主要选项的含义是：
depvar：被解释变量；varlist：所有解释变量；gmm(varlist,lag(# #) collapse)：varlist 为所有前定变量和内生变量，lag（# #）为限定滞后期，collapse 为压缩工具变量；iv(varlist)：严格外生的解释变量；nolevel：差分 GMM 或者系统 GMM；twostep：两步 GMM 估计；robust：稳健标准误；noconstant：不要常数项；small：小样本标识；or：是否采用正交变换。更多 options 和含义详见 help xtabond2。
为了阐述 AB 估计方法，我们可以用 xtabond2 重新估计模型，假设唯一的内生性存在的是涉及滞后因变量。回归命令为：xtabond2 n L(1/2).n L(0/1).w L(0/2).(k ys) yr*, gmm(L.n) iv(L(0/1).w L(0/2).(k ys) yr*) nolevel robust small
回归结果如下：
Dynamic panel-data estimation, one-step difference GMM
------------------------------------------------------------------
Group variable: id Number of obs = 611
Time variable : year Number of groups = 140
Number of instruments = 41 Obs per group: min = 4
F(0, 140) = . avg = 4.36
Prob > F = . max = 6
------------------------------------------------------------------
| Robust
n | Coef. Std. Err. t P>|t| [95% Conf. Interval]
--------+---------------------------------------------------------
n |
L1. | 0.6862 0.1469 4.67 0.000 0.3957 0.9767
L2. | -0.0854 0.0569 -1.50 0.136 -0.1979 0.0272
|
w |
--. | -0.6078 0.1811 -3.36 0.001 -0.9658 -0.2498
L1. | 0.3926 0.1707 2.30 0.023 0.0551 0.7301
|
k |
--. | 0.3568 0.0600 5.95 0.000 0.2383 0.4754
L1. | -0.0580 0.0744 -0.78 0.437 -0.2050 0.0890
L2. | -0.0199 0.0332 -0.60 0.549 -0.0857 0.0458
|
ys |
--. | 0.6085 0.1753 3.47 0.001 0.2619 0.9551
L1. | -0.7112 0.2355 -3.02 0.003 -1.1767 -0.2456
L2. | 0.1058 0.1435 0.74 0.462 -0.1779 0.3895
|
yr1978 | -0.0096 0.0105 -0.91 0.362 -0.0302 0.0111
yr1980 | 0.0125 0.0135 0.92 0.357 -0.0142 0.0391
yr1981 | -0.0213 0.0255 -0.84 0.405 -0.0718 0.0291
yr1982 | -0.0366 0.0261 -1.40 0.163 -0.0882 0.0150
yr1983 | -0.0309 0.0282 -1.10 0.275 -0.0866 0.0248
yr1984 | -0.0173 0.0296 -0.58 0.560 -0.0757 0.0412
------------------------------------------------------------------
Instruments for first differences equation
Standard
D.(w L.w k L.k L2.k ys L.ys L2.ys yr1976 yr1977 yr1978 yr1979 yr1980
yr1981 yr1982 yr1983 yr1984)
GMM-type (missing=0, separate instruments for each period unless collapsed)
L(1/8).L.n
------------------------------------------------------------------------------
Arellano-Bond test for AR(1) in first differences: z = -3.60 Pr > z = 0.000
Arellano-Bond test for AR(2) in first differences: z = -0.52 Pr > z = 0.606
------------------------------------------------------------------------------
Sargan test of overid. restrictions: chi2(25) = 67.59 Prob > chi2 = 0.000
(Not robust, but not weakened by many instruments.)
Hansen test of overid. restrictions: chi2(25) = 31.38 Prob > chi2 = 0.177
(Robust, but weakened by many instruments.)
Difference-in-Hansen tests of exogeneity of instrument subsets:
iv(w L.w k L.k L2.k ys L.ys L2.ys yr1976 yr1977 yr1978 yr1979 yr1980 yr1981 yr
> 1982 yr1983 yr1984)
Hansen test excluding group: chi2(11) = 12.01 Prob > chi2 = 0.363
Difference (null H = exogenous): chi2(14) = 19.37 Prob > chi2 = 0.151
如上结果表示，41 个工具变量已经被创建，其中 17 个对应的是 IV 型。此外，被解释变量滞后一期的系数为 0.686。该实例中有 25 个过度识别限制，如下面的第一列所示。hansen_df 代表 Hansen J 检验过度识别限制，该检验的 p 值 hansenp 所示。
-------------------------------------------------
(1) (2) (3)
All lags lags 2-5 lags 2-4
-------------------------------------------------
L.n 0.686*** 0.835** 1.107***
(4.67) (2.59) (3.94)
L2.n -0.085 0.262 0.231
(-1.50) (1.56) (1.32)
w -0.608*** -0.671*** -0.709***
(-3.36) (-3.18) (-3.26)
L.w 0.393** 0.436 0.608**
(2.30) (1.59) (2.32)
k 0.357*** 0.325*** 0.309***
(5.95) (4.95) (4.55)
L.k -0.058 -0.128 -0.211*
(-0.78) (-1.03) (-1.72)
L2.k -0.020 -0.170** -0.202***
(-0.60) (-2.60) (-3.07)
ys 0.609*** 0.640*** 0.698***
(3.47) (3.07) (3.45)
L.ys -0.711*** -0.776** -0.987***
(-3.02) (-2.14) (-2.85)
L2.ys 0.106 0.020 0.112
(0.74) (0.10) (0.56)
-------------------------------------------------
N 611 611 611
-------------------------------------------------
hansen_df 25 16 13
hansenp 0.177 0.676 0.714
-------------------------------------------------
t statistics in parentheses
* p<0.1, ** p<0.05, *** p<0.01
如上表所示，在调整滞后期后，Hansen J 检验的 p 值在限制滞后情况下也较大。另一方面，滞后因变量系数的估计似乎对滞后长度的选择相当敏感。
接着，我们采用正交变换，即 forward orthogonal deviations (FOD)再次进行参数估计。回归命令为：xtabond2 n L(1/2).n L(0/1).w L(0/2).(k ys) yr*, gmm(L.n) iv(L(0/1).w L(0/2).(k ys) yr*) nolevel robust small or
------------------------------------
(1) (2)
All lags FOD
------------------------------------
L.n 0.686*** 0.737***
(4.67) (5.14)
L2.n -0.085 -0.096
(-1.50) (-1.38)
w -0.608*** -0.563***
(-3.36) (-3.47)
L.w 0.393** 0.317**
(2.30) (1.99)
k 0.357*** 0.384***
(5.95) (6.85)
L.k -0.058 -0.083
(-0.78) (-1.18)
L2.k -0.020 -0.046
(-0.60) (-1.06)
ys 0.609*** 0.469***
(3.47) (2.72)
L.ys -0.711*** -0.625***
(-3.02) (-2.86)
L2.ys 0.106 0.042
(0.74) (0.30)
------------------------------------
N 611 611
------------------------------------
hansen_df 25 25
hansenp 0.177 0.170
------------------------------------
t statistics in parentheses
* p<0.1, ** p<0.05, *** p<0.01
相比于一阶变换，使用 FOD 变换时，结果似乎相当稳健，对大多数系数的估计稍微更精确。
最后，我们使用系统 GMM 方法再次进行参数估计。参考 Blundell 和 Bond，我们认为工资和资本是潜在的内生变量，给出如下的回归命令：xtabond2 n L.n L(0/1).(w k) yr*, gmm(L.(n w k)) iv(yr*, equation(level)) robust small
回归结果如下：
Dynamic panel-data estimation, one-step system GMM
---------------------------------------------------------------------
Group variable: id Number of obs = 891
Time variable : year Number of groups = 140
Number of instruments = 113 Obs per group: min = 6
F(12, 139) = 6359.46 avg = 6.36
Prob > F = 0.000 max = 8
---------------------------------------------------------------------
| Robust
n | Coef. Std. Err. t P>|t| [95% Conf. Interval]
--------+------------------------------------------------------------
n |
L1. | 0.9356 0.0266 35.21 0.000 0.8831 0.9881
|
w |
--. | -0.6310 0.1193 -5.29 0.000 -0.8668 -0.3951
L1. | 0.4826 0.1383 3.49 0.001 0.2091 0.7561
|
k |
--. | 0.4839 0.0544 8.89 0.000 0.3763 0.5915
L1. | -0.4244 0.0591 -7.18 0.000 -0.5412 -0.3076
|
yr1977 | -0.0214 0.0223 -0.96 0.339 -0.0655 0.0227
yr1978 | -0.0150 0.0098 -1.54 0.127 -0.0343 0.0043
yr1980 | -0.0147 0.0113 -1.31 0.192 -0.0370 0.0075
yr1981 | -0.0409 0.0144 -2.84 0.005 -0.0694 -0.0124
yr1982 | -0.0070 0.0159 -0.44 0.662 -0.0384 0.0245
yr1983 | 0.0065 0.0170 0.38 0.704 -0.0271 0.0400
yr1984 | 0.0027 0.0207 0.13 0.899 -0.0384 0.0437
_cons | 0.5495 0.1922 2.86 0.005 0.1694 0.9297
---------------------------------------------------------------------
Instruments for first differences equation
GMM-type (missing=0, separate instruments for each period unless collapsed)
L(1/8).(L.n L.w L.k)
Instruments for levels equation
Standard
yr1976 yr1977 yr1978 yr1979 yr1980 yr1981 yr1982 yr1983 yr1984
_cons
GMM-type (missing=0, separate instruments for each period unless collapsed)
D.(L.n L.w L.k)
------------------------------------------------------------------------------
Arellano-Bond test for AR(1) in first differences: z = -5.46 Pr > z = 0.000
Arellano-Bond test for AR(2) in first differences: z = -0.25 Pr > z = 0.804
------------------------------------------------------------------------------
Sargan test of overid. restrictions: chi2(100) = 157.54 Prob > chi2 = 0.000
(Not robust, but not weakened by many instruments.)
Hansen test of overid. restrictions: chi2(100) = 110.70 Prob > chi2 = 0.218
(Robust, but weakened by many instruments.)
Difference-in-Hansen tests of exogeneity of instrument subsets:
GMM instruments for levels
Hansen test excluding group: chi2(79) = 84.33 Prob > chi2 = 0.320
Difference (null H = exogenous): chi2(21) = 26.37 Prob > chi2 = 0.193
iv(yr1976 yr1977 yr1978 yr1979 yr1980 yr1981 yr1982 yr1983 yr1984, eq(level))
Hansen test excluding group: chi2(93) = 107.79 Prob > chi2 = 0.140
Difference (null H = exogenous): chi2(7) = 2.91 Prob > chi2 = 0.893
我们发现，被解释变量滞后一期的系数为 0.936，比 AB 估计要高得多，但处于 0.733-1.045 之间。其中，系统 GMM 创建了 113 个工具变量，在过度识别限制的检验中具有 100 度的自由度。注意，当工具变量较多时，可以使用collapse选项进行数量压缩。此外，在进行参数估计后，仍需进行上述提到的两项检验。$2
4. 参考资料
Boston College Course：Dynamic Panel Data estimators (Christopher F Baum), -Link- Anderson, T. W., C. Hsiao, 1982, Formulation and estimation of dynamic models using panel data, Journal of Econometrics, 18 (1): 47-82. -Link-, -PDF1- Arellano, M., S. Bond, 1991, Some tests of specification for panel data: Monte carlo evidence and an application to employment equations, The Review of Economic Studies, 58 (2): 277-297. -Link-, -PDF1- Blundell, R., S. Bond, 1998, Initial conditions and moment restrictions in dynamic panel data models, Journal of Econometrics, 87 (1): 115-143. -Link-, -PDF1-
5. 相关推文
Note：产生如下推文列表的 Stata 命令为：
lianxh 动态面板 xtabond 面板数据模型
安装最新版lianxh命令：
ssc install lianxh, replace
专题：面板数据 Stata：如何估计包含非时变变量的动态面板模型-xtseqreg 面板数据模型-xtdcce2：动态共同相关和截面相关 Stata：动态面板数据模型 OLS 估计的偏差 Stata 实操陷阱：动态面板数据模型 Stata：xtabond 命令中的 endog()选项设定问题 Stata: 面板数据模型一文读懂 Stata: 动态面板门槛模型 专题：空间计量 空间面板数据模型及 Stata 实现
New！ Stata 搜索神器：
lianxh和songblGIF 动图介绍
搜： 推文、数据分享、期刊论文、重现代码 ……
👉 安装：
. ssc install lianxh
. ssc install songbl
👉 使用：
. lianxh DID 倍分法
. songbl all
🍏 关于我们
连享会 ( www.lianxh.cn，推文列表) 由中山大学连玉君老师团队创办，定期分享实证分析经验。 直通车： 👉【百度一下：连享会】即可直达连享会主页。亦可进一步添加 「知乎」,「b 站」,「面板数据」,「公开课」 等关键词细化搜索。