嵌套刪失數(shù)據(jù)期望最大化的高斯混合聚類算法

余海燕 陳京京 邱 航 , 王 永 王若凡

高斯混合聚類[1?2]作為統(tǒng)計機器學習、模式識別和陣列數(shù)據(jù)分析等的重要模型,廣泛用于健康醫(yī)療[3?4]、故障診斷[5?6]等領(lǐng)域.然而,常因諸如截斷的數(shù)據(jù)、傳感器故障或傳輸錯誤等造成數(shù)據(jù)不完整問題[1],引起推斷偏差并使得聚類精度下降.例如在醫(yī)療決策智能支持中[7?8],需要依據(jù)患者的各項生理指標信息進行智能推理[9?10],然而由于記錄數(shù)據(jù)刪失或截斷等導致數(shù)據(jù)不完整,從而給數(shù)據(jù)分析帶來困難.在惡性淋巴瘤等疾病診斷[11]中,流式細胞儀記錄的數(shù)據(jù)因測量信號強度范圍有限而使得數(shù)據(jù)記錄在一個固定范圍內(nèi)(如0 到1 023 之間),如果測量值超出這一范圍,則該值將替換為最接近的值,小于0的值將被刪失記為0,大于1 023 的值將被刪失記為1023.類似的刪失數(shù)據(jù)還包括保險費理賠計算中,因一定數(shù)量免賠額的存在使得記錄成為刪失數(shù)據(jù)等.這類刪失數(shù)據(jù)處理不當會影響分析結(jié)果的可靠性,甚至使得聚類模型參數(shù)推斷出現(xiàn)較大偏差.又因這類數(shù)據(jù)的分布參數(shù)的精確估計,為處理變量或治療方案對觀察結(jié)果的因果效應分析[12]提供基礎(chǔ),甚至影響到后續(xù)的決策方案選擇.高斯聚類算法因能夠提供分布參數(shù)的估計,故而刪失數(shù)據(jù)的參數(shù)估計已成為高斯混合聚類的一個重要熱點問題.

刪失數(shù)據(jù)的處理方法?；谌笔?shù)據(jù)的處理機理.因數(shù)據(jù)缺失機制不同,處理方法也不盡相同.數(shù)據(jù)缺失可以分為隨機缺失(Missing at random,MAR)和非隨機缺失(Missing not at random,MNAR)兩大類[12].大多數(shù)傳統(tǒng)的缺失數(shù)據(jù)處理方法主要集中于使用樣本抽樣推斷、貝葉斯推斷和似然法推斷[13].其中貝葉斯推斷和似然法在實際數(shù)據(jù)中的應用更為普遍.當評估項目的長期性能數(shù)據(jù)隨機缺失且觀測數(shù)據(jù)也隨機缺失時,使用樣本抽樣估計數(shù)據(jù)集分布參數(shù)可以忽略缺失機制.當數(shù)據(jù)屬于隨機缺失且缺失機制參數(shù)不同于數(shù)據(jù)集分布參數(shù)時,使用貝葉斯推斷和似然法也可以忽略缺失機制.文獻[12]對非隨機缺失問題的探索,還包括不可忽略性無響應問題、不可忽略性缺失性問題,甚至被稱為有信息缺失的問題等.文獻[14]認為存在解決非隨機缺失的方法,但是通常難以檢驗,為此提出了懲罰驗證標準,通過懲罰未知參數(shù)過多的模型來防止模型過擬合.刪失數(shù)據(jù)作為一種非隨機性缺失數(shù)據(jù)[15?16],因其缺失機制(如刪失)的特殊性而不能直接使用一般的非隨機缺失方法直接計算[11].

刪失數(shù)據(jù)常包括右刪失和區(qū)間刪失等類型.對于右刪失數(shù)據(jù),文獻[17]基于一類廣義概率測度的誤差一致性,提供了適用于刪失數(shù)據(jù)的分類支持向量機并應用于刪失數(shù)據(jù)平均值、中位數(shù)、分位數(shù)的估計以及分類問題.針對區(qū)間刪失數(shù)據(jù),文獻[15]提出一種貝葉斯非參數(shù)化方法進行概率擬合.文獻[18]基于左截斷右刪失數(shù)據(jù)構(gòu)造了分位差的經(jīng)驗估計,并提出了分位數(shù)差的核光滑估計.針對刪失混合數(shù)據(jù),文獻[19]提出了一個加權(quán)最小二乘估計的一般族,并證明了現(xiàn)有的一致非參數(shù)方法屬于這個族,識別其估計量并分析其漸近性質(zhì).而在高斯混合聚類模型算法中,一般假設(shè)觀測值的特征向量對聚類有相同的權(quán)重[20].然而文獻[1]認為高斯混合聚類模型的每一個特征向量的權(quán)重并不一樣,提出競爭性懲罰期望最大化算法.該算法將特征選擇模型和高斯混合聚類模型結(jié)合在一起,使用馬爾科夫毯濾波器消除多余的特征項,找到最小的相關(guān)特征子集,同時確定高斯混合模型的混合成分個數(shù).文獻[21]提出了一種基于高斯混合聚類和模型平均的算法.對于缺失值,該方法將每一組成成分得出的估計值作為線性組合的概率估計權(quán)重,最終結(jié)果是混合成分的估計值的平均值.文獻[2]討論高斯混合聚類分析的過擬合問題.該文獻改變了以往認為不相關(guān)變量必須通過線性回歸方程依賴整個相關(guān)變量的做法,認為相關(guān)變量并不一定要解釋所有的不相關(guān)變量.該模型可以有效地提高聚類算法的性能且變量選擇的實現(xiàn)基于一個向后逐步算法.標準期望最大化(Expectation-maximization,EM)算法作為高斯混合模型中常用的缺失數(shù)據(jù)處理方法[22],更適用于處理隨機缺失數(shù)據(jù).本文在標準EM 的高斯混合聚類算法(EMGM)基礎(chǔ)上,提出了嵌套刪失數(shù)據(jù)期望最大化的高斯混合聚類算法(cenEMGM).

