混合函數(shù)型數(shù)據(jù)下Logistic回歸模型的懲罰估計(jì)

袁曉惠，金宛霖，曹儒雅

(長(zhǎng)春工業(yè)大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院， 吉林 長(zhǎng)春 130012)

0 引 言

函數(shù)型數(shù)據(jù)最早由Ramsay J O[1]提出，常以曲線或者圖像形式出現(xiàn)。由于其無限維的特性，傳統(tǒng)的數(shù)據(jù)處理方法具有局限性，需要探索新的方法和模型對(duì)函數(shù)型數(shù)據(jù)進(jìn)行分析，其中最經(jīng)典的是函數(shù)型線性回歸模型。Ramsay J O等[2]使用函數(shù)型線性模型對(duì)加拿大的溫度和降雨量進(jìn)行了實(shí)證分析。劉鋒等[3]利用低階基函數(shù)的線性組合研究了函數(shù)型線性模型在肉類光譜數(shù)據(jù)上的應(yīng)用。Lin Z等[4]提出一種稱為“fSCAD”(functional SCAD)的方法?；贐-樣條，結(jié)合粗糙度懲罰和稀疏懲罰得到回歸參數(shù)的最小二乘估計(jì)。徐夢(mèng)佳[5]構(gòu)造帶懲罰項(xiàng)的函數(shù)型多元線性回歸模型，并將該模型應(yīng)用到江西省的空氣質(zhì)量分析中。

對(duì)于協(xié)變量是函數(shù)，響應(yīng)變量是二值標(biāo)量的情形，經(jīng)典的統(tǒng)計(jì)模型為函數(shù)型Logistic回歸模型。Cardot H等[6]分析了標(biāo)量響應(yīng)和函數(shù)預(yù)測(cè)之間的聯(lián)系，在光滑性假設(shè)的條件下，利用樣條逼近的懲罰似然來估計(jì)模型中的系數(shù)函數(shù)。Liu Y等[7]提出一種稀疏的函數(shù)型Logistic模型來預(yù)測(cè)抑郁癥患者的治療緩解狀態(tài)，并將函數(shù)型數(shù)據(jù)分析應(yīng)用在重度抑郁這樣的疾病研究。王惠文等[8]研究含有函數(shù)型協(xié)變量的Logistic模型，結(jié)合B-樣條基函數(shù)，得到該模型回歸參數(shù)的極大似然估計(jì)。

對(duì)于響應(yīng)變量是標(biāo)量，協(xié)變量是混合函數(shù)型數(shù)據(jù)的情形，最常見統(tǒng)計(jì)的模型是部分函數(shù)型線性模型。程麗娟[9]對(duì)部分函數(shù)型線性模型在上證指數(shù)中的應(yīng)用進(jìn)行了研究。張雪[10]采用預(yù)平滑的方法得到部分函數(shù)型線性模型中系數(shù)函數(shù)的估計(jì)。丁輝[11]研究了部分函數(shù)型線性模型的局部稀疏估計(jì)，采用最小二乘法，結(jié)合粗糙度懲罰和fSCAD方法分別估計(jì)函數(shù)型系數(shù)和標(biāo)量型系數(shù)。

梳理文獻(xiàn)發(fā)現(xiàn)，盡管已有文獻(xiàn)給出了各種函數(shù)型回歸模型的研究和應(yīng)用，但對(duì)協(xié)變量同時(shí)包含函數(shù)型數(shù)據(jù)和數(shù)值型數(shù)據(jù)的Logistic模型研究還很少。文中針對(duì)響應(yīng)變量為二分類[12-13]，協(xié)變量為混合函數(shù)型數(shù)據(jù)的情形，構(gòu)造Logistic回歸模型，并求得回歸參數(shù)的懲罰似然估計(jì)，最后將該模型應(yīng)用于實(shí)例數(shù)據(jù)中。

1 函數(shù)型Logistic回歸模型介紹

數(shù)據(jù)類型為{Yi,Xi(t),Zi}，i=1,2,…,n,響應(yīng)變量Yi∈{0,1}，協(xié)變量為混合函數(shù)型數(shù)據(jù)的Logistic回歸模型為

P(Yi=1|Xi(t),Zi)=

(1)

i=1,2,…,n,

式中：μ----截距項(xiàng);

β(t)----系數(shù)函數(shù);

X1(t),X2(t),…,Xn(t)----定義在區(qū)間[0,T]上的函數(shù)型協(xié)變量;

Zi----p維向量,Zi=(Zi1,Zi2,…,Zip)T;

