二元選擇分位數(shù)回歸模型的貝葉斯估計(jì)方法及模擬研究

邸俊鵬，張曉峒

（1.上海社會(huì)科學(xué)院 數(shù)量經(jīng)濟(jì)研究中心，上海 200020；2.南開(kāi)大學(xué) 數(shù)量經(jīng)濟(jì)研究所，天津 300071）

0 引言

分位數(shù)回歸和貝葉斯估計(jì)屬于當(dāng)前計(jì)量經(jīng)濟(jì)學(xué)理論的前沿領(lǐng)域。分位數(shù)回歸作為一種不同于均值回歸、應(yīng)用更廣泛、提供信息更豐富的計(jì)量方法，雖然早在1978年就被Koenker等[1]提出，但關(guān)于它的理論研究和應(yīng)用研究方興未艾。貝葉斯估計(jì)尤其是采用馬爾科夫鏈蒙特卡羅（MCMC）方法，在小樣本性質(zhì)、假設(shè)檢驗(yàn)以及預(yù)測(cè)方面具有比傳統(tǒng)估計(jì)方法無(wú)可比擬的優(yōu)勢(shì)，因而受到越來(lái)越多的關(guān)注。Yu和Moyeed[2]在2001年首次將貝葉斯分析方法應(yīng)用于分位數(shù)回歸模型，提出貝葉斯分位數(shù)估計(jì)方法，并證實(shí)了該估計(jì)方法的有效性。

然而在貝葉斯分位數(shù)估計(jì)中大多數(shù)研究都是針對(duì)連續(xù)型因變量的，對(duì)離散型因變量，如二元變量的分位數(shù)回歸模型研究較少。目前離散型因變量模型代表性的文獻(xiàn)主要有：Manski（1975，1985）[3,4]定義了半?yún)?shù)的二元選擇分位數(shù)回歸估計(jì)量；Koenker和Hallock（2001）[5]主要致力于中位數(shù)上的二元選擇分位數(shù)回歸（Binary Quantile Regression）的相關(guān)研究；進(jìn)而，Kordas（2006）[6]將該研究擴(kuò)展到了各個(gè)分位數(shù)上，并闡明了二元數(shù)據(jù)的各分位數(shù)回歸可以提供解釋變量對(duì)因變量更豐富的影響。上述研究都是基于傳統(tǒng)的頻率學(xué)派方法的。Skouras（2003）[7]、Florios和Skouras（2008）[8]指出這些方法并不能保證得到目標(biāo)函數(shù)的全局最優(yōu)解。因此，采用頻率學(xué)派的方法求解二元選擇分位數(shù)回歸模型，無(wú)論在統(tǒng)計(jì)量的一致性方面還是在統(tǒng)計(jì)推斷方面都受到質(zhì)疑。

本文將針對(duì)二元選擇分位數(shù)回歸模型的貝葉斯估計(jì)方法進(jìn)行探索性研究。首先介紹基于ALD的貝葉斯二元選擇分位數(shù)回歸估計(jì)方法；進(jìn)而通過(guò)模擬實(shí)驗(yàn)，對(duì)不同先驗(yàn)設(shè)定和不同抽樣算法下二元選擇分位數(shù)回歸估計(jì)量性質(zhì)進(jìn)行比較研究；最后，比較頻率學(xué)派方法和貝葉斯估計(jì)方法對(duì)二元選擇分位數(shù)模型進(jìn)行估計(jì)時(shí)的不同表現(xiàn)。

1 二元選擇分位數(shù)回歸的貝葉斯估計(jì)方法

1.1 二元選擇分位數(shù)回歸

標(biāo)準(zhǔn)的二元選擇模型表達(dá)式為：

F（?）為累積分布函數(shù)。F（?）所采用的形式不同，二元選擇模型也相應(yīng)不同，常用的二元選擇模型如表1所示：

表1 常用的二元選擇模型