本文主要解決非隨機缺失下的刪失數(shù)據(jù)因利用率不高而導致聚類準確度不高的問題.本文的主要貢獻是:利用高斯混合模型聚類算法獨有的特性,在標準EM 算法的基礎(chǔ)上提出改進算法cenEMGM,并揭示了刪失率對模型算法的作用機制.將刪失數(shù)據(jù)和高斯混合模型聚類算法結(jié)合,更加準確地處理刪失數(shù)據(jù).通過調(diào)整刪失數(shù)據(jù)的分布函數(shù),使得刪失數(shù)據(jù)最大期望算法不斷更新均值、協(xié)方差和混合系數(shù)的估計值,從而使得聚類簇中心不斷接近真實的簇中心.cenEMGM 算法在標準EMGM 算法的基礎(chǔ)上進行改進,該方法更加靈活,對刪失和未刪失數(shù)據(jù)采取不同的處理方式.刪失數(shù)據(jù)EM 算法和高斯混合聚類相結(jié)合,使得該方法比原方法聚類效果更好,準確性更高.后續(xù)章節(jié)結(jié)構(gòu)如下:第1 節(jié)引入高斯混合聚類模型.第2 節(jié)論述刪失型缺失數(shù)據(jù)的相關(guān)概念.第3 節(jié)構(gòu)建高斯混合聚類的參數(shù)估計算法,包括標準EMGM 算法和cenEMGM 兩種算法,以及兩個模型校驗準則.第4 節(jié)使用數(shù)值實驗驗證算法.第5 節(jié)得出結(jié)論.

1 高斯混合聚類模型

對d維數(shù)據(jù)空間 Rd中,隨機變量y的觀察值為一個由n個樣本構(gòu)成的數(shù)據(jù)集,D={y1,y2,···,yn},其中yi為其第i個樣本.并將第j維數(shù)據(jù)記為y(j).假設(shè)樣本生成過程由包含K個成分的高斯混合分布確定.第k個成分fk的參數(shù)為 Θk=(πk,μk,Σk);其中,πk為其混合系數(shù),μk為均值,Σk為方差.全部參數(shù) Θ={Θ1,Θ2,···,ΘK}.y(j)為其第j維觀測值.對于y,定義高斯混合分布[20]如下:

其中,K為混合成分數(shù)量,且每個混合成分對應一個高斯分布 N(μk,Σk),相應的“混合系數(shù)”πk >0,.

樣本生成過程中,記π={π1,···,πK},首先根據(jù)π定義的先驗分布選擇高斯混合成分,且選擇第k個混合成分的概率為πk;然后,根據(jù)被選擇的混合成分的概率密度函數(shù)進行采樣,從而生成相應的樣本.