α----Zi的系數(shù),α=(α1,α2,…,αp)T。

以Logit為鏈接函數(shù)，則式(1)可表示為

logitP(Yi=1|Xi(t),Zi)=

i=1,2,…,n。

該模型描述了標(biāo)量型響應(yīng)變量和混合型協(xié)變量之間的關(guān)系，根據(jù)式(1)求得該模型的對(duì)數(shù)似然函數(shù)為

(2)

2 函數(shù)型Logistic回歸模型的懲罰似然估計(jì)

由于函數(shù)型數(shù)據(jù)是一種新型高維數(shù)據(jù)，需要對(duì)其進(jìn)行低維表示，常用的方法是選取基函數(shù)對(duì)其進(jìn)行基展開。在函數(shù)型數(shù)據(jù)分析中經(jīng)常使用的基函數(shù)有B-樣條基函數(shù)、傅里葉基函數(shù)、小波基函數(shù)等。文中選取B-樣條基函數(shù)對(duì)系數(shù)函數(shù)β(t)進(jìn)行基展開為

(3)

式中：B(t)----定義在區(qū)間[0,T]上的d次B-樣條基，B(t)=(B1(t),B2(t),…,BM+d(t))T;

b----對(duì)應(yīng)的系數(shù)向量,b=(b1,b2,…,bM+d)T。

在區(qū)間[0,T]上設(shè)置M+1個(gè)等距結(jié)點(diǎn)將區(qū)間分成M個(gè)等長(zhǎng)度的子區(qū)間,其中等距結(jié)點(diǎn)為0=t0

令U是一個(gè)n×(M+d)的矩陣，其中

且

U=(U1,U2,…,Un)T,i=1,2,…,n。

結(jié)合式(2)和式(3)，對(duì)數(shù)似然函數(shù)可改寫為

(4)

2.1 粗糙度懲罰

令V是一個(gè)(M+d)×(M+d)的矩陣，且

1≤i,j≤M+d，則粗糙度懲罰的懲罰函數(shù)可改寫為

γ‖Dmβ‖2=γbTVb。

(5)

2.2 函數(shù)型SCAD (fSCAD)

(6)

其中

令Wj是一個(gè)(M+d)×(M+d)的矩陣，具體形式為

2.3 結(jié)合粗糙度懲罰和稀疏懲罰的懲罰似然估計(jì)

采用Fan J等[15]提出局部平方近似(LQA)方法對(duì)式(6)中的稀疏懲罰項(xiàng)進(jìn)行近似。即給定初始值β(0)，有

G(β(0)),

(7)

其中

令

(8)

則有

(9)

nγbTVb-nbTW(0)b-

G(β(0))-nαTΣα。

可見β(t)的優(yōu)化問題等價(jià)于b的優(yōu)化問題。由于G(β(0))不依賴于β，因此其對(duì)b的優(yōu)化沒有影響，懲罰似然函數(shù)可改寫為

nγbTVb-nbTW(0)b-nαTΣα。

(10)

分別求得式(10)中b和α的一階導(dǎo)數(shù)和二階導(dǎo)數(shù)為：

2nγVb-2nW(0)b,

(11)

2nγV-2nW(0),

2nΣα,

2nΣ。

1)通過下式獲得b和α的初始估計(jì)值

nγbTVb-nαTΣα;

2)通過式(11)在每一次迭代中更新b和α，同時(shí)也對(duì)W(0)和Σ進(jìn)行迭代更新

b(m+1)=b(m)-

α(m+1)=α(m)-

3 數(shù)值模擬

通過模擬研究提出SLoS估計(jì)方法的數(shù)值性能，數(shù)據(jù)由以下模型產(chǎn)生

logitP(Yi=1|Xi(t),Zi)=μ+

(12)

i=1,2,…,n。

協(xié)變量函數(shù)Xi(t)基于等式

Xi(t)=∑aijBj(t)

來生成，其中aij服從標(biāo)準(zhǔn)正態(tài)分布，每個(gè)Bj(t)是101個(gè)具有98個(gè)等距結(jié)點(diǎn)的5階B-樣條基函數(shù)。對(duì)于數(shù)值型協(xié)變量，令Z1,Z2,…,Z5為來自獨(dú)立同分布的正態(tài)分布的樣本，回歸系數(shù)α=c(0,2,0,1,0)，且Xi(t)和Z1,Z2,…,Z5獨(dú)立。

在模擬中，令μ=0，考慮既包含非空子區(qū)間，又包含空子區(qū)間的β(t)，

