熊笛，何幼樺

（上海大學(xué)理學(xué)院，上海 200444）

半?yún)?shù)順序變量回歸模型

熊笛，何幼樺

（上海大學(xué)理學(xué)院，上海 200444）

在比例優(yōu)勢(shì)模型基礎(chǔ)上對(duì)順序變量回歸模型作更一般的推廣,建立了半?yún)?shù)順序變量回歸模型,構(gòu)造了模型中的線性和非線性部分的估計(jì)量,并證明了該估計(jì)量的弱相合性.通過數(shù)值模擬,考察了不同樣本容量下半?yún)?shù)順序變量回歸的判斷正確率和回歸函數(shù)的均方誤差.實(shí)驗(yàn)結(jié)果表明：半?yún)?shù)順序回歸模型在小樣本情況下仍具有較高精度,并且在實(shí)驗(yàn)點(diǎn)處的重復(fù)次數(shù)相對(duì)于觀察點(diǎn)個(gè)數(shù)對(duì)精度影響更大.通過對(duì)糧食預(yù)警實(shí)例的計(jì)算表明,半?yún)?shù)順序回歸模型較比例優(yōu)勢(shì)線性模型具有更好的外推效果.

比例優(yōu)勢(shì)模型；順序變量回歸；半?yún)?shù)回歸

順序變量是用于說明事物有序類別或者有序等級(jí)的一類以順序數(shù)據(jù)作為具體表現(xiàn)的變量，也是0，1二分類變量的擴(kuò)展，廣泛出現(xiàn)在各應(yīng)用領(lǐng)域的統(tǒng)計(jì)模型中.1959年Duncan［1］根據(jù)非相關(guān)選擇項(xiàng)的獨(dú)立性（independence from irrelevant alternation）特性首次提出了logit模型，該模型也是最早的離散因變量回歸模型.Cox［2］提出了0，1二分類logit模型，并對(duì)二分類變量回歸的線性模型形式進(jìn)行了詳細(xì)的分析.此后，自變量為數(shù)值變量、響應(yīng)變量為順序變量的多分類順序變量回歸問題被深入討論與研究.

設(shè)響應(yīng)變量Y為K個(gè)類別的順序變量，通過順序值1,2,…,K表征響應(yīng)變量所歸屬的類別或等級(jí).若對(duì)于d維解釋變量X，響應(yīng)變量Y屬于第j類的概率為pj（x）=P（Y=j|x=x），j=1,2,…,K，那么，響應(yīng)變量屬于第j類的累積概率可以表示為

可以按概率γj（x）和1-γj（x）把K個(gè)等級(jí)分成｛1,2,…,j｝和｛j+1,j+2,…,K｝兩類，在此基礎(chǔ)上以γj（x）/（1-γj（x））表示順序變量Y所屬級(jí)別或等級(jí)不大于j（Y 6 j）時(shí)的優(yōu)勢(shì)比.把｛1,2,…,j｝和｛j+1,j+2,…,K｝兩類視為一種兩個(gè)類別的數(shù)據(jù)形式，在Cox［2］的二分類logit模型基礎(chǔ)上得到了更一般的多類別模型：

或以線性模型形式

式中，θj=lgκj為第j個(gè)等級(jí)的基準(zhǔn)線.模型（1）為McCullagh［3］在1980年提出的比例優(yōu)勢(shì)模型（有序logit模型），它是二分類logit模型的擴(kuò)展.如果當(dāng)響應(yīng)變量只有兩類時(shí)，則比例優(yōu)勢(shì)模型就是一個(gè)線性logit模型.