在高斯混合聚類模型中,類似地存在K個簇,C={C1,C2,···,CK}.將yi是否被劃分到簇Ck中的隨機變量記為,簇指示變量∈{0,1}.當yi被劃分到簇Ck時,=1,意味著yi{由fk生成;否則=0.對于N個樣本總體,表示第k個 (k=1,2,···,K)高斯混合成分生成樣本y的指示變量值.因此,對于i=1,2,···,N,=1的概率對應于πk.根據(jù)貝葉斯定理,的后驗分布對應于

當高斯混合分布(1)已知時,高斯混合聚類將把樣本集D劃分為K個簇,樣本yi的簇標記λi.

可見,高斯混合聚類的本質(zhì)是采用概率模型(高斯分布)對原型進行刻畫,簇劃分則由原型對應后驗概率確定.因一個簇對應一個中心點,隸屬于每一個簇C的數(shù)據(jù)樣本將聚類在簇中心點附近.高斯混合聚類模型效果越好,所估計的簇中心點與實際簇的中心點之間距離將越小甚至重合.

2 刪失型數(shù)據(jù)缺失機制

2.1 數(shù)據(jù)缺失機制

依據(jù)文獻[12]將數(shù)據(jù)缺失機制分為四種類型,包括隨機缺失、完全隨機缺失、取決于未被觀測因素的缺失(可以通過未被觀察或記錄的數(shù)據(jù)進行預測的)以及和僅依賴于缺失值自身的缺失機制.后兩種缺失機制即為這里將定義的非隨機缺失.

在數(shù)據(jù)空間 Rd中,令A為一個實數(shù)集合,設(shè)為一個指示變量,表示y的元素在集合A中是否存在觀察值.若∈A,則=1,否則=0.這里yi不區(qū)分變量及其真實值,而將其觀測值記為.令作為yi中不存在缺失的部分,表示yi中存在缺失值的部分,那么.

定義1.如果對所有和參數(shù) Θ,

則缺失數(shù)據(jù)機制為隨機缺失.

定義2.如果對所有和參數(shù) Θ,

則缺失數(shù)據(jù)機制為非隨機缺失.

可見,對于隨機缺失數(shù)據(jù),其樣本數(shù)據(jù)及指示變量滿足交換性,而非隨機缺失數(shù)據(jù)不滿足這一性質(zhì)[12].當缺失數(shù)據(jù)是隨機缺失時,可直接使用標準EM 算法、多值插補、回歸等方法揭示缺失機制.下面引入一類非隨機性缺失數(shù)據(jù),即刪失數(shù)據(jù),并研究其缺失機制和參數(shù)估計方法.

2.2 刪失數(shù)據(jù)的似然函數(shù)

這里給出刪失數(shù)據(jù)的定義,并詳細闡述刪失數(shù)據(jù)的缺失機制和似然函數(shù).在數(shù)據(jù)空間 Rd中,[a,b]d為一個超矩陣[11],其中上邊界b=(b(1),···,b(d))T,下邊界a=（a(1),···,a(d)）T.

定義3.刪失數(shù)據(jù)(Censored data)是指yi的觀測值滿足分段函數(shù):

其中,a

換言之,yi中的缺失部分被分別賦予a或b對應維度上的元素值.為分析概率密度和估計參數(shù),假設(shè)的元素個數(shù)為J1,的元素個數(shù)為J2,且J1+J2=d.不妨進一步假設(shè),.對于刪失數(shù)據(jù),A=[a,b]d.為簡化,令δij=1?,當δij=1時,表示因刪失而存在缺失數(shù)據(jù),其對應觀測值被賦予邊界值;相應地,δij=0,表示不存在缺失數(shù)據(jù),即觀測值等同于真實值.y觀測值的樣本刪失率.對于一維數(shù)據(jù),刪失率pce=nce/n,其中nce是存在刪失的樣本數(shù).

根據(jù)刪失數(shù)據(jù)的定義,y1:n的部分真實值(如序數(shù)為n1+1,···,n的值)被修改.那么,其被修改后的數(shù)據(jù)(不存在缺失部分的值、和缺失部分的修改值)構(gòu)成新數(shù)據(jù)集,記為x1:n.對于?i,?j,有

與缺失數(shù)據(jù)機制對應,但因每一個樣本yi的刪失模式會不一樣,而使用im和io分別表示刪失和未刪失數(shù)據(jù)的坐標序號集,故分別指刪失部分的缺失值(缺失時的真實值)和刪失后的改寫值(簡稱刪失值),分別指原數(shù)據(jù)不存在缺失的部分與刪失型數(shù)據(jù)對應的部分值,盡管沒有刪失時它們值等同.那么.同時,.

為簡化,將y的數(shù)據(jù)空間劃分為{Yt|t=0,1,···,T},其中當,此時數(shù)據(jù)不存在刪失;而當,t>0 時,數(shù)據(jù)發(fā)生刪失.將刪失部分調(diào)整后的觀測值x的數(shù)據(jù)空間劃分為{Xt|t=1,···,T},注意,這里沒有涵蓋不存在刪失的部分,即x的數(shù)據(jù)空間劃分不涵蓋X0.對于yi∈Y0,觀測值xi的似然函數(shù)如下:

而對于yi缺失機制,有,ti >0,其似然函數(shù)如下:

并且關(guān)于f(xi)推導式(4)的右邊部分轉(zhuǎn)化為:

3 高斯混合聚類的參數(shù)估計

高斯混合聚類參數(shù)估計主要包括成分的期望、方差和對應的混合系數(shù).嵌套標準EM 的高斯混合聚類算法,這里簡記為EMGM.并將針對刪失數(shù)據(jù)所提出的改進算法,即嵌套刪失型數(shù)據(jù)期望最大化的高斯混合聚類算法,簡記為cenEMGM 算法.

3.1 基于高斯混合聚類的標準算法EMGM

對于獨立觀測變量集合y1:n,參數(shù)空間 Θ,第k個成分fk和簇指示變量,對數(shù)似然函數(shù)為:

其中,Θk=(πk,μk,Σk)表示第k個成分的參數(shù),(Σk)?1表示 Σk的倒數(shù),Const表示常數(shù),tr(·)表示矩陣的跡,In表示值全為1 的 1×n向量.

根據(jù)標準的期望最大化算法[23],其假設(shè)為數(shù)據(jù)存在隨機缺失.對于獨立觀測變量集合y1:n,Θ,Θold和 Θnew分別為參數(shù)空間,算法中更新前的參數(shù)及更新后的參數(shù).

算法第一步(步驟E):計算期望函數(shù)Q（Θ;Θold）=,步驟E 可以簡化為計算條件概率:

第二步(步驟M):尋找新的參數(shù)集 Θnew,使得Θnew=arg maxΘQ（Θ;Θold）.更新后的參數(shù)Θnew=,形成一個更新的閉環(huán)形式:

該算法不斷迭代E 步和M 步,直至收斂.以最后獲得的更新參數(shù)作為 Θ 的最優(yōu)估計值.

3.2 估計算法的有效性

引理1.通過最大似然估計方法獲得全數(shù)據(jù)的參數(shù),即求解全數(shù)據(jù)得分向量方程(yi,Θ)=0,得到

其中,IF(Θ)為全數(shù)據(jù)信息矩陣,IF(Θ)=E[SF(D,Θ)(D,Θ)].

通過正則漸近線性法(Regular and asymptotically linear,RAL)[24]獲得全數(shù)據(jù)的參數(shù)記為,即求解全數(shù)據(jù)得分向量方程=0.

引理2.對于RAL 方法估計的參數(shù),應滿足:

對于 arg maxΘQ（Θ;Θold）,根據(jù)全數(shù)據(jù)參數(shù)估計的引理,存在關(guān)于期望最大化算法估計刪失型缺失數(shù)據(jù)的定理.

定理1.令全數(shù)據(jù)D={y1,y2,···,yn},對應的刪失型缺失數(shù)據(jù),對缺失數(shù)據(jù)使用逐步更新的EMGM 算法估計參數(shù),可通過以下方程求解.

依據(jù)第3.2 節(jié)給定刪失數(shù)據(jù)及其似然函數(shù),cen-EMGM 算法首先計算完全對數(shù)似然函數(shù)的期望:

該式子可以由式(4)進一步推導出結(jié)果.

結(jié)合高斯混合分布定義(1),針對y(mi)的條件概率分布,,推導其條件分布期望.因為是正態(tài)密度函數(shù)且滿足

條件密度fk（y(mi)|x）是在Xc上的截尾正態(tài)密度函數(shù),那么計算關(guān)于Qc的充分統(tǒng)計量:

定理2.全數(shù)據(jù)D={y1,y2,···,yn},對應的刪失型缺失數(shù)據(jù),在給定缺失數(shù)據(jù)和RAL 估計參數(shù)下的原數(shù)據(jù)的得分向量,對缺失數(shù)據(jù)使用cenEMGM 算法估計參數(shù),滿足

其中,IF(Θ)為全數(shù)據(jù)信息矩陣,Θ 為數(shù)據(jù)的真實參數(shù),對于cenEMGM 算法.

證明.因cenEMGM 算法中刪失數(shù)據(jù)的對數(shù)似然函數(shù)期望為,那么其得分向量的期望

故而有

又因為

所以有

3.3 針對刪失數(shù)據(jù)的算法cenEMGM

根據(jù)定理2 獲得對數(shù)似然函數(shù)的期望Qc關(guān)于Θ最大化的解,即得到了Θ(t)=ar g maxΘQc（Θ;Θ(t?1)）的優(yōu)化解,.該算法的步驟t≥1,并且 Θ(0)表示初始值,可通過K-means 聚類方法獲得賦值.求解的高斯混合聚類的混合系數(shù)πk為:

同時,μk和 Σk關(guān)于 arg maxΘQc（Θ;Θ(t?1)）的優(yōu)化解分別為:

式(13)～(15)作為標準EM 算法式(6)～(8)針對刪失型缺失數(shù)據(jù)的改進.式(13)與(6)在形式上沒有變化,從理論上論證了刪失型算法cenEMGM與標準算法EMGM 在混合系數(shù)上一致.式(14)與(7)相比較發(fā)現(xiàn),在刪失數(shù)據(jù)算法cenEMGM 中,y1:n的刪失部分被條件均值代替.式(15)與(8)相比較發(fā)現(xiàn),刪失數(shù)據(jù)算法cenEMGM 的被樣本校正協(xié)方差所替代.標準算法EMGM 即為算法cenEMGM 處理不存在刪失數(shù)據(jù)時的特定情形.

3.4 模型檢驗準則

為了防止算法出現(xiàn)過擬合并計算估計值和真實值之間的距離,需要設(shè)定模型檢驗準則.這里引入信息散度(Kullback-Leibler divergence,KLD)和赤池弘次信息準則(Akaike＇s information criterion,AIC)[20,25].信息散度KLD 公式[25]為:

其中,p(y) 是y真實分布的概率密度函數(shù),q(y)是y估計分布的概率密度函數(shù).本文中y的概率密度函數(shù)由高斯混合分布(1)確定..在算法EMGM 中,p(y)由式(6)～(8)確定;在算法cenEMGM中,q(y)由式(13)～(15)確定.

對于AIC 準則,其值最小的模型即為最佳模型.假設(shè)模型的誤差服從獨立正態(tài)分布,AIC 可表示為:

其中,N(Θ)是模型算法參數(shù)的數(shù)量,d為D數(shù)據(jù)維度,K為高斯混合模型的成分數(shù)量,L(Θ)是參數(shù)集 Θ 的似然函數(shù).

3.5 cenEMGM 算法及分析

嵌套刪失型數(shù)據(jù)期望最大化的高斯混合聚類算法(cenEMGM)主要由高斯混合聚類和針對刪失數(shù)據(jù)的期望最大化算法構(gòu)成,如算法1 所示.第1)步初始化參數(shù),常使用k-means 算法.第2)～10)步,運行直至滿足停止條件,跳出循環(huán).其中第3)～4)步,cenEMGM 算法的E 步,計算后驗概率;第5)～9)步,cenEMGM 算法的M 步,計算新的模型參數(shù).第11)～13)步,劃分簇.算法流程的停止條件是,其中ε是一個小的正數(shù)(如1.0×10?6).其中,,k=1,2,···,K.cenEMGM 算法的計算復雜度(時間復雜度)受到樣本規(guī)模n和參數(shù)規(guī)模影響,其中d為D數(shù)據(jù)維度,K為高斯混合模型的成分數(shù)量.

