suest：跨模型比较与广义豪斯曼检验

???? 连享会主页：lianxh.cn

  连享会 · 名师讲坛

 
???? 空间计量 专题
⌚ 2020.12.10-13

???? 主讲：杨海生 (中山大学)；范巧 (兰州大学)

???? 课程主页：https://gitee.com/arlionn/SP

  连享会 · 计量专题

 
???? 因果推断-内生性 专题
⌚   2020.11.12-15
 
???? 主讲： 王存同 (中央财经大学)；司继春(上海对外经贸大学)
???? 课程主页：https://gitee.com/arlionn/YG

作者： 黄俊凯 (中国人民大学财政金融学院 )
邮箱： kopanswer@126.com

目录

温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。

在上一篇推文「SUR：似无相关估计是个啥东东？」 (微信版) 中，我们介绍了 SUR 模型的基本思想、估计方法和一些典型的应用。本文进一步介绍 SUR 模型的两个进阶用法：一是模型之间系数差异的检验；二是广义 Hausman 检验。

1. 理论回顾

sureg 命令的模型设定仅限于 GLS，无法进行标准误的调整，也无法使用其他常见的模型。本文我们介绍 suest 命令，它基于 “得分” 计算出系数的方差-协方差矩阵，因此并不局限于线性回归，也不局限于球形扰动项的严格前提。suest 命令默认 (对模型错误设定) 稳健的标准误，使用者也可以通过 vce(cluster clustervar) 指定 (允许组内相关的) 聚类稳健标准误。

suest 命令只能支持那些可以计算得分 (predict newvar, score ) 的回归，包括但不限于 logit、probit、heckman、poisson、zip 等基于极大似然估计的回归。但随后版本的 Stata 的 regress 命令也支持计算得分 (实际上是残差)，因此现在的 suest 也支持线性回归。

1.1 suest 的基准模型

通常，在相同的数据 (或有重叠的数据集) 上估计不同模型时，估计量会相关。此时，若要检验一个估计量之间的假设 ${\mathbb{H}}_{\mathbb{0}}$，就需要知道估计量的联合分布 (simultaneous distribution)。

考虑一个有 $n$ 个观测值和 $k$ 个模型 $\beta =({\beta }_1,{\beta }_2,\dots ,{\beta }_k)$ 的回归模型，对任意观测值 $i$ 的系数向量 ${\beta }_j$，我们称其对数似然函数的导数为 “得分”— $S_{ij}$：

对任意系数 ${\beta }_j$，其得分 $S_j$ 为：

极大似然估计量 ${\hat{\beta }}_j$ 的得分等于 $0$，有估计方程：

根据 White (1994)，在合适的正则条件 (suitable regularity conditions) 下，${\hat{\beta }}_j$ 是渐近正态分布的，其方差可以被如下 “三明治” 估计量一致的估计：

其中 $D_j$ 是得分向量 $S_j$ 的雅克比矩阵，在极大似然估计中 $D_j$ 可以被对数似然函数的负海塞矩阵 (或费雪信息矩阵 )一致的估计。如果模型以正确的方式设定，三明治的中间项将依概率收敛到 $D_j$，因此 $\mathbb{a}\mathbb{V}\mathbb{a}\mathbb{r}({\beta }_j)$ 可以被雅克比矩阵的逆矩阵 $D_j^{-1}$ 一致的估计。

Weesie (1999) 将得分 “堆栈” (stack) 的估计方程表示为：

在合适的正则条件下，$\beta$ 仍然是渐近联合正态分布的。此时雅克比矩阵 $D$ 与得分向量 $S$ 的关系为：

显然，$D$ 是由 $D_1,D2,\dots ,D_k$ 组成的分块对角矩阵，它的逆矩阵也是分块对角矩阵。我们可以将任意观测值 $i$ 在各个方程的得分级联 (concatenate) 在一起，得到：

从而得到 $\mathbb{a}\mathbb{V}\mathbb{a}\mathbb{r}(\beta )$，可以被如下 “三明治” 估计量一致的估计：

