面板数据的自适应Lasso分位回归方法研究

2014-05-12 10:23李子强田茂再罗幼喜
统计与信息论坛 2014年7期
关键词:分位先验贝叶斯

李子强,田茂再,罗幼喜

(1.湖北工业大学 理学院,湖北 武汉 430068;2.中国人民大学a.应用统计研究中心;b.统计学院,北京 100872)

面板数据的自适应Lasso分位回归方法研究

李子强1,田茂再2a,2b,罗幼喜1

(1.湖北工业大学 理学院,湖北 武汉 430068;2.中国人民大学a.应用统计研究中心;b.统计学院,北京 100872)

如何在对参数进行估计的同时自动选择重要解释变量,一直是面板数据分位回归模型中讨论的热点问题之一。通过构造一种含多重随机效应的贝叶斯分层分位回归模型,在假定固定效应系数先验服从一种新的条件Laplace分布的基础上,给出了模型参数估计的Gibbs抽样算法。考虑到不同重要程度的解释变量权重系数压缩程度应该不同,所构造的先验信息具有自适应性的特点,能够准确地对模型中重要解释变量进行自动选取,且设计的切片Gibbs抽样算法能够快速有效地解决模型中各个参数的后验均值估计问题。模拟结果显示,新方法在参数估计精确度和变量选择准确度上均优于现有文献的常用方法。通过对中国各地区多个宏观经济指标的面板数据进行建模分析,演示了新方法估计参数与挑选变量的能力。

面板数据;自适应Lasso;分位回归;切片Gibbs抽样

一、引 言

面板数据模型是当前学术界讨论最多的模型之一。传统的面板数据模型实际上是一种条件均值模型,即讨论在给定解释变量的条件下响应变量均值变化规律。这种模型的一个固有缺陷是只描述了响应变量的均值信息,其它信息则都忽略了。然而,数据的信息应该是全方位的,这种只对均值建模的方法有待改进。Koenker等提出的分位回归模型是对均值回归模型的一种有效改进,该模型可以在给定解释变量后对响应变量的任意分位点处进行建模,从而可以从多个层次刻画数据的分布信息[1]。同时,分位回归的参数估计是通过极小化加权残差绝对值之和得到,比传统均值回归模型下二次损失函数获得的最小二乘估计更为稳健[2]。

对于简单的线性模型,与分位回归方法相对应的参数点估计、区间估计、模型检验及预测已经有很多成熟的研究结果,但有关面板数据模型的分位回归方法研究文献还不多见。Koenker对固定效应的面板数据模型采用带Lasso惩罚的分位回归方法,通过对个体固定效应实施L1范数惩罚,该方法能够在各种偏态及厚尾分布下得到明显优于均值回归的估计,然而惩罚参数如何确定是该方法的一个难点[3];罗幼喜等也提出了3种新的固定效应面板数据分位回归方法,模拟显示,这些新方法在误差非正态分布情况下所得估计优于传统的最小二乘估计和极大似然估计,但新方法对解释变量在时间上进行了差分运算,当解释变量中包含有不随时间变化的协变量时,这些方法则无法使用[4];Tian等对含随机效应的面板数据模型提出了一种分层分位回归法,并利用EQ算法给出模型未知参数的估计,但该算法只针对误差呈正态分布而设计,限制了其应用范围[5]。以上文献均是直接从损失函数的角度考虑分位回归模型的建立及求解;Liu等利用非对称拉普拉斯分布与分位回归检验损失函数之间的关系,从分布的角度建立了含随机效应面板数据的条件分位回归模型,通过蒙特卡罗EM算法解决似然函数高维积分问题[6];Luo等则在似然函数的基础上考虑加入参数先验信息,从贝叶斯的角度解决面板数据的分位回归问题,模拟显示,贝叶斯分位回归法能有效地处理模型中随机效应参数[7];朱慧明等也考虑过将贝叶斯分位回归法应用于自回归模型,模拟和实证显示该方法能有效地揭示滞后变量对响应变量的位置、尺度和形状的影响[8]。