目前，比例優(yōu)勢(shì)模型成為順序變量回歸的主流方法之一.Pettitt［4］把比例優(yōu)勢(shì)模型應(yīng)用于生存數(shù)據(jù)研究中，并對(duì)比例優(yōu)勢(shì)模型的估計(jì)方法進(jìn)行了討論；Murphy等［5］研究了極大似然估計(jì)方法在比例優(yōu)勢(shì)模型上的運(yùn)用；Ibrahim等［6］對(duì)該模型在貝葉斯變量上的選擇方式進(jìn)行了分析；Lang［7］導(dǎo)出了順序回歸模型中混合連接函數(shù)的貝葉斯估計(jì)方法；Lam等［8］提出右刪失數(shù)據(jù)的比例優(yōu)勢(shì)模型的極大似然估計(jì)方法.國內(nèi)對(duì)該類模型的應(yīng)用問題也有較多的研究，如文獻(xiàn)［9-10］對(duì)二分類logit模型應(yīng)用的正確性進(jìn)行了探討，并將比例優(yōu)勢(shì)模型應(yīng)用于航空領(lǐng)域的加速壽命試驗(yàn)［11］、醫(yī)藥等［12］研究領(lǐng)域.

在很多實(shí)際問題中，自變量與因變量之間并不完全滿足線性關(guān)系，因此僅用線性回歸模型不能準(zhǔn)確地描述所討論的問題.在20世紀(jì)80年代中期，Engle等［13］提出了半?yún)?shù)（或稱為偏線性）模型：

式中，因變量u受到一些控制變量y∈Rp和x∈Rq以及隨機(jī)擾動(dòng)ε的影響，并且x對(duì)u的影響是線性的；f（·）為未知函數(shù)；β=（β1,β2,…,βq）T為未知參數(shù)；ε|（x,y）～（0,σ2）是隨機(jī)擾動(dòng)的，在很多情況下可以假設(shè)它是正態(tài)的.

因此，順序變量與影響因子之間的關(guān)系可以有更精細(xì)的描述.本研究考慮用一個(gè)連續(xù)非線性函數(shù)代替式（1）中的線性部分：

等式兩邊同時(shí)取對(duì)數(shù)，得到半?yún)?shù)順序變量回歸模型：

式中，θj=lgκj為第j個(gè)等級(jí)的基準(zhǔn)線，且θj＜θj+1，不失一般性，可設(shè)θ1=0.

1 半?yún)?shù)順序回歸模型的估計(jì)

令X為d維解釋變量，Y是順序響應(yīng)變量，Y=j（j=1,2,…,K）表示響應(yīng)變量歸屬于K個(gè)順序類別中的第j個(gè)類別.

針對(duì)樣本觀察值（xi,Yi）,Yi=Y（xi）∈｛1,2,…,K｝,i=1,2,…,n.記Rij=Rj（xi）= #｛xi|Y（xi）6 j｝，通過加權(quán)光滑方式得到γj（xi）的估計(jì)量

（2）ω（-x）=ω（x）；

那么，lg（γj（xi）/（1-γj（xi）））的估計(jì)則可以表示為則半?yún)?shù)順序回歸模型（5）的樣本模型為

對(duì)于擾動(dòng)項(xiàng)，假設(shè)在給定j的情況下，｛εij｝～（0,σ2）.

定理1 對(duì)于半?yún)?shù)順序回歸模型（5）中各θj的最佳線性無偏估計(jì)量為

證明 根據(jù)式（7），有

對(duì)于給定j，令ηj=εij-εi1～（0,Var（ηj）），根據(jù)最小二乘估計(jì)的基本結(jié)論，式（9）中各θj的最佳線性無偏估計(jì)為

注意到RijRi（j+1），且對(duì)于每一個(gè)j=1,2,…,K-1，至少存在一個(gè)i使得Rij＜Ri（j+1）（否則第j類和第j+1類可合并為一類），故

于是有

在根據(jù)定理1得到θj的估計(jì)量后，記

于是有關(guān)于f（x）的非參數(shù)回歸樣本模型：

將參數(shù)和非參數(shù)部分的估計(jì)代回式（5）得到

證明 首先，對(duì)于每一個(gè)j=1,2,…,K-1，根據(jù)大數(shù)定律，式（8）中的按概率收斂到θj，即

其次，根據(jù)局部線性回歸估計(jì)的殘差定理［14］：

得到

當(dāng)定理?xiàng)l件得到滿足時(shí)，有

