孫 興， 黃振生

(南京理工大學(xué) 理學(xué)院，南京 210094)

0 引 言

在醫(yī)療領(lǐng)域，糖尿病是備受關(guān)注的慢性病之一，也是全球嚴(yán)重的公共衛(wèi)生安全問題之一，它除了會(huì)給患者帶來痛苦之外，還會(huì)給家庭和社會(huì)帶來巨大的經(jīng)濟(jì)負(fù)擔(dān)和壓力。根據(jù)國際糖尿病聯(lián)盟(IDF)2017的調(diào)研結(jié)果顯示，全球共有約4.25億糖尿病患者，其中中國糖尿病患者占比超25%[1]。研究糖尿病與體質(zhì)測量數(shù)據(jù)之間的關(guān)系，可以更好地了解和預(yù)防糖尿病，因此具有較為重要的意義。研究的糖尿病數(shù)據(jù)集包含442個(gè)觀測樣本，其中響應(yīng)變量是糖尿病患者的定量測量數(shù)據(jù)，協(xié)變量分別為年齡(Age)、性別(Sex)、體質(zhì)指數(shù)(BMI)、平均血壓(BP)和6種血清指標(biāo)測量數(shù)據(jù)，分別以符號(hào)TC、LDL、HDL、TCH、LTG和GLU表示。

部分變系數(shù)單指標(biāo)模型(PVCSIMs)是一類重要的半?yún)?shù)模型，它不僅保留了非參數(shù)模型的特點(diǎn)，還能避免“維數(shù)災(zāi)禍”問題，因此是統(tǒng)計(jì)分析的重要工具，尤其在處理高維數(shù)據(jù)時(shí)十分有優(yōu)勢，模型的結(jié)構(gòu)如下：

Y=θT(U)Z+g(βTX)+ε

(1)

可以看出，該模型具有一般性，它不僅兼有變系數(shù)模型和單指標(biāo)模型的特點(diǎn)，還包含了許多其他重要的半?yún)?shù)模型作為特例：對于變系數(shù)部分，函數(shù)θ(·)表示的是Z和U的相互作用，如果假設(shè)θ(·)是常數(shù)向量，那么模型就可以看作是部分線性單指標(biāo)模型，若進(jìn)一步令系數(shù)函數(shù)向量θ(·)的維數(shù)q等于1，那么就可以得到單指標(biāo)模型。對于單指標(biāo)部分，若令聯(lián)系函數(shù)g(βTX)≡βTX，那么模型就變成了部分線性變系數(shù)模型，進(jìn)一步取參數(shù)β的維數(shù)p為1，模型就退化成了變系數(shù)模型。

