井曉丹,朱永忠,井 睿
(1.河海大學(xué) 理學(xué)院,南京 211100;2.西北農(nóng)林科技大學(xué) 經(jīng)濟(jì)管理學(xué)院,陜西 咸陽(yáng) 712100)
模型選擇是建模的基本問(wèn)題之一,對(duì)模型進(jìn)行選擇的目的是選擇輸入變量的一個(gè)最小子集,并使最終模型具有較好的預(yù)測(cè)性能。在最優(yōu)模型選擇過(guò)程中,要充分考慮計(jì)算的復(fù)雜程度和模型的預(yù)測(cè)性能。
到目前為止,許多科研工作者對(duì)貝葉斯模型選擇進(jìn)行了大量的研究,得出了很多模型選擇方法[1]。在模型選擇時(shí),常用的方法是根據(jù)預(yù)測(cè)性能對(duì)模型進(jìn)行選擇。交叉驗(yàn)證法[2]和廣義信息準(zhǔn)則[3]都可以得到模型預(yù)測(cè)性能的一個(gè)漸近無(wú)偏估計(jì)[4],但不是真正無(wú)偏估計(jì)模型的泛化能力。當(dāng)給定的數(shù)據(jù)集較小時(shí),這兩個(gè)估計(jì)中含有隨機(jī)誤差項(xiàng),較大的方差可能會(huì)導(dǎo)致過(guò)擬合現(xiàn)象,在模型選擇過(guò)程中導(dǎo)致選擇非最優(yōu)模型以及在所選模型的性能估計(jì)中產(chǎn)生偏倚[5-7]。模型選擇的另一種方法是構(gòu)建一個(gè)完整的參考模型[8-10],這種方法給出的預(yù)測(cè)模型幾乎與參考模型一致。如果完整的模型過(guò)于復(fù)雜或觀察所有變量的成本太高,那么該模型可以通過(guò)投影方法穩(wěn)健地簡(jiǎn)化,將全模型的信息投影到子模型上。投影方法考慮到形成全模型時(shí)模型的不確定性[11],將參考模型的后驗(yàn)信息投影到候選模型,從而使候選模型中的方差減小,降低模型中的過(guò)擬合現(xiàn)象。在預(yù)測(cè)性能方面,投影方法似乎優(yōu)于其他的模型選擇方法,但對(duì)于不同大小的數(shù)據(jù)集,投影方法對(duì)模型中自變量數(shù)目的選擇[12]存在不確定性。鑒于此,本文在對(duì)模型進(jìn)行選擇時(shí)將k折交叉驗(yàn)證與投影預(yù)測(cè)法相結(jié)合,數(shù)值實(shí)驗(yàn)結(jié)果表明,改進(jìn)后的方法能選出具有良好預(yù)測(cè)性能的最簡(jiǎn)模型。
統(tǒng)計(jì)學(xué)家們通常用效用函數(shù)描述模型的質(zhì)量,一種常用的衡量候選模型M預(yù)測(cè)性能好壞的效用函數(shù)[13]是:
其中,y?表示未來(lái)預(yù)測(cè)值,表示訓(xùn)練數(shù)據(jù)集。
由于未來(lái)預(yù)測(cè)值是未知的,導(dǎo)致效用函數(shù)u(M,?)無(wú)法求值。因此通常將式(1)的計(jì)算用下式代替:
(1)留一交叉驗(yàn)證
式(2)可以通過(guò)使用已獲得的數(shù)據(jù)集D代替真實(shí)數(shù)據(jù)進(jìn)行計(jì)算,生成概率分布pt(?)的一個(gè)估計(jì)。但使用與擬合模型相同的數(shù)據(jù)集D,會(huì)導(dǎo)致泛化能力有一個(gè)樂(lè)觀的估計(jì)。為了解決這一問(wèn)題,1974年Geisser[14]等提出留一交叉驗(yàn)證法,即將給定的數(shù)據(jù)集分成兩部分,一部分作為訓(xùn)練集,用于擬合模型;另一部分作為測(cè)試集,用于檢測(cè)模型的預(yù)測(cè)性能。其基本思想是:每次從數(shù)據(jù)集D=
中取一個(gè)樣本(xi,yi)作為測(cè)試集,用剩下的n-1 個(gè) 樣 本(xn,yn)}組成訓(xùn)練集。將此步驟重復(fù)n次,依次取遍所有的n個(gè)樣本分別作為測(cè)試集,最后用n個(gè)測(cè)試誤差的均值作為泛化誤差的估計(jì)。
由于模型需要被擬合n次,導(dǎo)致留一交叉驗(yàn)證法的計(jì)算量可能很大。如果n很大,不僅擬合次數(shù)很多,也會(huì)使得擬合單個(gè)模型的過(guò)程變得特別慢,那么,這種方法將非常耗時(shí)。
(2)k折交叉驗(yàn)證
1979年,Geisser和Eddy[2]提出了k折交叉驗(yàn)證。k折交叉驗(yàn)證是一種常用的樣本重用方法,它與留一交叉驗(yàn)證
其中,pt(?)是概率分布,下標(biāo)t代表此概率分布是由真實(shí)數(shù)據(jù)計(jì)算得到的。這個(gè)表達(dá)式粗略地描述了候選模型M的預(yù)測(cè)性能。
法類似,其本質(zhì)區(qū)別在于數(shù)據(jù)的切分方式不同。
k折交叉驗(yàn)證的基本思想是:將數(shù)據(jù)集等分成k個(gè)子集I1,…,Ik,每次從中取一個(gè)數(shù)據(jù)集作為測(cè)試集,而其余k-1個(gè)數(shù)據(jù)集作為訓(xùn)練集,重復(fù)此步驟k次。平均k次的結(jié)果作為泛化誤差的估計(jì):
其中,Is(i)表示將第i個(gè)樣本集作為測(cè)試集,DIs(i)表示在訓(xùn)練集中將第i個(gè)樣本集刪除。
k折交叉驗(yàn)證中,較小的k會(huì)使方差變大,導(dǎo)致模型預(yù)測(cè)效果變差。統(tǒng)計(jì)學(xué)家們通過(guò)利用大量數(shù)據(jù)集,使用不同方法進(jìn)行大量的實(shí)驗(yàn),表明10折是獲得最好誤差估計(jì)的恰當(dāng)選擇。但這也不是絕對(duì)的,爭(zhēng)議仍然存在。令k=10其優(yōu)勢(shì)在于:第一,方便計(jì)算。使用10折交叉驗(yàn)證只需要把模型擬合十次,可行性更高;第二,考慮到偏差-方差的權(quán)衡問(wèn)題,令k=10時(shí),會(huì)產(chǎn)生一個(gè)中等程度的偏差,使得測(cè)試誤差的估計(jì)不會(huì)有過(guò)大的偏差或方差[15]。本文中將使用k=10進(jìn)行討論。
廣義信息準(zhǔn)則可以估計(jì)模型的預(yù)測(cè)性能,包括參數(shù)的不確定性,也可以用于異常模型。其估計(jì)模型泛化能力的公式為:
其中,V是函數(shù)方差,由下式給出:
其中,·||·是布爾運(yùn)算中的“或”運(yùn)算,θ是模型參數(shù)。
參考模型M*和候選模型M之間的差異定義如下:
這里,兩個(gè)期望都采用了后驗(yàn)概率p(θ|D,M)。
模型選擇的另一種方法是建立一個(gè)參考模型M*,這種方法稱為參考模型法[8-10]。從預(yù)測(cè)角度看,用所有的變量形成參考模型將會(huì)得到更好的預(yù)測(cè)效果。參考模型法被認(rèn)為是對(duì)未來(lái)預(yù)測(cè)值的最好擬合。投影預(yù)測(cè)方法是參考模型方法的一種表現(xiàn)形式。
投影預(yù)測(cè)方法的基本思想是:將參考模型M*的后驗(yàn)信息投影到候選模型M上,因此候選模型的預(yù)測(cè)分布與參考模型的分布非常接近。候選模型的參數(shù)是由擬合的參考模型決定的,而不是由數(shù)據(jù)決定。
給出參考模型的參數(shù)θ*,在候選模型M的參數(shù)空間上,投影參數(shù)θ⊥定義如下:
計(jì)算KL散度,并根據(jù)式(6)分別計(jì)算出投影參數(shù)那么式(7)近似等于:
Dupuis和Robert[16]提出了模型選擇的一種方法:
其中,M0代表空模型,它與參考模型的差異是最大的,從變量選擇角度看,M0是自由變量模型,M*代表參考模型,M代表候選模型。Dupuis和Robert指出,在模型選擇過(guò)程中,要選擇具有較好解釋能力的最簡(jiǎn)模型,但他們沒(méi)有討論模型預(yù)測(cè)性能這一閾值的影響。一般情況下,對(duì)子模型的預(yù)測(cè)性能來(lái)說(shuō),?(M)是一個(gè)不可靠的指標(biāo)。這是因?yàn)橥ǔG闆r下,參考模型M*與真實(shí)數(shù)據(jù)生成的模型Mt是不同的,且M和M*之間存在差異,但M,M*卻與Mt有相同的預(yù)測(cè)性能。
投影方法考慮到形成全模型時(shí)的不確定性,能最好地保留全模型的預(yù)測(cè)性能。但在一般情況下,參考模型和子模型之間預(yù)測(cè)性能的差異是一個(gè)不可靠的指標(biāo)。這個(gè)屬性使得投影方法不容易確定最終模型的大小。
在模型空間上有如下分布:
其預(yù)測(cè)通過(guò)貝葉斯模型平均(BMA)獲得:
嚴(yán)格來(lái)說(shuō),BMA假設(shè)候選模型是真實(shí)數(shù)據(jù)生成的模型。此外,在變量選擇的背景下,BMA已被Raftery[17]證明在理論上和實(shí)踐中都具有良好的預(yù)測(cè)性能。
投影方法的預(yù)測(cè)效果雖然很好,但它決定最終模型中自變量數(shù)目的效果并不理想。針對(duì)這一點(diǎn),本文在對(duì)模型進(jìn)行選擇時(shí)將k折交叉驗(yàn)證與投影預(yù)測(cè)法相結(jié)合對(duì)其進(jìn)行改進(jìn)。
其基本思想為:使用投影方法進(jìn)行變量搜索,將顯著性變量作為訓(xùn)練集;交叉驗(yàn)證時(shí),訓(xùn)練集用來(lái)確定模型的大小,測(cè)試集對(duì)最終選擇的模型在獨(dú)立數(shù)據(jù)集上的性能進(jìn)行預(yù)測(cè)。
選擇中的過(guò)擬合現(xiàn)象和選擇的偏倚取決于模型預(yù)測(cè)性能估計(jì)中的方差和用于比較的候選模型數(shù)量。對(duì)于后者,與在變量組合中用逐步回歸法選擇模型相比,改進(jìn)后的方法中,交叉驗(yàn)證法通過(guò)比較p+1個(gè)模型(變量已經(jīng)排序),降低了模型的比較數(shù),當(dāng)僅用來(lái)決定模型大小時(shí),它產(chǎn)生的過(guò)擬合是相當(dāng)小的。
對(duì)模型中自變量進(jìn)行選擇時(shí),應(yīng)盡可能地使自變量個(gè)數(shù)達(dá)到最少,而且模型具有較高的預(yù)測(cè)精度。改進(jìn)后的方法允許在模型的預(yù)測(cè)性能和模型中自變量個(gè)數(shù)之間進(jìn)行取舍,在損失較少預(yù)測(cè)精度的基礎(chǔ)上簡(jiǎn)化模型。
對(duì)于臨界值U和α,變量選擇時(shí)應(yīng)滿足的條件為:
其中,ΔMLPD(m)表示交叉驗(yàn)證時(shí)訓(xùn)練集中m變量的樣本外預(yù)測(cè),U表示為了減少模型中自變量數(shù)量而損失的預(yù)測(cè)精度,α是置信水平。
方法改進(jìn)后,k折交叉驗(yàn)證令搜索重復(fù)進(jìn)行k次,然后對(duì)預(yù)測(cè)模型進(jìn)行擬合,需要注意的是:在這個(gè)過(guò)程中參考模型也被擬合了k次,每次也是用相同測(cè)試集中的數(shù)據(jù)對(duì)模型的預(yù)測(cè)性能進(jìn)行估計(jì)。因此可以比較預(yù)測(cè)模型和參考模型在獨(dú)立數(shù)據(jù)集上的預(yù)測(cè)性能,并進(jìn)行估計(jì),進(jìn)而判斷選擇的模型是否是最優(yōu)模型。
改進(jìn)后的方法優(yōu)點(diǎn)在于:第一,降低了模型的復(fù)雜程度,能夠保證模型中所含自變量的準(zhǔn)確性;第二,考慮到形成全模型時(shí)的不確定性,能夠最好地保留全模型的預(yù)測(cè)性能,準(zhǔn)確給出模型預(yù)測(cè)性能的一個(gè)漸近無(wú)偏估計(jì)。
2.2.1 模型
下文將討論線性模型中的回歸問(wèn)題。對(duì)于回歸,這里采用標(biāo)準(zhǔn)的高斯模型:
其中,x是輸入變量的p維向量,ω是對(duì)應(yīng)的權(quán)重,σ2是噪聲的方差。因?yàn)槌瑓?shù)τ2是共軛先驗(yàn)的,因此這個(gè)模型在大部分情況下都可以求解。對(duì)于回歸模型(13),通過(guò)在輸入向量x=(x0,x1,…,xp)和相應(yīng)的權(quán)向量ω=(ω0,ω1,…,ωp)中增加常數(shù)項(xiàng),將截距項(xiàng)計(jì)入。參數(shù)ατ,βτ,ασ,βσ的值將在數(shù)據(jù)集的描述中一起給出。
由于考慮變量的選擇問(wèn)題,因此子模型中輸入變量的數(shù)目不同,使得x和ω的維數(shù)是變化的。在預(yù)測(cè)性能的比較中,參考模型M*使用BMA[16]。
2.2.2 數(shù)據(jù)仿真
本文將進(jìn)行兩個(gè)仿真實(shí)驗(yàn),分別采用不同評(píng)價(jià)指標(biāo)比較 LOO-CV、CV-10、WAIC、BMA-Proj、IMP-Method的預(yù)測(cè)性能和變量選擇,并由此說(shuō)明方法改進(jìn)后的優(yōu)越性。為了在模型選擇方法之間產(chǎn)生有效的比較,在每次實(shí)驗(yàn)中,各種模型選擇方法都在相同的訓(xùn)練集和測(cè)試集對(duì)上進(jìn)行。
首先,做一個(gè)模擬的變量選擇實(shí)驗(yàn),用以說(shuō)明不同方法之間的差異。數(shù)據(jù)分布如下:
仿真實(shí)驗(yàn)中,將變量總數(shù)設(shè)為p=100。對(duì)所有變量進(jìn)行分組,每組五個(gè)變量。每個(gè)變量xj的均值為0,方差為1,且與同一組中其他變量相關(guān),相關(guān)系數(shù)設(shè)為ρ,但與其他組的變量均不相關(guān)。將前三組變量的權(quán)重設(shè)為(ω1∶5,ω6∶10,ω11∶15) =(ξ,0.5ξ,0.25ξ),其余變量的權(quán)重設(shè)為0。常量ξ調(diào)整數(shù)據(jù)集的信噪比。也就是說(shuō),在數(shù)據(jù)集中有15個(gè)相關(guān)變量,85個(gè)無(wú)關(guān)變量。將訓(xùn)練集的大小設(shè)為n=100,300,500,令相關(guān)系數(shù)ρ=0,0.5,0.9,為了得到不同水平相關(guān)系數(shù)ρ的比較結(jié)果,將ξ的值設(shè)為σ2Var[y]。在回歸模型(13)中,令參數(shù)ασ=0.5,使用BMA作為參考模型M[18]。*
對(duì)每個(gè) (n,ρ)的組合,模型選擇方法見(jiàn)表1。LOO-CV、CV-10、WAIC、BMA-Proj中用逐步回歸法對(duì)變量進(jìn)行搜索,即從空模型開(kāi)始,使每一步選入的變量都最大限度地增加了模型的預(yù)測(cè)性能。
表1 數(shù)據(jù)仿真中使用的模型選擇方法
將預(yù)測(cè)模型在數(shù)據(jù)量為n?=1000且與訓(xùn)練集數(shù)據(jù)不相關(guān)的測(cè)試集上進(jìn)行測(cè)試。對(duì)模型的預(yù)測(cè)性能進(jìn)行評(píng)價(jià)時(shí),使用對(duì)數(shù)預(yù)測(cè)密度的均值(MLPD):
為了更直觀地比較不同方法的預(yù)測(cè)性能,用下式表示性能差異:
其中,M是選擇的子模型,M*是參考模型(BMA)。若ΔMLPD(M)=0,則說(shuō)明所用方法與BMA有相同的預(yù)測(cè)性能;出現(xiàn)負(fù)值說(shuō)明所用方法的預(yù)測(cè)性能比BMA糟糕;相對(duì)應(yīng)的,正值表示所用方法的預(yù)測(cè)性能更好。由于BMA已經(jīng)被證明具有良好的預(yù)測(cè)性能[11,17,19],因此對(duì)模型進(jìn)行選擇的目的是希望找到與BMA預(yù)測(cè)性能接近的最簡(jiǎn)模型。
圖1模型選擇方法的比較
圖1左列中正方形、圓形以及三角形對(duì)應(yīng)的橫坐標(biāo)表示模型中自變量的平均數(shù),灰色線表示數(shù)據(jù)集中相關(guān)變量的數(shù)目是15、右列中正方形、圓形以及三角形對(duì)應(yīng)的橫坐標(biāo)表示不同模型選擇方法的預(yù)測(cè)性能與BMA預(yù)測(cè)性能的差異大小,三條垂直的點(diǎn)線代表空模型的性能(僅指截距項(xiàng)),正方形、圓形以及三角形代表變量間相關(guān)水平的強(qiáng)弱(相關(guān)程度見(jiàn)圖例)。
從圖1可以看出:仿真實(shí)驗(yàn)中使用的所有方法的預(yù)測(cè)性能都不如BMA方法的預(yù)測(cè)性能好。從預(yù)測(cè)的角度來(lái)看,模型平均通常會(huì)產(chǎn)生最好的預(yù)測(cè)效果。因此,對(duì)模型進(jìn)行選擇的主要目的是簡(jiǎn)化模型,并且基本上不會(huì)影響預(yù)測(cè)的準(zhǔn)確性,而不是試圖通過(guò)考慮模型的不確定性來(lái)提高預(yù)測(cè)精度。
其次,對(duì)于小數(shù)據(jù)集(n=100)的模型選擇,LOO-CV、CV-10、WAIC選出的模型預(yù)測(cè)性能比較差,甚至比沒(méi)有變量的模型(右列中三條垂直的點(diǎn)線)預(yù)測(cè)性能更糟。這是因?yàn)槿狈?shù)據(jù),預(yù)測(cè)性能估計(jì)中的高方差導(dǎo)致選擇了過(guò)擬合的模型。過(guò)度擬合是一個(gè)潛在的問(wèn)題,阻礙模型的選擇,尤其是當(dāng)數(shù)據(jù)集較小時(shí),會(huì)導(dǎo)致估計(jì)中的高方差,高方差對(duì)模型預(yù)測(cè)精度的提高是不利的。此外,從圖1可以看出:變量間的高相關(guān)性在一定程度上改善了預(yù)測(cè)效果。上述三種方法僅對(duì)大數(shù)據(jù)集(n=500)的預(yù)測(cè)性能較好。相比之下,投影預(yù)測(cè)法的預(yù)測(cè)效果明顯更好,特別是對(duì)于小的數(shù)據(jù)集(n=100),其預(yù)測(cè)性能接近BMA。但從左列可以看出:在不同大小的數(shù)據(jù)集上,投影方法選擇的模型中自變量個(gè)數(shù)的變化比較大,在小數(shù)據(jù)集上模型復(fù)雜程度較大。而方法改進(jìn)后,模型中自變量個(gè)數(shù)趨于穩(wěn)定,在較小數(shù)據(jù)集上模型的復(fù)雜程度變小,且模型的預(yù)測(cè)性能仍接近BMA的預(yù)測(cè)性能。
采用逐步回歸法對(duì)變量進(jìn)行排序,然后重新進(jìn)行實(shí)驗(yàn),并將模型的CV效用和測(cè)試效用作為評(píng)價(jià)模型性能的指標(biāo)。
圖2 CV效用和測(cè)試效用
圖2中虛曲線表示CV效用,代表模型選擇中的過(guò)擬合現(xiàn)象;實(shí)曲線表示變量排序后的測(cè)試效用。兩條曲線之間的距離表示選擇引起的偏倚。垂直虛線表示模型中自變量的平均數(shù)。
令ρ=0.5,用逐步回歸法對(duì)變量進(jìn)行選擇。CV效用(10折)用選出的顯著性變量進(jìn)行計(jì)算,測(cè)試效用在數(shù)據(jù)集中余下的其他數(shù)據(jù)上進(jìn)行計(jì)算。
由圖2可知,前三種方法中,所選模型預(yù)測(cè)性能的變化是非常大的。在不同大小的數(shù)據(jù)集上進(jìn)行變量選擇時(shí),模型選擇中出現(xiàn)過(guò)擬合現(xiàn)象,而且當(dāng)訓(xùn)練集變大時(shí),選擇過(guò)程中的過(guò)擬合現(xiàn)象減?。粡目漳P烷_(kāi)始,在預(yù)測(cè)模型中一次添加一個(gè)變量,模型的CV效用比較高,但測(cè)試效用卻很糟糕。也就是說(shuō),用逐步回歸法對(duì)變量進(jìn)行選擇時(shí),模型的預(yù)測(cè)性能被高估。由此可知,經(jīng)過(guò)變量選擇后,CV效用是一個(gè)樂(lè)觀的估計(jì)。但是對(duì)空模型和包含所有變量的模型而言,CV效用和測(cè)試效用幾乎相同,因?yàn)檫@些模型不涉及任何選擇。
與LOO-CV、CV-10、WAIC相比,投影預(yù)測(cè)方法中測(cè)試效用的變化程度小,且不容易產(chǎn)生過(guò)擬合現(xiàn)象。但從圖2可以看出:投影方法對(duì)模型中自變量個(gè)數(shù)的選取并不理想。
方法改進(jìn)后,模型的復(fù)雜程度降低,且仍具有良好的預(yù)測(cè)性能。當(dāng)增加多個(gè)變量時(shí),子模型會(huì)越來(lái)越接近參考模型(BMA),從而避免了預(yù)測(cè)精度下降的問(wèn)題。
表2 改進(jìn)前后模型中自變量個(gè)數(shù)及CV值的比較
由表2可知,方法改進(jìn)后模型的大小趨于穩(wěn)定,但預(yù)測(cè)精度略有下降;樣本量越大的數(shù)據(jù)集確定的模型越合適。
通過(guò)數(shù)據(jù)仿真可以看出:改進(jìn)后的方法產(chǎn)生了較好的預(yù)測(cè)效果,它簡(jiǎn)化了模型,而且不會(huì)失去太多的預(yù)測(cè)精度。在尋找預(yù)測(cè)性能良好的最簡(jiǎn)模型時(shí)似乎是最穩(wěn)健的方法。
本文將k折交叉驗(yàn)證與投影預(yù)測(cè)法相結(jié)合對(duì)模型進(jìn)行選擇,通過(guò)仿真實(shí)驗(yàn)說(shuō)明了使用這種方法對(duì)模型進(jìn)行選擇時(shí),能選出具有良好預(yù)測(cè)性能的最簡(jiǎn)模型。但改進(jìn)后的方法也存在一些問(wèn)題,方法改進(jìn)后不僅增加了計(jì)算量,而且預(yù)測(cè)精度略有下降。針對(duì)出現(xiàn)的不足,還需要進(jìn)行更深入的研究。