基于估計穩定性的變系數模型選擇

呂曉玲，劉擷芯，戴秀紅

（中國人民大學a.應用統計研究中心；b.數據挖掘中心；c.統計學院，北京100872）

0　引言

變系數模型最初由Hastie和Tibshirani（1993）[1]提出，是一類非常重要的非參數回歸模型。它考慮了指示變量與協變量之間的交互效應，與常規的線性模型相比具有更強的適應性和解釋性。它在計量經濟、生物統計、社會科學等多個領域中都有著廣泛的應用，已成為處理多元非參數、半參數回歸問題的有力工具[2-4]。

在對實際問題進行回歸建模時，為了減小可能存在的模型誤差，研究者在初始建模時往往會引入很多可能與被解釋變量相關的協變量。但為了提高模型的預測精度、增強模型的可解釋性，研究者需要判別對因變量具有顯著影響的重要變量。因此，變量選擇已成為當今統計分析中一個重要的研究課題。各種各樣的正則化估計方法應運而生，也即在傳統損失函數的基礎上加入懲罰函數，從而實現變量選擇和參數估計的同時進行。

與其他正則化估計方法一樣，在對變系數模型的正則化估計中，調節參數的選擇至關重要。常用的選參準則包括交叉驗證（Cross-validation,CV）、貝葉斯信息準則（Bayesian Information Criterion,BIC）、赤池信息量準則（Akaike Information Criterion,AIC）等，這些方法在判別顯著變量、提高模型預測準確性等多個方面都非常有效。不過，這些方法所確定的模型或多或少都不具有穩定性，尤其是在高維數據設定下。對此，Lim和Yu（2013）[5]針對線性模型的LASSO問題中正則化參數的選擇提出了ESCV（estimation stability cross validation）方法，這一方法有效彌補了以往選參方法在高維數據分析中不穩定這一不足。因此，本文將ESCV方法作為一種選參準則引入到變系數模型的正則化估計中，以期提高變系數模型的穩定性。

1　變系數模型及KLASSO估計

其中ei∈R1是隨機噪聲，滿足E(ei|Xi,Zi)=0，系數向量β(z)={β1(z),...,βd(z)}T∈Rd是未知的，并且是Z的光滑函數。其真實值可給定為β0(z)={β01(z),...,β0d(z)}T∈Rd。不失一般性的假定存在整數d0≤d，對于任意的j≤d，有0＜E{(Z)}＜∞，但對于任意的j＞d，0i0E{(Z)}=0，簡單來說，就是假定前d個預測變量與響i0應變量是真實相關的，其余的是不相關的。

Wang和Xia（2013）[4]提出的KLASSO（Kernel Lasso）估計，是一種將流行的核光滑方法與加罰估計結合起來的估計方法，其基本思想是將一個典型的收縮方法即LASSO算法的局部連續核估計應用于變系數模型，KLASSO估計方法如下：

對于任意的指標變量Zi∈[0,1]，β(z)可以通過最小化下面的局部加權最小二乘函數來估計：

對于B0={β0(Z1),…,β0(Zn)}T∈Rn×d，可通過最小化如下全局最小二乘函數來估計：

注意到，在模型假設下，矩陣B0的最后(d-d0)列應該都是0，因此變量選擇就等價于在矩陣B0中辨別出稀疏列。借用Yuan和Lin（2006）[6]提出的Group LASSO的方法來判別稀疏列，提出下面的加罰估計：

其中，bj是B的第j列，‖.‖表示常用的歐幾里德范數。Wang和Xia（2013）[4]使用局部二階近似算法得到上述估計的解，并證明了KLASSO方法有很好的理論性質。

上述方法涉及到調節參數（核函數K的窗寬h以及懲罰函數的λj,1≤j≤d）的選擇問題。文中第一步使用了留一交叉驗證方法選取h，然后簡使用BIC準則選取λ0。

