gsreg：自动模型设定和变量筛选

发布时间：2021-02-09 阅读 14939

Stata连享会主页 || 视频 || 推文 || 知乎

温馨提示： 定期清理浏览器缓存，可以获得最佳浏览体验。

New！ lianxh 命令发布了：
随时搜索连享会推文、Stata 资源，安装命令如下：
. ssc install lianxh
使用详情参见帮助文件 (有惊喜)：
. help lianxh

✌ 课程详情： https://gitee.com/arlionn/Course | lianxh.cn

✨ 课程主页： https://gitee.com/arlionn/Course

作者：肖泽林 (中山大学)
邮箱：xiaozlin@mail2.sysu.edu.cn

编者按：本文摘译自下文，特此致谢！

Source：Gluzmann P, Panigo D. Global search regression: A new automatic model-selection technique for cross-section, time-series, and panel-data regressions[J]. The Stata Journal, 2015, 15(2): 325-349. -PDF-

1. 背景介绍

计量经济学者通常面临全局优化问题，即从大量可选的计量经济模型中识别真正的数据生成过程 (DGP)。但是大多数模型自动选择技术 (AMST) 无法保证模型选择中的 “全局最优”，即根据不同搜索参数和搜索起点，将获得不同的最终结果。

一些新的 AMST (如 PC-GETS 或 RETINA) 通过多路径-多样本向后和向前查找的方法避免这一问题。但由于未探索的还原路径、顺序测试的数据大小-速率权衡和累积类型 I 错误，他们仍然无法保证全局最优。

近年来，计算机的计算能力呈指数增长，使计量经济模型选择问题中的详尽搜索得以实现 (如 vselect 命令)。然而，他们都无法提供高精度 (当使用样本外选择标准或需要假设检验时)，以及没有为敏感性分析提供详尽的结果。

本推文将要介绍的 GSREG (全局搜索回归) 是一种较新的横截面、时间序列和面板数据回归模型自动选择技术。与其他穷举搜索算法 (如 VSELECT) 一样，GSREG 避免了标准后向和前向搜索方法的特征路径依赖缺陷。其优势如下：

当使用样本外选择标准时，可以保证最优性；
允许对每个备选方案进行残差测试；
保留使用样本内选择标准时应用 VSELECT 的 "leaps-and-bound" 子算法的快捷思路；
为每个备选模型提供结果统计的完整信息数据集。

2. 其他常用回归模型自动选择技术

根据不同模型选择技术的特点，在过往研究分类的基础上，Gluzmann 和 Panigo (2015) 总结出以下 “概念树”：

2.1 PC-GETS

根据 Hoover (2006) 的定义，GETS 方法从尽可能广泛的一般情况开始，然后在可能的限制空间内搜索，以找到最简洁的情况。

在早期的程序 (例如 STEPWISE) 中，变量输入 (或删除) 的顺序和候选协变量的数量决定了最终的结果。但随着 PC-GETS 代码的出现，小样本限制、路径依赖和累积类型 I 错误等问题已经被部分克服。

虽然 PC-GETS 足够灵活，可以应用于许多其他方面，但其主要还是用于大型模型选择问题，这些问题目标通常是最大限度地提高样本内拟合优度和残差稳健性。

2.2 RETINA

另一种自动启发式模型选择技术是输入网络方法的相关转换算法 (RETINA)，由 Perez-Amaral (2003) 开发，后被 Marinucci (2008) 拓展的 GAUSS 软件包。

RETINA 与 PC-GETS 在方法和目标上有所相同，RETINA 算法是向前查找 (从具体到一般)，其主要目标是预测。在小样本模型选择问题中，RETINA 特别有用，因为解释变量的条件均值函数形式是未知的。同时，凭借他在样本外良好的特性，迅速吸引了大量用户。

然而，上述 PC-GETS 的许多缺陷也适用于 RETINA。 RETINA 仍然不能保证模型选择的最优性，即无论是路径依赖还是累积类型 I 错误的问题都不能完全消除。

2.3 vselect

VSELECT 是 Lindsey 和 Sheter (2010) 开发的一种用于模型选择的 Stata 代码。他包括 Stata 内置逐步算法的一个稍微不同的版本 (其中的 Wald 显着性检验被一组可选信息标准如 $R^{2}$ 、AIC、BIC 等取代)，以及 Furnival 和 Wilson (1974) 的 “Leaps & Bound” 代码 (在线性回归中执行 “有效详尽” 的样本内搜索模型选择) 改编版本。