因?yàn)榫€性概率模型不能保證條件概率的預(yù)測(cè)值在0和1之間，即使加以約束，其預(yù)測(cè)結(jié)果也往往與現(xiàn)實(shí)不符，故應(yīng)用較少。Logit模型假設(shè)條件概率分布的累積分布函數(shù)；Probit模型假設(shè)條件概率為標(biāo)準(zhǔn)正態(tài)分布的累積分布函數(shù)；而對(duì)于互補(bǔ)對(duì)數(shù)模型，假設(shè)條件概率為極值分布的累積分布函數(shù)。Probit曲線和logit曲線都是在概率為0.5處存在拐點(diǎn)，但logit曲線在兩個(gè)分布的尾部要比Probit曲線厚。

其中，Q（?）和F（?）分別表示潛變量的條件分位數(shù)函數(shù)和條件分布函數(shù)，βτ為在τ分位數(shù)下的參數(shù)向量。因?yàn)闈撟兞渴怯^測(cè)不到的，所以不能直接用模型（2）進(jìn)行估計(jì)。

為了估計(jì)離散選擇的分位數(shù)回歸模型，Powell（1986）[9]指出可以借助分位數(shù)回歸的同變性（Equivariance），這是分位數(shù)回歸的一大優(yōu)勢(shì)。假設(shè)h（?）是因變量yi的變換（transition）函數(shù)，是定義在實(shí)數(shù)空間上的非遞減函數(shù)，那么對(duì)于任意的隨機(jī)變量y，則有：，即隨機(jī)變量先變換再進(jìn)行分位數(shù)估計(jì)與先進(jìn)行分位數(shù)估計(jì)再變換是等價(jià)的。因此隨機(jī)變量y的單調(diào)變換不影響分位數(shù)估計(jì)結(jié)果。在離散選擇模型中，潛變量與因變量yi之間滿足這種單調(diào)變換關(guān)系，因此潛變量與解釋變量的分位數(shù)估計(jì)結(jié)果可以通過(guò)可觀測(cè)的因變量與解釋變量的分位數(shù)估計(jì)得到。

1.2 頻率學(xué)派二元選擇分位數(shù)回歸方法存在的問(wèn)題

頻率學(xué)派在估計(jì)二元選擇分位數(shù)模型時(shí)，首先要對(duì)其進(jìn)行極大得分估計(jì)（maximum score estimate）。這個(gè)方法幾乎對(duì)誤差分布不作過(guò)多假設(shè)，只要求誤差以自變量為條件的中位數(shù)為0。因此，極大得分估計(jì)不知道解釋變量和誤差項(xiàng)的函數(shù)形式，因而它也適應(yīng)存在異方差的情形。Manski（1975）[3]最早只是關(guān)注中位數(shù)上的估計(jì)，后來(lái)Manski（1985）[4]將其擴(kuò)展到更一般的分位數(shù)情形。極大得分估計(jì)如下：

其中sgn（）為符號(hào)函數(shù)，ρτ（?）為損失函數(shù)。

Kim和Pollard（1990）[10]指出極大得分估計(jì)量收斂速度低而且漸近分布復(fù)雜，而極限分布的復(fù)雜性限制了它在統(tǒng)計(jì)推斷方面的應(yīng)用，同時(shí)也使得目標(biāo)函數(shù)不存在漸近的一階條件。Delgado等（2001）[11]采用子采樣（subsampling）的方法克服上述問(wèn)題。他們從理論上證明了子采樣對(duì)于極大得分估計(jì)是有效的，并給出了模擬證據(jù)。但該方法的一個(gè)主要缺點(diǎn)是計(jì)算量大，因而它只適用于解決小樣本和維度低的情形。該方法起初用于中位數(shù)估計(jì)，Kordas（2006）[6]將其擴(kuò)展到一般的分位數(shù)情形，但是這個(gè)平滑的估計(jì)量對(duì)誤差分布假設(shè)過(guò)于嚴(yán)格。模擬實(shí)驗(yàn)表明，即使樣本容量很大時(shí)，用正態(tài)分布近似也是不準(zhǔn)確的；即使采用自舉也很難得到估計(jì)量的標(biāo)準(zhǔn)誤差[12]。Skouras（2003）[7]、Florios和Skouras（2008）[8]則側(cè)重于就對(duì)目標(biāo)函數(shù)（3）的優(yōu)化，但他們指出沒(méi)有哪個(gè)算法可以保證能得到全局最優(yōu)解，即使是采用改進(jìn)的算法，如混合整數(shù)規(guī)劃（mixed integer programs）。總之，采用頻率學(xué)派的方法求解二元選擇分位數(shù)回歸，無(wú)論在統(tǒng)計(jì)量的一致性方面還是在統(tǒng)計(jì)推斷方面都受到質(zhì)疑，即使提出了改進(jìn)的方法，仍存在某些缺陷。