那么當(dāng)X=x時(shí)，響應(yīng)變量Y屬于第j個(gè)類別的隸屬概率的估計(jì)為

定義記

則稱

本研究采用判斷正確率（correct rate，CR）和均方誤差（mean squared error，MSE）兩個(gè)指標(biāo)作為估計(jì)結(jié)果優(yōu)良性的評(píng)判標(biāo)準(zhǔn).

2 數(shù)值模擬

設(shè)定一個(gè)函數(shù)作為原始模型進(jìn)行數(shù)值模擬，隨機(jī)產(chǎn)生一系列實(shí)驗(yàn)樣本點(diǎn)，然后利用這些樣本數(shù)據(jù)對(duì)半?yún)?shù)順序回歸進(jìn)行估計(jì).

重復(fù)N次如下試驗(yàn)：

步驟1 隨機(jī)產(chǎn)生n個(gè)解釋變量值x1,x2,…,xn，在每個(gè)xi處對(duì)Y重復(fù)m次觀察，共產(chǎn)生m×n組樣本；

步驟2 利用半?yún)?shù)順序回歸模型中的式（12），計(jì)算

根據(jù)上述步驟得到的結(jié)果計(jì)算出半?yún)?shù)順序回歸模型判斷正確率和均方誤差兩個(gè)指標(biāo)：

（1）固定實(shí)驗(yàn)觀察點(diǎn)個(gè)數(shù)n=30，改變?cè)诿總€(gè)xi處對(duì)Y重復(fù)觀察次數(shù)m=1,2,5，將得到的判斷正確率CR（x）和均方誤差MSE（x）的數(shù)值進(jìn)行比較，結(jié)果如圖1所示.

（2）固定每個(gè)xi處重復(fù)次數(shù)m=1，改變實(shí)驗(yàn)觀察點(diǎn)個(gè)數(shù)n=30,100，將得到的判斷正確率CR（x）和均方誤差MSE（x）數(shù)值進(jìn)行比較（見圖2）.

（3）固定實(shí)驗(yàn)樣本容量m×n=60，討論n=60，m=1，n=30,m=2以及n=20，m=3這3種情況，將得到的判斷正確率CR和均方誤差MSE數(shù)值進(jìn)行比較（見圖3）.

實(shí)驗(yàn)結(jié)果表明：

（1）當(dāng)實(shí)驗(yàn)觀察點(diǎn)個(gè)數(shù)n不變時(shí)，每個(gè)xi處重復(fù)次數(shù)（m）越多，判斷正確率就越高，回歸函數(shù)的均方誤差越??；

（2）當(dāng)每個(gè)xi處重復(fù)次數(shù)m不變時(shí)，實(shí)驗(yàn)觀察點(diǎn)個(gè)數(shù)（n）越多，判斷正確率越高，回歸函數(shù)的均方誤差越??；

（3）當(dāng)實(shí)驗(yàn)樣本容量m×n不變時(shí)，在每個(gè)xi處重復(fù)次數(shù)（m）對(duì)判斷正確率和回歸函數(shù)的均方誤差兩個(gè)指標(biāo)的影響比觀察點(diǎn)個(gè)數(shù)n對(duì)它們的影響相對(duì)更大，即在樣本容量相同的情況下，實(shí)驗(yàn)點(diǎn)xi處重復(fù)次數(shù)（m）越多，判斷正確率就越高，回歸函數(shù)的均方誤差越??；

（4）當(dāng)x靠近樣本集邊界時(shí)，判斷正確率和回歸函數(shù)的均方誤差兩個(gè)指標(biāo)均不如x位于樣本集內(nèi)部時(shí)的情形，此時(shí)判斷正確率相對(duì)更低，回歸函數(shù)的均方誤差相對(duì)較大.

圖1 判斷正確率CR和均方誤差MSE（n=30,m=1,2,5）Fig.1 CR and MSE（n=30,m=1,2,5）

圖2 判斷正確率CR和均方誤差MSE（m=1,n=30,100）Fig.2 CR and MSE（m=1,n=30,100）