近年來，關(guān)于模型式(1)的研究成果已經(jīng)十分豐富。為了研究化學(xué)污染物水平與每天因呼吸系統(tǒng)疾病住院的總?cè)藬?shù)的關(guān)系，以及溫度和相對濕度對入院人數(shù)的影響，Wong等[2]首次提出PVCSIMs，他們結(jié)合二元局部線性方法、平均方法和一步回?cái)M技術(shù)(One-step Back-fitting Technique)得到函數(shù)和參數(shù)的有效估計(jì)?；谒麄兊难芯?，Huang和Zhang[3]進(jìn)一步利用廣義似然法(GLR)解決了模型中變系數(shù)部分的檢驗(yàn)問題。Li和Zhang[4]通過將系數(shù)函數(shù)和聯(lián)系函數(shù)樣條化，提出了模型的懲罰樣條估計(jì)方法，該方法可以同時(shí)得到未知參數(shù)和函數(shù)的估計(jì)值。Wang和Xue[5]的研究指出，Wong等[2]用二元局部線性光滑進(jìn)行估計(jì)可能會(huì)導(dǎo)致估計(jì)量不相合，因此，提出一種較為穩(wěn)定的逐步估計(jì)法對模型進(jìn)行估計(jì)，其基本思想是：假設(shè)參數(shù)已知，將模型轉(zhuǎn)換成變系數(shù)模型，并利用Nadaraya-Watson核估計(jì)方法和局部線性回歸方法逐步得到系數(shù)函數(shù)和聯(lián)系函數(shù)的初步估計(jì)，然后根據(jù)這些初始估計(jì)量計(jì)算未知參數(shù)的估計(jì)值，文章還討論了估計(jì)量的漸近性質(zhì)，并且建立了逐點(diǎn)置信區(qū)間和置信域。Huang[6]通過經(jīng)驗(yàn)似然方法研究了單指標(biāo)參數(shù)的極大似然估計(jì)，并且利用截面經(jīng)驗(yàn)似然方法構(gòu)造了各參數(shù)分量的置信區(qū)間。Huang[7]等結(jié)合SCAD(Smoothly Clipped Absolute Deviation)懲罰和逐步估計(jì)法研究了指標(biāo)參數(shù)β的變量選擇問題，在一定的正則化條件下，還構(gòu)建了估計(jì)量的大樣本性質(zhì)。最近，受到圖像數(shù)據(jù)分析的啟發(fā)，Li等[8]討論了函數(shù)型數(shù)據(jù)下PVCSIMs的估計(jì)問題，利用局部線性方法逐步迭代得到了系數(shù)函數(shù)、聯(lián)系函數(shù)、指標(biāo)參數(shù)以及方差函數(shù)的估計(jì)值，并且證明提出的方法相較于Wong等[2]和Wang和Xue[5]中的方法更加穩(wěn)定。

雖然半?yún)?shù)模型可以用來解決大部分回歸擬合問題，然而在實(shí)際應(yīng)用中，由于操作人員的失誤和測量工具不精確等問題，收集到數(shù)據(jù)中往往帶有測量誤差，因此進(jìn)一步研究半?yún)?shù)誤差模型是很有必要的。扭曲測量誤差(Distorted Measurement Error)作為誤差的常見形式，是近年來學(xué)者們研究的熱點(diǎn)問題。Sentürk 和 Müller[9]提出協(xié)變量調(diào)整回歸(Covariate-adjusted Regression, CAR)模型，即假設(shè)響應(yīng)變量和協(xié)變量都含有扭曲測量誤差，通過將線性誤差模型轉(zhuǎn)換成變系數(shù)模型，并結(jié)合分箱法給出了對未知參數(shù)的估計(jì)。Cui等[10]提出一種非參數(shù)誤差回歸模型的一般估計(jì)方法，該方法通過核估計(jì)得到誤差函數(shù)的估計(jì)量并以此計(jì)算受污染變量的校正值，然后根據(jù)校正后的變量對目標(biāo)參數(shù)進(jìn)行估計(jì)。Delaigle[11]等進(jìn)一步討論了在不同假設(shè)條件下非參數(shù)誤差回歸模型的估計(jì)問題，在弱化對未知變量或扭曲函數(shù)的假設(shè)后，提出了更一般的估計(jì)方法，并且建立了相應(yīng)估計(jì)量的漸近性質(zhì)。此外，Qian和Huang[12]研究了含扭曲測量誤差的部分非線性變系數(shù)模型的統(tǒng)計(jì)推斷問題。Dai和Huang[13]將協(xié)變量含扭曲測量誤差的情形推廣到部分非線性變系數(shù)模型。

對于糖尿病數(shù)據(jù)集，根據(jù)數(shù)據(jù)本身的特點(diǎn)，利用單指標(biāo)部分含扭曲測量誤差的PVCSIMs模型進(jìn)行擬合。該模型具有復(fù)雜的結(jié)構(gòu)，可以靈活地?cái)M合變量之間的關(guān)系。進(jìn)一步考慮了測量誤差的存在，使得模型更加符合實(shí)際情形，可以更好地挖掘變量之間潛在的聯(lián)系。模型的估計(jì)方法主要參考Cui等[10]和Huang[6]中的思想。

1 模型建立

含有扭曲測量誤差的部分變系數(shù)單指標(biāo)模型具有如下形式：

(2)

(1) ‖β0‖=1，其中‖·‖表示Euclid模；