2　基于估計穩定性的新的變量選擇方法

估計穩定性對于一個合理的估計過程來說是一個必要的性質，如果隨著樣本的不同，估計的值變動相當大，那么這個估計是沒有意義的[5，7]。當用L2誤差來度量不同樣本間的差異時，估計穩定性顯然與方差相關，然而在統計學上人們傾向于用穩定性而非變異性來形容不同條件和環境對所估計模型的影響，這就是說穩定性是一個比方差或者變異性更廣泛的概念。現有文獻中研究變系數模型變量選擇和估計方法的文章很多，但是研究變系數模型穩定性的文章卻很少。然而模型穩定性對于任何模型來說都是重要的，尤其是在數據采集技術及數據存儲技術日益強大的今天，人們經常可以收集到非常多的變量和樣本數據，數據往往呈現海量或高維的形態。在分析這些大數據和高維數據時，統計方法的不穩定性出現得更為普遍。

在對有限樣本且無模型假定的數據建模時，交叉驗證（CV）是建模常用方法，它依賴數據重抽樣來評估候選模型的預測誤差。具體做法是：在給定的建模樣本中，拿出大部分樣本作為訓練集建立模型，留小部分樣本作為測試集，用訓練集所建立的模型對預測集進行預測，并求出測試樣本的預測誤差，記錄它們的誤差平方和，這個過程一直進行，直到所有的樣本都被作為測試集測試了一次而且僅被測試一次時，選出預測誤差平方和最小的模型作為最終模型。交叉驗證的目的是為了得到可靠穩定的模型，然而，數據重抽樣會引發模型的不穩定性，尤其是數據為大數據或者高維數據時。在正則化估計如LASSO估計方法中經常用CV方法來選擇調節參數，然而CV通常會導致模型不穩定，從而不利于可靠性解釋。Lim和Yu（2013）[5]提出了ESCV方法，即將數據可信度需求加入到交叉驗證中，ESCV是一個基于估計穩定性ES（Estimation stability）并將其與CV結合起來的一種無需模型假定的變量選擇方法。

在變系數模型KLASSO估計實際計算中，需要選擇合適的調節參數，調整參數的選擇在加罰估計的變量選擇過程中起著極其重要的作用。當調整參數λ=0時,所有的變量都被選進模型；當調整參數λ=∞時,那么模型中不含有任何變量。λ起到了控制模型復雜度的作用。λ取值越大，得到的模型越簡單。反之，λ取值越小，得到的模型越復雜。大的λ給出的估計的方差比較小,而小的λ對應的模型偏差會比較小。因此，在KLASSO估計中調節參數λ的確定對模型的穩定性有重要影響，估計方差和模型偏差之間一個好的平衡就需要選出一個比較理想的λ，如何選擇一個合適的調整參數使得模型在預測性和解釋性上都能達到一個理想結果就成為人們所關心的重要問題。

Wang和Xia（2009）[4]提出的變系數模型的KLASSO估計中，確定收縮參數λ，是采用BIC最小準則，BIC雖然易于計算，但其有效性依賴于模型假定，而且它是漸近性結果，因此在樣本量有限的情況下，BIC模擬結果表現不一定很好，且BIC在統計性能上是不穩定的[8]，當數據是高維數據時，即樣本量n小于變量維度p時，Lim和Yu（2013）[5]將ESCV、CV（cross validation）和BIC應用與Lasso方法，并對這三種方法所估計的模型的穩定性進行比較，結果表明ESCV方法在多個指標上面都表現較好。基于此，本文將ESCV作為一種選參方法引入到變系數模型加罰估計中，以期提高變系數模型在傳統變量選擇方法BIC下的模型穩定性，挑選λ的準則是選擇具有局部最小標準化方差的[λ]，即就是要使ES(λ)值最小。

