兩種概率密度導(dǎo)數(shù)估計(jì)方法的比較研究*

王文武， 孔德茹

(曲阜師范大學(xué)統(tǒng)計(jì)與數(shù)據(jù)科學(xué)學(xué)院,273165,山東省曲阜市)

0 引 言

在探索樣本總體信息時(shí)，概率分布函數(shù)及密度函數(shù)為統(tǒng)計(jì)研究的基本問題. 相應(yīng)地，其估計(jì)理論已經(jīng)得到較好的完善，但密度函數(shù)的導(dǎo)數(shù)卻鮮少被關(guān)注. 雖然密度導(dǎo)數(shù)沒有受到足夠的重視，但是其作用不容小覷. 例如，一階導(dǎo)函數(shù)為密度函數(shù)的斜率，反映了密度函數(shù)的波動(dòng)幅度；二階導(dǎo)函數(shù)為密度函數(shù)的曲率，反映了密度函數(shù)的彎曲程度以及光滑程度. 密度函數(shù)的導(dǎo)數(shù)不僅直觀反映了概率密度的圖像走勢(shì)，而且被廣泛應(yīng)用于經(jīng)濟(jì)學(xué)與統(tǒng)計(jì)學(xué)等領(lǐng)域：在經(jīng)濟(jì)學(xué)中，常需要GDP(國(guó)民生產(chǎn)總值)進(jìn)行建模并做出導(dǎo)數(shù)估計(jì)，利益的最值問題也需要導(dǎo)數(shù)或者偏導(dǎo)數(shù)的計(jì)算[1]；在統(tǒng)計(jì)學(xué)中，概率密度函數(shù)的導(dǎo)數(shù)也發(fā)揮了重要作用，例如，導(dǎo)數(shù)作為重要的一環(huán)，出現(xiàn)在最速下降及牛頓等算法中.

因此，密度函數(shù)導(dǎo)數(shù)需要得到關(guān)注并且估計(jì)方法也應(yīng)被深入研究. 從密度導(dǎo)數(shù)本身出發(fā)，H?rdle等提出了平均導(dǎo)數(shù)估計(jì)并給出了實(shí)例應(yīng)用[2]；次年，又提出了一種廣義交叉驗(yàn)證，通過利用差商的核平滑來估計(jì)一階導(dǎo)數(shù)[3]. 2009年，Hansen在《Lecture Notes on Nonparametrics》書中總結(jié)了密度函數(shù)導(dǎo)數(shù)的簡(jiǎn)單估計(jì)方法——根據(jù)導(dǎo)數(shù)定義對(duì)密度函數(shù)的導(dǎo)數(shù)進(jìn)行核估計(jì)[4]. 2013年，Chacón等提出了多變量數(shù)據(jù)下的核密度導(dǎo)數(shù)估計(jì)帶寬選擇器，優(yōu)化了核方法在導(dǎo)數(shù)估計(jì)中的應(yīng)用[5]. 近年來，回歸函數(shù)的導(dǎo)數(shù)估計(jì)受到了重視，這為密度函數(shù)導(dǎo)數(shù)估計(jì)提供了新的思路.Charnigo等在2011年提出通過GCp標(biāo)準(zhǔn)選擇合適的調(diào)整參數(shù)從而提高導(dǎo)數(shù)估計(jì)的效果，并且該方法可用于任何非參數(shù)回歸方法[6]. 2013年，Brabanter等闡述了經(jīng)驗(yàn)一階導(dǎo)數(shù)在局部多項(xiàng)式回歸框架中的使用，并推導(dǎo)了偏差和方差[7]. 2016年，Wang等將局部加權(quán)最小二乘方法以及差分方法應(yīng)用到導(dǎo)數(shù)估計(jì)中，得到了更小的均方誤差[8]. 2019年，Cattaneo等采用局部多項(xiàng)式回歸估計(jì)各階導(dǎo)數(shù)，給出了較為全面的軟件包，并應(yīng)用到斷點(diǎn)回歸的臨界點(diǎn)檢驗(yàn)中[9，10].

