吳 笛,劉 文
(武漢理工大學(xué)理學(xué)院,湖北武漢 430070)
概率密度估計(jì)既是傳統(tǒng)的概率論與數(shù)理統(tǒng)計(jì)的重點(diǎn),也是統(tǒng)計(jì)學(xué)習(xí)理論的重要研究內(nèi)容[1]。在解決統(tǒng)計(jì)學(xué)習(xí)問題的傳統(tǒng)模式中,模式識(shí)別和回歸估計(jì)都建立在密度估計(jì)的基礎(chǔ)之上[2]。且概率密度估計(jì)在實(shí)際中有廣泛的應(yīng)用,如電子器件壽命估計(jì)和排隊(duì)論等。但在實(shí)際應(yīng)用中很多時(shí)候并不知道概率密度的分布,這時(shí)可根據(jù)樣本點(diǎn)進(jìn)行回歸分析得到實(shí)際概率密度的一個(gè)近似估計(jì)。目前的概率密度估計(jì)方法主要分為參數(shù)估計(jì)和非參數(shù)估計(jì)兩大類。參數(shù)估計(jì)方法具有較大的局限性,其前提是已知數(shù)據(jù)密度符合某種分布,但問題在于如何確定密度函數(shù)中的參數(shù),這種方法強(qiáng)烈依賴于前提假設(shè),一旦假設(shè)錯(cuò)誤,估計(jì)值就無法較好地反映真實(shí)值。非參數(shù)估計(jì)具有更廣的應(yīng)用范圍,并且已經(jīng)得到了廣泛研究,提出了核估計(jì)和近鄰估計(jì)等方法。然而上述方法需要用到所有的訓(xùn)練樣本以估計(jì)出概率密度,當(dāng)訓(xùn)練樣本非常多時(shí),計(jì)算量非常大并不實(shí)用[3]。長期以來,人們希望尋求一種既可減少計(jì)算量又可保持估計(jì)精度的求解方法。正則化理論[4]的引入較好地解決了這一問題,利用懲罰項(xiàng)來得到偏移-方差平衡,防止概率密度過度擬合的發(fā)生。筆者將概率密度的估計(jì)轉(zhuǎn)化成一類線性算子方程問題的求解,并根據(jù)算子方程核矩陣奇異值的性質(zhì),構(gòu)建了概率密度估計(jì)的TSVD正則化方法[5-7],同時(shí)將其與估計(jì)概率密度的K+線性Bregman迭代正則化方法[8-10]進(jìn)行比較分析。
根據(jù)概率密度的定義,密度估計(jì)是一個(gè)求解第一類Fredholm積分方程的問題,如:
對式(1),概率密度函數(shù)p(x)須同時(shí)滿足如下條件:
在實(shí)際應(yīng)用中,往往只考慮區(qū)間[a,b]上的概率密度估計(jì)問題,對于已知的獨(dú)立同分布數(shù)據(jù)x=(x1,x2,…,xM),此時(shí)的經(jīng)驗(yàn)分布函數(shù) FM(x)可由x來構(gòu)造,即:
其中指示函數(shù)θ(x)為:
筆者的目的是將概率密度估計(jì)問題轉(zhuǎn)化為線性算子方程,然后利用截?cái)嗥娈愔捣纸?truncated singular value decomposition,TSVD)方法進(jìn)行求解。根據(jù)式(2)~式(4)的定義,可將式(1)離散成如下的線性算子方程:
其中,核矩陣K為0-1矩陣。在實(shí)際問題中右端項(xiàng)F(x)通常包含噪聲,給右端項(xiàng)F(x)一個(gè)很小的擾動(dòng)都會(huì)使得式(5)的解產(chǎn)生巨大的震蕩現(xiàn)象,無法對概率密度函數(shù)p(x)進(jìn)行準(zhǔn)確估計(jì)。這時(shí)式(1)和式(5)就變?yōu)橐粋€(gè)不適定問題。因此,筆者將采用計(jì)算量較小的TSVD正則化方法估計(jì)概率密度函數(shù)p(x)。
假設(shè)實(shí)際觀察中式(5)右端項(xiàng)為Fδ(x),則:
其中,δ為誤差水平。
根據(jù)核矩陣K基于指示函數(shù)θ(x)的定義,對核矩陣K進(jìn)行奇異值分解(singular value decomposition,SVD),以便分析其性質(zhì),構(gòu)建求解式(5)的正則化方法。以100×100大小的核矩陣K為例,對應(yīng)奇異值的分布情況如圖1所示。
圖1 核矩陣K(100×100)奇異值的分布情況
由圖1可知,核矩陣K的能量主要集中在少數(shù)幾個(gè)數(shù)值較大的奇異值上。KIRSCH指出不適定問題解的不穩(wěn)定性的根源在于緊算子的奇異值有趨于零的性質(zhì),由此引入TSVD正則化方法來減弱或?yàn)V掉奇異值趨于零的性質(zhì)對解的穩(wěn)定性影響,并對核矩陣K進(jìn)行奇異值分解得到:
其中,奇異值 si=S(i,i)(i=1,2,…,M),且ui和vi分別為酉矩陣U和V的行向量,滿足以下性質(zhì):
根據(jù)式(7),此時(shí)式(5)可轉(zhuǎn)化為:
利用行向量ui和vi的性質(zhì),根據(jù)式(9)和式(10)可求得概率密度函數(shù)為:
式(11)中的奇異值si隨著i增加逐漸趨于0,而1/si將趨于無窮大。此時(shí)式(5)的解不連續(xù)依賴于右端的數(shù)據(jù),當(dāng)右端數(shù)據(jù)F(x)存在誤差時(shí),方程的解與真解間存在較大的誤差,無法得到準(zhǔn)確的概率密度函數(shù)p(x)。因此,有必要在保證概率密度函數(shù)估計(jì)值滿足ε的基礎(chǔ)上,對奇異值si進(jìn)行截?cái)鄟肀WC方程解的穩(wěn)定性,則此時(shí)的概率密度函數(shù)估計(jì)值為:
其中,k(1≤k≤M)為正則化參數(shù),以控制奇異值si的截?cái)唷?/p>
研究重點(diǎn)在于根據(jù)概率密度的定義,將概率密度的估計(jì)轉(zhuǎn)化成第一類Fredholm算子方程的求解問題,然后利用能保證正則解誤差具有漸進(jìn)最優(yōu)階的TSVD正則化方法進(jìn)行求解,以體現(xiàn)TSVD正則化方法在處理概率密度估計(jì)問題中的優(yōu)越性。對于TSVD正則化方法的誤差估計(jì)、正則化參數(shù)先驗(yàn)選取及其收斂性和收斂階分析,黃小為等已作了深入研究。針對式(5),其給出了TSVD正則化方法在概率密度估計(jì)問題中的收斂定理1。
定理1 設(shè)式(5)的精確解p+(x)=(K*·K)vz∈R(K*K)v,z∈X(p)且‖z‖2≤E,則對于TSVD正則化方法,若選取正則化參數(shù)k=k(δ)使得成立,則有:
其中:K*為K的伴隨算子;為概率密度的估計(jì)值;E∈R為常數(shù);X(p)為實(shí)Hilbert空間,且滿足p(x)∈X(p)。該定理說明,TSVD正則化方法可保證概率密度估計(jì)值的誤差具有漸進(jìn)最優(yōu)階,且數(shù)值計(jì)算簡單易行,其具體證明過程及相關(guān)性質(zhì)分析請參見文獻(xiàn)[5-6]。
由式(13)產(chǎn)生100個(gè)服從正態(tài)分布的隨機(jī)樣本:
其中,均值分別為μ1=0和μ2=2,標(biāo)準(zhǔn)差分別為 σ1=1和 σ2,變量 x的取值范圍為[-5,5]。
為了加快求解式(14)的收斂速度,在此利用K+線性Bregman迭代正則化方法對其進(jìn)行求解,具體迭代過程如下:
其中,K+=KT(KKT)-1,參數(shù) α >0,μ >0,并定義向量閾值算子為:
根據(jù)前文對TSVD正則化方法的介紹,得到求解概率密度函數(shù)的TSVD方法的迭代過程如下:
迭代循環(huán):
假設(shè)右端數(shù)據(jù)F(x)分別受均值為0,方差為0.001、0.010和0.050的高斯白噪聲的干擾,并利用線性Bregman迭代正則化方法和TSVD正則化方法分別對式(5)進(jìn)行求解。其中線性Bregman迭代正則化方法中的參數(shù)分別為α=0.8和μ=0.005,則實(shí)驗(yàn)結(jié)果如圖2和表1所示。
圖2 實(shí)驗(yàn)結(jié)果比較分析(噪聲方差為0.001)
由圖2可知,當(dāng)式(5)的右端項(xiàng)F(x)受噪聲干擾時(shí),TSVD正則化方法可以較好地估計(jì)概率密度函數(shù)pT(x),而K+線性Bregman迭代正則化方法易受噪聲的干擾,估計(jì)的概率密度函數(shù)pB(x)存在較為明顯的波動(dòng),不能較好地逼近真實(shí)的概率密度函數(shù)p(x)。由表1可知,當(dāng)噪聲方差分別為0.001、0.010和0.050時(shí),基于 TSVD的概率密度估計(jì)要好于線性Bregman迭代正則化方法,顯示出更強(qiáng)的穩(wěn)定性,且更易于編程實(shí)現(xiàn)。
表1 線性Bregman與TSVD的實(shí)驗(yàn)結(jié)果比較
其中,假設(shè) f∈IRN:=IR{1,2,…,N},對于 1≤p <∞,則lp范數(shù)可定義為:
筆者在分析核矩陣K的奇異值性質(zhì)的基礎(chǔ)上,構(gòu)建了概率密度估計(jì)的TSVD正則化方法,將概率密度估計(jì)問題轉(zhuǎn)化成l1正則化問題,并構(gòu)建了求解該最優(yōu)化問題的K+線性Bregman迭代正則化方法。由實(shí)驗(yàn)結(jié)果可以看出,用TSVD估計(jì)的概率密度可以較好地逼近真實(shí)的概率密度,較線性Bregman迭代正則化方法表現(xiàn)出更好的估計(jì)效果,且更易于編程實(shí)現(xiàn),在不同噪聲水平下表現(xiàn)出更強(qiáng)的穩(wěn)定性。
[1]VLADIMIR N V.統(tǒng)計(jì)學(xué)習(xí)理論的本質(zhì)[M].張學(xué)工,譯.北京:清華大學(xué)出版社,2000:12-98.
[2]VLADIMIR N V.統(tǒng)計(jì)學(xué)習(xí)理論[M].許建華,張學(xué)工,譯.北京:電子工業(yè)出版社,2004:87-125.
[3]FUKUNAGA K,HAYES R R.The reduced Parzen classifier[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,1989,11(4):423-425.
[4]ENGL W,HANKE M,NEUBAUER A.Regularization of inverse problems[M].Dordrecht:Kluwer Acdemic Publishers,1996:19-76.
[5]黃小為,吳傳生,朱華平.求解不適定問題的TSVD正則化方法[J].武漢理工大學(xué)學(xué)報(bào),2005,27(2):90-92.
[6]黃小為,吳傳生,李卓球.TSVD正則化方法的參數(shù)選取及數(shù)值計(jì)算[J].華中師范大學(xué)學(xué)報(bào):自然科學(xué)版,2006,40(2):154-157.
[7]BOUHAMIDI A,JBILOU K,REICHEL L,et al.An extrapolated TSVD method for linear discrete ill-posed problems with Kronecker structure[J].Linear Algebra and its Applications,2011(7):1677-1688.
[8]CAI J F,OSHER S,SHEN Z W.Linearized Bregman iterations for compressed sensing[J].Mathematics of Computation,2009(78):1515-1536.
[9]張慧,成禮智.A-線性Bregman迭代算法[J].計(jì)算數(shù)學(xué),2010,32(1):97-104.
[10]YIN W T,OSHER S,GOLDFARB D,et al.Bregman iterative algorithms for l1-minimization with applications to compressed sensing[J].SIAM Journal on Imaging Sciences,2008,1(1):143-168.