本文的分析比較中，選用模型預測均方誤差（MSE）、模型大小（MS）以及顯著性變量個數（NOSV）及其百分比（PSV）四個方面來度量模型的穩定性。均方誤差是度量模型穩定性的首要標準，模型預測能力不好，則模型不可靠。模型大小即所選變量的個數，在高維數據分析中，需要控制模型復雜度，若所選變量過多，模型太過復雜，模型的穩定性就可能得不到保證。在眾多變量中對模型有顯著性影響的自變量對模型穩定性有重要影響，顯著性變量個數及其百分比是指挑選多個變量的情況下，對模型有顯著影響的自變量個數及其占所選全部變量的比例。

3　模擬研究

3．1　正態分布下變系數模型模擬

本文的次模擬是模擬實際中常見的數據分布形式，即自變量服從或近似服從正態分布的情形，模擬所采用的模型如下：

其中假定X=(Xi1,Xi2,Xi3)服從正態分布N(0,1)，ei服從正態分布N(0,0.8)，σe=1.2，并設定不顯著變量(Xi4,…,Xip)服從正態分布N(0,0.8)，全部變量之間的協示變量Z服從

i均勻分布U(0,1)。

在自變量來自正態分布的變系數模型中，本文用KLASSO方法進行參數估計和變量選擇，在估計過程中，調節參數分別選用ESCV準則和傳統的BIC準則進行確定。為了比較在不同變量維度下BIC和ESCV方法進行變量選擇對模型穩定性的影響，本文設定總變量數p∈{10,30,60,70,90,100}。

在模型樣本量n=50不變，變量數p不斷增大的情況下，將每個模型隨機模擬100次，結果如表1所示。

表1　正態分布下模擬結果

從表1中可以看出，在樣本量n=50保持不變而總變量數p變化時，兩種變量選擇方法的均方誤差（MSE）都隨著變量總數的增加而增大。當變量維度p小于樣本量n，即當p為10和30時，ESCV方法估計的預測誤差、變量個數以及顯著性百分比都不如BIC方法，但在高維數據情形下，即當變量維度p大于樣本量n時，ESCV方法的預測誤差、變量個數以及顯著性變量百分比優于BIC方法且這種優勢隨著變量維度p的增大越發明顯。

當樣本量n=50，變量維度p=70時，ESCV的100次模擬平均預測誤差為2.29，BIC的100次模擬平均預測誤差為2.80，ESCV所選模型的MSE小于BIC所選模型且較BIC所選模型的MSE降低了18.21%，同時ESCV的100次模擬所選變量個數平均為16.92，BIC的100次模擬所選變量個數平均為42.88，ESCV所選模型的變量個數不到BIC所選模型變量個數的一半，ESCV方法較BIC大大縮減了模型變量維度，在顯著性變量占所選變量百分比上，ESCV所選模型的顯著性變量百分比為14.36%，是BIC所選模型的兩倍。當p=100時，ESCV所選模型的MSE較BIC所選模型降低約20%，ESCV所挑選的變量個數僅占全部變量數的16.05%，而BIC所選變量個數占全部變量數的66.29%，ESCV所選變量個數大約是BIC所選變量個數的四分之一，在顯著性變量百分比上，ESCV所選模型的顯著性變量百分比為14.70%，是BIC所選模型的三倍。由上述分析可知，當變量維度p大于樣本量n時，ESCV方法在模型穩定性上的表現優于BIC，且在樣本量不變的情況下，隨著變量維度p的增加優勢越發明顯。

3．2　稀疏情況下變系數模型模擬

Lim和Yu（2013）[5]給出在高維稀疏情況下，ESCV方法較BIC估計方法更能突顯模型的穩定性優勢，因此，改變變系數模型中自變量的分布，將自變量分布稀疏化來探討在自變量稀疏情況下，ESCV方法與傳統變量選擇方法BIC在模型穩定性上的不同表現。