然而,上述方法均不能对模型中自变量进行选择,但在实际的经济问题中,人们在建立模型之前经常会面临较多解释变量,且对哪个解释变量最终应该留在模型中没有太多信息。如果将一些不重要的噪声变量包含在模型之中,不仅会影响其它重要解释变量估计的准确性,也会使模型可解释性和预测准确性降低。Park等在研究完全贝叶斯分层模型时提出了一种新的贝叶斯Lasso方法,通过假定回归系数有条件Laplace先验信息给出了参数估计的Gibbs抽样算法,这一工作使得一些正则化的惩罚方法都能够纳入到贝叶斯的框架中来,通过特殊的先验信息对回归系数进行压缩,该方法能够在估计参数的同时对模型中自变量进行选择[9-10]。Alhamzawi等将贝叶斯 Lasso方法引入到面板数据分位回归模型中来,使得在估计分位回归系数的同时能够对模型中重要解释变量进行自动选择[11-12]。但是,上述研究中均假设回归系数先验分布所依赖的条件参数对所有解释变量都是相同的,也即对所有分量压缩程度一样,正如Zou所指出,这样得到的回归系数估计将不是无偏估计[13]。为了改进这一缺陷,本文拟构造一种自适应的贝叶斯Lasso分位回归方法,即假定回归系数的每个分量先验分布都依赖不同的条件参数,从而对不同的解释变量施加不同的惩罚权重,这不仅能够改进回归系数估计偏差,而且能够自动压缩模型中非重要解释变量回归系数为0,达到变量选择的目的。虽然面临需要估计更多参数的困境,但本文通过对Laplace分布的分解和引进辅助变量构造的切片Gibbsl抽样算法能够快速有效地解决这一问题[14]。

二、模型及方法

(一)面板数据的贝叶斯分位回归模型

定义1 考虑含多重随机效应的面板数据模型,定义给定τ时的条件分位回归函数如下:

(二)非对称Laplace分布分解与自适应先验信息的选取

显然,给定适当的先验信息后,上述模型(4)即可以通过一般的 MCMC方法进行求解。然而,考虑到非对称Laplace分布没有共轭先验,这将为MCMC算法的估计带来极大的计算负担,为此给出非对称Laplace分布的一个重要分解:

利用引理1,ALD分布可以表示为正态和指数两个常见分布的混合,这为后面建立未知参数的Gibbs抽样算法带来了极大方便。关于先验信息,选取的方法很多,其中共轭先验信息选取法由于其计算推导简洁应用最为广泛。对于随机效应通常假定αi|φ~N(0,φ),φ~IG(k0,w0);对于尺度参数σ~IG(c0,d0),其 中 IG(a,b)表 示 参 数 为 a,b 的 逆Gamma分布。

对于参数β,如果按照通常共轭先验信息的选取方法则为正态分布,但这一先验分布无法起到变量选择的作用。Alhamzawi等将Laplace先验引入到贝叶斯分位回归模型中来,使得在估计分位回归系数的同时能够对模型中重要解释变量进行自动选择,改进了正态分布先验的缺陷。需要指出的是,虽然他们提出的先验能够对解释变量系数进行压缩起到变量选择的作用,但其所依赖的条件参数λ对β的所有分量都是相同的,也即对所有分量压缩程度一样,这显然会限制了β变化的灵活性,与实际中不同的解释变量应该有不同的权重也不符。为了改进这一缺陷,本文在其基础上提出一种自适应的β先验信息分布假设:

从而可以获得与自适应Lasso惩罚相对应的贝叶斯自适应Lasso分位回归方法。另外值得一提的是,式(6)中先验条件依赖于σ也很重要,正如Park等所指出,它能够保证后验密度是单峰的。如果没有条件依赖于σ,则后验密度可能不再是单峰的,缺少单峰性质不仅会使得后面的Gibbs抽样算法收敛速度变慢,而且所得估计结果意义不大[9]。

(三)参数估计的切片Gibbs抽样算法构造

上述提出的自适应β先验虽然具有自动的对模型中解释变量进行选择且考虑到了不同的解释变量应该对应不同的权重两大优点,但同时也带来了两个难题:一是条件Laplace分布不是共轭先验,从而其后验密度的推导及计算均较为复杂;二是每个解释变量都引进了不同的压缩参数,待估计的未知参数增多,从而加大了计算量。下面将通过引进辅助变量构造一种易于实施的切片Gibbs抽样算法来解决推导及计算上带来的困难[14]。