1.3 基于ALD的貝葉斯二元選擇分位數(shù)回歸估計(jì)方法

τ為關(guān)注的分位數(shù)，如τ=0.5，則是二元選擇中位數(shù)回歸。同樣地：

這里，F(xiàn)y*（?）設(shè)定為非對(duì)稱拉普拉斯變量y*的累積分布函數(shù)。

π（β）為回歸系數(shù)的先驗(yàn)分布，I（?）為指示函數(shù)。這個(gè)后驗(yàn)分布不是我們所熟悉的分布，不可能對(duì)其進(jìn)行直接抽樣，而MCMC方法可以解決這個(gè)問(wèn)題。以β為條件的y*的全條件后驗(yàn)分布是可知的：

當(dāng)yi=1時(shí)，

當(dāng)yi=0 時(shí)，

由（5）可知，以y*、τ和觀測(cè)值為條件的β的后驗(yàn)分布是：

這個(gè)全條件后驗(yàn)分布不同于潛變量y*的全條件后驗(yàn)分布，不存在解析形式，因此采用MCMC方法，如Gibbs抽樣、M-H抽樣來(lái)獲取后驗(yàn)分布。具體而言，給定觀測(cè)值、參數(shù)的先驗(yàn)分布和感興趣的分位點(diǎn)，借助式（4）和（5）給出的條件分布，采用Gibbs或者M(jìn)-H抽樣，可以得到聯(lián)合后驗(yàn)分布（6）。以極大似然估計(jì)值為初值，并刪除開(kāi)始的部分估計(jì)值，最后得到一個(gè)抽樣值序列，從而可以方便地得到參數(shù)的點(diǎn)估計(jì)和置信區(qū)間。這也是貝葉斯方法相對(duì)于頻率方法的一大優(yōu)勢(shì)。

關(guān)于先驗(yàn)分布，Yu和Moyeed（2001）[2]指出，回歸參數(shù)的先驗(yàn)分布π（β）可以是任意的，即使是一個(gè)不合適的均勻分布，得到的后驗(yàn)分布也是合適的。此外，由式（6）可知，后驗(yàn)分布是由誤差服從非對(duì)稱拉普拉斯分布的假設(shè)所決定的，這也表明模型參數(shù)的估計(jì)量和統(tǒng)計(jì)推斷是受這個(gè)假設(shè)影響的。然而，在下文的模擬中可以看出，即使誤差不服從這個(gè)假設(shè)分布，相關(guān)結(jié)論也是相當(dāng)穩(wěn)定的。

1.4 模擬實(shí)驗(yàn)

下面通過(guò)一個(gè)模式實(shí)驗(yàn)展示采用貝葉斯方法分析二元選擇分位數(shù)回歸的過(guò)程。在二元離散選擇模型中，如果存在異方差，貝葉斯分位數(shù)回歸方法可以有效地捕捉和反應(yīng)這種異方差。參數(shù)先驗(yàn)采用模糊的標(biāo)準(zhǔn)正態(tài)分布，β～N（0,10），以此來(lái)減弱對(duì)后驗(yàn)分布的影響。MCMC抽樣值的收斂性由從邊際分布中得到的抽樣值的時(shí)間序列是否平穩(wěn)來(lái)判定。

首先通過(guò)下面的異方差模型，生成n=200數(shù)據(jù)：

在上述數(shù)據(jù)生成過(guò)程（7）的基礎(chǔ)上，定義離散二元變量yi：當(dāng)時(shí)yi=0，當(dāng)時(shí)yi=1。采用貝葉斯二元選擇分位數(shù)回歸方法對(duì)其進(jìn)行估計(jì)。由于為潛變量，是不可觀測(cè)的，所以本文使用的數(shù)據(jù)是二元離散變量yi和解釋變量x，來(lái)估計(jì)潛變量中隱含的異方差。在0.05、0.25、0.5、0.75、0.95分位數(shù)上的抽樣值軌跡圖和核密度圖如圖2所示。