在估計(jì)概率密度導(dǎo)數(shù)時(shí)，盡管核方法和局部多項(xiàng)式回歸方法經(jīng)常被使用，但在應(yīng)用中并不清楚哪種方法更有優(yōu)勢(shì). 故本文從理論和模擬兩個(gè)角度進(jìn)行比較研究. 首先，對(duì)這兩種非參數(shù)估計(jì)方法的理論知識(shí)進(jìn)行了梳理總結(jié)并做出了合理的比較. 其次，選擇均勻分布、正態(tài)分布、卡方分布，在不同樣本量下對(duì)概率密度導(dǎo)數(shù)進(jìn)行數(shù)值模擬并根據(jù)經(jīng)驗(yàn)積分均方誤差比較估計(jì)效果. 為了更充分地對(duì)比估計(jì)性質(zhì)，每種方法均通過兩個(gè)軟件包實(shí)現(xiàn). 通過觀察兩種方法在4個(gè)R包下均方誤差值的變化以及估計(jì)結(jié)果圖，最終證實(shí)局部多項(xiàng)式回歸方法相較于核方法具有更加優(yōu)良的性質(zhì).

1 兩種估計(jì)方法

在總體信息模糊的情況下，常采用非參數(shù)方法來達(dá)到目的. 常見的密度函數(shù)估計(jì)方法有直方圖法、核方法、局部多項(xiàng)式回歸法、樣條估計(jì)法以及最近鄰估計(jì)法等，其中核方法和局部多項(xiàng)式回歸法最為常用[11]. 本部分將簡(jiǎn)單介紹核方法和局部多項(xiàng)式回歸方法的一般理論.

1.1 核方法

首先采用核方法估計(jì)概率密度函數(shù)，即核密度估計(jì). 核密度估計(jì)通過核函數(shù)來估計(jì)未知的概率密度，其基本思想是按照其余觀測(cè)值距估計(jì)點(diǎn)的遠(yuǎn)近對(duì)該估計(jì)點(diǎn)賦予不同的權(quán)重：距估計(jì)點(diǎn)較近的觀測(cè)值包含該估計(jì)點(diǎn)較多的信息，所以賦予更大的權(quán)重；距估計(jì)點(diǎn)較遠(yuǎn)的觀測(cè)值包含信息較少，故賦予較小的權(quán)重甚至0權(quán)重[10].

其中k(·)為核函數(shù).通過核密度估計(jì)的最終形式可以看出核密度估計(jì)結(jié)果實(shí)際上是對(duì)協(xié)變量x附近的觀測(cè)值進(jìn)行加權(quán)平均，核函數(shù)為權(quán)重函數(shù). 另外，核函數(shù)需要滿足非負(fù)性、對(duì)稱性、正則性等條件. 在有關(guān)核的應(yīng)用中，均勻核函數(shù)、Bieweight核函數(shù)以及Gaussian核函數(shù)較為常見.

在核密度估計(jì)的基礎(chǔ)上由導(dǎo)數(shù)定義，對(duì)估計(jì)得到的密度函數(shù)求r階導(dǎo)[4]，即

1.2 局部多項(xiàng)式回歸方法

局部多項(xiàng)式回歸方法作為非參數(shù)統(tǒng)計(jì)的一種重要估計(jì)方法，應(yīng)用范圍極其廣泛. 其核心思想是在局部區(qū)域用某一多項(xiàng)式近似目標(biāo)函數(shù)，從而得到該函數(shù)及其各階導(dǎo)數(shù)的估計(jì)[12]. 估計(jì)過程中用到的多項(xiàng)式通常是采用泰勒展開得到的. 因此，當(dāng)估計(jì)同一函數(shù)時(shí)，采用不同階數(shù)的多項(xiàng)式將得到不同的估計(jì)效果[13].