从式(10)可以看到,通过引入辅助变量S,将(β,S)的联合先验变为了正态与指数分布的混合,假定η2l~IG(a0,b0),则构造了一种易于实施的切片Gibbs抽样算法,该算法中所有未知参数均来自常见分布,从而可利用专门的Gibbs抽样软件WinBUGS来实现后验样本的抽取,而且速度快、操作方便,很容易从抽取的后验样本获得各待估参数的点估计和区间估计以及响应变量预测值。WinBUGS还提供了各种诊断和监视抽取样本是否收敛的图形和方法。

(四)切片Gibbs抽样算法中各参数条件后验密度的推导

实事上,与李翰芳等提出的贝叶斯Lasso方法相比,此处切片Gibbs抽样算法中条件后验密度首先发生改变的是参数η2l、sl,l=1,2,…,k,模型的似然函数可重写为:

9.用每次上步生成的数据值代入下一步生成新的数据,重复2~8步直至收敛。

三、蒙特卡罗模拟

Luo等在假设固定效应系数先验为正态分布的情况下对贝叶斯分位回归估计(BQR)、混合数据普通最小二乘估计(LS)、考虑固定效应的最小二乘估计(LSFE)、限制极大似然估计(REML)、混合数据的分位回归估计(QR)、个体效应惩罚的分位回归估计(PQR)进行了比较,结果表明BQR方法均优于其它方法,特别是当误差分布为非正态分布时,分位回归估计明显优于传统的均值估计[3]。李翰芳等在BQR的基础上改进固定效应系数先验分布假设,提出了与Lasso方法等价的贝叶斯Lasso分位回归估计(BLQR),使得在建立面板数据分位回归模型的同时能够对模型中解释变量进行自动选择,通过蒙特卡罗模拟比较发现,BLQR较BQR有更强的排除无关解释变量的能力[12]。考虑到BLQR方法对所有解释变量的压缩权重均相同的缺陷,本文提出了一种贝叶斯自适应Lasso分位回归法(BALQR)进行改进,即对不同的解释变量实施不同的压缩权重。下面将通过蒙特卡罗模拟来比较BQR、BLQR以及本文提出的BALQR在固定效应系数估计、随机效应方差分量估计及排除模型中无关解释变量的能力。另外,比Luo等的研究更进一步的是,前两者均只考虑了含单个随机截距的面板数据模型,本文则将其推广至含多重随机效应的面板数据模型[12]。

下面利用如下含多重随机效应的面板数据模型生成数据:

取β0=0,N=30,T=10。对于模型中固定效应部分,首先利用标准正态分布随机生成8个解释变量,X′it= (Xit1,Xit2,…,Xit8),且任意两个解释变量Xl与Xk之间的相关系数为ρ|l-k|,设ρ=0.5。对于解释变量的系数,分别取两组值进行模拟:(1)β=(β1,β2,…,β8)′= (2,4,6,8,0,0,0,0)′; (2)β=(β1,β2,…,β8)′= (5,5,5,5,0,0,0,0)′,其中两组中均只有前4个解释变量与响应变量相关,后4个解释变量与响应变量无关,以便测试BALQR法变量选择的能力。另外,第(1)组中所有相关解释变量权重系数均不同,而第(2)组中所有相关解释变量权重系数均相同,以便比较权重系数对BLQR与BALQR估计的影响。对于模型中随机效应部分,取Z′it= (Xit1,Xit2,Xit3,Xit4)′,αi~ N4(0,φ2I4),φ=2,即每个与响应变量相关的解释变量都受到了一个个体随机效应的干扰;对于模型中特异误差εit,考虑到其来自正态分布N(0,22),每个模拟重复500次。

在模拟过程中,BQR法中取β有弱先验信息βi~ N(0,106),i=0,1,…,8,另外,3种方法中其它参数也均假定有弱先验信息:σ~IG(10-6,10-6),φ~ IG(10-6,10-6),η2~ IG(10-6,10-6),η2l~IG(10-6,10-6),l=1,2,…,k。3种方法的 Gibbs抽样均迭代10 000次,舍弃前面5 000个点,取后面的5 000个点用于计算各参数的估计值。

