基于非對稱拉普拉斯分布的混合分位數回歸參數估計

張發趕,何幼樺

(上海大學理學院,上海 200444)

傳統線性回歸要求誤差項服從正態分布,然而現實生活中許多數據尤其是金融經濟數據大多數是尖峰厚尾的,此時模型的估計結果將不具有穩健性.1978 年Bassett 等[1]提出分位數回歸模型.

1986 年,Powell[2]解決了分位數回歸中存在缺失數據的問題.1998 年,Thompson 等[3]簡單介紹了貝葉斯方法.2001 年,Yu 等[4]進一步完善了貝葉斯方法,完整地將貝葉斯框架引入分位數回歸模型中,并使用非對稱拉普拉斯分布來描述誤差項.2010 年,Taddy 等[5]解決了非參數貝葉斯分位數回歸問題.

混合回歸模型是研究多個子聚類混合的統計模型.Goldfeld 等[6]首次提出混合回歸模型.2014 年,Yao 等[7]和Song 等[8]分別使用t 分布以及Laplace 分布給出新的混合回歸模型,相比于傳統模型而言,該類模型具有更好的穩健性.Park[9]在混合回歸模型中率先將均值和方差綜合起來考慮,提出了聯合均值方差模型.2016 年,Wu 等[10]提出了混合分位數回歸模型.2017 年,詹金龍等[11]利用混合Laplace 分布提出了聯合位置和尺度參數的回歸模型.

上述研究基本不涉及尺度參數的結構,詹金龍等[11]雖然考慮了尺度參數的結構,但只限于拉普拉斯分布.在實際問題中往往出現不僅類別不同,每個聚類本身也存在差異的情況,例如男女身高的區別,不同學歷收入的區別等,基于此本工作提出非對稱拉普拉斯分布下的混合分位數回歸模型,并同時考慮位置參數和尺度參數的回歸問題.

1 混合分位數回歸模型

設隨機向量y 服從非對稱拉普拉斯分布(asymmetric Laplace distribution,ALDp),并且由m 個子聚類混合而成,概率密度函數如下:

對每一個ALDp(yi;μj,),概率密度為

式中:p 為非對稱參數,0

若隨機變量x 服從ALDp(x;μ,σ2),則有P(x<μ)=p,P(x>μ)=1 ?p,即位置參數μ就是分布ALDp(x;μ,σ2)的p 分位數,所以在上述假設下估計模型的參數與求y 的p 分位數是等價的.

在很多情況下,數據本身存在異方差性,這使得傳統的混合分位數回歸模型的估計存在偏差,針對這種情況,本工作同時對位置參數和尺度參數進行回歸,提出了基于ALDp的混合分位數回歸模型,

式中:xi={xi1,xi2,···,xir}T和hi={hi1,hi2,···,hiq}T是回歸方程的解釋變量;相應的{yi}ni=1是獨立的被解釋變量;βj={βj1,βj2,···,βjr}T是第j 個子聚類中維數為r×1 的位置模型的未知參數;γj={γj1,γj2,···,γjq}T是第j 個子聚類中維數為q×1 的尺度模型的未知參數.βj,γj,πj可以與p 有關,為符號簡潔起見,βj即為βj(p),γj即為γj(p),πj即為πj(p).直接求解上述模型是比較困難的,本工作擬采用期望最大化(expectation maximization,EM)算法對上述參數進行估計.

2 基于非對稱拉普拉斯分布的EM 算法

EM 算法[12]是含有隱變量的概率模型參數的極大似然估計方法,本工作利用EM 算法對參數進行估計.

假定子聚類數m 是固定和已知的,引入隱變量zij對混合比例進行刻畫,若zij屬于m 個子聚類的第j 類,則zij等于1;若zij不屬于m 個子聚類的第j 類,而屬于其余m ?1 類中的某一類,則zij等于0.

式(2)在完全數據下關于參數Θ=(β1,γ1,π1,···,βm,γm,πm)T的對數似然函數可以寫成

EM 算法是一種迭代算法,流程主要分為兩個步驟:E 步和M 步.E 步是計算對數似然函數的期望;M 步是尋找能使E 步產生的似然期望最大化的參數值;重復執行E 步和M 步,直至參數Θ 收斂.具體操作如下.

步驟一 給定參數迭代初始值

步驟二(E 步) 利用第k 次迭代得到的結果估計Θ(k),計算,

式中:

步驟三(M 步) 用Q(Θ;y,x,h,Θ(k))對Θ 求最大值,將得到

這里采用牛頓法對上述問題進行求解,由于目標函數含有絕對值無法直接求導,因此下文推導過程中對應部分使用差商代替求導,這種做法也方便之后進行數值計算.

令θ=(β,γ),Q(θ)=Q(β,γ),

設計如下迭代過程:

為了計算方便,

步驟四 重復E 步和M 步,直至參數Θ 收斂.

3 Monte Carlo 數值模擬

利用數值模擬對本工作提出的非對稱拉普拉斯分布的EM算法進行驗證,以分析樣本量和p 對估計效果的影響.

首先,產生服從ALDp分布的n 個隨機數,若要產生來自分布F(x)的隨機數,需先產生U(0,1)的隨機數u,然后計算F?1(u).具體計算步驟如下:(1) 從U(0,1)隨機產生u;(2) 計算x=F?1(u),其中F?1(u)=inf{x:F(x)≥u}.

首先求ALDp的累計概率密度函數

式中:

計算反函數F?1(u),

式中:p 為非對稱參數;μ為位置參數,σ 為尺度參數.

通過對不同分類數的模型進行模擬,結果表明樣本量和分位數對參數估計精度的影響關系與兩分類問題相同.此處僅以m=2 為例,考慮如下混合分位數回歸模型:

