基于R語言的互補雙對數模型分析

戴建國 楊劍紅

摘 要： 針對分類數據模型的誤差不滿足正態分布，Logistic分布等常用分布，而是滿足極值分布這類特殊情況時，采用互補雙對數模型進行分析。以R語言作為分析工具，并將其與幾種常用模型進行對比分析，結果表明互補雙對數模型優于其他模型。

關鍵詞： 互補雙對數模型； R語言； 有序分類

中圖分類號： O 212.1 文獻標志碼： A 文章編號： 1671-2153（2017）04-0087-03

1 問題提出

大數據時代，數據挖掘與數據分析過程中遇到的數據類型主要有連續型，離散型，以及混合型，其涉及的方法眾多，如機器學習，深度學習等。當然，對于離散型數據，回歸模型也是比較流行的統計分析方法之一，而且在多數情況下，因變量有可能是有序二分類或多分類的離散型變量，在醫療衛生統計或社會統計尤為常見，如對某事的評價可能分為不好、好、非常好這樣三個等級。對于這類數據常用的模型有Logistic回歸模型（包括累積Logistic，鄰近Logistic等模型），Probit回歸模型，但它們分別是假設模型誤差分布是Logistic分布和正態分布的，而有些情況誤差并不滿足這樣的假設，其誤差分布是非對稱的極值分布或稱Gumbel 分布[1]。如某年某河面最高水位這樣一個極值，這時需要尋求一個更為有效的模型，此時采用互補雙對數線性模型與其他幾種模型進行對比分析。并且近年來R語言已成為當下最流行做數據分析和挖掘的工具之一，它不僅有豐富的函數包，而且能與Hadoop，Python等軟件結合，使得成為重要的數據分析工具。下面基于R語言用互補雙對數模型對一組實際數據進行對比分析。

2 互補雙對數線性模型[1]

該模型的一般形式為

log{-log[1-P（Y≤j）]}=αj+β'X， （1）

從而可得求概率的表達式為

式中：αj為解釋變量不能解釋的部分；β'為參數構成的向量。當只有一個自變量，并且當X=xi（i=1，2）時，有性質：

P（Y>j|x1）=P（Y>j|x2）exp[β（x1-x2）]。

對于各自變量的回歸系數通常采用極大似然法來估計。由于互補雙對數函數是嚴格凹的，海森矩陣便是負的，因此通過牛頓—拉夫森方法[2]得到參數估計，對參數進行檢驗，模型擬合優度檢驗[3]，這些都可以由R語言的summary（）函數得到。

3 實例分析與R程序

下面實例是研究膚色和性別對壽命的影響，其中壽命作為因變量分為5個類別（0～20歲、20～40歲、40～50歲、50～60歲、>60歲），分別記為1，2，3，4，5；性別和膚色作為兩個自變量，記為X，Z。對于自變量X包含男性、女性兩個類別，分別用0，1表示；自變量Z包含白人、黑人兩個類別，同樣用0，1表示。且每行的總數均為1000，數據見表1。由表1可以看出，在壽命的給類別上分布是很不平衡的，如比較壽命分別對應0～20和>60的兩列，下面給出幾種模型的對比分析。

3.1 建立數據文件

首先將表1中數據的每一列構成因子，然后用cbind（）函數將其按列合并，最后用as.data.frame（）將合并后的數據轉化成數據框，并命名為w，程序如下：

>gender<-c（0， 0， 1， 1）； race<-c（0， 1， 0， 1）；y1<-c（13， 26， 9， 18）；y2<-c（28， 49， 13， 24）；y3<-c（32， 56， 19， 37）；>y4<-c（122， 201， 80， 129）； y5<-c（805， 668， 879， 792）

>w<-as.data.frame（cbind（gender，race，y1，y2，y3，y4，y5））

3.2 安裝和載入VGAM包

在數據文件建立后，需要安裝和載入所需要的函數包，以下用一種簡單有效的方法載入和安裝所需的VGAM包。

>install.package（“VGAM”）； library（VGAM）

3.3 模型的構建及預測

在完成前面的準備工作后，使用VGAM包中的vglm（）函數對表中的數據用上述模型和其他兩種常用模型進行分析及預測。首先構建互補雙對數模型，模型的表達式已在式（1）給出。以下代碼中設置link = cloglog表示互補雙對數，將構建的模型記為CLL_fit，可用summary（CLL_fit）得到模型參數的估計和檢驗，并用predict（）所構建模型對數據進行預測，用round（）結果保留一位小數，預測結果見表1中括號內的數據。相關代碼如下：

>CLL_fit<-vglm（cbind（y1，y2，y3，y4，y5）～gender+race，family = cumulative（link = cloglog，parallel = T），data = w）

>summary（CLL_fit）

>fitted<-predict（CLL_fit，type = "response"）；round（fitted*1000，1）

3.3.1 鄰近logit回歸模型

該模型的一般形式為