算法1.嵌套刪失型數(shù)據(jù)期望最大化的高斯混合聚類算法cenEMGM

cenEMGM 算法的核心步驟主要基于式(13)～(15).與之對應的標準EMGM 算法,其核心是式(6)～(8).cenEMGM 算法是針對刪失型缺失數(shù)據(jù)的改進算法,先根據(jù)新均值向量計算新樣本規(guī)模,然后計算新混合系數(shù).因為樣本規(guī)模改變,所以樣本方差、刪失率、觀測數(shù)據(jù)均值等參數(shù)同步做出改變.針對刪失數(shù)據(jù)修改的這些內(nèi)容,使cenEMGM 算法更靈活,更能適應含有刪失數(shù)據(jù)的高斯混合聚類.

高斯混合分布中,πk是選擇第k個混合成分的概率,由式(8)和式(13)可以看出,樣本刪失率間接地通過樣本容量影響著πk,所以pce對πk產(chǎn)生影響.數(shù)據(jù)質(zhì)量可以衡量采樣機制產(chǎn)生的選擇偏差程度[26],其不僅和估計準確度有關(guān),更是與刪失率有關(guān).為了提高模型的準確性,可以根據(jù)刪失率調(diào)整并確定樣本規(guī)模n.關(guān)于樣本規(guī)模在實驗設(shè)計中已有討論[27].這里給出樣本方差未知時刪失率pce與樣本規(guī)模n的結(jié)論.根據(jù)統(tǒng)計推斷理論,檢驗水準α時,預測能力 (1?β)表示,當所考慮的總體與原假設(shè)H0確有差別時,按照檢驗水準α能夠發(fā)現(xiàn)拒絕它的概率.總體方差未知時,在刪失數(shù)據(jù)缺失率為pce的情況下,估計樣本容量大小如下:n0=,其中δ表示估計精度(即允許誤差),為數(shù)據(jù)分布中的真實缺失率,t為檢驗統(tǒng)計量.對于一定規(guī)模的同一數(shù)據(jù)集,隨著樣本刪失率pce上升,參數(shù)估計模型的估計能力下降,導致準確性也降低.因此,數(shù)據(jù)分析中要求樣本容量不小于n0.隨著數(shù)據(jù)感知和收集成本下降,數(shù)據(jù)可得性變高,統(tǒng)計機器學習模型使用的數(shù)據(jù)規(guī)模選取常會超過模型的測試能力要求,且通常會考慮數(shù)據(jù)的缺失機制[12].