(2)E|ψ(U)|=1,E|φr(U)|=1；

(3) 模型誤差ε的方差是有限的。

其中，r=1,2…,p。假設(shè)條件(1)保證了單指標(biāo)參數(shù)的唯一性，假設(shè)條件(2)保證了誤差函數(shù)的可識(shí)別性，假設(shè)條件(3)保證了估計(jì)量的漸近性質(zhì)。

接下來，介紹模型(2)的估計(jì)方法。首先估計(jì)誤差函數(shù)，根據(jù)Cui等[10]和Delaigle等[11]中提出的方法，結(jié)合假設(shè)條件(2)可以得到：

因此可以用N-W核估計(jì)方法可以得到誤差函數(shù)的估計(jì)式分別為

根據(jù)校正后的變量，可以得到模型(1)的一個(gè)近似形式

(3)

模型(3)的估計(jì)過程主要參考文獻(xiàn)[6]中的思想和方法，簡要表述過程如下：

令B={β∈Rp:‖β‖=1}，可以推出β0∈B。由于假設(shè)條件(1)的存在，目標(biāo)函數(shù)

β∈B，在β0處的一階導(dǎo)數(shù)不存在，因此考慮使用Zhu和Xue[14]中的“去一分量”法得到β0的有效估計(jì)。不妨假設(shè)β的第r個(gè)分量βr>0，定義

β(r)=(β1,β2,…,βr-1,βr+1,…βp)T

則有

接著，可以計(jì)算出β關(guān)于β(r)的Jacobian矩陣為

不難推出，{ηi(β(r)),i=1,2,…,n}是相互獨(dú)立的且E(ηi(β(r)))=0。因此β(r)的經(jīng)驗(yàn)對數(shù)似然比函數(shù)可以定義為

(4)

接下來，只要求解函數(shù)的初始估計(jì)值即可。對于函數(shù)的初始估計(jì)，考慮利用局部線性光滑方法分步求解系數(shù)函數(shù)和聯(lián)系函數(shù)，具體的估計(jì)步驟如下：首先假設(shè)參數(shù)β0已知，對模型(3)的兩邊求條件期望，經(jīng)過簡單的計(jì)算可以得到：

其中

且有

參考Einmahl和Mason[15]的方法，應(yīng)用N-W核估計(jì)求解兩個(gè)二元函數(shù)的估計(jì)表達(dá)式。選擇核函數(shù)K1(t,u)=K(t)·K(u)和帶寬h2(n)→0，ωj(·,·)(j=1,2)的估計(jì)量分別為

接著，令a=(a1,…,aq)T和b=(b1,…,bq)T,在u的鄰域內(nèi)有

通過最小化加權(quán)平方和

可以得到系數(shù)函數(shù)θ(·)及其一階導(dǎo)數(shù)θ′(·)的初始估計(jì)量為

其中κ(u)=diag(Kh2(U1-u),Kh2(U2-u),…,Kh2(Un-n))以及

然后，帶入系數(shù)函數(shù)估計(jì)量到式(3)中，通過局部線性方法，類似地，可以得到聯(lián)系函數(shù)g(·)及其一階導(dǎo)數(shù)g′(·)的初始估計(jì)方程估計(jì)分別為

以及

2 模型擬合

糖尿病數(shù)據(jù)集共包含442個(gè)樣本觀測值和11個(gè)變量，數(shù)據(jù)集來自Efron等[16]。Zhang等[17]利用含測量誤差的部分線性單指標(biāo)模型研究了該數(shù)據(jù)集，并且證明變量Sex和6種血清的化驗(yàn)數(shù)據(jù)與變量BMI可能存在非線性關(guān)系，受到該實(shí)驗(yàn)結(jié)果的啟發(fā)，考慮用部分變系數(shù)單指標(biāo)模型重新分析該數(shù)據(jù)集，觀察變量Sex和6種血清的化驗(yàn)數(shù)據(jù)變量的系數(shù)是否受到BMI的影響。此外，考慮到病人的BMI可能與他的年齡和血壓存在關(guān)系，因此假設(shè)變量Age和BP受到BMI的污染，具體的變量意義見表1。實(shí)驗(yàn)前，對各個(gè)變量進(jìn)行了標(biāo)準(zhǔn)化處理。