首先从表1中4个相关解释变量在模型中权重系数均不相同的情况来看,显然无论是在中位点处还是极端分位点处,本文提出的BALQR法MSE值均是最小的,也即其整体估计效果最优;从对取值不为0的4个权重系数估计准确性来看,BALQR和BQR法相当,且偏差及标准差都明显小于BLQR法,这是因为BLQR法和普通的Lasso方法一样,对模型中所有解释变量的压缩程度都是一样的,而本模拟中各个相关解释变量在模型中权重系数均不相同,这种采取同样策略的压缩显然会给估计带来较大的偏差。而本文提出的具有自适应性的BALQR法则克服了该缺点,对每个权重系数都引进一个单独可以灵活变动的压缩参数,从而使得估计精确度大大提高;从对4个无关解释变量的选择来看,BLQR和BALQR明显优于BQR法,两方法对无关解释变量的权重系数估计都与0非常接近,也即这两种方法都具有变量选择的功能,从而验证了Laplace先验比正态先验有更高的概率使得变量落在0周围的性质;从对随机效应方差分量φ的估计来看,3种方法估计精度几乎相当,即他们在处理随机效应干扰方面的能力一致,这也是因为理论上3种方法对于随机效应的先验假设均是一致的。另外,从中位点和极端分位点的估计对比来看,3种方法在中位点处的估计精度较极端分位点要高。

表1 3种方法在相关解释变量系数不等时的估计结果比较表

从表中4个相关解释变量在模型中权重系数都相同的情况来看,BLQR法在两个分位点处对固定效应系数总体估计精度均是最高的,BALQR法次之,BQR法最差。这一点其实不难理解,因为BLQR法将所有解释变量都同等对待,而本模拟中真实参数的设置即所有相关解释变量权重系数相同,这正好满足其理论假设,从而获得比较优良的估计。虽然本文提出的BLAQR法在这种情况下精度逊于BLQR法,但不难看到它们之间的差距并不明显,且BLAQR与BLQR法一样,在对无关“噪声”解释变量的排除能力上都优于BQR法。

表2 3种方法在相关解释变量系数相等时的估计结果比较表

总结以上模拟分析结果可知:1.3种贝叶斯分位回归法均能够有效处理面板数据模型中的多重随机效应,即能够有效解决个体样本之间相关对模型估计带来的影响这一问题,并且对解释变量间存在的相关性并不敏感;2.本文提出的BALQR法与李翰芳等提出的BLQR法一样,不仅能有效地处理模型中随机效应带来的影响,而且能够在估计固定效应系数的同时对模型中解释变量进行自动选择,即通过对非重要解释变量系数实施更大程度的压缩,使之与0较为接近,从而排除其对建模时带来的不利影响和干扰,也正是由于这一点,使得BLQR法和BALQR法能够在模型的整体估计效果上明显优于BQR法;3.相比BLQR法,本文提出的BALQR法能够对模型中不同的解释变量实施不同的惩罚权重系数,从而获得更为精确的估计,即使对于具有相同重要程度的解释变量,其也能够获得与BLQR法相当的估计结果。另外,由于本文构造的一种易于实施的切片Gibbs抽样算法可在专门的抽样软件WinBUGS中实现,所以虽然待估参数较BQR法和BLQR法多,但在计算时间消耗上并无明显差别。

四、真实数据分析

为了便于与李翰芳等研究结果比较,下面以1998—2009年中国各个省市地区的宏观经济指标面板数据为例,利用本文提出的贝叶斯自适应Lasso分位回归方法对近些年来各宏观解释指标与GDP关系进行建模分析,探讨各指标对GDP的贡献程度。根据GDP的主要内在和外在影响因素,并考虑到中国目前的经济结构,初步选取了以下8个宏观解释指标:总固定资产投资额 (Finvest)、城镇居民全年平均消费性支出(Consume)、进出口总额(Imexport)、财 政 支 出 (Finac)、外 商 直 接 投 资(FDI)、就业总人数 (Employ)、R&D 经费支出(R&D)、能源消耗 (Energy)。数据来源于《中国统计年鉴 (1999—2010)》,其中由于西藏缺少数据较多故删去。为了降低各时间序列的不稳定性并使模型各解释变量系数之间有可比性,分别对GDP和所有宏观指标取对数后并标准化,为方便起见,经变换后各指标的记号不变。