模擬模型假定與第一種情況相同，只是假定自變量X來自于服從均勻分布的隨機稀疏矩陣sprand(ss,p,d),其中ss為樣本量，p為總變量個數，d為非零元素分布密度的大小，設為0.4。

模型在樣本量n=50保持不變，總變量數p∈{10,30,60,70,90,100}的情況下，對每組模型設定都進行100次模擬，在KLASSO估計方法下，分別選用ESCV選參準則和BIC選參準則來選擇調節參數，結果如表2所示。

表2　稀疏情況下模擬結果

從表2可以看出，自變量來自于稀疏分布情形時，在樣本量p小于n的情況下，即p=10或30時，在均方誤差MSE上，ESCV方法所選模型對應的均方誤差較BIC方法要小，在模型大小以及顯著性變量百分比上，ESCV方法對應的模型較BIC要大，但隨著變量維度的增加，當變量維度p大于樣本量n，即p∈{60,70,90,100}時，ESCV篩選變量的優勢大大增強，均方誤差也越來越小。當p=100時，ESCV方法平均100次模擬的均方誤差為1.06，BIC方法對應的均方誤差為1.18，ESCV方法的對應模型的MSE較BIC減少了10.17%，且ESCV方法為模型所挑選的變量個數平均為4.09，大大少于BIC方法所選變量個數33.77，ESCV所選變量個數大約是BIC方法所選變量個數的四分之一，在顯著性變量百分比上，ESCV所選顯著性變量占全部所選變量的50.86%，而BIC方法所選顯著性變量占所選全部變量的7.26%，即BIC方法較ESCV方法更多的選擇了不顯著變量。

與自變量來自于正態分布下的變系數模型相比，稀疏分布下兩種方法擬合的模型在模型均方誤差、模型大小、顯著性變量百分比上的表現都更好，但與BIC方法相比，ESCV方法對模型穩定性的影響更為顯著，即ESCV方法下的模型均方誤差、模型大小減小幅度更大，顯著性變量百分比增大幅度則更多，如稀疏情形下，BIC方法對應的模型大小是正態分布下模型大小的二分之一，而ESCV方法對應的模型大小是正態分布下的四分之一。且稀疏分布情形下，低維時，ESCV方法對應的模型均方誤差小于BIC，這與自變量來自于正態分布且數據為低維情況下ESCV預測誤差大于BIC相反。高維時，不論自變量來自哪種分布，ESCV在均方誤差、模型大小、顯著性變量百分比上都優于BIC。在高維稀疏數據分析中，ESCV較BIC方法的穩定性優勢更加明顯。

4　實例分析

4．1　Boston Housing數據分析

本文將分析Boston住房數據，該數據是在1970年Boston地區收集的506個人口普查區的房價信息。本文沿用Fan和Huang（2005）[9]的變量設定，將MEDV（業主自用房子的中位數，以$1000為單位）作為響應變量Y，LSTAT（地區較低地位人群占總體的百分比）作為指示變量Z，數據集中的其他七個變量作為自變量INT（截距，X1），CRIM（鎮上人均犯罪率，X2），RM（每座房子的平均房間數，X3），PTRATIO（鎮上學生-老師人數比，X4），NOX（氮氧化物濃度，X5），TAX（房間全價值物業稅率，以$10000為單位，X6）和AGE（業主占用房子建造早于1940年的比例，X7）。

將全部506個樣本單元隨機分成十份，每次選取其中九份做訓練集，另一份為預測集。在應用模型之前，需要將自變量X和指示變量LSTAT標準化處理。由抽取的訓練集數據建立變系數模型，并在KLASSO方法下分別選用ESCV與BIC兩種方法進行變量選擇，記錄10次抽樣數據擬合中每個變量被選入模型的次數，以及每次模型預測集的均方誤差MSE，結果如表3所示。