在介紹局部多項(xiàng)式回歸方法估計(jì)概率密度導(dǎo)數(shù)之前，作出兩個(gè)一般假定[9].

假定1假定x1,…,xn是來自分布函數(shù)為F(x)的一組隨機(jī)觀測(cè)值，概率密度函數(shù)為f(x)，其中，協(xié)變量x的定義域?yàn)閇xlower,xupper](xupper≠xlower)，分布函數(shù)F(x)至少p+1階可導(dǎo)并且F(r+1)(x)=f(r)(x),r≤p.

假定2為了優(yōu)化估計(jì)效果的表達(dá)式，本文要求核函數(shù)k(·)除滿足一般性質(zhì)外，還需要滿足該條件：無論核函數(shù)為哪種形式，均只在[-1,1]上有意義，其余位置點(diǎn)為0.

注對(duì)于假定1的條件，幾乎所有數(shù)據(jù)都可以滿足.由于距離x越近的觀測(cè)點(diǎn)包含的信息越多，給出的權(quán)重也應(yīng)越大[14]，故假定2把核函數(shù)的取值范圍限制在了[-1,1]上，即在[x-h,x+h]上核函數(shù)非零.對(duì)于高斯核函數(shù)，雖不能滿足該條件，但在[-1,1]之外的區(qū)域，核函數(shù)取值接近于0.因此可認(rèn)為高斯核函數(shù)滿足假定2中的限制.此外，當(dāng)帶寬接近于0時(shí)，假定2將核函數(shù)限制在估計(jì)點(diǎn)附近，從而可優(yōu)化邊界估計(jì)效果.

目標(biāo)函數(shù)F(x)為x的分布函數(shù)，yi為F(x)的樣本值，則有局部最小二乘估計(jì)

1.3 兩種方法的估計(jì)效果

在衡量估計(jì)效果時(shí)，常采用積分均方誤差MISE(Mean Integrated Squared Error),

從上式可以看出均方誤差既包含了方差也包含了偏差，能夠很好地度量估計(jì)性質(zhì)：積分均方誤差越小說明估計(jì)性質(zhì)越好.

通過變量代換和r次分部積分可以得到核方法下導(dǎo)數(shù)估計(jì)的偏差和方差[3]，從而積分均方誤差為

通過使積分均方誤差達(dá)到最小，得到密度函數(shù)各階導(dǎo)數(shù)的最優(yōu)帶寬為

h=Cr,v(k,f)n-1/(1+2r+2v)，Cr,v(k,f)=R(f(r+v))-1/(1+2r+2v)Ar,v(k)，

通過上述結(jié)果可以看出，密度函數(shù)導(dǎo)數(shù)核估計(jì)使用的最優(yōu)帶寬為O(n-1/(1+2r+2v)),而密度函數(shù)使用的最優(yōu)帶寬為O(n-1/(1+2v)).這是因?yàn)楣饣烙?jì)的最優(yōu)值取決于估計(jì)對(duì)象和分析目的，故估計(jì)的均方誤差會(huì)發(fā)生變化，最優(yōu)帶寬也會(huì)隨之發(fā)生變化[16].

到得草原中來，綠色王國(guó)，花草天堂，但在漫無涯際之中，此番我的觀光焦點(diǎn)，不全在草原山崗河流等大景致，反而破天荒地“微觀化”了。

引理1假定密度函數(shù)f(x)至少r+v+1階可導(dǎo)，核函數(shù)滿足k(s)(∞)=0,k(s)(-∞)=0，s=0,1,2,…,r,h→0，當(dāng)n→∞時(shí)，nh→∞，則核方法的漸近偏差與漸近方差分別為

引理2滿足假定1和假定2，2≤r≤p，h→0，當(dāng)n→∞時(shí)，nh→∞且nh2p+1=O(1)，則局部多項(xiàng)式回歸方法的漸近偏差與漸近方差分別為