李翰芳等研究发现,如果直接将所有数据混合建立简单的线性回归模型,虽然可以得到较高的模型拟合优度,但有众多解释变量系数为负无法通过显著性检验,进一步的研究还发现,各个宏观经济变量之间的相关性比较强,存在着严重的共线性,因而通常的最小二乘估计已经失效。李翰芳等考虑了带个体随机效应的面板数据模型并利用贝叶斯Lasso分位回归法来解决解释变量之间的多重共线性问题[12]。考虑到中国各地区经济发展极不平衡,各个指标均存在较大的地区效应,建立如下更一般的多重随机效应面板数据模型:

该模型是一个既含随机截距也含随机斜率的面板数据模型,假设随机效应向量α= (αi0,αi1,…,αi8)′~ N9(0,φ2I9)。下面利用本文提出的贝叶斯自适应分位回归法对上述模型进行估计。取τ=0.1,0.5,0.9分别计算3个分位点处的估计结果,其中各个参数先验与蒙特卡罗模拟中一样,均取弱先验信息,Gibbs抽样共迭代40 000次,舍弃前面的20 000个点,取后面的20 000个点用于计算各参数的估计值,估计结果如表3。

表3 贝叶斯自适应Lasso在3个分位点处的估计结果表

表3的结果表明,各个宏观指标对GDP的贡献度都不一样,且在高、中、低分位点处也都不尽相同。首先与李翰芳等给出的中位点处后验均值估计结果对比来看,虽然各个指标系数估计差别不大,但从第二部分的模拟分析来看,当各个解释指标所占权重系数不同时,本文提出的自适应Lasso方法精度更高一些[12]。外商直接投资和R&D经费支出对GDP贡献度最小,其次是进出口总额,特别是在低分位点处,0均包含在这3个指标的95%置信区间内,即相比其它指标,这3个指标为非重要指标。

从权重系数较大的几个指标来看,在低分位点τ=0.1处,从大到小依次为就业总人数、能源消耗、财政支出、城镇居民全年平均消费性支出、总固定资产投资额;在中位点τ=0.5处,从大到小依次为就业总人数、能源消耗、城镇居民全年平均消费性支出、财政支出、总固定资产投资额;在高分位点τ=0.9处,从大到小依次为城镇居民全年平均消费性支出、就业总人数、能源消耗、财政支出、总固定资产投资额。不难看到,就业总人数和能源消耗两个指标无论在哪个分位点处权重系数都排在前3名之内,一方面说明扩大就业始终是拉动GDP增长的最关键因素,另一方面也能看出目前中国经济增长对能源消耗的依赖度还很高。从各个分位点指标重要程度排名变化情况来看,城镇居民全年平均消费性支出变化最大,从低分位点处的排名第4跃升致高分位点处的排名第1,充分说明了消费对拉动GDP高速增长的巨大潜力。

从各个指标在高中低3个分位点处的权重系数大小来看,随着分位点增加,总固定资产投资额和能源消耗系数是逐渐变小的,而其它指标均是逐渐增大的,也即要使得GDP能够长期持续增长,则需要降低其对固定资产投资额和能源消耗的依赖度,这也给当前经济结构调整及转型提供了重要的启示。

[1] Koenker R,Bassett G.Regression Quantiles[J].Econometrica,1978,46(1).

[2] 陈建宝,丁军军.分位数回归技术综述[J].统计与信息论坛,2008,23(3).

[3] Koenker R.Quantile Regression for Longitudinal Data[J].Journal of Multivariate Analysis,2004,91(1).

[4] 罗幼喜,田茂再.面板数据的分位回归方法及其模拟研究[J].统计研究,2010(10).

[5] Tian M Z,Chen G M.Hierarchical Linear Regression Models for Conditional Quantiles[J].Science in China Series A:Mathematics,2006,49(12).

[6] Liu Y,Bottai M.Mixed-Effects Models for Conditional Quantiles with Longitudinal Data[J].The International Journal of Biostatistics,2009(5).