根據式(11)產生模擬數據,yi服從混合ALDp分布,xi服從U(1,2),hi服從U(1,2),彼此相互獨立.考慮分位數對模型參數的影響,即β1=50p,σ1=exp(?p)/5,β2=30p +2,σ2=2p2+0.1,其中p 表示非對稱參數.給定混合比例π1=0.3,π2=0.7,取樣本量n=100,200,400,600,800,1 000,p=0.1,0.2,···,0.9,重復模擬300 次.利用均方誤差(mean square error,MSE)來觀察估計效果,結果如表1～4 所示.

表1 參數β1 在不同樣本量下的MSETable 1 MSE of β1 under different sample sizes

表2 參數σ1 在不同樣本量下的MSETable 2 MSE of σ1 under different sample sizes

通過計算參數在不同樣本量下的均方誤差,可以發現模型參數估計的整體MSE 較小,參數估計效果良好.值得注意的是,在小樣本下,參數的估計相對于大樣本而言精確度會低一些,原因之一是在小樣本下高低兩側分位數的數據不足可能會導致部分參數估計的結果存在偏差.但是隨著樣本量的增加,所有參數在各個分位數下的MSE 均明顯降低,說明隨著樣本量的逐漸增加,模型的估計效果越來越好,并且隨著n 的增大,分位數對估計結果的影響也越來越小,模型在各個分位數上的估計穩定性越來越高.顯然地,當樣本量足夠大時,高低兩側分位數也已經擁有了足夠多的數據.

表3 參數β2 在不同樣本量下的MSETable 3 MSE of β2 under different sample sizes

表4 參數σ2 在不同樣本量下的MSETable 4 MSE of σ2 under different sample sizes

4 實證分析

4.1 數 據

下面將利用模型(2)對波士頓房價數據進行分析,該數據來源于UCI(University of California Irvine)數據庫.波士頓房價數據于1978 年開始統計,包含了波士頓地區房屋的一系列特征(離散變量已由該數據庫處理完成),例如犯罪率、一氧化氮濃度、城鎮師生比例等,每個特征有506 個樣本點.首先對特征進行篩選,通過相關性分析,排除對房價沒有影響的變量,進而篩選出其中最重要的兩個變量,分別為每棟住宅房間數和該地區房東屬于低收入階層的比例,對應散點如圖1 所示.

4.2 模型建立與結果分析

由圖1 可以看出,房價的波動與每棟住宅房間數(x1)和該地區房東低收入階層比例(x2)都有關,利用該數據集構建混合分位數回歸模型(12),將數據(房屋)類型分為兩類,對這兩類數據進行研究發現房價與房間數以及房價與房東屬于低收入階層比例的依賴關系有很大不同.

圖1 各解釋變量與被解釋變量散點圖Fig.1 Scatter plot of each explanatory variable and explained variable

利用EM 算法對式(12)中的參數進行估計,得到不同分位數下的位置參數和尺度參數估計結果如表5 所示.

表5 不同分位數下模型位置參數和尺度參數以及混合比例的估計Table 5 Estimation of model position parameters,scale parameters and mixing ratio under different quantiles

進一步分析表5 的實際意義,考察不同分位數下解釋變量x1和x2對不同類型房屋的房價邊際影響如圖2 所示.

由圖2(a)可知,對于第二類房屋,β21>0,表明該類房屋的房價和房間數的多少呈正相關,并且房間數對房價的邊際影響明顯高于第一類房屋,但是這種邊際影響隨著房價的升高逐漸減小.對于第一類房屋,β11在p>0.3 時接近于0,說明此時房間數的多少對房價的影響不大.

圖2 不同分位數下x1,x2 對房價的邊際影響Fig.2 Marginal influence of x1 and x2 on housing prices at different quantiles

對于圖2(b),β12<0,β22<0,表明在房價的所有分位數點上,該地區房東屬于低收入階層的比例對房價的邊際影響是負向的,并且整體而言這種影響會隨著房價的升高而逐漸增大,同時第二類房屋的負向影響高于第一類房屋.綜上所述,可以將上述房屋分為兩類,高檔小區和普通小區.對高檔小區（第二類）而言,房間數的增加會明顯提高房價,可見每個房間的價值都很高,比如市中心的房屋.對于普通小區而言,房間數的多少對房價影響不大,每個房間的價值較低,比如郊區的房屋.同時,對于高檔小區而言,隨著房價的升高,房間數對房價的邊際影響逐漸減小,這表明當房價足夠高時,人們會更多地關心居住地的品質,而不僅僅考慮房間的大小.

另一方面,波士頓房東屬于低收入階層的比例會對房價造成負面影響,并且對高檔小區的整體負面影響高于普通小區,特別是當房價處于極高的位置時,二者的差距會被進一步放大.這是顯然的,高檔小區(如別墅)的定位是高收入人群,該區域低收入購房者的比例越高,對此類房屋價格的提高越是不利的,事實證明富人更喜歡和富人居住.

本工作提出的模型相對于傳統混合分位數回歸模型而言,主要優點在于不僅考慮了解釋變量對位置的影響,而且考慮了解釋變量對尺度的影響,因此可以在實際運用中對模型的各個部分進行更有效的控制.在實證分析中,使用每棟住宅房間數和該地區房東屬于低收入階層的比例這兩個解釋變量,對波士頓的房屋價格進行分析,發現每棟住宅房間數對高檔小區的正向邊際影響高于普通小區,低收入階層的比例對高檔小區的負向影響高于普通小區.隨著房價逐漸升高,每棟住宅房間數對房屋價格的邊際影響最終會有所降低,而低收入階層的比例對房屋價格的負向影響會逐漸增大.