因?yàn)闈M足林德伯格條件[17]，所以由林德伯格-萊維中心極限定理可得：核方法和局部多項(xiàng)式回歸方法下估計(jì)的概率密度導(dǎo)數(shù)服從漸近正態(tài)分布，關(guān)于更加細(xì)致的證明詳見參考文獻(xiàn)[4，9].

核方法估計(jì)密度導(dǎo)數(shù)的漸近偏差與漸近方差的證明過程詳見文獻(xiàn)[4]，局部多項(xiàng)式回歸方法估計(jì)性質(zhì)的證明在Cattaneo和Jansson (2019)中有詳盡的敘述[9]，這里僅對(duì)結(jié)果進(jìn)行分析. 從式子上看，核方法與局部多項(xiàng)式回歸方法下的估計(jì)偏差與方差并不能合理比較，因此在第2節(jié)中，通過數(shù)值模擬結(jié)合實(shí)際比較二者的估計(jì)性質(zhì).

2 數(shù)值模擬

本節(jié)通過數(shù)值模擬的方式，以密度函數(shù)的一階導(dǎo)數(shù)估計(jì)為例，分別采用核方法與局部多項(xiàng)式回歸方法估計(jì)，對(duì)比同一分布下2種方法的經(jīng)驗(yàn)均方誤差，并且直觀呈現(xiàn)密度導(dǎo)數(shù)估計(jì)結(jié)果[13].

2.1 均勻分布

從圖1可以看出，在樣本量一致的情況下，無論使用哪個(gè)包估計(jì)，局部多項(xiàng)式回歸方法的估計(jì)積分均方誤差始終小于核估計(jì)的誤差，故對(duì)導(dǎo)數(shù)估計(jì)，局部多項(xiàng)式回歸方法優(yōu)于核方法. 橫向比較趨勢(shì)，局部多項(xiàng)式回歸方法隨著樣本量增加，均方誤差存在大幅度或小幅度的減少；而核方法的兩個(gè)估計(jì)結(jié)果卻越來越差. 在不同樣本量下，lpdensity包的經(jīng)驗(yàn)均方誤差值最小，且隨著樣本量增加，該估計(jì)的均方誤差越來越小，即大樣本量下，會(huì)有更小的均方誤差. 此外，通過圖2，可以看出在4個(gè)軟件包中，lpdensity的估計(jì)結(jié)果最接近于真實(shí)值0. 雖然局部多項(xiàng)式回歸方法在邊界處的估計(jì)仍存在進(jìn)步空間，但相較于核方法在-1和1附近的估計(jì)結(jié)果，局部多項(xiàng)式估計(jì)的表現(xiàn)還是明顯更優(yōu)的. 局部多項(xiàng)式對(duì)內(nèi)點(diǎn)的估計(jì)效果也優(yōu)于核估計(jì).

圖1 均勻分布導(dǎo)數(shù)估計(jì)的均方誤差

圖2 均勻分布導(dǎo)數(shù)估計(jì)結(jié)果

2.2 正態(tài)分布

采用核估計(jì)方法和局部多項(xiàng)式估計(jì)方法，選擇均值為0，方差為1的正態(tài)分布進(jìn)行100次數(shù)值模擬，給出樣本量為2 000時(shí)，呈現(xiàn)兩種方法中最優(yōu)的估計(jì)結(jié)果見圖3和圖4.

圖3 正態(tài)分布導(dǎo)數(shù)估計(jì)的均方誤差