[7] Luo Y X,Lian H,Tian M Z.Bayesian Quantile Regression for Longitudinal Data Models[J].Journal of Statistical Computation and Simulation,2012,82(11).

[8] 朱慧明,王彦红,曾惠芳.基于逆跳 MCMC的贝叶斯分位自回归模型研究[J].统计与信息论坛,2010,25(1).

[9] Park T,Casella G.The Bayesian Lasso[J].Journal of the American Statistical Association,2008,103(482).

[10]Hans C.Bayesian Lasso Regression[J].Biometrika,2009,96(4).

[11]Alhamzawi R,Yu K.Bayesian Lasso-mixed Quantile Regression[J].Journal of Statistical Computation and Simulation,2014(4).

[12]李翰芳,罗幼喜,田茂再.面板数据的贝叶斯Lasso分位回归方法[J].数量经济技术经济研究,2013(2).

[13]Zou H.The Adaptive Lasso and Its Oracle Properties[J].Journal of the American Statistical Association,2006,101(476).

[14]Neal R.Slice Sampling(with Discussion)[J].Annals of Statistics,2003,31(3).

[15]Kozumi H,Kobayashi G.Gibbs Sampling Methods for the Bayesian Quantile Regression[R].Technical Report,Kobe University,2009.

[16]Andrews D F,Mallows C L.Scale Mixtures of Normal Distributions[J].Journal of the Royal Statistical Society,Ser.B,1974,36(1).LI Zi-qiang1,TIAN Mao-zai2a,2b,LUO You-xi1

Study on Adaptive Lasso Quantile Regression for Panel Data Models

(1.School of Science,Hubei University of Technology,Wuhan 430068,China;a.Center for Applied Statistics;b.School of Statistics,2.Renmin University of China,Beijing 100872,China)

How to do parameter estimation and variable selection simultaneously is a hot issue in the study of quantile regression for panel data models.On the base of the assumption that the fixed effect coefficients are subject to a novel conditional Laplace prior,the paper constructs a hierarchical Bayesian quantile regression model and gives the Gibbs sample algorithm for the unknown parameter estimation.In consideration of different explain variables should have different shrinkage degree,the proposed prior has the property of adaptivity,which could select the important explain variables in the model automatically.Furthermore,the slice Gibbs sample algorithm that the paper proposed is able to estimate the posteriori mean estimation of unknown parameter quickly and efficiently.Monte Carlo simulation study indicates that the proposed method is obviously superior to the existing methods in literatures on the accuracy of parameter estimation and variable selection.Finally,the paper gives a research of modeling the panel data including several macroeconomic indicators of our country and demonstrates the new method's capability of estimating parameters and doing variable selection.

panel data;adaptive Lasso;quantile regression;slice Gibbs sampler

O212∶F064.1

A

1007-3116(2014)07-0003-08

2014-02-28

国家自然科学基金项目《基于当代分位回归与鞍点逼近方法的复杂数据分析》(11271368);教育部人文社会科学青年基金项目《面板数据的分位回归方法及其变量选择问题研究》(10XNL018);湖北省教育厅人文社科项目《面板数据的分位回归方法及其应用研究》(2012G078);湖北工业大学博士科研启动基金《高维复杂纵向数据的分位回归建模研究》(BSQD13050)

李子强,男,湖北武穴人,副教授,研究方向:应用数理统计;

田茂再,男,湖南凤凰人,理学博士,教授,博士生导师,研究方向:数理统计;

罗幼喜,男,湖北红安人,经济学博士,副教授,研究方向:数理统计及计量经济建模。

(责任编辑:崔国平)

猜你喜欢
分位先验贝叶斯
1980~2019年夏季青藏高原中东部极端日降水分布特征
BOP2试验设计方法的先验敏感性分析研究*
最清晰的加仓信号:估值
一类低先验信息度的先验分布选择研究
基于贝叶斯解释回应被告人讲述的故事
基于动态贝叶斯估计的疲劳驾驶识别研究
贵州近50 a来夏季不同等级强降水事件特征研究
基于自适应块组割先验的噪声图像超分辨率重建
基于“业绩与薪酬双对标”的国有企业负责人薪酬研究
基于互信息的贝叶斯网络结构学习