4 數(shù)值實驗分析

這里使用人工數(shù)值實驗與真實數(shù)據(jù)分析,驗證方法的有效性.

4.1 人工數(shù)值實驗分析

實驗從預設(shè)分布生成數(shù)據(jù)集,并對數(shù)據(jù)進行刪失處理.在刪失數(shù)據(jù)上,分別采用嵌套標準EM 的高斯混合聚類算法EMGM 和嵌套刪失型數(shù)據(jù)cenEM 的高斯混合聚類算法cenEMGM 進行實驗分析.實驗結(jié)果通過聚類的真實參數(shù)與估計參數(shù)比較、KL 散度等統(tǒng)計指標進行比較分析.

為在多變量上比較算法,這里設(shè)計兩個含有三個成分的二元高斯混合模型的實驗.在兩個實驗中,實驗數(shù)據(jù)集D S-a 的觀測值 (Y1,Y2)被設(shè)置在[10,50]×[5,45]的矩形窗中,用于右刪失型數(shù)據(jù)和雙邊刪失型數(shù)據(jù)在EMGM 算法和cenEMGM 算法上的實驗;實驗數(shù)據(jù)集DS-b 的觀測值 (Y1,Y2)被設(shè)置在 [?20,60]×[?10,60] 的矩形窗中,用于左刪失型數(shù)據(jù)和雙邊刪失型數(shù)據(jù)在EMGM 算法和cenEMGM算法上的實驗.右(左)刪失型缺失是指在變量值域范圍內(nèi),設(shè)定了觀測值上(下)界,且大(小)于該上(下)界的其他值被賦予該上(下)界值,但并無給定的下(上)界.雙邊刪失型缺失是指在變量值域范圍內(nèi),同時設(shè)定了觀測值上界和下界值,大于該上界的其他值被賦予該上界值,且小于該下界的其他值被賦予該下界值.這里生成的兩組數(shù)據(jù)分別采用了兩種刪失機制,并非只討論一組數(shù)據(jù)的左刪失、右刪失及雙邊刪失,以便體現(xiàn)刪失數(shù)據(jù)邊界的多樣性和實驗的可重復性.

在實驗數(shù)據(jù)集DS-a 中,三個分量的中心都在對應的矩形窗內(nèi),參數(shù)設(shè)置如下:成分權(quán)重為π=(0.25,0.40,0.35);均值為μ1=(23.50,23.50),μ2=(33.50,23.50),μ3=(40.50,40.50);方差中,成分1 與成分2 在兩個變量之間不存在相關(guān)性:

成分3 的兩個變量之間存在相關(guān)性:

在實驗數(shù)據(jù)集DS-b 中,雖然三個成分的中心都在對應的矩形窗內(nèi),但有兩個成分的中心落在了下界之外.參數(shù)設(shè)置如下:成分權(quán)重和方差分別與實驗數(shù)據(jù)集DS-a 對應一致.但它們的均值分別為μ1=(?3.50,23.50),μ2=(33.50,?3.50),μ3=(40.50,40.50).