從表3中可以看出，對10次擬合進行平均，ESCV的均方誤差MSE與BIC方法對應模型的MSE相當，但在10次擬合中，ESCV方法傾向于選擇前三個變量，變量X1，X2，X3被選中的次數分別為10、9、10，即ESCV方法能夠穩定的選擇出前三個自變量，而BIC方法在每次模擬時所選變量的個數以及傾向于選擇哪些變量都不穩定，即抽樣數據不同時，BIC不能保證所估計模型的穩定性。

為了進一步理解抽樣數據分析中的不穩定性，給出第6次抽樣下，BIC方法對應模型的估計系數圖，如圖1（上）；ESCV方法所對應模型的估計系數圖，如圖1（下）。

從圖1（上）中可以看出，在BIC方法下，用此次抽樣所得的455個數據進行變量選擇，所選取的變量數為1，即僅第一個變量INT顯著不為0，其他六個變量的全部估計為0。用此次數據建立模型并對余下的51個數據進行預測，所對應的平均預測誤差為0.5115；圖1（下）給出了相同樣本數據下，用ESCV方法進行變量選擇，所選取的變量個數為3，即INT、RM、CRIM這3個變量顯著不為0，用所建模型對余下的51個數據進行預測，所對應的平均預測誤差為0.4082。對比兩個圖可以看出在一次抽樣數據擬合中，對同一個自變量的估計，例如自變量INT，ESCV方法對該變量估計的波動程度要顯著小于BIC方法，即ESCV方法估計的變量系數更為穩定。

表3　Boston住房數據結果

圖1 BIC（上）和ESCV（下）方法下自變量估計系數變動情況

從Boston housing的數據分析中可以看到，當所抽取的樣本數據發生變化時，BIC方法選擇的變量個數就會隨之發生大的變動，而ESCV卻能穩定地選出對因變量有重要影響的自變量，且對所選變量的系數估計也更穩定。

4．2　新浪新聞數據分析

利用爬蟲技術在新浪新聞網站獲取2013年7月1日至9月30日（共12周）財經（標簽為+1）和健康（標簽為-1）兩類新聞文檔。共390篇，兩類的比例為1：1。利用分詞軟件將原始文本數據轉化為文檔詞頻矩陣。隨機選取150篇文檔為訓練集，剩余240篇為測試集。利用LASSO初步篩選出162個關鍵詞。指示變量z為時間并以周為時間單位。即假定關鍵詞對新聞類別的影響與其出現的時間有關。

在變系數模型KLASSO估計中，分別用ESCV和BIC兩種方法確定調節參數，用所確定模型的分類準確率來衡量模型預測誤差，用所選擇變量的個數來確定模型大小。將數據進行10次抽樣，并對每次抽樣所得數據進行擬合，每次擬合模型所選變量數，以及模型分類準確率所得結果如表4所示：

表4　新聞數據分類準確率和模型大小

從表4中可以看出，在樣本量n=150小于變量維度p=162時，盡管10次所抽取的樣本數據不一樣，但用ESCV方法選擇出的變量所建模型的分類準確率始終要大于BIC方法的分類準確率，且在模型大小上，BIC方法傾向于選擇幾乎所有的變量，而ESCV傾向于選擇固定的20多個變量，ESCV挑選出的變量個數始終顯著小于BIC方法下的變量個數，在此次新聞數據分析中，ESCV變量選擇方法有明顯的降維作用，在分類準確率上較傳統BIC方法也有優勢，采用ESCV變量選擇方法可以顯著提高模型的穩定性。

為了更好地理解ESCV選擇變量的穩定性要優于BIC，給出某次抽樣下，ESCV選擇變量的情況，如在第四次抽樣中，ESCV在162個關鍵詞中，選擇了22個對分類有重要影響的關鍵詞，即選擇的變量個數為22。分析這22個關鍵詞，大致可以分為三大類，第一類是明顯跟經濟類相關的詞，如：美元、經濟、投資、下降、漲、中國、企業、發展、部門；第二類是跟健康類相關詞，如：醫院、性、肌肉、成分、疾病、健康、破壞、食物、效果；第三類是對分類沒有很重要作用，但在兩類文章中都會出現的詞，如：發布、公布、好。