其中，对角元 $\mathbb{a}\mathbb{V}\mathbb{a}\mathbb{r}({\hat{\beta }}_j)$ 的一致估计正如我们之前讨论的样子：

而 “堆栈” 的好处就在于，我们还能方便的得到非对角元 $\mathbb{a}\mathbb{C}\mathbb{o}\mathbb{v}({\hat{\beta }}_{j_1},{\hat{\beta }}_{j_2})$ 的一致估计 (White, 1982 )：

该表达式也是 Rogers (1993) 提出的聚类修正的三明治估计量的一个应用。我们考虑这样一个简单的情形：数据集中前一半的观测值用 logit 模型，后一半用 OLS 模型，并且前一半和后一半之间存在一对一的聚类。如果两个模型之间没有共同的参数，“堆栈” 模型的得分在它不属于的那一半中等于 $0$。而在 Rogers 模型中，我们必须在聚类中加入观测值的得分，实际上就是将得分级联到聚类层面。在上文中，我们正是将得分聚类到观测值层面。

1.2 suest 的假设检验

在估计出 $\beta$ 的渐近联合分布后，我们就可以进行假设检验。考虑一个最简单的情形，检验两个系数 ${\beta }_1$ 和 ${\beta }_2$ 是否相等，即 ${\mathbb{H}}_{\mathbb{0}}:b_1-b_2=0$。我们可以计算出 ${\beta }_1-{\beta }_2$ 的渐近方差：

因此可以方便的为 ${\mathbb{H}}_{\mathbb{0}}$ 构造 Wald 检验统计量：

统计量 $W$ 服从自由度为 $\mathbf{r}\mathbf{a}\mathbf{n}\mathbf{k}(\mathbb{a}\mathbb{V}\mathbb{a}\mathbb{r}({\beta }_1-{\beta }_2))$ 的 ${\chi }^2$ 分布。如果 $W$ 大于临界值，则拒绝原假设。

豪斯曼检验常被用于判断估计量的一致性和有效性。豪斯曼 (1978) 证明了当 ${\beta }_1$ 是有效的估计量 (最小的渐近协方差矩阵)时，有：

所以有：

而由于 ${\beta }_1$ 是有效的估计量，故 $\mathbb{a}\mathbb{V}\mathbb{a}\mathbb{r}({\beta }_1)-\mathbb{a}\mathbb{V}\mathbb{a}\mathbb{r}({\beta }_2)$ 是渐近正定矩阵。但对有限样本，或者说对全部的实际应用场景，$\mathbb{V}\mathbb{a}\mathbb{r}({\beta }_1-{\beta }_2)$ 的正定性是无法保证的。并且，由于豪斯曼检验的定义模糊，常常使实证研究者陷入困惑。而 Weesie (1999) 的研究则利用 “得分” 统计量构造豪斯曼检验中的 ${\chi }^2$ 统计量。正如第一节的推导所展示的那样，这种方法不再要求 ${\beta }_1$ 是有效估计量，就能保证计算出来的 $\mathbb{V}\mathbb{a}\mathbb{r}({\beta }_1-{\beta }_2)$ 始终是正定的，从而改进了豪斯曼检验，该方法又被称为 “广义豪斯曼检验” (generalized hausman test)。

2. suest 实例

2.1 系数差异检验 (不同模型)

logit 模型和 probit 模型的系数有可比性吗？我们用 suest 命令来一探究竟。

sysuse auto, clear

qui logit foreign price rep78
est sto m1
qui probit foreign price rep78
est sto m2

suest m1 m2
estat vce, corr

通过 estat vce, corr 命令，我们可以查看估计量的相关系数矩阵。

Correlation matrix of coefficients of suest model

             | m1_for~n                     | m2_for~n                     
        e(V) |    price     rep78     _cons |    price     rep78     _cons 
-------------+------------------------------+------------------------------
m1_foreign   |                              |                              
       price |   1.0000                     |                              
       rep78 |   0.0660    1.0000           |                              
       _cons |  -0.4275   -0.9164    1.0000 |                              