式中：c為因變量的類別數；X為解釋變量構成的向量；πj為類別j在給定解釋變量下對應的概率；βj'為對應的回歸系數，該模型的假設是模型誤差服從Logistic分布。同樣使用VGAM包中的vglm（）函數對表中的數據進行分析及預測，此時需設置family=acat（），將構建后的模型記為AL_fit，由predict（）命令得到預測結果見表2。代碼如下：

>AL_fit<-vglm（cbind（y1，y2，y3，y4，y5）～gender+race，family=acat（reverse=TRUE，parallel=FALSE）， data=w）

>fitted<-predict（AL_fit，type ="response"）；round（fitted*1000，1）

3.3.2 累積Probit模型

模型的一般形式為

P（Y≤j）=？椎（αj+β'X） （j=1，…，c-1）， （3）

式中：c為因變量的類別數；X為解釋變量構成的向量；β'為對應的回歸系數；？椎為正態分布的分布函數；P（Y≤j）為累積概率，該模型的假設是模型誤差服從正態分布。可同互補雙對數模型一樣用代碼構建模型，只需將連接函數設為probit，即link=probit，同樣，得到的擬合值如表2所示，代碼如下：

>fit.cprobit<-vglm（cbind（y1，y2，y3，y4，y5）～gender+race，family=cumulative（link=probit，parallel=TRUE）， data=w）

>fitted<-predict（fit.cprobit，type = "response"）；round（fitted*1000，1）

3.4 模型擬合評價

本文分別用BIC和AIC準則[4-5]對上述模型進行評價，對于具體的公式及證明在此不討論，有興趣者可參考相關文獻。當模型的BIC和AIC值越小，說明模型擬合的越好。下面分別給出代碼計算各模型的AIC和BIC值，結果如表3所示，代碼如下：

>AIC（CLL_fit）；BIC（CLL_fit）；AIC（AL_fit）；BIC（AL_fit）；AIC（fit.cprobit）；BIC（fit.cprobit）

綜合上面代碼得到表2和表3的結果，比較發現：表1中模型（1）的擬合數據比表2中另外兩種模型的擬合值更接近真實值，尤其是在年齡大于60的那一列，而且比較表3中三個模型的AIC和BIC值可知模型（1）的值明顯小于另外兩種模型，因此，表明對于這樣一類特殊數據，互補雙對數模型具有明顯的優勢。

由互補雙對數模型擬合的輸出結果，可以得到4個回歸方程分別為

log{-log[1-P（Y≤1）]}=-4.21-0.54X+0.61Z，

log{-log[1-P（Y≤2）]}=-3.19-0.54X+0.61Z，

log{-log[1-P（Y≤3）]}=-2.58-0.54X+0.61Z，

log{-log[1-P（Y≤4）]}=-1.52-0.54X+0.61Z。

由互補雙對數模型給出的性質有P（Y>j|X=1，Z）=P（Y>j|X=0，Z）exp（0.54），即在給定膚色的情況下，男性壽命長于固定時間的比例是女性比例的exp（0.54）=1.71次方；在給定性別時，黑人壽命長于某一固定時間的比例是白人比例的exp（0.61）=1.84次方，并且可以得到對于壽命長于某一固定時間，黑種男人的比例是白種女人比例的3.55次方。

4 結束語

對于互補雙對數線性模型主要用于誤差分布是極值分布的情況，如上例所示，其可以對自變量進行選擇，也可以估計概率的大小，從而得到擬合值，當因變量是二分變量時，可以用于預測分類，在文獻[6-7]中給出了相關實例分析，說明了在誤差分布非對稱時，互補雙對數模型體現一定的優勢，而且R語言能夠很有效，簡潔完成分析過程，從而也說明R語言是一款很好的統計分析工具。

參考文獻：

[1] AGRESTI，ALAN. Analysis Of Ordinal Categorical Data，Second Edition[M]. Canada：Wiley，2010.

[2] BOLARINWA，OLUBUSOYE O E. On the sensitivity of probit and complementary log-log models to violated tolerance[J]. 2013，12：351-360.

[3] AGRESTI，ALAN. An Introduction to Categorical Data Analysis[M]. Wiley，New York，Duxbury Press，1996.

[4] BURNHAM K P，ANDERSON D R. Multimodel Inference understanding AIC and BIC in model selection[J]. Sociological Methods & Research，2004，33（33）：261-304.

[5] KUHA J. AIC and BIC comparisons of assumptions and performance[J]. Sociological Methods & Research，2004，33（2）：188-229.

[6] 肖珉，周宗放. 基于互補雙對數模型的企業集團財務困境預測研究[C]. 中國災害防御協會風險分析專業委員會第四屆年會論文集，2010：976-981.

[7] 史芳芳，楊海珍，徐昭，等. 基于互補雙對數模型的國際資本流入風險影響因素研究[J].數學的實踐與認識，2016，46（10）：91-96.