veselect 命令通过 best 选项提供了 “Leaps & Bound” 子算法选择模型。该算法将所有可能的模型组织成树结构并进行扫描，跳过那些绝对不是最优的。树中的每个节点对应于两组预测器 (predictor)。预测器列表是基于原始回归中所有预测器 t 统计量的自动排序创建的。当算法检查一个节点时，他将该对预测器列表的回归与已经进行的每个预测器的最优回归进行比较。根据结果，该算法可以跳过该节点的所有或部分后代，最终算法能够在找到最优预测器列表并只检查所有可能回归的一小部分之后就完成。

VSELECT-BEST 克服了启发式方法最重要的弱点，并成为中等大小模型选择问题的最佳方案，其样本内解释和非样本外预测是主要目标。但此算法还存在一些问题：

不适用于样本外模型选择问题；
大型模型选择问题变得不可行和非常耗时；
由于其设计不是为了稳健性分析，因此他只保留每个子集大小的最佳模型结果。

3. gsreg 命令介绍

3.1 gsreg 的步骤

gsreg 命令有两个主要步骤：

第一步，他将创建一组集合，每个集合包含一个可能的解释变量集。因此，全套集合包含所有可能的候选解释变量组合。
第二步，命令将对先前创建的每个集合执行回归。

第一阶段