-------------+------------------------------+------------------------------
m2_foreign   |                              |                              
       price |   0.9957    0.0710   -0.4223 |   1.0000                     
       rep78 |   0.0543    0.9996   -0.9090 |   0.0607    1.0000           
       _cons |  -0.4371   -0.9151    0.9990 |  -0.4357   -0.9084    1.0000 

接下来，我们检验 logit 模型和 probit 模型的系数是否相等。

test [m1_foreign]price = [m2_foreign]price, notest
test [m1_foreign]rep78 = [m2_foreign]rep78, acc

 (1)  [m1_foreign]price - [m2_foreign]price = 0
 (2)  [m1_foreign]rep78 - [m2_foreign]rep78 = 0

             chi2(2) =   13.07
         Prob > chi2 =   0.0015

可以看出，${\chi }^2$ 值高达 $13.07$，$p$ 值低于 $0.01$，logit 模型和 probit 模型对相同样本的回归结果显著不同。

然而，logit 模型和 probit 模型是基于不同的量纲，标准逻辑分布的标准差为 $\sqrt{{\pi }^2/3}$，而标准正态分布的标准差为 $1$。因此，有理由认为 logit 模型的系数和 probit 模型的系数之间成比例，这就必需用非线性形式的 Wald 检验。原假设为：

非线性的 Wald 检验对模型设定相当敏感，我们也可以构造一个更加 “线性化” 的原假设：

我们用 testnl 命令执行非线性 Wald 检验。

testnl [m1_foreign]price/[m1_foreign]rep78 = [m2_foreign]price/[m2_foreign]rep78

      chi2(1) =    0.85
  Prob > chi2 =    0.3567

testnl [m1_foreign]price*[m2_foreign]rep78 = [m2_foreign]price*[m1_foreign]rep78

      chi2(1) =    0.83
  Prob > chi2 =    0.3630

此时，${\chi }^2$ 值始终不大于 $1$，$p$ 值始终大于 $0.1$，说明 logit 模型和 probit 模型的系数之间存在等比例关系。

2.2 系数差异检验 (不同样本)

事实上，suest 用于系数差异检验在「Stata: 如何检验分组回归后的组间系数差异？」 (知乎版) 一文已经介绍过。

eststo clear
sysuse auto, clear

reg price rep78 headroom trunk weight length if foreign == 0
est sto m1
reg price rep78 headroom trunk weight length if foreign == 1
est sto m2

suest m1 m2

test [m1_mean = m2_mean]
test [m1_mean = m2_mean], constant
test [m1_mean]rep78 = [m2_mean]rep78
test [m1_mean]length = [m2_mean]length

上面的代码表明，我们既可以用 test [m1_mean = m2_mean] 检验两组的全部系数是否相等。检验结果如下，两组系数在 $5\mathrm{\%}$ 的水平上存在显著差异：

 (1)  [m1_mean]rep78 - [m2_mean]rep78 = 0
 (2)  [m1_mean]headroom - [m2_mean]headroom = 0
 (3)  [m1_mean]trunk - [m2_mean]trunk = 0
 (4)  [m1_mean]weight - [m2_mean]weight = 0
 (5)  [m1_mean]length - [m2_mean]length = 0

             chi2(5) =   13.15
         Prob > chi2 =    0.0220

如果还要求包括常数项，你还可以加入 constant 选项，此时命令变成 test [m1_mean = m2_mean], constant，检验结果如下，在加入常数项后两组系数的差异更加显著，${\chi }^2$ 值高达 $61$：

 (1)  [m1_mean]rep78 - [m2_mean]rep78 = 0
 (2)  [m1_mean]headroom - [m2_mean]headroom = 0
 (3)  [m1_mean]trunk - [m2_mean]trunk = 0
 (4)  [m1_mean]weight - [m2_mean]weight = 0
 (5)  [m1_mean]length - [m2_mean]length = 0
 (6)  [m1_mean]_cons - [m2_mean]_cons = 0

             chi2(6) =   61.00
         Prob > chi2 =    0.0000