將上述假設(shè)代入模型(12)中，考慮了三種不同樣本容量N=100,450,1 000下的系數(shù)函數(shù)β(t)的估計(jì)。繪制樣本量N=1 000時(shí),β(t)的SLoS估計(jì)曲線與β(t)的真實(shí)曲線,如圖1所示。

圖1 β(t)的SLoS估計(jì)曲線

β(t)的估計(jì)效果由空子區(qū)間和非空子區(qū)間上的積分平方誤差的大小來決定，分別定義如下:

式中：l0----β(t)的空子區(qū)間長(zhǎng)度;

l1----β(t)的非空子區(qū)間長(zhǎng)度。

表1 β(t)的ISE0和ISE1

根據(jù)表1發(fā)現(xiàn)，β(t)的空子區(qū)間ISE0和非空子區(qū)間ISE1都有隨著樣本量的增大而減小的趨勢(shì)。結(jié)合圖1可知，SLoS方法可以較準(zhǔn)確識(shí)別空子區(qū)間和非空子區(qū)間，是一種比較有效的估計(jì)方法。

α估計(jì)值的表現(xiàn)通過標(biāo)準(zhǔn)差、偏差和均方誤差體現(xiàn)。將α在樣本量為N=100,450,1 000時(shí)的標(biāo)準(zhǔn)差、偏差和均方誤差，分別匯總見表2～表4。

表2 α的標(biāo)準(zhǔn)差(×10-3)

表3 α的偏差(×10-4)

表4 α的均方偏差(×10-6)

由上述表中可以明顯看出,隨著樣本量的增大，α估計(jì)的標(biāo)準(zhǔn)差、偏差和均方誤差都在減小。綜上可知，函數(shù)型Logistic模型對(duì)系數(shù)函數(shù)β(t)和回歸系數(shù)α都給出了良好的估計(jì)，證明該模型是有效的。

4 實(shí)證分析

數(shù)據(jù)來源于http://lib.stat.cmu.edu/datasets/tecator的Tecator數(shù)據(jù)集。該數(shù)據(jù)集包含以百分比為單位的215個(gè)碎肉樣本的脂肪含量、水含量和蛋白質(zhì)含量。光譜數(shù)據(jù)的每個(gè)樣本記錄了在波長(zhǎng)為[850 nm,1 050 nm]上的100個(gè)觀測(cè)通道下的光譜吸光記錄。將該數(shù)據(jù)集中脂肪含量高于28%的54塊碎肉判定為具有一定肥胖特征的樣本，脂肪含量低于28%的161塊碎肉判定為不具有肥胖特征的樣本。從兩類樣本中各隨機(jī)抽取30條光譜曲線如圖2所示。

圖2 兩類樣本的部分光譜曲線

首先對(duì)數(shù)據(jù)進(jìn)行預(yù)處理，然后使用函數(shù)型Logistic模型進(jìn)行實(shí)例研究

logitP(Yi=1|Xi(t),Zi)=μ+

(13)

式中:Yi=1----樣本具有一定的肥胖特征;

Xi(t)----光譜數(shù)據(jù);

Z1----水的含量;

Z2----蛋白質(zhì)的含量。

在對(duì)該實(shí)例進(jìn)行分析時(shí)，主要探究光譜數(shù)據(jù)的各段波長(zhǎng)范圍是否對(duì)判定肥胖特征有顯著影響。在對(duì)碎肉樣本是否具有肥胖特征的研究中使用提出的SLoS估計(jì)方法分析上述問題，得到β(t)的估計(jì)曲線如圖3所示。

圖3 β(t)的估計(jì)曲線

由圖3可以發(fā)現(xiàn)，波長(zhǎng)在[960 nm,980 nm]時(shí)，光譜數(shù)據(jù)對(duì)肥胖特征的判定沒有影響。

5 結(jié) 語

面向協(xié)變量為函數(shù)型和數(shù)值型混合數(shù)據(jù)的二分類問題，提出一種基于Logit變換的函數(shù)型 Logistic回歸模型。通過數(shù)值模擬驗(yàn)證了該方法的有效性，并將該模型應(yīng)用于Tecator數(shù)據(jù)集進(jìn)行實(shí)例分析。實(shí)證結(jié)果表明，波長(zhǎng)范圍大約在[960 nm,980 nm]時(shí)，光譜吸收率對(duì)肥胖特征的判定沒有影響。文中僅考慮了函數(shù)型協(xié)變量為單一變量的回歸模型,針對(duì)多元函數(shù)型協(xié)變量的函數(shù)型Logistic回歸模型可作為后續(xù)研究。