橫向比較圖3，當(dāng)樣本量逐漸增大，4個(gè)軟件包下的估計(jì)均方誤差均存在減小的趨勢(shì)，即局部多項(xiàng)式回歸方法和核方法的估計(jì)效果越來越好. ks的估計(jì)均方誤差隨著樣本量增加越來越接近0，在圖4中的估計(jì)結(jié)果也十分優(yōu)秀. 可見對(duì)正態(tài)分布來說，合適的核方法可以得出較好的估計(jì). 每種樣本量下，局部多項(xiàng)式估計(jì)的兩個(gè)包的MISE雖然不是最優(yōu)，但也不是最差. 另外，由于KernSmooth包locpoly函數(shù)的帶寬需要提前設(shè)定，這可能導(dǎo)致了經(jīng)驗(yàn)均方誤差有小幅增加. 通過圖4可以看出在正態(tài)分布的密度導(dǎo)數(shù)估計(jì)中，ks包中的kdde函數(shù)在邊界附近和內(nèi)點(diǎn)的擬合均明顯優(yōu)于lpdensity的估計(jì). 雖然lpdensity的表現(xiàn)并不最優(yōu)，但得到的估計(jì)曲線較為光滑，趨勢(shì)也接近真實(shí)概率密度導(dǎo)函數(shù).

圖4 正態(tài)分布導(dǎo)數(shù)估計(jì)結(jié)果

2.3 卡方分布

以卡方分布作為偏態(tài)分布的代表，對(duì)χ2(2)進(jìn)行100次數(shù)值模擬，兩種方法下的經(jīng)驗(yàn)均方誤差結(jié)果和樣本量為2 000的估計(jì)結(jié)果與真實(shí)值，見圖5和圖6.

圖5 卡方分布導(dǎo)數(shù)估計(jì)的均方誤差

從以上結(jié)果可以清晰地看出偏態(tài)分布下，局部多項(xiàng)式估計(jì)依然優(yōu)于核估計(jì). 圖中的數(shù)值趨勢(shì)與均勻分布類似：無論在哪個(gè)包下，局部多項(xiàng)式回歸方法的均方誤差均小于核方法的均方誤差；隨著樣本量增加，MISElp表現(xiàn)越來越好，估計(jì)結(jié)果越來越接近真實(shí)值，而MISEk卻越來越大. 通過圖6，可以直觀看出核估計(jì)在邊界處的估計(jì)結(jié)果與真實(shí)值的偏差較大，但局部多項(xiàng)式回歸方法在0附近明顯有更好的表現(xiàn)效果. 在內(nèi)點(diǎn)，核估計(jì)出現(xiàn)了欠光滑，但局部多項(xiàng)式估計(jì)幾乎完全擬合真實(shí)密度導(dǎo)數(shù)，估計(jì)性質(zhì)十分優(yōu)秀.

圖6 卡方分布導(dǎo)數(shù)估計(jì)結(jié)果

2.4 小結(jié)

對(duì)比3個(gè)分布下的估計(jì)結(jié)果圖，lpdensity對(duì)卡方分布的密度導(dǎo)數(shù)估計(jì)最優(yōu)，其次是均勻分布，正態(tài)分布. 雖然在正態(tài)分布下，ks的估計(jì)在邊界點(diǎn)和內(nèi)點(diǎn)更勝一籌，但局部多項(xiàng)式估計(jì)方法在這3種分布下一直表現(xiàn)良好，說明該方法不僅有邊界自適應(yīng)性，還具有穩(wěn)健性. 核方法在均勻分布和卡方分布的估計(jì)均方誤差較大，這說明核方法不具有穩(wěn)健性. 當(dāng)分布未知時(shí)，核方法并不是最優(yōu)選擇. 正態(tài)分布下，ks包估計(jì)效果優(yōu)于lpdensity. 但在其他兩個(gè)分布下，ks的估計(jì)效果并不好，而且與事實(shí)相悖：隨著樣本量增加，均方誤差卻存在或大或小的增幅. 這可能與軟件包的內(nèi)置函數(shù)有關(guān)，不同R包的最優(yōu)適應(yīng)范圍不同，為了更好地適應(yīng)某類數(shù)據(jù)需要做出適當(dāng)調(diào)整. 另外，ks包可以用于1至6維數(shù)據(jù)，在多維數(shù)據(jù)中也許會(huì)具有一定的穩(wěn)健性.