suest 也支持单变量比较，以变量 rep78 为例，使用的命令为 test [m1_mean]rep78 = [m2_mean]rep78，下方的结果表明，$p$ 值为 $0.1929$，说明分组检验中的变量 rep78 的系数没有显著差异。

 (1)  [m1_mean]rep78 - [m2_mean]rep78 = 0

           chi2(  1) =    1.70
         Prob > chi2 =    0.1929

直播回放：实证研究设计 (2.4小时, 88元)

???? https://lianxh.duanshu.com

2.3 广义豪斯曼检验

suest 是一个估计后命令，suest 命令后既可以使用 test、testnl、lincom、nlcom 等命令进行假设检验，也可以进行聚类稳健的标准误调整或调查设计效果。suest 适用于不同模型组成的方程组，比如 OLS 模型和 probit 模型组成的方程组。suest 也适用于由对不同样本的估计组成的方程组，比如有交集的数据或无交集的数据。suest 还可以与 vce(robust) 和 vce(cluster clustvar) 选项连用。此外，尽管是估计后命令，suest 的结果还可以用 eststo 命令保存。

在下面这个例子中，我们将借助「Stata-多元 Logit 模型详解 (mlogit)」 (知乎版) 所需的无关方案独立性检验 (Independence of Irrelevant Alternatives, IIA) 比较豪斯曼检验和广义豪斯曼检验。在该例中，被解释变量 insurance 是一个分类变量 (categorical variable )：等于 $1$ 时，代表消费者选择了事后赔偿金 (Indemnity) 的方案；等于 $2$ 时，代表消费者选择了事前预付款 (Prepaid) 的方案；等于 $3$ 时，代表消费者选择了不保险的方案，其中 indemnity 是多项 logit 回归中的基准方案。

eststo clear
*数据下载地址：https://gitee.com/arlionn/data/blob/master/data01/sysdsn4.dta
use sysdsn4.dta, clear
tab insure

     insure | Freq.     Percent     Cum.
------------+---------------------------
  Indemnity |   294       47.73    47.73
    Prepaid |   277       44.97    92.69
   Uninsure |    45        7.31   100.00
------------+---------------------------
      Total |   616      100.00

在无关方案独立性 (IIA) 假设下，当我们从方案集合 (Alternatives) 中剔除某个方案时，回归系数不发生系统性的改变。因此，我们分别对全样本、排除 uninsure 方案的子样本和排除 prepaid 方案的子样本进行多项 logit 回归，用 est sto 命令保存，并用 est tab 命令展示三个样本的回归结果。

qui mlogit insure age male
est sto m1, title(all three insurance forms)
qui mlogit insure age male if insure != "Uninsure":insure
est sto m2, title(insure != "Uninsure":insure)
qui mlogit insure age male if insure != "Prepaid":insure
est sto m3, title(insure != "Prepaid":insure)

est tab m1 m2 m3, b(%8.3f) star stats(N ll) keep(Prepaid: Uninsure:)

--------------------------------------------------------
    Variable |     m1            m2            m3       
-------------+------------------------------------------
Prepaid      |
         age |   -0.010        -0.010                   
        male |    0.510**       0.514**                 
       _cons |    0.263         0.268                   
-------------+------------------------------------------
Uninsure     |
         age |   -0.005                      -0.004     
        male |    0.475                       0.459     
       _cons |   -1.757***                   -1.802***  
-------------+------------------------------------------
Statistics   |                                          
           N |      615           570           338     
          ll | -551.328      -390.486      -131.768     
--------------------------------------------------------
                legend: * p<0.05; ** p<0.01; *** p<0.001

获得回归结果后，我们可以进行豪斯曼检验：

hausman m2 m1, alleqs constant
hausman m3 m1, alleqs constant

全样本与排除 uninsure 方案的子样本之间的豪斯曼检验结果如下表所示，可以看出模型 m1 和 m2 的系数 (包括常数项的系数) 无差异 (${\chi }^2=0.08$，$p=0.994$ )。

         ---- Coefficients ----
       |    (b)          (B)          (b-B)     sqrt(diag(V_b-V_B))
       |     m2           m1       Difference          S.E.
-------+------------------------------------------------------------
   age | -.0101521    -.0100251     -.0001269               .
  male |  .5144003     .5095747      .0048256        .0123338
 _cons |  .2678043     .2633838      .0044205               .