比較使用二元選擇貝葉斯分位數(shù)回歸（圖2）和采用潛變量的貝葉斯分位數(shù)回歸得到的各分位數(shù)下估計(jì)值（圖1），見(jiàn)表2。

表2 貝葉斯二元選擇分位數(shù)估計(jì)量與真值的比較

圖2貝葉斯二元選擇分位數(shù)回歸系數(shù)的抽樣軌跡圖和核密度圖

表2中，β的真值與圖1中擬合直線的斜率相對(duì)應(yīng)。從表2可以看出，通過(guò)可觀測(cè)的變量采用貝葉斯二元選擇分位數(shù)回歸方法得到的估計(jì)量雖然與真值在數(shù)值上有偏差，但是它仍可以反應(yīng)各個(gè)分位數(shù)上潛變量與自變量的相關(guān)關(guān)系，比如在低分位數(shù)上負(fù)相關(guān)，且隨著分位數(shù)的提高，相關(guān)系數(shù)逐漸增大。在下文中比較不同抽樣方法對(duì)二元選擇分位數(shù)回歸估計(jì)結(jié)果時(shí)，采用某一分位數(shù)下的值占所有分位數(shù)下值之和的比例，來(lái)討論哪種抽樣方法更能反映各個(gè)分位數(shù)下變量之間的關(guān)系。

2 不同先驗(yàn)設(shè)定下二元選擇分位數(shù)回歸估計(jì)量性質(zhì)的比較

針對(duì)一般化的形式：

誤差項(xiàng)存在異方差。為比較不同先驗(yàn)分布條件下估計(jì)量的差異，本文設(shè)定回歸系數(shù)的先驗(yàn)分布為正態(tài)分布，而且隨著方差的減小，信息由弱到強(qiáng)。形式分別如下：

先驗(yàn)設(shè)定1：β～N（0,100）

先驗(yàn)設(shè)定2：β～N（0,10）

先驗(yàn)設(shè)定3：β～N（0,1）σ為參數(shù)化的尺度變量，設(shè)定其先驗(yàn)分布σ～χ2（3）。

根據(jù)上述數(shù)據(jù)生成過(guò)程，結(jié)合不同的先驗(yàn)信息，運(yùn)用Gibbs抽樣算法對(duì)模型參數(shù)進(jìn)行估計(jì)，考察不同的先驗(yàn)分布對(duì)估計(jì)結(jié)果的影響。

表3 不同先驗(yàn)分布下二元選擇分位數(shù)回歸估計(jì)量分布的數(shù)字特征

由表3可知，與連續(xù)數(shù)據(jù)貝葉斯分位數(shù)回歸不同，二元選擇模型的貝葉斯分位數(shù)回歸參數(shù)在各個(gè)分位點(diǎn)下的偏誤和標(biāo)準(zhǔn)差均不受先驗(yàn)分布的影響。因此，在對(duì)二元選擇模型進(jìn)行貝葉斯分位數(shù)回歸時(shí)，不必考慮先驗(yàn)分布的選取，即先驗(yàn)可以是無(wú)信息先驗(yàn)。

3 不同抽樣方法對(duì)二元選擇分位數(shù)回歸估計(jì)量的影響

在貝葉斯分析中，Gibbs抽樣和M-H抽樣是目前較流行的抽樣算法。而不同的抽樣方法施行不同的算法和抽樣規(guī)則，本文將比較采用Gibbs抽樣和M-H抽樣對(duì)二元選擇分位數(shù)回歸估計(jì)量的影響。

數(shù)據(jù)生成過(guò)程為：

在Gibbs抽樣過(guò)程中，參數(shù)β的先驗(yàn)分布為，將尺度變量σ參數(shù)化，并設(shè)定其先驗(yàn)分布為χ2（3）。在M-H抽樣中，參數(shù)的先驗(yàn)分布與Gibbs設(shè)定相同，建議分布為高斯分布。每種抽樣算法下共迭代6000次，為消除初值的影響，去掉前2000個(gè)抽樣值。通過(guò)抽樣值時(shí)序圖和自相關(guān)圖可以對(duì)兩種抽樣方法進(jìn)行初步的比較。圖3是為樣本容量為200時(shí)，各分位數(shù)下Gibbs和M-H抽樣值軌跡圖和自相關(guān)圖。左列為Gibbs抽樣結(jié)果，右列為M-H抽樣結(jié)果；在每個(gè)分位數(shù)下，第一行為4000個(gè)抽樣值的時(shí)序軌跡圖，第二行為自相關(guān)圖。

