郭媛媛，楊雪梅，孫志華,3?

(1 中國(guó)科學(xué)院大學(xué)數(shù)學(xué)科學(xué)學(xué)院， 北京 100049； 2 華北電力大學(xué)數(shù)理學(xué)院， 北京 102206；3 中國(guó)科學(xué)院大數(shù)據(jù)挖掘與知識(shí)管理重點(diǎn)實(shí)驗(yàn)室， 北京 100190)(2019年5月14日收稿； 2019年10月9日收修改稿)

分位回歸模型具有穩(wěn)健的特點(diǎn)，并且能夠?qū)憫?yīng)變量的分布做出精細(xì)的描述， 因此獲得很多學(xué)者的關(guān)注。 目前， 分位回歸方法已經(jīng)成為分析數(shù)據(jù)的一個(gè)非常重要的工具，廣泛地應(yīng)用在金融、醫(yī)療、生物研究等領(lǐng)域。 分位回歸模型的研究和應(yīng)用可參考文獻(xiàn)[1-4]及相關(guān)文獻(xiàn)。 參數(shù)分位回歸模型是在實(shí)際中應(yīng)用廣泛的一類(lèi)回歸模型[5]。但很多時(shí)候無(wú)法對(duì)參數(shù)分位回歸模型進(jìn)行正確的設(shè)定?；谡`定的參數(shù)分位回歸的統(tǒng)計(jì)推斷結(jié)果經(jīng)常是不可信的， 文獻(xiàn)[6-7] 闡述了參數(shù)分位回歸需要進(jìn)行檢驗(yàn)的必要性。非參數(shù)分位回歸不存在誤定的風(fēng)險(xiǎn)，但當(dāng)樣本量比較小且協(xié)變量比較多時(shí)，非參數(shù)分位回歸方法可能會(huì)受到維數(shù)禍根的問(wèn)題的困擾。

對(duì)分位回歸，構(gòu)建目標(biāo)函數(shù)時(shí)用到的損失函數(shù)ρτ(r)=τr-I{r<0}具有不光滑的特點(diǎn), 從而使得求解目標(biāo)函數(shù)的最小值比較困難，且可能出現(xiàn)多個(gè)最小值點(diǎn)的情況，參見(jiàn)文獻(xiàn)[5,8-11]。一種解決上面問(wèn)題的方法是將分位回歸模型的求解問(wèn)題轉(zhuǎn)化為線性規(guī)劃問(wèn)題，再利用單純形法或內(nèi)點(diǎn)法進(jìn)行計(jì)算。不管是單純形法，還是內(nèi)點(diǎn)法，運(yùn)算效率都不能令人滿意。2000 年，Hunte和Lange[12]提出一種新的用于求解分位回歸問(wèn)題的算法，即MM算法。MM算法概念簡(jiǎn)單，易于執(zhí)行，且數(shù)值穩(wěn)定，比內(nèi)點(diǎn)法擁有更強(qiáng)的數(shù)值計(jì)算能力。文獻(xiàn)[11]對(duì)4種求解分位回歸問(wèn)題的算法，即內(nèi)點(diǎn)法、MM算法、坐標(biāo)下降法和ADMM算法進(jìn)行比較研究，驗(yàn)證了MM 算法具有數(shù)值穩(wěn)定和計(jì)算效率高的特點(diǎn)。

單指標(biāo)分位回歸模型具有降維的效果，同時(shí)保持了非參數(shù)分位回歸的穩(wěn)健性，其估計(jì)問(wèn)題的研究吸引了很多研究者的興趣。文獻(xiàn)[13-14] 提出基于兩步迭代的估計(jì)方法，文獻(xiàn)[15-16]進(jìn)一步提出不需要迭代的估計(jì)方法，文獻(xiàn)[17-18]又提出基于貝葉斯方法的估計(jì)方法，文獻(xiàn)[19-21]探討單指標(biāo)分位回歸模型的變量選擇以及加權(quán)復(fù)合單指標(biāo)分位回歸模型的估計(jì)。然而這些文獻(xiàn)所提的估計(jì)方法都基于內(nèi)點(diǎn)法來(lái)實(shí)現(xiàn)，內(nèi)點(diǎn)法在計(jì)算分位回歸模型時(shí)，計(jì)算效率低、耗時(shí)久，尤其在樣本量較大的情況下，這種缺點(diǎn)更為明顯。MM 算法在求解分位回歸模型的估計(jì)時(shí)比較高效和便捷，這在文獻(xiàn)[22-24]中均有體現(xiàn)，但是沒(méi)有文獻(xiàn)研究單指標(biāo)分位回歸模型的MM算法，故本文研究單指標(biāo)分位回歸模型估計(jì)的MM算法。

我們借鑒文獻(xiàn)[12] 的方法，對(duì)單指標(biāo)分位回歸模型的每一步迭代程序中目標(biāo)函數(shù)構(gòu)建其替代函數(shù)，從而將復(fù)雜的優(yōu)化問(wèn)題簡(jiǎn)單化。然后，基于優(yōu)化函數(shù)再進(jìn)行求解計(jì)算得到估計(jì)值。我們構(gòu)建的優(yōu)化函數(shù)是光滑的，并能夠保證每次迭代目標(biāo)函數(shù)是下降的。數(shù)值模擬和實(shí)例分析結(jié)果表明基于MM 算法的估計(jì)程序具有較好的穩(wěn)定性，能夠得到比較準(zhǔn)確的估計(jì)結(jié)果，并且相較于傳統(tǒng)的內(nèi)點(diǎn)算法具有更強(qiáng)的數(shù)值計(jì)算能力，用時(shí)更短。

1 單指標(biāo)分位回歸模型的估計(jì)介紹

對(duì)于給定的分位數(shù)τ∈(0,1)，在給定x的條件下，響應(yīng)變量y的τ分位數(shù)θτ(x)與協(xié)變量x之間的關(guān)系如下：

θτ(x)=g(xTγ),

其中x∈d是d維協(xié)變量，g(·)表示未知的一元聯(lián)系函數(shù)。另外γ=(γ1,…,γd)T為未知的單指標(biāo)向量，滿足‖γ‖=1且γ1>0,‖·‖表示Euclidean范數(shù)。這個(gè)約束條件是為了模型的可識(shí)別性[25]，已廣泛應(yīng)用在有關(guān)單指標(biāo)模型的文獻(xiàn)中。

本文采用局部線性方法對(duì)γ和g(·)進(jìn)行估計(jì)，詳細(xì)內(nèi)容可參考文獻(xiàn)[13]，具體算法如下：

(1)

(2)

step 4重復(fù)step 1和step 2直至收斂。

(3)

2 估計(jì)程序的MM算法

2.1 MM算法介紹

下面介紹MM算法的基本思想。假設(shè)需要最小化的目標(biāo)函數(shù)為L(zhǎng)(θ):p→，θk為第k步的迭代值。MM算法每次迭代分兩步來(lái)進(jìn)行。首先，構(gòu)造目標(biāo)函數(shù)的優(yōu)化函數(shù)Q(θ|θk):p×p→ 滿足

Q(θk|θk)=L(θk),

Q(θ|θk)≥L(θ)?θ.

(4)

然后, 對(duì)優(yōu)化函數(shù)Q(θ|θk)進(jìn)行最小化, 得到下一步的迭代值θk+1，則有

Q(θk+1|θk)≤Q(θk|θk).

(5)

綜合式(4)和式(5)，可知L(θk+1)≤L(θk)。這種下降趨勢(shì)保證了MM算法具有顯著的數(shù)值穩(wěn)定性。

觀察目標(biāo)函數(shù)式(1)～式(3)，可以發(fā)現(xiàn)它們均為非光滑函數(shù)， 因而不易得到最優(yōu)解，故借用文獻(xiàn)[12] 提出的MM 算法的思想來(lái)處理這個(gè)問(wèn)題。其主要的處理方式如下：首先給ρτ(r)加一擾動(dòng)ε，得到其近似函數(shù)

2.2 非參數(shù)部分的估計(jì)

首先，將式(1)中的目標(biāo)函數(shù)修正為

(6)

(7)

由此，可以將MM算法總結(jié)為如下步驟：

1) 選擇迭代初始值θ0和一個(gè)較小的正常數(shù)ε，置k=0；

3) 令k=k+1，判斷是否滿足收斂準(zhǔn)則，若滿足收斂準(zhǔn)則，即： 當(dāng)

時(shí)，可終止迭代，其中δ是預(yù)先取定的足夠小的數(shù)。否則返回2)繼續(xù)迭代，直到滿足收斂準(zhǔn)則。

2.3 單指標(biāo)向量的估計(jì)

首先定義式(2)中目標(biāo)函數(shù)的近似函數(shù)：

在γk處的優(yōu)化函數(shù)可以構(gòu)建為

(8)

最后，可將基于MM算法的單指標(biāo)模型的估計(jì)總結(jié)為如下步驟：

1)參考第1節(jié)step 1所提供的方法，得到γ的初始估計(jì)；

4)重復(fù)2)、3)步驟，直至收斂。

2.4 聯(lián)系函數(shù)g(·)的估計(jì)

對(duì)于聯(lián)系函數(shù)的估計(jì)，式(6)中目標(biāo)函數(shù)式的近似函數(shù)可定義為

3 數(shù)值模擬

3.1 模擬1

借鑒文獻(xiàn)[13]模擬1的模型設(shè)置，考慮模型

在不同的分位數(shù)下，計(jì)算上述評(píng)價(jià)指標(biāo)值，所得結(jié)果列于表1和表2，可以看出，用MM算法計(jì)算單指標(biāo)分位回歸模型，無(wú)論是單指標(biāo)向量的估計(jì)還是聯(lián)系函數(shù)的估計(jì)，都有比較好的結(jié)果，且在較小和較大的分位數(shù)下依然有良好的表現(xiàn)，這表明本文所提出的計(jì)算方法是有效的。將該估計(jì)結(jié)果與文獻(xiàn)[13] 中關(guān)于該模型的模擬結(jié)果進(jìn)行對(duì)比，可以看出，用MM 算法計(jì)算單指標(biāo)分位回歸模型，與內(nèi)點(diǎn)法相比，估計(jì)的精度是相似的，估計(jì)系數(shù)的偏差大小都在10-3～10-2，估計(jì)系數(shù)的樣本標(biāo)準(zhǔn)差數(shù)量級(jí)均為10-2，但在計(jì)算效率上，本文所提出的方法大大優(yōu)于內(nèi)點(diǎn)法，這將在模擬3 中展示和說(shuō)明。

表1 模擬1中不同的τ下，估計(jì)的Bias、SE、MSE、CI、CPTable 1 The Bias、SE、MSE、CI、CP of under different choices of τ in simulation 1

表2 模擬1中不同τ的選擇下，的ASE、AAE的Mean、SETable 2 Outcomes of for the models under different choices of τ in simulation 1

3.2 模擬2

借鑒文獻(xiàn)[13]模擬2的模型設(shè)置，考慮模型

驅(qū)動(dòng)層介于系統(tǒng)層和硬件設(shè)備之間，提供應(yīng)用程序訪問(wèn)硬件設(shè)備資源的接口，同時(shí)也為存儲(chǔ)軟件提供了基礎(chǔ)環(huán)境和接口。linux系統(tǒng)將設(shè)備分為3類(lèi)：字符設(shè)備、塊設(shè)備、網(wǎng)絡(luò)設(shè)備?？紤]到安全存儲(chǔ)的硬件加解密設(shè)備與塊設(shè)備、網(wǎng)絡(luò)設(shè)備特性的差異，而與面向流的字符設(shè)備類(lèi)似，驅(qū)動(dòng)層在采用字符設(shè)備的基礎(chǔ)上提供相應(yīng)的庫(kù)文件供應(yīng)用程序?qū)崿F(xiàn)加解密功能[15]。

3.3 與內(nèi)點(diǎn)法比較

由前兩個(gè)模擬可以發(fā)現(xiàn)，MM算法在單指標(biāo)分位回歸估計(jì)的計(jì)算問(wèn)題中表現(xiàn)良好，接下來(lái)比較MM算法與內(nèi)點(diǎn)法的計(jì)算效率，考慮如下3個(gè)模型：

表3 模擬2中不同的τ下，估計(jì)的Bias、SE、MSE、CI、CP Table 3 The Bias、SE、MSE、CI、CP of under different choices of τ in simulation 2

表4 模擬2中不同τ的選擇下，的ASE、AAE的Mean、SE Table 4 Outcomes of for the models under different choices of τ in simulation 2

表5 n=100， 模型1、2、3的估計(jì)結(jié)果比較 Table 5 Estimation comparison among models 1， 2， and 3 with n=100

可以發(fā)現(xiàn)MM算法所用的時(shí)間遠(yuǎn)遠(yuǎn)少于內(nèi)點(diǎn)法，且隨著樣本量的增大，這種計(jì)算效率上的優(yōu)勢(shì)更加明顯。這是由于用內(nèi)點(diǎn)法解決分位回歸問(wèn)題，是將目標(biāo)函數(shù)及約束條件轉(zhuǎn)化為線性規(guī)劃問(wèn)題，再用內(nèi)點(diǎn)法來(lái)求解該問(wèn)題，但轉(zhuǎn)化之后的線性規(guī)劃問(wèn)題，協(xié)變量維數(shù)與樣本量的大小有著正相關(guān)的關(guān)系，這種方法極大地增加了算法的計(jì)算量與所用時(shí)間。本文第3.2節(jié) 中的問(wèn)題轉(zhuǎn)化為線性規(guī)劃后，協(xié)變量的維數(shù)為2n+p，第3.3 節(jié) 中的問(wèn)題轉(zhuǎn)化為線性規(guī)劃后，維數(shù)為2n2+p，具體轉(zhuǎn)化方法及維數(shù)的增加量可見(jiàn)文獻(xiàn)[11]。而MM 算法只需對(duì)p維矩陣做運(yùn)算，故兩種方法的計(jì)算效率隨樣本量的增加會(huì)產(chǎn)生越來(lái)越大的差距。

表6 n=200，模型1、2、3的估計(jì)結(jié)果比較Table 6 Estimation comparison among models 1， 2， and 3 with n=200

4 實(shí)例分析

θτ(MEDV|RM,TAX,PTRATIO,LSTAT)=

g(γ1RM+γ2log(TAX)+γ3PTRATIO)+

γ4log(LSTAT).

用本文提出的方法對(duì)該問(wèn)題進(jìn)行估計(jì)，計(jì)算在不同分位數(shù)下系數(shù)的估計(jì)值，并采用bootstrap方法估計(jì)標(biāo)準(zhǔn)差，方法如下，具體細(xì)節(jié)可參考文獻(xiàn)[13,29]。

重復(fù)模擬100次計(jì)算標(biāo)準(zhǔn)差，所得結(jié)果列于表7。從表7可以發(fā)現(xiàn)，RM的系數(shù)在不同的分位數(shù)下皆為正，這表明每棟房屋的房間數(shù)量越多，房?jī)r(jià)就越高且收入越多的家庭更加在意每棟房屋的房間數(shù)量；log(TAX)的系數(shù)為負(fù)且隨分位數(shù)逐漸變大，這表明不動(dòng)產(chǎn)的稅率越高，房?jī)r(jià)越低且收入較低的家庭更加在意不動(dòng)產(chǎn)稅率的大?。籔TRATIO 的系數(shù)為負(fù)且隨分位數(shù)變化較小，這表明學(xué)生與教師的比例越大，即教師資源越匱乏，房?jī)r(jià)越低且低收入家庭與高收入家庭對(duì)教育的重視程度是同樣大的；log(LSTAT) 系數(shù)為負(fù)且隨分位數(shù)逐漸變小，這表明一個(gè)地區(qū)低收入人群所占的百分比越高，房?jī)r(jià)越低且收入較高的家庭更加在意一個(gè)地區(qū)的低收入人群比例。

表7 波士頓數(shù)據(jù)集在單指標(biāo)分位回歸模型下的系數(shù)估計(jì)及標(biāo)準(zhǔn)差估計(jì)Table 7 Coefficient estimation and standard deviationestimation of Boston data set under the single-indexquantile regression mode

圖1 聯(lián)系函數(shù)g(u)及在不同的τ下,的估計(jì)Fig.1 The link function g(u) and the estimation of under different choices of τ

5 結(jié)論

本文研究單指標(biāo)分位回歸模型估計(jì)方法的MM算法。相比于內(nèi)點(diǎn)法，MM算法極大地縮短了計(jì)算時(shí)間，提高了運(yùn)算效率。此外，本文給出單指標(biāo)分位回歸模型在MM 算法下的參數(shù)估計(jì)公式，在每次迭代過(guò)程中，將協(xié)變量與響應(yīng)變量的觀測(cè)值直接代入公式，即可得到參數(shù)的估計(jì)值，避免了每次迭代都要優(yōu)化目標(biāo)函數(shù)的麻煩。