圖3 判斷正確率CR和均方誤差MSE（n=60,m=1和n=30,m=2以及n=20,m=3）Fig.3 CR and MSE（n=60,m=1 and n=30,m=2 and n=20,m=3）

3 半?yún)?shù)順序變量回歸模型的應(yīng)用

作為一個(gè)應(yīng)用實(shí)例，對(duì)糧食預(yù)警問題建立一個(gè)半?yún)?shù)順序變量回歸模型，將影響糧食價(jià)格波動(dòng)的警源作為解釋變量，對(duì)糧食警情等級(jí)進(jìn)行預(yù)報(bào).

將1978—2012年的糧食價(jià)格作為研究對(duì)象，取其價(jià)格相對(duì)變動(dòng)作為糧食價(jià)格警情的指標(biāo).以當(dāng)年糧食播種面積增長率、當(dāng)年糧食畝產(chǎn)增長率、當(dāng)年受災(zāi)面積增長率作為影響當(dāng)年糧食價(jià)格波動(dòng)的警源（樣本數(shù)據(jù)見附錄）.

對(duì)于警情則采用多數(shù)原則，即把計(jì)算得到的糧食波動(dòng)率數(shù)值從小到大排列，從第一個(gè)數(shù)據(jù)開始，將占總數(shù)2/3的數(shù)據(jù)作為安全警限，即無警警限，依次在剩下的波動(dòng)率數(shù)據(jù)中劃分輕警、中警、重警、巨警，根據(jù)實(shí)際劃分情況，本研究將余下4個(gè)警限按照等距劃分并依次將這5個(gè)警級(jí)命名為警級(jí)1，2，3，4，5.因此依據(jù)附錄中糧食價(jià)格波動(dòng)情況，結(jié)合多數(shù)原則將糧食價(jià)格警度進(jìn)行劃分（見表1）［15］.

表1 多數(shù)原則下的糧食價(jià)格警度警限Table 1 Grain price warning degree under principle of majority

首先取1978—2012年數(shù)據(jù)作為樣本點(diǎn)用半?yún)?shù)順序回歸模型進(jìn)行內(nèi)插檢驗(yàn)，其判斷正確率為100%.再以1978—2007年數(shù)據(jù)為訓(xùn)練樣本，用比例優(yōu)勢(shì)線性模型和半?yún)?shù)順序回歸模型對(duì)2008—2012年糧食警級(jí)進(jìn)行外推，判斷正確率分別為60%和100%（見表2和3）.

表2 順序回歸線性模型外推糧食價(jià)格警級(jí)結(jié)果Table 2 Extrapolation of grain price warning degree using ordinal regression linear model

表3 半?yún)?shù)順序回歸模型外推糧食價(jià)格警級(jí)結(jié)果Table 3 Extrapolation of grain price warning degree using semi-parametric ordinal regression model

4 結(jié)束語

本研究所建立的半?yún)?shù)順序變量回歸模型是在傳統(tǒng)的線性順序變量回歸模型基礎(chǔ)上考慮了非線性部分，擴(kuò)展了模型的實(shí)際應(yīng)用范圍.同時(shí)，本研究通過半?yún)?shù)順序回歸模型對(duì)糧食價(jià)格進(jìn)行了預(yù)警，從預(yù)警結(jié)果來看半?yún)?shù)順序回歸模型具有很好的預(yù)測(cè)效果.后續(xù)工作將從以下兩個(gè)方向進(jìn)行：①對(duì)于γj（x）估計(jì)的改進(jìn).當(dāng)在每一個(gè)x處重復(fù)觀察一次或很少時(shí)，所采用的估計(jì)方法（6）會(huì)有較大的誤差，這個(gè)誤差直接影響了模型估計(jì)的最終效果.②研究模型（4）基準(zhǔn)量κj的一般化問題.如假設(shè)κj依賴于其他外生變量或與解釋變量X有一定相關(guān)性，則整個(gè)估計(jì)方法會(huì)有較大的改變，模型的適用范圍可以更大.

附錄