圖3各分位數(shù)下Gibbs和M-H抽樣值軌跡圖和自相關(guān)圖

從0.05、0.25和0.5分位數(shù)的時(shí)序軌跡圖可以看出，兩種抽樣方法得到的馬爾科夫鏈都收斂，但Gibbs抽樣值的分布呈左偏形態(tài)，比較而言M-H抽樣值的分布較為對(duì)稱。從自相關(guān)圖角度看，Gibbs抽樣值的自相關(guān)度下降快于M-H抽樣值的自相關(guān)度。如在0.05分位數(shù)下，Gibbs抽樣值的自相關(guān)度在滯后10期時(shí)基本趨于零，而M-H抽樣值的自相關(guān)度在滯后35期才趨于零。在高分位數(shù)0.75和0.95下，抽樣結(jié)果表明Gibbs抽樣明顯優(yōu)于M-H抽樣：Gibbs抽樣值的時(shí)序圖集中在某一個(gè)穩(wěn)定值上下波動(dòng)，且自相關(guān)圖快速衰減，基本在5期滯后衰減為零；而M-H抽樣值從時(shí)序圖可以看出馬爾科夫鏈?zhǔn)諗啃圆?，存在自相關(guān)。此外不衰減的自相關(guān)圖也說(shuō)明了這一點(diǎn)。因此，從抽樣的時(shí)序圖和自相關(guān)圖可以直觀地初步判斷，在高分位下Gibbs抽樣明顯優(yōu)于M-H抽樣，而在低分位數(shù)下需要進(jìn)一步考察。下面對(duì)抽樣得到的后驗(yàn)分布的統(tǒng)計(jì)特征進(jìn)行進(jìn)一步的分析。

依照上述數(shù)據(jù)生成過(guò)程分別生成25、75、100、200、500、800容量的樣本，并在每個(gè)樣本容量下分別進(jìn)行Gibbs和M-H抽樣，抽樣結(jié)果見(jiàn)表4。

表4 不同抽樣算法下二元選擇分位數(shù)回歸結(jié)果比較

通過(guò)以上分析，得出相關(guān)結(jié)論：①關(guān)于標(biāo)準(zhǔn)差：在給定樣本容量下，各個(gè)分位數(shù)上，Gibbs抽樣得到的后驗(yàn)分布的標(biāo)準(zhǔn)差均小于M-H的標(biāo)準(zhǔn)差。這表明在貝葉斯二元選擇分位數(shù)回歸中，采用Gibbs抽樣得到的估計(jì)量精度更高，統(tǒng)計(jì)推斷更準(zhǔn)確。②關(guān)于偏誤：在兩種抽樣算法下，抽樣后驗(yàn)分布的均值對(duì)真值的偏誤都比較大，其結(jié)果與表3類似（表中未列出），而且Gibbs的偏誤比M-H的偏誤更大。但需要強(qiáng)調(diào)的是，在二元離散選擇模型中，每個(gè)分位數(shù)上估計(jì)量偏誤的絕對(duì)值并不重要，重要的是各個(gè)分位數(shù)上估計(jì)量的相對(duì)值，因此在表4中本文給出了各個(gè)分位數(shù)下估計(jì)量的偏誤占所有分位數(shù)上估計(jì)量偏誤之和的比重（特定樣本容量下，第一行數(shù)值），同時(shí)計(jì)算真值在各個(gè)分位數(shù)下占所有分位數(shù)上真值之和的比重（特定樣本容量下，第二行數(shù)值）。如果采用某抽樣方法時(shí)各個(gè)分位數(shù)上參數(shù)估計(jì)量的偏誤占總偏誤的比重更接近于對(duì)真值所對(duì)應(yīng)的比重，則該抽樣方法更優(yōu)。如樣本容量為100下，Gibbs抽樣得到的β在0.05、0.25、0.5、0.75、0.95分位數(shù)上的偏誤占總偏誤的比重分別為-0.270、0.012、0.151、0.232、0.542；而M-H 抽樣得到的β在上述分位數(shù)上的偏誤占總偏誤的比重分別為-0.074、0.010、0.046、0.144、0.874；真值在各分位數(shù)上的偏誤占總偏誤的比重分別為-0.187、0.025、0.202、0.397、0.563，顯然Gibbs抽樣下偏誤的相對(duì)值與真值偏誤的相對(duì)值更接近。由此可知，貝葉斯二元選擇分位數(shù)回歸采用Gibbs抽樣方法得到的估計(jì)量更能描述潛變量與自變量關(guān)系的全貌。