--------------------------------------------------------------------
                b = consistent under Ho and Ha; obtained from mlogit
 B = inconsistent under Ha, efficient under Ho; obtained from mlogit

    Test:  Ho:  difference in coefficients not systematic

                  chi2(3) = (b-B)'[(V_b-V_B)^(-1)](b-B)
                          =        0.08
                Prob>chi2 =      0.9944
                (V_b-V_B is not positive definite)

全样本与排除 prepaid 方案的子样本之间的豪斯曼检验结果如下表所示，可以看出，模型 m1 和 m3 的系数也无差异，但 ${\chi }^2$ 却令人困惑的为负，Stata 也提醒我们「see suest for a generalized test」。

         ---- Coefficients ----
       |    (b)          (B)         (b-B)     sqrt(diag(V_b-V_B))
       |     m3           m1      Difference          S.E.
-------+-----------------------------------------------------------
   age | -.0041055    -.0051925      .001087        .0021355
  male |  .4591074     .4748547    -.0157473               .
 _cons | -1.801774    -1.756843    -.0449311        .1333421
-------------------------------------------------------------------
               b = consistent under Ho and Ha; obtained from mlogit
B = inconsistent under Ha, efficient under Ho; obtained from mlogit

Test:  Ho:  difference in coefficients not systematic

       chi2(3) = (b-B)'[(V_b-V_B)^(-1)](b-B)
               =    -0.18    chi2<0 ==> model fitted on these
                             data fails to meet the asymptotic
                             assumptions of the Hausman test;
                             see suest for a generalized test

接下来，我们首先用 suest 储存三个模型的回归结果，并分别用 test 命令对 m1 和 m2、m1 和 m3、以及 m1 m2 和 m3 的系数进行广义豪斯曼检验。可以看出，全部 ${\chi }^2$ 值的 $p$ 值最小也有 $0.6845$，因此该被解释变量 insure 满足 IIA 假设。

suest m1 m2, noomitted
test [m1_Prepaid = m2_Prepaid], cons

             chi2(3) =    0.89
         Prob > chi2 =    0.8266

suest m1 m3, noomitted
test [m1_Uninsure = m3_Uninsure], cons

             chi2(3) =    1.49
         Prob > chi2 =    0.6845

suest m1 m2 m3, noomitted
test [m1_Prepaid = m2_Prepaid], cons notest
test [m1_Uninsure = m3_Uninsure], cons acc

             chi2(6) =    1.95
         Prob > chi2 =    0.9240

如果观测值是聚类在城市层面的 (cityid)，经典的「McFadden-Hausman」 检验不再合适，因为没有一个估计量是有效的。但我们可以通过选项 vce(cluster cityid) 来轻松地校正豪斯曼统计需要的统计量，具体如下：

suest m1 m2, vce(cluster cityid) noomitted
test [m1_Prepaid = m2_Prepaid], cons

 (1)  [m1_Prepaid]age - [m2_Prepaid]age = 0
 (2)  [m1_Prepaid]male - [m2_Prepaid]male = 0
 (3)  [m1_Prepaid]_cons - [m2_Prepaid]_cons = 0

             chi2(3) =    0.79
         Prob > chi2 =    0.8529

统计结果再次表明，我们的模型满足 IIA 假定。

温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。

扩展命令

其他似不相关估计的扩展命令

• nlsur 命令用可行广义非线性最小二乘法 (FGNLS) 拟合非线性方程组。该命令提供了灵活的语法，甚至可以用该命令快速编写函数计算器程序。
• xtsur 命令对似不相关的非平衡面板数据集进行单向随机效应估计，最早由 Biorn (2004) 提出。
• bireprob 命令拟合一个二值随机效应 Probit 模型，该模型考虑了随机效应误差项和异质性冲击的相关性，并最大化模拟的似然估计。该命令在使用前需要将数据声明为面板数据 xtset，并要求安装 mdraws 命令。
• sudcd命令用于联合估计一个离散选择模型和一个久期 (生存) 模型。其基本原理类似于 SUR 模型，但允许使用者估计服从指数分布、威布尔分布或对数正态分布的似不相关离散选择-久期方程组。
• rifsureg 命令用再中心化影响函数的估计方法执行似不相关的分位数回归，这是一种重要的无条件分位数回归模型。