在此次抽樣中，BIC方法所選變量個數為160，幾乎所有由LASSO初步篩選出的關鍵詞都被引入到模型中。因此BIC所建模型比ESCV更復雜。在兩種選參方法下，“健康”一詞都被選中，但兩種方法對其重要性的估計不同，下面本文給出兩種方法下，關鍵詞“健康”在12周中估計系數的變化對比圖，如圖2所示。

圖2ESCV和BIC方法下關鍵詞“健康”估計系數變化圖

從圖2中可以看出，在12周中，“健康”一詞在分類上始終是有重要作用的詞匯，但在BIC方法下估計出的系數值即關鍵詞重要性波動較大，而在ESCV方法下估計系數值變動較小，在前4周幾乎沒有變動，第4周略有下降，但在后4周中又開始固定不變，系數值基本保持在0.4995的水平上，ESCV比BIC方法對該詞的估計更為穩定。

從新浪新聞數據分析可以看出，在高維數據情況下，ESCV在模型預測、變量選擇上較BIC表現得更好，即在變系數模型KLASSO估計下，選用ESCV準則比選用BIC準則進行變量選擇所確立模型的穩定性更強。

5　總結

本文是基于Yu和Lim（2013）[5]提出的ESCV方法以及對模型穩定性的度量標準，將ESCV方法引入到變系數模型加罰估計中，以期提高變系數模型的穩定性。本次研究主要是基于Wang和Xia（2009）[4]提出的變系數模型KLASSO估計，在KLASSO估計實際計算中，分別應用ESCV方法與BIC方法進行調節參數的選擇，并對比不同選參方法對模型穩定性的影響，而模型穩定性主要從模型預測誤差、模型大小和顯著性變量百分比上來進行比較。

本文雖然找到了一種能夠提高變系數模型穩定性的方法，但同樣存在很多問題：首先對于模型穩定性統計學上還沒有給出標準的定義，本文只能直觀地從預測誤差、選擇變量等方面來衡量模型是否穩定；其次變系數模型有很多估計方法，而此次研究僅限于KLASSO估計，在其他變系數模型估計方法下，ESCV方法是否能夠比BIC方法表現好還有待進一步研究；最后數據在低維情形時，多次抽樣情況下，雖然ESCV方法的變量選擇穩定性要優于BIC，但是BIC方法的平均預測誤差要小于ESCV，ESCV方法可能存在總是漏選某個重要變量的情況。

參考文獻：

[1]Hastie T,Tibshirani R.Varying Coefficient Models[J].Journal of Royal Statistical Society:Series B,1993,(55).

[2]Fan J,Zhang W.Statistical Estimation in Varying Coefficient Models[J].Journal of the American Statistical Association,1999,(27).

[3]Chiang C,Rice J A,Wu C O.Smoothing Spline Estimation for Varying Coefficient Models With Repeatedly Measured Dependent Variables[J].Journal of American Statistical Association,2001,(96).

[4]Wang H,Xia Y.Shrinkage Estimation of the Varying Coefficient Model[J].Journal of the American Statistical Association,2009,(104).

[5]Lim C,Yu B.Estimation Stability With Cross Validation(ESCV)[J].arXiv,2013,(1303).

[6]Yuan M,Lin Y.Model Selection and Estimation in Regression With Grouped Variables[J].Journal of the Royal Statistical Society:Series B,2006,(68).

[7]Yu B.Stability[J].Bernoulli,2013,19(4).

[8]Breiman L.Heuristics of Instability and Stabilization in Model Selection[J].Annals of Statistics,1996,(24).

[9]Fan J,Huang T.Profile Likelihood Inferences on Semiparametric Varying Coefficient Partially Linear Models[J].Bernoulli,2005,(11).