表1 糖尿病數(shù)據(jù)集變量及其含義Table 1 The variables and meanings of diabetes data set

在非參數(shù)回歸方法中，帶寬對估計(jì)精度較大，因此使用合適的帶寬選擇方法十分重要。由于對于誤差函數(shù)的估計(jì)方法比較簡單，因此可以采用基于經(jīng)驗(yàn)的拇指規(guī)則(Rule of Thumb)來選擇帶寬。選取h1=n-1/3SE(V),其中

而對于系數(shù)函數(shù)向量θ(·)和聯(lián)系函數(shù)g(·)的估計(jì)比較復(fù)雜，因此使用交叉驗(yàn)證法(Cross-validation)選擇最優(yōu)帶寬hcv，并令h2=hcv,h3=hcv(nlogn)-1/20。

根據(jù)第2節(jié)的模型估計(jì)方法，本次實(shí)驗(yàn)的具體算法流程如下：

Step2 選擇符合假設(shè)(1)的初始參數(shù)βini；

使用模型(2)對糖尿病數(shù)據(jù)集進(jìn)行擬合，變量選擇如表1所示。首先得到的是誤差函數(shù)ψ(·)和φr(·)(r=1,2)的估計(jì)曲線，見圖1和圖2，其中圖2(a)和圖2(b)分別代表函數(shù)φ1(·)和φ2(·)的估計(jì)曲線。從估計(jì)曲線的趨勢可以看出，各條估計(jì)曲線都不是水平的，這說明變量(Y,X1,X2)與混淆因子V存在非線性的關(guān)系，這驗(yàn)證了之前的實(shí)驗(yàn)假設(shè)：糖尿病患者定量測量數(shù)據(jù)，Age和BP受到了混淆因子BMI的污染。

圖1 ψ(·)的估計(jì)值Fig. 1 The estimation of ψ(·)

(a)φ1(·) (b)φ2(·)

(a)Z1 (b)Z2 (c)Z3

(d)Z4 (e)Z5 (f)Z6 (g)Z7

圖4 聯(lián)系函數(shù)g(·)的估計(jì)曲線Fig. 4 The estimation curve of contact function g(·)

3 結(jié)束語

對于高維數(shù)據(jù)分析，相較于傳統(tǒng)的參數(shù)模型和非參數(shù)模型，半?yún)?shù)模型不僅可以較好地?cái)M合數(shù)據(jù)，還可以避免“維數(shù)災(zāi)禍”問題，因此廣泛地應(yīng)用醫(yī)藥和經(jīng)濟(jì)等領(lǐng)域。此外，在實(shí)際應(yīng)用中，由于外界因素的干擾，很難避免測量誤差的產(chǎn)生，如果忽略誤差的影響，就有可能導(dǎo)致模型擬合產(chǎn)生偏差，因此研究帶有測量誤差的模型應(yīng)用問題是比較重要的。這里，針對糖尿病數(shù)據(jù)集，應(yīng)用部分變系數(shù)單指標(biāo)模型進(jìn)行擬合，并假設(shè)變量Y和X受到混淆因子BMI的乘積污染。觀察實(shí)驗(yàn)結(jié)果發(fā)現(xiàn)，6種血清測量數(shù)據(jù)和變量Sex的系數(shù)會(huì)隨著BMI的變化而變化，并且對比帶有測量誤差和不含測量誤差兩種情形的結(jié)果發(fā)現(xiàn)，糖尿病人定量測量值、Age和BP均受到BMI的污染。這些結(jié)果說明選擇單指標(biāo)部分帶有測量誤差的部分變系數(shù)單指標(biāo)模型對該數(shù)據(jù)集進(jìn)行擬合是合理的，并且相較于不含測量誤差的半?yún)?shù)模型，可以更好地挖掘數(shù)據(jù)中的信息。