联立方程的扩展命令

• reg3 是三阶段最小二乘法 3SLS 的命令。3SLS 最早由 Zellner and Theil (1992) 提出，是联立方程的完全信息估计方法，虽然有偏但比 2SLS 更有效。reg3 的功能十分广泛，除了支持 3SLS，还提供支持似不相关估计的选项 sure、多元回归mvreg、支持逐方程 OLS 的选项 ols 和支持 2SLS 的选项 2sls，使用者可以方便的比较不同方法的异同。估计后命令可以参考 help reg3 postestimation。

参考文献和资料

• 连享会推文，Stata: 如何检验分组回归后的组间系数差异？-Link-
• 连享会推文，Stata: 多元 Logit 模型详解 (mlogit) -Link-
• 连享会推文，SUR：似无相关估计是个啥东东？
• Biørn, Erik. "Regression systems for unbalanced panel data: a stepwise maximum likelihood procedure." Journal of Econometrics 122.2 (2004): 281-291. -Link-
• Boehmke, Frederick J. "The influence of unobserved factors on position timing and content in the NAFTA vote." Political Analysis (2006): 421-438. -Link-
• Hausman, Jerry A. "Specification tests in econometrics." Econometrica: Journal of the econometric society (1978): 1251-1271. -Link-
• Hausman, Jerry, and Daniel McFadden. "Specification tests for the multinomial logit model." Econometrica: Journal of the Econometric Society (1984): 1219-1240. -Link-
• Rogers, William. "Regression standard errors in clustered samples." Stata technical bulletin 3.13 (1994). -Link-
• Weesie, J. 1999. sg121: Seemingly unrelated estimation and the cluster-adjusted sandwich estimator. Stata Technical Bulletin 52: 34–47. Reprinted in Stata Technical Bulletin Reprints, vol. 9, pp. 231–248. College Station, TX: Stata Press. -Link-
• White, Halbert. "Maximum likelihood estimation of misspecified models." Econometrica: Journal of the Econometric Society (1982): 1-25. -Link-
• Zellner, Arnold, and Henri Theil. "Three-stage least squares: simultaneous estimation of simultaneous equations." Henri Theil’s Contributions to Economics and Econometrics. Springer, Dordrecht, 1992. 147-178. -Link-

温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。

???? 因果推断-内生性 专题   ⌚ 2020.11.12-15
???? 主讲： 王存同 (中央财经大学)；司继春(上海对外经贸大学)
???? 课程主页：https://gitee.com/arlionn/YG

???? 空间计量 专题   ⌚ 2020.12.10-13
???? 主讲：杨海生 (中山大学)；范巧 (兰州大学)
???? 课程主页：https://gitee.com/arlionn/SP

???? ???? ???? ????
连享会主页： ???? www.lianxh.cn
直播视频： lianxh.duanshu.com
 

免费公开课：

• 直击面板数据模型：https://gitee.com/arlionn/PanelData - 连玉君，时长：1小时40分钟
• Stata 33 讲：https://gitee.com/arlionn/stata101 - 连玉君, 每讲 15 分钟.
• 部分直播课课程资料下载 ???? https://gitee.com/arlionn/Live (PPT，dofiles等)

温馨提示： 文中链接在微信中无法生效。请点击底部「阅读原文」。

关于我们

• ???? 连享会 ( 主页：lianxh.cn ) 由中山大学连玉君老师团队创办，定期分享实证分析经验。
• ???? 直达连享会：【百度一下：连享会】即可直达连享会主页。亦可进一步添加 主页，知乎，面板数据，研究设计 等关键词细化搜索。

????  连享会小程序：扫一扫，看推文，看视频……

???? 扫码加入连享会微信群，提问交流更方便

???? 连享会学习群-常见问题解答汇总：
???? https://gitee.com/arlionn/WD