其他
不可不知的门槛回归模型代码合集
计量经济学服务中心专辑汇总！计量百科·资源·干货：
Stata |Python |Matlab |Eviews |R Geoda |ArcGis |GeodaSpace |SPSS 一文读懂 |数据资源 |回归方法 |网络爬虫 门限回归 |工具变量 |内生性 |空间计量 因果推断 |合成控制法 |倾向匹配得分 |断点回归 |双重差分 面板数据 | 动态面板数据
不可不知的门槛回归模型代码合集
进行回归分析，一般需要研究系数的估计值是否稳定。很多经济变量都存在结构突变问题，使用普通回归的做法就是确定结构突变点，进行分段回归。这就像我们高中学习的分段函数。但是对于大样本、面板数据如何寻找结构突变点。所以本文在此讲解面板门限回归的问题，门限回归也适用于时间序列（文章后面将介绍stata15.0新命令进行时间序列的门限回归）。
门限效应，是指当一个经济参数达到特定的数值后，引起另外一个经济参数发生突然转向其它发展形式的现象（结构突变）。作为原因现象的临界值称为门限值。例如，成果和时间存在非线性关系，但是在每个阶段是线性关系。有些人将这样的模型称为门槛模型，或者门限模型。如果模型的研究对象包含多个个体多个年度，那么就是门限面板模型。
一、history
常见模型如下：门槛回归模型（threshold regression，也称门限回归）：
汉森（Bruce E. Hansen）在门限回归模型上做出了很多贡献。Hansen于1996年在《Econometrica》上发表文章《Inference when a nuisance parameter is not identified under the null hypothesis》，提出了时间序列门限自回归模型（TAR）的估计和检验。之后，他在门限模型上连续追踪，发表了几篇经典文章，尤其是1999年的《Threshold effects in non-dynamic panels: Estimation, testing and inference》（Hansen (1999) 首次介绍了具有个体效应的面板门限模型的计量分析方法, 该方法以残差平方和最小化为条件确定门限值, 并检验门限值的显著性, 克服了主观设定结构突变点的偏误。具体思路是:选定某一变量作为门限变量, 根据搜寻到的门限值将回归模型区分为多个区间, 每个区间的回归方程表达不同, 根据门限划分的区间将其他样本值进行归类, 回归后比较不同区间系数的变化。），2000年的《Sample splitting and threshold estimation》和2004年与他人合作的《Instrumental Variable Estimation of a Threshold Model》。
在这些文章中，Hansen介绍了包含个体固定效应的静态平衡面板数据门限回归模型，阐述了计量分析方法。方法方面，首先要通过减去时间均值方程，消除个体固定效应，然后再利用OLS（最小二乘法）进行系数估计。如果样本数量有限，那么可以使用自举法（Bootstrap）重复抽取样本，提高门限效应的显著性检验效率。在Hansen（1999）的模型中，解释变量中不能包含内生解释变量，无法扩展应用领域。Caner和Hansen在2004年解决了这个问题。他们研究了带有内生变量和一个外生门限变量的面板门限模型。与静态面板数据门限回归模型有所不同，在含有内生解释变量的面板数据门限回归模型中，需要利用简化型对内生变量进行一定的处理，然后用2SLS（两阶段最小二乘法）或者GMM（广义矩估计）对参数进行估计。
二．显著性检验
门槛回归模型显著性检验的目的是，检验以门檻值划分的两组样本其模型估计参数是否显著不同。
因此，不存在门槛值的零假设为：Ho：两个系数相同。同时构造LM统计量：
其中，So是在零假设下的残差平方和。由于LM统计量并不服从标准的分布。因此， Hansen（2000）提出了通过“自举法”（ Bootstrap）来获得渐进分布的想法，进而得出相应的概率p值，也称为 Bootstrap P值。
这种方法的基本思想是：在解释变量和门槛值给定的前提下，模拟（ Simulate）产生一组因变量序列，并使其满足N（0，e2），其中e是式（4）的残差项。每得到一个自抽样样本，就可以计算出一个模拟的エM统计量。将这一过程重复1000次。Hansen（1996）认为模拟产生的LM统计量大于式（6）的次数占总模拟次数的百分比就是“自举法”估计得到的P值。这里的Bootstrap P值类似于普通计量方法得出的相伴概率P值。例如，当 Bootstrap P值小于0.01时，表示在1 %的显著性水平下通过了LM检验，以此类推。
三．置信区间
以上的检验过程为只有一个门槛值的检验过程，为了能确定是否存在两个门槛值或者是更多的门槛值，我们应当检验是否存在两个门槛值，拒绝意味着至少存在一个门槛值。我们可以假设己经估计的第一个门槛值，然后开始寻找第二个门槛值。在确定有两个门槛值后，再寻找第三个门槛值，方法都和前面的一样，直至我们不能拒绝零假设。
1、门槛回归xthreg命令
// 双重门槛回归模型
use hansen1999, clear
tabstat i q1 c1 d1, stat(min p25 p50 p75 max) col(stat)
* 单一门槛回归模型single threshold model
xthreg i q1 q2 q3 d1 qd1, rx(c1) qx(d1) thnum(1) grid(400) trim(0.01) bs(300)
* 三重门槛回归模型
xthreg i q1 q2 q3 d1 qd1, rx(c1) qx(d1) thnum(3) grid(400) ///
trim(0.01 0.01 0.05) bs(0 300 300) thgiven nobslog noreg
* 绘制LR图
_matplot e(LR21), columns(1 2) yline(7.35, lpattern(dash)) ///
connect(direct) msize(small) mlabp(0) mlabs(zero) ///
ytitle("LR Statistics") xtitle("First Threshold") recast(line) name(LR21) nodraw
_matplot e(LR22), columns(1 2) yline(7.35, lpattern(dash)) ///
connect(direct) msize(small) mlabp(0) mlabs(zero) ///
ytitle("LR Statistics") xtitle("Second Threshold") recast(line) name(LR22) nodraw
graph combine LR21 LR22, cols(1)
// estimate triple thrshold model directly
xthreg i q1 q2 q3 qd1 d1, rx(c1) qx(d1) thnum(3) grid(400) ///
trim(0.01 0.01 0.05) bs(300 300 300)
结果为：
. *========================================
. * 高级计量经济学
. *========================================
.
.
. * 计量经济学服务中心
.
.
. *-------------------------------------------------------------------------------
. * 参考资料：
. * 《初级计量经济学及Stata应用：Stata从入门到进阶》
. * 《高级计量经济学及Stata应用：Stata回归分析与应用》
. * 《量化社会科学方法》
. * 《社会科学因果推断》
. * 《面板数据计量分析方法》
. * 《时间序列计量分析方法》
. * 《高级计量经济学及Eviews应用》
. * 《R、Python、Mtalab初高级教程》
. * 《空间计量入门：空间计量在Geoda、GeodaSpace中的应用》
. * 《零基础|轻松搞定空间计量：空间计量及GeoDa、Stata应用》
. * 《空间计量第二部：空间计量及Matlab应用课程》
. * 《空间计量第三部：空间计量及Stata应用课程》
. * 《空间计量第四部：《空间计量及ArcGis应用课程》
. * 《空间计量第五部：空间计量经济学》
. * 《空间计量第六部：《空间计量及Python应用》
. * 《空间计量第七部：《空间计量及R应用》
. * 《空间计量第八部：《高级空间计量经济学》
. *-------------------------------------------------------------------------------
. // 双重门槛回归模型
.
. use hansen1999, clear
(The Value and Performance of U.S.Corporations (B.H.Hall & R.E.Hall, 1993))
. tabstat i q1 c1 d1, stat(min p25 p50 p75 max) col(stat)
Variable | Min p25 p50 p75 Max
-------------+--------------------------------------------------
i | .00051 .04877 .075855 .11295 1.65509
q1 | .02119 .37089 .674795 1.30773 111.796
c1 | -.93734 .12401 .21585 .32084 8.70715
d1 | 0 .08906 .205505 .3198 4.67335
----------------------------------------------------------------
.
. * 单一门槛回归模型single threshold model
. xthreg i q1 q2 q3 d1 qd1, rx(c1) qx(d1) thnum(1) grid(400) trim(0.01) bs(300)
Estimating the threshold parameters: 1st ...... Done
Boostrap for single threshold
.................................................. + 50
.................................................. + 100
.................................................. + 150
.................................................. + 200
.................................................. + 250
.................................................. + 300
Threshold estimator (level = 95):
-----------------------------------------------------
model | Threshold Lower Upper
-----------+-----------------------------------------
Th-1 | 0.0154 0.0141 0.0167
-----------------------------------------------------
Threshold effect test (bootstrap = 300):
-------------------------------------------------------------------------------
Threshold | RSS MSE Fstat Prob Crit10 Crit5 Crit1
-----------+-------------------------------------------------------------------
Single | 17.7818 0.0023 35.20 0.0000 13.0421 15.6729 21.1520
-------------------------------------------------------------------------------
Fixed-effects (within) regression Number of obs = 7910
Group variable: id Number of groups = 565
R-sq: Within = 0.0951 Obs per group: min = 14
Between = 0.0692 avg = 14.0
Overall = 0.0660 max = 14
F(7,7338) = 110.21
corr(u_i, Xb) = -0.3972 Prob > F = 0.0000
------------------------------------------------------------------------------
i | Coefficient Std. err. t P>|t| [95% conf. interval]
-------------+----------------------------------------------------------------
q1 | .0105555 .0008917 11.84 0.000 .0088075 .0123035
q2 | -.0202872 .0025602 -7.92 0.000 -.025306 -.0152683
q3 | .0010785 .0001952 5.53 0.000 .0006959 .0014612
d1 | -.0229482 .0042381 -5.41 0.000 -.031256 -.0146403
qd1 | .0007392 .0014278 0.52 0.605 -.0020597 .0035381
|
_cat#c.c1 |
0 | .0552454 .0053343 10.36 0.000 .0447885 .0657022
1 | .0862498 .0052022 16.58 0.000 .076052 .0964476
|
_cons | .0628165 .0016957 37.05 0.000 .0594925 .0661405
-------------+----------------------------------------------------------------
sigma_u | .03980548
sigma_e | .04922656
rho | .39535508 (fraction of variance due to u_i)
------------------------------------------------------------------------------
F test that all u_i=0: F(564, 7338) = 6.90 Prob > F = 0.0000
.
. * 三重门槛回归模型
. xthreg i q1 q2 q3 d1 qd1, rx(c1) qx(d1) thnum(3) grid(400) ///
> trim(0.01 0.01 0.05) bs(0 300 300) thgiven nobslog noreg
Estimating the threshold parameters: 2nd ...... 3rd ...... Done
Boostrapping for threshold effect test: 2nd ...... 3rd ...... Done
Threshold estimator (level = 95):
-----------------------------------------------------
model | Threshold Lower Upper
-----------+-----------------------------------------
Th-1 | 0.0154 0.0141 0.0167
Th-21 | 0.0154 0.0141 0.0167
Th-22 | 0.5418 0.5268 0.5473
Th-3 | 0.4778 0.4755 0.4823
-----------------------------------------------------
Threshold effect test (bootstrap = 0 300 300):
-------------------------------------------------------------------------------
Threshold | RSS MSE Fstat Prob Crit10 Crit5 Crit1
-----------+-------------------------------------------------------------------
Single | 17.7818 0.0023 35.20 0.0000 13.0421 15.6729 21.1520
Double | 17.7258 0.0022 24.97 0.0133 12.7394 15.1986 25.3522
Triple | 17.7119 0.0022 6.20 0.5700 15.1892 19.3399 25.4014
-------------------------------------------------------------------------------
.
. * 绘制LR图
.
. _matplot e(LR21), columns(1 2) yline(7.35, lpattern(dash)) ///
> connect(direct) msize(small) mlabp(0) mlabs(zero) ///
> ytitle("LR Statistics") xtitle("First Threshold") recast(line) name(LR21) nodraw
.
. _matplot e(LR22), columns(1 2) yline(7.35, lpattern(dash)) ///
> connect(direct) msize(small) mlabp(0) mlabs(zero) ///
> ytitle("LR Statistics") xtitle("Second Threshold") recast(line) name(LR22) nodraw
.
. graph combine LR21 LR22, cols(1)
2、门槛回归xtthres操作与应用
. use hansen1999, clear // 调入 Hansen99 数据
(The Value and Performance of U.S.Corporations (B.H.Hall & R.E.Hall, 1993))
. *-Table 1: 描述性统计分析
. tabstat i q1 c1 d1, s(min p25 p50 p75 max) ///
> format(%6.3f) c(s)
variable | min p25 p50 p75 max
-------------+--------------------------------------------------
i | 0.001 0.049 0.076 0.113 1.655
q1 | 0.021 0.371 0.675 1.308 111.796
c1 | -0.937 0.124 0.216 0.321 8.707
d1 | 0.000 0.089 0.206 0.320 4.673
----------------------------------------------------------------
.
. *-门槛回归估计
.
end of do-file
. do "C:\Users\admin\AppData\Local\Temp\STD67b8_000000.tmp"
. xtthres i q1 q2 q3 d1 qd1, th(d1) d(c1) ///
> min(120) bs1(3) bs2(3) bs3(2)
Begin Time:
+------------------------------------------+
| ---单一门槛面板模型--- |
+------------------------------------------+
第一个门槛估计值(gamma1)： 0.016
Fixed-effects (within) regression Number of obs = 7910
Group variable: id Number of groups = 565
R-sq: within = 0.0951 Obs per group: min = 14
between = 0.0692 avg = 14.0
overall = 0.0660 max = 14
F(7,7338) = 110.22
corr(u_i, Xb) = -0.3971 Prob > F = 0.0000
------------------------------------------------------------------------------
i | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
q1 | .0105533 .0008917 11.84 0.000 .0088053 .0123012
q2 | -.0202819 .0025603 -7.92 0.000 -.0253008 -.015263
q3 | .0010782 .0001952 5.52 0.000 .0006955 .0014609
d1 | -.0229513 .0042381 -5.42 0.000 -.0312592 -.0146435
qd1 | .0007396 .0014278 0.52 0.604 -.0020592 .0035385
c1_1 | .0552464 .0053324 10.36 0.000 .0447933 .0656995
c1_2 | .0862637 .0052019 16.58 0.000 .0760665 .0964608
_cons | .0628205 .0016955 37.05 0.000 .0594969 .0661441
-------------+----------------------------------------------------------------
sigma_u | .0398016
sigma_e | .0492263
rho | .39531091 (fraction of variance due to u_i)
------------------------------------------------------------------------------
F test that all u_i=0: F(564, 7338) = 6.89 Prob > F = 0.0000
Note: c1_1: c1*I(d1<.016)
c1_2: c1*I(d1>=.016)
STATA 自抽样中，请等待 ... ...
+------------------------------------------+
| ---双重门槛面板模型--- |
+------------------------------------------+
---搜索第二个门槛值---
第二个门槛估计值(gamma2)： 0.777
---重新搜索第一个门槛值---
更新后的第一个门槛估计值(gamma1)： 0.013
Fixed-effects (within) regression Number of obs = 7910
Group variable: id Number of groups = 338
R-sq: within = 0.0689 Obs per group: min = 14
between = 0.2050 avg = 23.4
overall = 0.1003 max = 70
F(8,7564) = 69.91
corr(u_i, Xb) = 0.0551 Prob > F = 0.0000
------------------------------------------------------------------------------
i | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
q1 | .0086772 .0008421 10.30 0.000 .0070264 .010328
q2 | -.0204894 .002588 -7.92 0.000 -.0255626 -.0154161
q3 | .001235 .0002005 6.16 0.000 .000842 .0016279
d1 | .0054882 .0041086 1.34 0.182 -.0025657 .0135422
qd1 | -.0016447 .0015209 -1.08 0.280 -.004626 .0013367
c1_1 | .0092101 .0037334 2.47 0.014 .0018916 .0165286
c1_2 | .0698294 .0046948 14.87 0.000 .0606262 .0790326
c1_3 | .0014582 .0170336 0.09 0.932 -.0319323 .0348488
_cons | .0646393 .0015652 41.30 0.000 .061571 .0677075
-------------+----------------------------------------------------------------
sigma_u | .03031803
sigma_e | .05380741
rho | .2409758 (fraction of variance due to u_i)
------------------------------------------------------------------------------
F test that all u_i=0: F(337, 7564) = 5.37 Prob > F = 0.0000
Note: c1_1: c1*I(d1<.013)
c1_2: c1*I(.013<=d1<.777)
c1_3: c1*I(d1>=.777)
STATA 自抽样中，请等待 ... ...
+------------------------------------------+
| ---三重门槛面板模型--- |
+------------------------------------------+
第三个门槛估计值(gamma3)： 0.539
Fixed-effects (within) regression Number of obs = 7910
Group variable: id Number of groups = 365
R-sq: within = 0.0978 Obs per group: min = 14
between = 0.0570 avg = 21.7
overall = 0.0759 max = 84
F(8,7537) = 102.17
corr(u_i, Xb) = -0.2914 Prob > F = 0.0000
------------------------------------------------------------------------------
i | Coef. Std. Err. t P>|t| [95% Conf. Interval]
-------------+----------------------------------------------------------------
q1 | .0077292 .0008106 9.54 0.000 .0061402 .0093183
q2 | -.0165064 .0025753 -6.41 0.000 -.0215547 -.011458
q3 | .0008924 .0002021 4.42 0.000 .0004962 .0012885
d1 | .0072797 .0041604 1.75 0.080 -.0008758 .0154352
qd1 | -.0015131 .0014997 -1.01 0.313 -.0044529 .0014267
c1_1 | .0530393 .0048158 11.01 0.000 .0435991 .0624796
c1_2 | .0883785 .0047951 18.43 0.000 .0789789 .0977782
c1_3 | 0 (omitted)
c1_4 | .0112092 .0169942 0.66 0.510 -.0221042 .0445226
_cons | .0590864 .0016366 36.10 0.000 .0558781 .0622946
-------------+----------------------------------------------------------------
sigma_u | .0348355
sigma_e | .0532425
rho | .29976019 (fraction of variance due to u_i)
------------------------------------------------------------------------------
F test that all u_i=0: F(364, 7537) = 5.60 Prob > F = 0.0000
Note: c1_1: c1*I(d1<.013)
c1_2: c1*I(.013<=d1<.777)
c1_3: c1*I(.777<=d1<.777)
c1_4: c1*I(d1>=.777)
STATA 自抽样中，请等待 ... ...
+------------------------------------------------+
| ---门槛估计值和置信区间--- |
+------------------------------------------------+
-------------------------------------------------------------------------------
门槛估计值 95% 置信区间
-----------------------------------------------------------------------------
单一门槛模型(g1) 0.016 [ 0.013 , 0.017 ]
-----------------------------------------------------------------------------
双重门槛模型：
Ito1(g1) 0.777 [ 0.002 , 0.777 ]
Ito2(g2) 0.013 [ 0.013 , 0.013 ]
-----------------------------------------------------------------------------
三重门槛模型(g3)： 0.539 [ 0.280 , 0.572 ]
-------------------------------------------------------------------------------
Note: g# denotes gamma#, the estimated threshold values, #=1,2,3
+---------------------------------------------+
| ---门槛效果自抽样检验--- |
+---------------------------------------------+
-------------------------------------------------------------------------------
临界值
------------------------------------------------------------------
模型 F值 P值 BS次数 1% 5% 10%
-------------------------------------------------------------------------------
单一门槛 35.342*** 0.000 3 2.305 2.305 2.305
双重门槛 24.634*** 0.000 3 1.568 1.568 1.568
三重门槛 0.000*** 0.000 2 0.000 0.000 0.000
-------------------------------------------------------------------------------
+------------------------------------------------+
| ---门槛模型系数估计结果--- |
+------------------------------------------------+
--------------------------------------------------------------------
(1) (2) (3)
Single Double Triple
--------------------------------------------------------------------
L.q 0.00912*** 0.00868*** 0.00773***
(10.82) (10.30) (9.54)
L.q^2/100 -0.0213*** -0.0205*** -0.0165***
(-8.26) (-7.92) (-6.41)
L.q^3/1000 0.00127*** 0.00123*** 0.000892***
(6.37) (6.16) (4.42)
L.d 0.000862 0.00549 0.00728*
(0.22) (1.34) (1.75)
L.q*L.d -0.00227 -0.00164 -0.00151
(-1.49) (-1.08) (-1.01)
c1_1 0.00831** 0.00921** 0.0530***
(2.25) (2.47) (11.01)
c1_2 0.0655*** 0.0698*** 0.0884***
(14.14) (14.87) (18.43)
c1_3 0.00146 0
(0.09) (.)
c1_4 0.0112
(0.66)
Constant 0.0663*** 0.0646*** 0.0591***
(43.78) (41.30) (36.10)
--------------------------------------------------------------------
r2_w 0.0660 0.0689 0.0978
r2_b 0.210 0.205 0.0570
r2_o 0.0982 0.100 0.0759
N 7910 7910 7910
--------------------------------------------------------------------
t statistics in parentheses
* p<0.1, ** p<0.05, *** p<0.01
你可以输入或点击如下命令查看各个门槛的置信区间图：
Single Threshold: xttr_graph
Double Threshold(1st round): xttr_graph, m(22)
Double Threshold(2ed round): xttr_graph, m(21)
Triple Threshold: xttr_graph, m(3)
For details, see: xttr_graph
Over Time: Time used: 140.398s
.
end of do-file