表4中各項(xiàng)波動(dòng)率、增長率是根據(jù)歷年的統(tǒng)計(jì)年鑒（http：//data.stats.gov.cn/workspace/ index？m=hgnd）計(jì)算得到的.

表4 1978—2012年糧食數(shù)據(jù)表Table 4 1978—2012 grain's data%

［1］DUNCAN L R.Individual choice behavior：a theoretical analysis［M］.New York：John Wiley& Sons,1959.

［2］COx D R.The analysis of multivariate binary data［J］.Royal Statistical Society,1972,21（2）：113-120.

［3］MCCULLAGH P.Regression models for ordinal data［J］.Journal of the Royal Statistical Society, 1980,42（2）：109-142.

［4］PETTITT A N.Inference for the linear model using a likelihood based on ranks［J］.Journal of the Royal Statistical Society,1982,44（2）：234-243.

［5］MURPHY S A,ROSSINI A J.Maximum likelihood estimation in the proportional odds model［J］. Journal of the American Statistical Association,1997,92（439）：968-976.

［6］IBRAHIM J G,CHEN M H,MACEACHERN S N.Bayesian variable selection for proportional hazards models［J］.The Canadian Journal of Statists,1999,27（4）：701-717.

［7］LANG J B.Bayesian ordinal and binary regression models with a parametric family of mixture links［J］.Computational Statistics&Data Analysis,1999,31（1）：59-87.

［8］LAM K F,LEUNG T L.Marginal likelihood estimation for proportional odds models with right censored data［J］.Lifetime Data Analysis,2001,7（1）：39-54.

［9］馮國雙,陳景武,周春蓮.logistic回歸應(yīng)用中容易忽視的幾個(gè)問題［J］.中華流行病學(xué)雜志,2004, 25（6）：544-545.

［10］趙宇東,劉嶸,劉延齡.多元logistic回歸的共線性分析［J］.中國衛(wèi)生統(tǒng)計(jì),2001,17（5）：259-261.

［11］黃婷婷,姜同敏.基于比例危險(xiǎn)-比例優(yōu)勢(shì)模型的加速壽命試驗(yàn)設(shè)計(jì)［J］.北京航空航天大學(xué)學(xué)報(bào), 2010,36（5）：570-579.

［12］唐俐玲,翟曉紅.累積比數(shù)logit模型在有序資料中的正確應(yīng)用［J］.徐州醫(yī)學(xué)院學(xué)報(bào),2010,30（9）：577-579.

［13］ENGLE R F,GRANGER C W J,RICE J,et al.Semiparametric estimates of the relationship between weather and electricity sales［J］.Journal of the American Statistical Association,1986, 81（394）：310-320.

［14］RUPPERT D,WAND M P.Multivariate locally weighted least squares regression［J］.The Annals of Statistics,1994,22（3）：1346-1370.

［15］吳璇.中國糧食價(jià)格預(yù)警系統(tǒng)研究［D］.北京：中國農(nóng)業(yè)大學(xué),2003.

Semi-parametric ordinal variable regression model

XIONG Di，HE Youhua
（College of Sciences，Shanghai University，Shanghai 200444，China）

Based on a proportional odds model,the ordinal variable regression model is generalized,a semi-parametric ordinal regression model is established,and consistency of the estimators both in linear and nonlinear parts are proved in this paper.Simulation is conducted to analyze the correct rate and mean square error in the semi-parametric ordinal variable regression model with different sample sizes.The result shows that the semi-parametric ordinal regression model has high accuracy even with small samples.Compared to the number of observation points,the repeat number of experimental points has greater influence on accuracy.Calculation of the grain price warning problem shows that the semi-parametric ordinal regression model provides better extrapolation results than the proportional odds model.

proportional odds model；ordinal variable regression；semi-parametric regression

O 212

A

1007-2861（2016）04-0477-09

10.3969/j.issn.1007-2861.2014.04.010

2014-11-21

國家自然科學(xué)基金資助項(xiàng)目（11371242）

何幼樺（1960—），男，教授，博士，研究方向?yàn)楦怕式y(tǒng)計(jì).E-mail：heyouhua@shu.edu.cn