根据用户指定的原始变量集合和基于候选变量的其他协变量确定候选变量的总集 Lvc；
若未指定 ncomb(#1,#2) 选项，则通过不重复选取候选变量所有可能的组合创建集合 SL。否则，SL 是通过不重复选取数量从 #1 到 #2 的候选变量得到的。此时， $S L = ⟮ L_{i n t 1}, \dots, L_{i n t 2} ⟯$ ，其中每个 $L_{i}$ 是候选变量集 Lvc 的特定子集；
若指定 square 选项，将对上一点的每个 $L_{i}$ 中创建平方项集合 $S q L_{i}$ ，并添加到集合 SL 中。指定cubic同理；
若指定 interactions 选项，将为每个 $L_{i}$ 创建一个 $I n t L_{i}$ 集合，其中不仅包含 $L_{i}$ ，还包括所有 $L_{i}$ 变量之间可能的交互组合，然后将 $I n t L_{i}$ 集合添加到 SL 中；
通过指定 fixinteractions 选项，可以为每个 $L_{i}$ 创建一个新 $F i n t L_{i}$ 集合，他不仅包括所有 $L_{i}$ 变量，而且还包括 $L_{i}$ 和 fixvar 变量之间所有可能的交互项；
若指定 schange 选项，将从 SL 创建一个新集合 SC 以测试是否存在结构变化，其中包括 SL 变量与用户定义的结构变化变量 (如虚拟变量) 的交互项。

第二阶段

在第二阶段，gsreg 执行每个 SL 集合的回归，并在一个 Stata.dbf 文件中保存系数和不同的 (默认和用户定义的) 统计数据。下图总结了 gsreg 程序的步骤流程：

3.2 gsreg 语法结构

gsreg 命令是 Gluzmann 和 Panigo (2015) 根据以上步骤所编写的 Stata 新命令。默认情况下，gsreg 执行将替代 OLS 回归，寻找最佳的被解释变量数据生成过程，迭代指定的解释变量之间的所有可能组合。回归结果将自动保存在名为 gsreg.dta 的文件中。

gsreg 命令属于外部命令，安装方法如下：

ssc install gsreg, replace

gsreg命令的语法结构如下：

gsreg depvar varlist_ocand [weight] [if] [in] [,          ///
      ncomb(#1,#2) samesample vselect dlags(numlist)      ///
      ilags(numlist) lags(numlist) fixvar(fixvarlist)     ///
      schange(varschange) interactions squares cubic      ///
      fixinteractions outsample(#) cmdest(commandname)    ///
      cmdoptions(commandoptions) cmdstat(commandstats)    ///
      cmdiveq(varlist_end=valist_inst) aicbic hettest     ///
      hettest_o(hettestmoptions) archlm                   ///
      archlm_o(archlmoptions) bgodfrey                    ///
      bgodfrey_o(bgodfreyoptions) durbinalt               ///
      durbinalt_o(durbinaltoptions) dwatson sktest        ///
      sktest_o(sktestoptions) swilk swilk_o(swilkoptions) ///
      sfrancia testpass(#) resultsdta(newbasename)        ///
      replace double nocount compact nindex(lcimplist)    ///
      mindex(lcimplist) best(#) backup(#) part(#1, #2)]

下面以使用 gsreg 命令的例子简单演示其使用方法，具体命令如下：

sysuse auto.dta, clear
gsreg mpg weight foreign

在这种情况下，有两个候选协变量 (weight 和 foreign)，gsreg 将从中不重复地执行所有可能的组合：

reg mpg weight
reg mpg foreign
reg mpg weight foreign

在以上三种解释变量组合中，gsreg 命令对三个回归结果的调整 $R^{2}$ 自动排序，并报告调整 $R^{2}$ 最高的回归结果。

----------------------------------------------------
Total Number of Estimations: 3
----------------------------------------------------
Computing combinations...
Preparing regression list...
Doing regressions...
Estimation number 1 of 3
Estimation number 2 of 3
Estimation number 3 of 3
Saving results...
file gsreg.dta saved
----------------------------------------------------
Best estimation in terms of r_sqr_a 
Estimation number 3
----------------------------------------------------

  Source |       SS           df       MS      Number of obs   =        74
---------+----------------------------------   F(2, 71)        =     69.75
   Model |   1619.2877         2  809.643849   Prob > F        =    0.0000
Residual |  824.171761        71   11.608053   R-squared       =    0.6627
---------+----------------------------------   Adj R-squared   =    0.6532
   Total |  2443.45946        73  33.4720474   Root MSE        =    3.4071
--------------------------------------------------------------------------
     mpg |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
---------+----------------------------------------------------------------
  weight |     -0.007      0.001   -10.34   0.000       -0.008      -0.005
 foreign |     -1.650      1.076    -1.53   0.130       -3.796       0.495
   _cons |     41.680      2.166    19.25   0.000       37.362      45.998
--------------------------------------------------------------------------

除返回最优估计的结果外，gsreg 命令还自动将三次回归的相关数据保存在 gsreg.dta 文件中。通过输入以下命令，查看数据结构：

. use gsreg.dta, clear
. list, noobs

+----------------------------------------------------------------------------------------------------------------------------+
| order       v_1_b       v_1_t       v_2_b       v_2_t      v_3_b      v_3_t   obs   nvar    r_sqr_a    rmse_in      nindex |
|----------------------------------------------------------------------------------------------------------------------------|
|     3   -.0065879   -10.34022   -1.650029   -1.533493    41.6797   19.24673    74      3   .6532016   3.407059    .5884358 |
|     1   -.0060087   -11.60251           .           .   39.44028   24.43631    74      2   .6466914    3.43889    .5661933 |
|     2           .           .    4.945804    3.630848   19.82692   26.69507    74      2   .1430221    5.35582   -1.154629 |
+----------------------------------------------------------------------------------------------------------------------------+

表中的 order 表示不同组合回归的 ID，v_1_b、v_2_b、v_3_b 表示对应解释变量的估计系数 (最后一个为截距系数)，v_1_t、v_2_t、v_3_t 表示对应系数的 t 统计量，obs 表示观察值个数，nvar 为包含截距系数在内的待估系数个数，r_sqr_a 为调整 $R^{2}$ ，rmse_in 为样本内均方根误差，nindex 表示归一化估计的线性组合指数。

3.3 具体选项

常规选项

ncomb(#1,#2)：指定要包含在过程中的协变量最小和最大数目。gsreg 将在个数范围从 #1 到 #2 的候选变量之间执行所有可能的组合，例如：

gsreg depvar ocand1 ocand2 ocand3, ncomb(1,2)

*以上命令限制解释变量个数为1到2个，即不重复地进行以下6个回归
regress depvar ocand1
regress depvar ocand2
regress depvar ocand3
regress depvar ocand1 ocand2
regress depvar ocand1 ocand3
regress depvar ocand2 ocand3

*而不进行以下同时包含三个解释变量的回归
regress depvar ocand1 ocand2 ocand3

backup(#)：在执行 gsreg 期间创建 # 个备份文件，每个备份文件将包含回归总数的大约 1/# 个回归结果；
part(#1,#2)：将所有解释变量组合的回归大约等分为 #2 组，运行其中特定的第 #1 个分组；
samesample：对同一观察样本 (定义为最大公共样本) 进行所有回归；
nocount：从屏幕上隐藏正在估计的回归数量。如果未指定此选项，gsreg 将显示正在估计的回归序号和要估计的回归总数；
vselect：运行 vselect命令 (best 选项) ，根据最佳样本内 RMSE 获得最优的 n 个解释变量。使用 vselect 时，唯一可以指定的附加选项是 fixvar 选项。

gsreg depvar ocand1 ocand2 ocand3, vselect

*将执行以下命令
vselect depvar ocand1 ocand2 ocand3, best
*并且只运行能够根据样本内RMSE、AIC和BIC获得第1、2、3个协变量系数（加上截距）的最佳模型所需要的最小数量回归。

滞后结构选项

dlags(numlist)：允许在候选协变量中包含因变量滞后项。使用此选项时必须指定 tsset。以下为使用选项 dlags 的一个例子：

gsreg depvar ocand1 ocand2, dlags(1/4)
*在这个例子中，因变量的1到4阶滞后将添加到候选解释变量
*即包含6个候选解释变量：ocand1 ocand2 L1.depvar L2.depvar L3.depvar L4.depvar

ilags(numlist)：允许包含原始自变量的滞后项。用法和选项 dlags 相同；
lags(numlist)：允许同时包含原始因变量和自变量的滞后项。用法和选项 dlags 相同。选项 lags 不得与 dlags 或 ilags 一起使用。

固定变量选项

fixvar(varlist_fix)：允许用户指定必须包含在所有回归中的协变量子集。

gsreg depvar ocand1 ocand2, ncomb(1) fixvar(fixvar1 fixvar2)
*以上命令限制了解释变量必须包含 fixvar1 和 fixvar2
*即仅执行以下两次回归
regress depvar ocand1 fixvar1 fixvar2
regress depvar ocand2 fixvar1 fixvar2

转换和交互选项

schange(varname_sc)：验证引入虚拟变量与原始解释变量 (以及截距项) 交互项的结构变化。虚拟变量 varname_sc 不能包含在原始候选解释变量中；

gsreg depvar ocand1 ocand2, ncomb(1) schange(varschange)
*schange的选项引入了虚拟变量varschange
*上式等同于执行以下6次回归并比较
regress depvar ocand1
regress depvar ocand1 c.ocand1#varschange
regress depvar ocand1 c.ocand1#varschange varschange
regress depvar ocand2
regress depvar ocand2 c.ocand2#varschange
regress depvar ocand2 c.ocand2#varschange varschange

interactions：将评估在候选解释变量中所有可能的交互项组合。当指定了 dlags、ilags 或 lags，滞后项也会同样被交互。当与 schange 选项一起使用时，只有当这些交互项包含在原始解释变量范围中时，才会使用他们与虚拟变量交互项的结构变化；
squares：将候选解释变量 (如果指定了 lags，也包括滞后项) 的平方项添加为新的候选解释变量。但只有在回归式中存在一次项时才会引入该变量的平方项；
cubic：将候选解释变量 (如果指定了 lags，也包括滞后项) 的立方项添加为新的候选解释变量。但只有在回归式中同时存在一次项和平方项时才会引入该变量的立
fixinteractions：类似选项 interactions，但他只包括原始候选解释变量 varlist_ocand (以及滞后项，如果指定 lags) 和 varlist_fix 中的每个固定变量之间的所有可能交互，需要和选项 fixvar(varlist_fix) 一起使用。

时间序列和面板数据预测选项

outsample(#)：用于时间序列和面板模型。他把样本一分为二，第一个子样本用于回归，第二个子样本用于评估预测精度。outsample(#) 留下最后的 # 个时期进行预测。指定此选项时，gsreg 将计算并存储时期 1 到 N-# 之间的样本内均方根误差，以及时期 N-# 到 N 之间的样本外均方根误差。使用此选项时必须指定 tsset。

回归命令选项

cmdest(commandname)：允许选择要使用的回归命令。如果未指定该选项，commandname 默认值为 regression。此选项允许使用 regress、xtreg、probit、logit、areg、qreg和plreg，也接受任何服从 regress 语法的回归命令；
cmdoptions(commandoptions)：允许为每个回归添加受该 commandname 支持的附加选项；
cmdstat(commandstats)：允许 gsreg 保存由特定回归命令commandname 保存为标量 e() 的其他回归统计数据信息；
cmdiveq(varlist_end = varlist_inst)：一个服务于 IV 估计的特殊选项。当估计命令为 ivregress 时，包含内生变量 varlist_end 和工具变量 varlist_inst。

估计后续选项

aicbic：在每次估计后执行 estat ic，得到 AIC (赤池信息准则) 和 BIC (贝叶斯信息准则)；
hettest：在每次回归后执行 estat hettest，并保存 p 值；
hettest_o(hettestmoptions)：允许添加在 hettest 中的选项；
archlm：在每次回归后执行 estat archlm，并保存 p 值。使用此选项时必须指定 tsset；
archlm_o(archlmoptions)：允许添加在 archlm 中的选项。

序列自相关检验

bgodfrey：在每次回归后执行estat bgodfrey，并保存 p 值。使用此选项时必须指定 tsset；
bgodfrey_o(bgodfreyoptions)：允许添加在bgodfrey中的选项；
durbinalt：在每次回归后执行estat durbinalt，并保存 p 值。使用此选项时必须指定 tsset；
durbinalt_o(durbinaltoptions)：允许添加在durbinalt中的选项；
dwatson：在每次回归后执行estat dwatson，并保存 Durbin-Watson 统计量。使用此选项时必须指定 tsset。

残差的正态性检验

sktest：在每次回归后执行 sktest，并保存正态分布的偏度和峰度联合概率的 p 值。使用此选项时必须指定 tsset；
sktest_o(sktestoptions)：允许添加在 sktest 中的选项；
swilk：在每次回归后执行 swilk，并保存 Shapiro-Wilk 正态性检验的 p 值。使用此选项时必须指定 tsset；
swilk_o(swilkoptions)：允许添加在swilk中的选项；
sfrancia：在每次回归后执行sfrancia，并保存 Shapiro-Francia 正态性检验的 p 值。使用此选项时必须指定 tsset；
testpass(#)：允许通过只保存那些完成用户指定的所有残差测试的回归 (在 # 的显着性水平下) 来减少结果数据库的大小。

结果输出选项

resultsdta(newbasename)：允许在 newbasename 中定义结果输出文件名称。默认情况下名称将是 gsreg.dta；
replace：如果 (同名) 结果文件已经在当前工作目录中创建，则替换；
double：强制结果以双精度格式创建和保存；
compact：通过删除所有系数和 t 统计量来减少结果数据库的大小；
nindex(lcimplist)：允许指定标准化精度的索引 nindex。默认情况下 nindex 将基于调整 $R^{2}$ 。lcimplist 允许创建多项排序标准 (例如 gr_sqr_a、rmse-in、rmse-out、AIC、BIC 等)，但必须有用户定义的实数权重，如 nindex(0.3 r_sqr_a -0.3 aic -0.4 bic)；
mindex(lcimplist)：必须和选项 best(#) 同时使用。与选项 nindex 语法相同，但是其参数标准化是通过使用最优的 #+1 个回归的平均值获得的。因此，mindex 随着每一个新回归的加入而更新，只有最好的 # 个回归结果被保存。同时使用 mindex 和 best 选项可以大幅减少结果数据库的大小，从而使更大规模的模型选择问题变得可行。

4. Stata 实操

gsreg 具有很多用途，本节将简单介绍三个不同的应用。在此不会详细探讨选项的使用，感兴趣的朋友可以 help gsreg 找到深入的解释和选项使用示例。

4.1 示例 A：模型选择与残差检验

在第一个例子中，我们将了解 gsreg 如何在样本内拟合优度方面获得最佳模型，前提是回归的残差满足一些理想的性质。

假设我们希望通过使用两个协变量 $x$ 和 $z$ 的一些组合，获得最好的模型来解释 $y$ 。假设以下数据生成过程 (DGP)：

y_{t} = β_{0} + β_{1 t} x_{1 t} + β_{2 t} z_{1 t} + u_{t}

t = 1, \dots, 1000; β_{0} = 1

β_{1 t} = {\begin{cases} 0.9 & t < 800 \\ 0 & t \geq 800 \end{cases}, β_{2 t} = {\begin{cases} 0 & t < 600 \\ 0.6 & t \geq 600 \end{cases}

z \sim U (0, 1); u \sim N (0, 1)

? ? ? ? \sim {\begin{cases} ? ? ? ? (0, 2) & t < 600, \\ ? ? ? ? (0, 2.4) & t \geq 600 \end{cases}

以上数据生成过程 Stata 命令如下：

clear
set obs 1000
set seed 500
gen t=_n
tsset t
gen z=runiform()
gen u=rnormal(0,1)
gen beta_0=1
gen x=runiform(0,2) in 1/599
replace x=runiform(0,2.4) in 600/1000
gen beta_1=0.9 in 1/799
replace beta_1=0 in 800/1000
gen beta_2=0 in 1/599
replace beta_2=0.6 in 600/1000
gen y=beta_0+beta_1*x+beta_2*z+u

通过构造，不难发现变量 x 具有比 z 更高的解释力，但会产生异方差问题。接下来，我们将使用 gsreg 来估计所有可能的组合。因为有两个候选协变量，所以只有 3 个可能的模型。对于每个回归，将生成和保存关于：

AIC 和 BIC 的信息 (通过 aicbic 选项)；
标准异方差检验的 p 值 (通过 hettest 选项)；
在屏幕上显示基于调整后 $R^{2}$ 、AIC 和 BIC 的多项式归一化 nindex 的最佳回归。

gsreg y x z, resultsdta(res1) replace       ///
      hettest aicbic                        ///
      nindex(0.4 r_sqr_a -0.3 aic -0.3 bic)

----------------------------------------------------
Total Number of Estimations: 3
----------------------------------------------------
Computing combinations...
Preparing regression list...
Doing regressions...
Estimation number 1 of 3
Estimation number 2 of 3
Estimation number 3 of 3
Saving results...
file res1.dta saved
----------------------------------------------------
Best estimation in terms of 0.4 r_sqr_a -0.3 aic -0.3 bic 
Estimation number 3
----------------------------------------------------

  Source |       SS           df       MS      Number of obs   =     1,000
---------+----------------------------------   F(2, 997)       =     85.26
   Model |  194.460346         2  97.2301729   Prob > F        =    0.0000
Residual |  1137.03902       997   1.1404604   R-squared       =    0.1460
---------+----------------------------------   Adj R-squared   =    0.1443
   Total |  1331.49937       999   1.3328322   Root MSE        =    1.0679
--------------------------------------------------------------------------
       y |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
---------+----------------------------------------------------------------
       x |      0.688      0.053    12.92   0.000        0.584       0.793
       z |      0.242      0.119     2.03   0.043        0.008       0.476
   _cons |      1.014      0.089    11.40   0.000        0.840       1.189
--------------------------------------------------------------------------

由上可知，最佳模型同时包括变量 x 和 z。然而，res1.dta 文件显示了一些有趣的结果：

use res1.dta, clear
sort order 
gen model = "x"
replace mode = "z" if order == 2
replace model = "x 与 z" if order == 3
order order model
list order model r_sqr_a aic bic hettest nindex, noobs

+------------------------------------------------------------------------+
| order    model    r_sqr_a        aic        bic    hettest      nindex |
|------------------------------------------------------------------------|
|     1        x   .1416655    2974.42   2984.235   .0000115    .5719652 |
|     2        z   .0020743   3125.105    3134.92   .0606163   -1.154568 |
|     3   x 与 z   .1443331   2972.305   2987.028   .0000316    .5826026 |
+------------------------------------------------------------------------+

可以看出，只有 x 作为解释变量，BIC 是最好的，而在调整 $R^{2}$ 和 AIC 方面，最好的模型都是同时使用 x 和 z 作为解释变量。然而，这两种模型都未能满足残差同方差要求，hettest 的 p 值均低于0.01。

另一方面，z 模型 (模型 2) 虽然在任何选择准则下都是最差的，但他是唯一不能拒绝同方差原假设的模型。

当用户关注估计稳健性时，残差的要求变得至关重要，gsreg 提供了一个比其他模型选择命令更好的方案，能够实现在含白噪声残差中的模型中找到最优模型的目标。

4.2 示例 B：样本外预测

如果样本外精度是模型选择的主要关注点，那么完整的搜索方法 (如 gsreg) 是必不可少的。gsreg 能够同时保证样本内以及样本外模型选择的最优性，减少用户对多元关系中结构断裂的担忧。

为了说明这一点，假设我们希望得到基于两个协变量 $x$ 和 $z$ 的最佳模型 (根据一些样本外标准) 来解释 $y$ 。假设有 100 个时间序列观测值 (使用最后 20 个进行样本外模型评估)，并假设以下数据生成过程 (DGP)：

y_{t} = β_{0} + β_{1 t} x_{1 t} + β_{2 t} z_{1 t} + u_{t}

t = 1, \dots, 100; β_{0} = 1; β_{1} = 1

β_{2 t} = {\begin{cases} 1 & t \leq 70 \\ 0 & t > 70 \end{cases}

x, z, u \sim N (0, 1)

以上数据生成过程 Stata 命令如下：

clear
set obs 100
set seed 500
gen t=_n
tsset t
gen u=rnormal()
gen x=rnormal()
gen z=rnormal()
gen beta_0=1
gen beta_1=1
gen beta_2=1 in 1/70
replace beta_2=0 in 70/100
gen y=beta_0+beta_1*x+beta_2*z+u

通过构造，这两个协变量都具有较高的样本内解释能力，但对于样本外评估目标来说， $z$ 变得不显著。如果结构变化是未知的，并且我们不使用 gsreg 来评估预测精度，那么最好的模型表示将明显包括 $x$ 和 $z$ 作为解释变量。

相反，关注潜在结构中断风险影响的用户将利用一些数据库子样本来检查参数稳定性 (例如最后 20 次观测值)，并使用gsreg检查每个备选模型的解释力和预测精度。对于本例，最简单的命令是：

gsreg y x z, outsample(20) replace

默认情况下，gsreg 结果保存为 gsreg.dta，并根据调整 $R^{2}$ 在屏幕上显示的最佳模型，为同时包含 $x$ 和 $z$ 作为解释变量的模型。结果如下所示：

----------------------------------------------------
Total Number of Estimations: 3
----------------------------------------------------
Computing combinations...
Preparing regression list...
Doing regressions...
Estimation number 1 of 3
Estimation number 2 of 3
Estimation number 3 of 3
Saving results...
file gsreg.dta saved
----------------------------------------------------
Best estimation in terms of r_sqr_a 
Estimation number 3
----------------------------------------------------

      Source |       SS           df       MS      Number of obs   =        80
-------------+----------------------------------   F(2, 77)        =     71.02
       Model |  164.161773         2  82.0808864   Prob > F        =    0.0000
    Residual |  88.9914596        77  1.15573324   R-squared       =    0.6485
-------------+----------------------------------   Adj R-squared   =    0.6393
       Total |  253.153232        79   3.2044713   Root MSE        =    1.0751

------------------------------------------------------------------------------
           y |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
           x |      1.050      0.112     9.35   0.000        0.827       1.274
           z |      0.842      0.128     6.58   0.000        0.587       1.097
       _cons |      0.837      0.124     6.75   0.000        0.590       1.084
------------------------------------------------------------------------------

默认情况下，他还将计算样本内和样本外的均方根误差。

use gsreg.dta, clear
sort order 
gen model = "x"
replace mode = "z" if order == 2
replace model = "x 与 z" if order == 3
order order model
list order model r_sqr_a-nindex, noobs

+-------------------------------------------------------------+
| order    model    r_sqr_a    rmse_in   rmse_out      nindex |
|-------------------------------------------------------------|
|     1        x   .4437445   1.335105   .5969427     .013723 |
|     2        z   .2399295   1.560649   .8685482   -1.006791 |
|     3   x 与 z   .6393373    1.07505   .6877004    .9930679 |
+-------------------------------------------------------------+

可以看出，样本内标准 (调整 $R^{2}$ 和样本内均方根误差) 的最佳模型 (模型 3，以 x 和 z 为解释变量) 在样本外均方根误差准则下并不是最优的。相反，模型 1 (仅包括 x 作为解释变量) 的样本内性能相对较差，但样本外均方根误差却最小。通过分别将 rmse_out 或 rmse_in 作为排序依据，用户能够详尽地交叉比较模型的最优性。

4.3 示例 C：参数稳定性分析

本例为参数稳定性分析 (跨越不同的控制变量模型结构)。通过生成一个包含所有回归备选方案详尽信息的数据库，GSERG 是参数稳定性分析的唯一工具。

在该例中，将根据 Gluzmann 和 Guzman (2011) 的 crisis_fr 数据集 (包含关于 1973-2005 年期间 89 个国家的金融危机、金融改革和一系列控制变量的信息) 来评估多个控制变量结构下参数的稳定性。

作者对未来 5 年金融危机的概率进行了一个混合数据线性回归 (针对拉丁美洲国家、亚洲新兴国家和转型经济国家) ，并对金融改革指数 ifr 及其最近的变化 d_ifr 进行了分析。

*数据下载地址：https://gitee.com/arlionn/data/blob/master/data01/crisis_fr.dta

use crisis_fr.dta, clear
reg fc5 ifr d_ifr if EA_LA_TR==1

  Source |       SS           df       M      Number of obs   =        928
---------+----------------------------------  F(2, 925)       =      13.97
   Model |   5.08852674        2  2.54426337  Prob> F         =     0.0000
Residual |   168.410396      925  .182065293  R-squared       =     0.0293
---------+----------------------------------  Adj R-squared   =     0.0272
   Total |   173.498922      927  .187161729  Root MSE        =     .42669
--------------------------------------------------------------------------
     fc5 |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
---------+----------------------------------------------------------------
     ifr |   -.00711    .0025616    -2.78   0.006    -.0121372   -.0020829
   d_ifr |  .0549173    .0110602     4.97   0.000     .0332113    .0766233
   _cons |  .2793652    .0255525    10.93   0.000     .2292176    .3295128
--------------------------------------------------------------------------

可以看出，fc5 和 ifr 之间显著负相关，以及 d_ifr 的显著正相关。在 Gluzmann 和 Guzman（2011）的文章中，还另外确定了23个理论上相关的控制变量。与以前的例子不同，我们不会在这里根据样本内或样本外信息准则使用 gsreg 来获得最佳模型，但会使用整套结果来评估 ifr 和 d_ifr 回归系数和 t 统计量分布。

为了达到这一目的，我们使用 crisis_fr 数据集中的现有信息，以 ifr 和 d_ifr 作为固定变量，运行所有可能的回归，并限制 GSREG 为每种备选方案使用三个控制变量，命令如下所示：

gsreg fc5 v1-v23 if EA_LA_TR ==1, ncomb(3)  ///
      fixvar(ifr d_ifr) replace nocount

*以下为显示器上返回的回归总数
Total Number of Estimations: 1771
file gsreg.dta saved

命令的主要结果可以很容易地通过下面的核密度图来描述：

从上图可以看出，ifr 系数分布很大一部分集中在 0 左右，而 d_ifr 系数则几乎完全分布在正值上。

5. 参考文献

Gluzmann P, Panigo D. Global search regression: A new automatic model-selection technique for cross-section, time-series, and panel-data regressions[J]. The Stata Journal, 2015, 15(2): 325-349. -PDF-
Lindsey C, Sheather S. Variable selection in linear regression[J]. The Stata Journal, 2010, 10(4): 650-669. -PDF-
Furnival G M, Wilson R W. Regressions by leaps and bounds[J]. Technometrics, 2000, 42(1): 69-79. -PDF-
Gluzmann P A, Guzman M M. Reformas financieras e inestabilidad financiera[J]. 2011. -PDF-

6. 相关推文

Note：产生如下推文列表的命令为：
lianxh AIC 筛选挑选, m
安装最新版 lianxh 命令：
ssc install lianxh, replace

专题	嘉宾	直播/回看视频
⭐ 最新专题		因果推断, 空间计量，寒暑假班等
⭕ 数据清洗系列	游万海	直播, 88 元，已上线
研究设计	连玉君	我的特斯拉-实证研究设计，-幻灯片-
面板模型	连玉君	动态面板模型，-幻灯片-
面板模型	连玉君	直击面板数据模型 [免费公开课，2小时]