在每種情形下繪制1 000 個數(shù)據(jù)點后,根據(jù)刪失型缺失的預設(shè)邊界,邊界外的所有數(shù)據(jù)都刪失.在DS-a 中,針對右刪失缺失型數(shù)據(jù),其上界值設(shè)為43.5,表明刪失類型的(超)矩形窗為[10,43.5]×[5,43.5],其中10 和5 為小于其觀測值最小值的一個數(shù),來源于觀測值的矩形窗下界,并不表示刪失數(shù)據(jù)的下界,并觀察到約862 個數(shù)據(jù)點未刪失,并使用EMGM 算法和cenEMGM 算法進行實驗,如圖1 所示;若其還存在左刪失,如將其下界值設(shè)為15,形成雙邊刪失型缺失數(shù)據(jù),表明刪失類型的(超)矩形窗為 [15,43.5]×[15,43.5],約818 個數(shù)據(jù)點未刪失,如圖2 所示.類似地,在DS-b 中,針對左刪失缺失型數(shù)據(jù),其下界值設(shè)為0,表明刪失類型的(超)矩形窗為 [0,60]×[0,60],其中60 為大于其觀測值最大值的一個數(shù),來源于觀測值的矩形窗上界,并不表示刪失數(shù)據(jù)的上界,約484 個數(shù)據(jù)點未刪失,如圖2 所示;若其還存在右刪失,例如其上界值設(shè)為40,形成雙邊刪失型缺失數(shù)據(jù),表明刪失類型的(超)矩形窗為 [0,40]×[0,40],約241 左右的數(shù)據(jù)點未刪失,如圖3 所示.圖中小十字表示刪失后的數(shù)據(jù)點,‘o’ 和實心橢圓是每個成分在算法估計后的聚類中心和距離為1 的等高曲線.其距離使用成對馬氏(Mahalanobis)距離計算.‘+’ 和虛線橢圓表示高斯混合模型成分的真實聚類中心和等高曲線.

圖1 在數(shù)據(jù)集DS-a 右刪失上的兩種算法比較Fig.1 Comparison of the two algorithms on the dataset DS-a with right censoring

圖2 在數(shù)據(jù)集DS-a 雙邊刪失上的兩種算法比較Fig.2 Comparison of the two algorithms on the dataset DS-a with double-side censoring

圖3 在數(shù)據(jù)集DS-b 左刪失上的兩種算法比較Fig.3 Comparison of the two algorithms on the dataset DS-b with left censoring

圖1 顯示EMGM 算法和cenEMGM 算法在二維合成數(shù)據(jù)DS-a 右刪失上的實驗結(jié)果.EMGM 算法在該數(shù)據(jù)集上的結(jié)果(圖1 (a))顯示,‘o’ 和實心橢圓所表示的估計的聚類中心和距離為1 的等高曲線與 ‘+’ 和虛線橢圓表示高斯混合模型成分的真實聚類中心和等高曲線之間存在顯著差異.而cenEMGM 算法在該數(shù)據(jù)集上的結(jié)果(圖1(b))顯示,cenEMGM 算法估計的聚類中心和等高曲線與真實聚類中心和等高曲線之間的差異明顯減小,其結(jié)果明顯優(yōu)于EMGM 算法.

圖2 顯示EMGM 算法和cenEMGM 算法在二維合成數(shù)據(jù)DS-a 雙邊刪失上的實驗結(jié)果.EMGM算法在該數(shù)據(jù)集上的結(jié)果(圖2 (a))顯示,聚類中心和距離為1 的等高曲線比EMGM 算法(圖1 (a))明顯更接近于真實值.因為這里除了存在右刪失外,還存在左刪失.盡管缺失率更高,但觀測到的數(shù)據(jù)(未刪失部分)的均值更接近真實值.同時可見,cenEMGM 算法估計(圖2 (b))的聚類中心和真實聚類中心之間的差異也明顯更小,其結(jié)果進一步表明cenEMGM 算法在處理刪失數(shù)據(jù)聚類問題上明顯優(yōu)于EMGM 算法.

圖3 顯示EMGM 算法和cenEMGM 算法在二維合成數(shù)據(jù)DS-b 左刪失上的實驗結(jié)果.‘+’ 表示高斯混合模型成分的真實聚類中心,其中兩個已落在了值域的下界之外.EMGM 算法的結(jié)果(圖3 (a))顯示,其估計的聚類中心(‘o’)和等高曲線(實心橢圓)沒有超出值域的下界,表明估計值與對應的真實值之間存在顯著差異.而cenEMGM 算法的估計結(jié)果(圖3 (b))顯示,其估計的聚類中心和等高曲線與真實值之間的差異明顯更小.對于圖3 (b)圖中靠近Y2坐標軸的成分,盡管其估計值與真實值之間尚存在一些差異,但這一差異與EMGM 算法所表現(xiàn)出的差異已經(jīng)小很多,且另外兩個成分的估計值與真實值之間幾乎無差異,因此這些結(jié)果進一步表明cenEMGM 算法在這類數(shù)據(jù)聚類上更優(yōu)于EMGM 算法.