4 不同估計(jì)方法對(duì)二元選擇分位數(shù)回歸估計(jì)量的影響

針對(duì)二元選擇分位數(shù)回歸頻率學(xué)派的主要方法是基于二次抽樣標(biāo)準(zhǔn)差的二元選擇分位數(shù)回歸方法[4]和基于漸近標(biāo)準(zhǔn)差的二元選擇分位數(shù)回歸方法[6]。下面比較了不同樣本容量下貝葉斯二元選擇分位數(shù)回歸方法Bayes（LAD）與頻率學(xué)派這兩種方法的參數(shù)統(tǒng)計(jì)性質(zhì)。

數(shù)據(jù)生成過(guò)程為：

差）

在不同樣本容量下，比較兩種傳統(tǒng)二元分位數(shù)回歸于貝葉斯二元分位數(shù)回歸得到的估計(jì)結(jié)果，實(shí)驗(yàn)結(jié)果見(jiàn)表5。

表5 不同估計(jì)方法下二元選擇分位數(shù)回歸估計(jì)量比較

主要結(jié)論為：①隨著樣本容量的增加，三種估計(jì)方法下估計(jì)量的偏誤、均方誤以及置信區(qū)間都減小。這表明可采用的數(shù)據(jù)越多，參數(shù)的估計(jì)量越趨向于真值，而且估計(jì)的不確定性越小。②貝葉斯二元選擇分位數(shù)回歸（Bayes（LAD）估計(jì)量比頻率學(xué)派得到的估計(jì)量擁有更小的偏誤、更小的均方誤和更精準(zhǔn)的置信區(qū)間。③隨著樣本容量的增加，頻率方法和貝葉斯方法得到的二元選擇分位數(shù)回歸估計(jì)量的差異會(huì)減??；BRQ和sBRQ在大樣本下才能進(jìn)行可靠的統(tǒng)計(jì)推斷；同時(shí)這也說(shuō)明，在小樣本情形下，貝葉斯方法用于二元選擇分位數(shù)回歸估計(jì)效果更好。上述結(jié)果也驗(yàn)證了Benoit等（2012）[13]、Abrevaya和 Huang（2005）[12]、Kottas和Krnjajic（2009）[14]得出的結(jié)論。

5 總結(jié)

本文對(duì)二元選擇貝葉斯分位數(shù)回歸方法進(jìn)行了研究，模擬結(jié)果表明：二元選擇模型的貝葉斯分位數(shù)回歸參數(shù)在各個(gè)分位點(diǎn)下的偏誤和標(biāo)準(zhǔn)差均不受先驗(yàn)分布的影響，因此在對(duì)二元選擇模型進(jìn)行貝葉斯分位數(shù)回歸時(shí)，不必考慮先驗(yàn)分布的選取。在進(jìn)行二元選擇分位數(shù)模型的貝葉斯估計(jì)時(shí)，與M-H抽樣相比，采用Gibbs抽樣得到的估計(jì)量精度更高，統(tǒng)計(jì)推斷更準(zhǔn)確，更能描述潛變量與自變量關(guān)系的全貌，而且在高分位下Gibbs抽樣的優(yōu)勢(shì)更明顯。貝葉斯二元選擇分位數(shù)回歸估計(jì)量比頻率學(xué)派得到的估計(jì)量擁有更小的偏誤、更小的均方誤和更精準(zhǔn)的置信區(qū)間，尤其是在小樣本情形下，采用貝葉斯方法對(duì)二元選擇分位數(shù)回歸模型進(jìn)行估計(jì)效果更好。