通過觀察數(shù)值模擬得到的3個(gè)均方誤差圖以及估計(jì)結(jié)果對(duì)比圖，可以得到以下兩個(gè)結(jié)論.

第一，由于正態(tài)分布的取值范圍是全體實(shí)數(shù)，局部多項(xiàng)式回歸方法在邊界點(diǎn)擬合效果并沒有完美體現(xiàn). 但另外兩個(gè)分布存在明確的邊界，并且相對(duì)于核估計(jì)，局部多項(xiàng)式估計(jì)更接近真實(shí)值. 即在一般情況下，局部多項(xiàng)式估計(jì)具有邊界自適應(yīng)性.

第二，雖然局部多項(xiàng)式回歸方法在正態(tài)分布密度導(dǎo)數(shù)下的估計(jì)性質(zhì)不是最優(yōu)，但這并不能否定該方法在另外兩個(gè)分布下的良好表現(xiàn). 即在大多分布下，相較于核方法，局部多項(xiàng)式估計(jì)方法更穩(wěn)健、高效.

3 討 論

本文基于核方法和局部多項(xiàng)式回歸方法，得到概率密度函數(shù)的導(dǎo)數(shù)估計(jì)，并通過數(shù)值模擬對(duì)比兩種方法下的估計(jì)性質(zhì)，最終得到結(jié)論：樣本量一定時(shí)，相較于核方法，無論分布為均勻、正態(tài)還是偏態(tài)，局部多項(xiàng)式回歸方法的估計(jì)性質(zhì)均更加優(yōu)良. 局部多項(xiàng)式回歸方法無需任何復(fù)雜的數(shù)據(jù)預(yù)處理就可以減輕邊界偏倚并且具有邊界自適應(yīng)性、估計(jì)性質(zhì)穩(wěn)健等優(yōu)點(diǎn). 此外，該方法的理論體系愈見成熟，各階導(dǎo)數(shù)估計(jì)偏差和估計(jì)方差的計(jì)算和證明也都有著系統(tǒng)的闡述. 為了進(jìn)一步對(duì)比核方法和局部多項(xiàng)式回歸方法，在第二部分以均勻分布、正態(tài)分布和卡方分布為例，在不同樣本量下采用兩種方法估計(jì)概率密度函數(shù)的一階導(dǎo)數(shù)，分析對(duì)比了均方誤差和各個(gè)點(diǎn)上的估計(jì)結(jié)果. 從數(shù)值與圖形兩個(gè)方面，驗(yàn)證了局部多項(xiàng)式回歸方法不僅在內(nèi)點(diǎn)的估計(jì)性質(zhì)優(yōu)越，在邊界上的估計(jì)效果也十分良好，整體估計(jì)效果明顯優(yōu)于核估計(jì). 在估計(jì)密度函數(shù)的導(dǎo)數(shù)時(shí)，該方法會(huì)更加穩(wěn)健、高效、簡(jiǎn)單. 此外，局部多項(xiàng)式估計(jì)不僅可以應(yīng)用在操縱檢驗(yàn)中[10]，還可以應(yīng)用到極小化算法中[18]，這會(huì)使檢驗(yàn)結(jié)果與算法應(yīng)用更加精準(zhǔn)有效.

在導(dǎo)數(shù)估計(jì)中，本文僅討論了一元局部多項(xiàng)式回歸模型和核方法的估計(jì)性質(zhì)，未比較兩種方法在多元情形下的效果，這一比較會(huì)更加復(fù)雜. 另外，因?yàn)閹掃x擇對(duì)最終估計(jì)結(jié)果具有顯著影響，所以這將會(huì)是另一個(gè)有意義的研究課題. 在統(tǒng)計(jì)推斷方面，關(guān)于局部多項(xiàng)式回歸的置信帶也是一個(gè)不錯(cuò)的探索方向.