圖4 顯示EMGM 算法和cenEMGM 算法在二維合成數(shù)據(jù)DS-b 雙邊刪失上的實驗結(jié)果.三個成分的聚類中心真實值(‘+’)都在下界或上界之外.EMGM 算法在該數(shù)據(jù)集上的結(jié)果(圖4 (a))顯示,三個成分的估計的聚類中心和距離為1 的等高曲線與真實值之間都存在顯著差異.與此相反,cenEMGM算法在該數(shù)據(jù)集上的結(jié)果(圖4 (b))顯示,其估計值也可以位于上下界之外,更接近真實聚類中心和等高曲線,即估計值與真實值之間的差異明顯變小.結(jié)果表明cenEMGM 算法在處理這類刪失數(shù)據(jù)聚類時明顯優(yōu)于EMGM 算法.

圖4 在數(shù)據(jù)集DS-b 雙邊刪失上的兩種算法比較Fig.4 Comparison of the two algorithms on the dataset DS-b with double-side censoring

此外,進行 100 次重復實驗,記錄多次實驗結(jié)果在KLD 值與AIC 值上的平均值和方差.實驗合成數(shù)據(jù)集真實分布和估計分布之間的KLD 值見表1,對于參數(shù)估計的兩種算法AIC 值比較見表2.結(jié)果表明,對于兩種算法在同一數(shù)據(jù)集上的表現(xiàn),不論是KLD 值還是AIC 值,cenEMGM 算法的值都小于對應EMGM 算法的值,說明在刪失型缺失數(shù)據(jù)參數(shù)估計上cenEMGM 算法優(yōu)于EMGM 算法.對于同一算法在不同數(shù)據(jù)集上的表現(xiàn),因雙邊刪失比對應的單邊刪失因缺失而擁有更少的樣本數(shù)據(jù),雙邊刪失的AIC 值小于對應的單側(cè)刪失的AIC 值.

表1 實驗合成數(shù)據(jù)集真實分布和估計分布之間的KLD 值Table 1 Kullback-Leibler divergence (KLD)between the true densities and the estimated densities of the synthetic data set

表2 實驗合成數(shù)據(jù)集參數(shù)估計的兩種算法AIC 比較Table 2 AIC comparison of the two estimation algorithms on the synthetic data set

4.2 真實數(shù)據(jù)分析

數(shù)據(jù)來源于某大型醫(yī)院信息系統(tǒng)中的臨床數(shù)據(jù)[4].這些數(shù)據(jù)樣本包括554 個相關(guān)屬性,其中有106 個建檔屬性、23 個檢驗數(shù)據(jù)屬性、157 個來自實驗室信息系統(tǒng)的試驗結(jié)果屬性以及268 個電子健康檔案中病案首頁的屬性.根據(jù)醫(yī)學領(lǐng)域?qū)＜乙庖姾臀墨I進行屬性篩選,經(jīng)過數(shù)據(jù)清理后所得數(shù)據(jù)集包括50 個屬性,具體包括年齡、婚齡、孕婦體重指數(shù)、紅細胞計數(shù)、谷氨酰轉(zhuǎn)肽酶、空腹血糖水平值等屬性.根據(jù)驗證的目的,這里所使用的數(shù)據(jù)集為原臨床數(shù)據(jù)集中提取的包含4 個屬性的數(shù)據(jù).這些屬性具體為關(guān)于孕婦在篩檢妊娠期糖尿病過程中的血糖水平值和醫(yī)生給出的診斷結(jié)果,即是否患有妊娠期糖尿病.其中包括關(guān)于血糖水平值的3 個屬性分別為口服糖耐量試驗中的空腹血糖水平值(Fasting blood sugar level,FBSL)、1 小時血糖水平值(1h-blood sugar level,1h-BSL)和2 小時后的血糖水平值.根據(jù)國際妊娠合并糖尿病研究組織建議,妊娠期糖尿病的診斷標準為[4],空腹血糖水平值高于5.1 mmol/L、1 小時血糖水平值高于10 mmol/L 和2 小時血糖水平值高于8.5 mmol/L,滿足以上三項中的任一項即診斷為患有妊娠期糖尿病,數(shù)據(jù)記錄聚類為患病簇,否則為正常簇.在電子病歷記錄與數(shù)據(jù)聯(lián)結(jié)整合中,小于等于10 mmol/L 的血糖水平值記錄為原始測量值,而高于10 mmol/L 的空腹血糖水平值和1 小時血糖水平值的數(shù)據(jù)被記錄為“>10mmol/L”型刪失型數(shù)據(jù).雖然這些刪失型數(shù)據(jù)能夠為診斷結(jié)果提供直接的臨床證據(jù),但是這些數(shù)據(jù)的刪失對于進一步探索關(guān)于妊娠期糖尿病的風險因子,以及這些因子對血糖水平值影響的因果關(guān)系研究構(gòu)成困難.又因妊娠期糖尿病的主要治療方案包括膳食改變、增加鍛煉甚至胰島素等的藥物治療[28],但這些治療方案對以血糖水平值作為結(jié)果的影響作用大小是有差異的.為后續(xù)研究這些影響作用,在使用這些刪失型的血糖水平值數(shù)據(jù)時,需要對這些數(shù)據(jù)的分布參數(shù)進行較為精確的估計.本文的聚類算法正是針對這些刪失型數(shù)據(jù)提供分布參數(shù)的估計.

從原數(shù)據(jù)中選擇了917 例數(shù)據(jù)進行數(shù)值計算,其中756 例樣本屬于正常簇,161 例樣本屬于患病簇.在917 例樣本數(shù)據(jù)中,以空腹血糖水平值和1小時血糖水平值進行分析,發(fā)現(xiàn)78 例樣本數(shù)據(jù)屬于刪失型數(shù)據(jù),主要存在于1 小時血糖水平值上.對這一數(shù)據(jù)集,分別采用EMGM 算法和cenEMGM算法進行高斯混合聚類,結(jié)果如圖5 所示.

圖5 在血糖測試數(shù)據(jù)右刪失上兩種算法比較Fig.5 Comparison of the two algorithms on the dataset of blood sugar tests with right-side censoring

圖5 顯示了EMGM 算法和cenEMGM 算法在刪失型血糖水平值數(shù)據(jù)上的聚類結(jié)果.橫坐標為空腹血糖水平值,縱坐標為1 小時血糖水平值,其樣本數(shù)據(jù)關(guān)于“>10 mmol/L”刪失.真實數(shù)據(jù)中一個成分的聚類中心真實值(‘+’)在樣本數(shù)據(jù)所展示的范圍內(nèi),為(4.96,7.16);另一個成分的聚類中心真實值(‘+’)在樣本數(shù)據(jù)的上界之外,為(6.09,11.16),即中心值在1 小時血糖水平值上“>10 mmol/L”.圖5 (a)顯示EMGM 算法在該數(shù)據(jù)集上存在一個成分的估計聚類中心和距離為1 的等高曲線與真實值之間存在顯著差異,即估計值所在的聚類中心在1 小時血糖水平值以下,而真實值所在的聚類中心在1 小時血糖水平值以上.不同的是,圖5 (b)顯示cenEMGM 算法在該數(shù)據(jù)集上的估計值也可以位于上界之外,使得其更接近真實聚類中心,說明估計值與真實值之間的差異明顯變小.在模型檢驗準則上,對于這一真實數(shù)據(jù)集,EMGM 算法在真實分布與估計分布之間的KLD 值(12.7)高于cenEMGM算法的KLD 值(9.1),同時后者的AIC 值(4 263)低于前者的AIC 值(4 366).因此,這些結(jié)果說明cenEMGM 算法在處理真實的刪失數(shù)據(jù)聚類時優(yōu)于EMGM 算法.

此外,為進一步驗證方法的有效性,對于真實數(shù)據(jù)調(diào)整刪失率進行拓展,動態(tài)改變刪失率而進行計算,并對聚類中心、AIC 與KLD 值進行定量對比,如表3 所示.

表3 真實數(shù)據(jù)及其拓展數(shù)據(jù)的兩種算法比較Table 3 Comparison of the two algorithms with the real data and its extended data

表3 結(jié)果表明,當右側(cè)刪失率從8.51%增加到11.67%時,兩種算法的聚類中心估計值與真實值(4.96,7.16)和(6.09,11.16)之間的差異增大,KLD 值與AIC 值減小.cenEMGM 算法的KLD 值與AIC 值比EMGM 算法的對應值小,說明其在處理刪失數(shù)據(jù)聚類時仍然優(yōu)于EMGM 算法.當將數(shù)據(jù)拓展為雙邊刪失型數(shù)據(jù)時,即在右邊刪失的基礎(chǔ)上增加左邊刪失6.54 %,總體上刪失15.05 %時,兩種算法的聚類中心估計值與真實值之間的差異進一步增大,且KLD 值增大而AIC 值減小.總體上,隨著刪失率的增加,算法處理的能力在一定程度上逐漸減弱,但是cenEMGM 算法的聚類中心估計值與真實值相對更接近,且KLD 值與AIC 值比EMGM算法的對應值更小,進一步說明其通過聚類在處理刪失數(shù)據(jù)的參數(shù)估計時仍然優(yōu)于EMGM 算法.

5 結(jié)論

刪失型數(shù)據(jù)處理特別是在機器學習或數(shù)據(jù)挖掘等數(shù)據(jù)處理中,作為工程實踐和管理中數(shù)據(jù)處理的焦點問題.由于刪失數(shù)據(jù)處理的知識有限性,需要根據(jù)刪失模式制定合適的算法模型.盡管當前數(shù)據(jù)智能處理所面臨的數(shù)據(jù)規(guī)模較大,但選取高價值的實驗數(shù)據(jù)或稀有事件等所面臨的刪失數(shù)據(jù)處理仍然顯得較為重要.然而,現(xiàn)有的缺失數(shù)據(jù)處理問題主要集中在隨機缺失,對非隨機缺失下的刪失型數(shù)據(jù)研究不深,因此本文根據(jù)估計算法的有效性理論,針對刪失數(shù)據(jù)期望最大化的高斯混合聚類算法(cenEMGM),通過關(guān)于得分向量期望的方程得出算法估計的最優(yōu)參數(shù).與嵌套標準EM 的高斯混合聚類算法(EMGM)相比,本方法根據(jù)刪失數(shù)據(jù)的指示變量調(diào)整樣本似然函數(shù),進而改進參數(shù)估計的期望最大化算法,使得高斯混合聚類模型參數(shù)估計準確性更高,AIC 信息準則值更小,聚類效果更好.并通過數(shù)值實驗論證了本方法相對于EMGM 算法的優(yōu)越性.更多類型數(shù)據(jù)中的刪失型缺失機制(模式)識別、不同刪失情形下多種算法有效性分析及其高斯混合聚類算法拓展是下一步工作重點.