呂書龍,劉文麗,梁飛豹,薛美玉
(福州大學 數學與計算機科學學院,福州 350108)
參數估計和假設檢驗是統計中最常見的兩類統計推斷問題,在點估計方面,矩估計法和極大似然估計法是基本方法,而區間估計和假設檢驗則通常是通過構造樞軸量[1-2]得以實現。隨著非參數理論和方法的深入研究,各種基于非參數理論的統計推斷方法也相繼出現,如秩方法[2-3]、Bootstrap法[4-6]、隨機加權法[7]等。關于分布參數的參數型和非參數型的統計推斷方法的研究很多,成果也很豐碩,文獻[8-16]中列舉了近幾年的一些研究。大量的研究成果極大地促進了經典的統計推斷理論和方法的發展。
經典的參數統計推斷問題的背景可描述為:假設總體X~F(x,θ)(或f(x,θ)),且F(x,θ)(或f(x,θ))的函數形式已知,但θ未知,θ∈Θ,其中Θ為θ的參數空間。已知該總體的一個樣本X1,X2,…,Xn及其觀測值x1,x2,…,xn,對于未知參數θ,經典的參數推斷問題為:① 關于θ的點估計和置信度為1-α的區間估計;② 關于θ的假設檢驗問題。
現有的參數型或非參數型的很多方法都可以解決上述問題,其中有些方法都已成為教學和實際應用中的經典。在講授“應用概率統計”和“統計計算”課程[1,17]的參數估計和假設檢驗內容時,為了激發學生對統計方法的探索和對R統計軟件的應用水平,以密度函數、核估計、Bootstrap法3個關鍵詞,要求構造一套基于直觀統計理論和隨機模擬的,有一定新意的統計推斷方法。這個問題來自是教學過程中的一個突發的想法,本意是希望通過問題式、探究式教學來促進學生的統計計算思維。遺憾的是,這個綜合了統計建模、實驗設計、程序設計和隨機模擬的問題沒有在學生群體中得到突破,反而成了本文研究的一個起點。希望通過此文對類似課程的問題式、探究式、實驗式教學提供一種參考,以便提升教學效果。
核密度是對總體密度函數f(x,θ)的一種估計實現,經驗分布函數是對總體分布函數F(x,θ)的一種估計實現,不妨將這兩個估計統一稱為擬合分布,而F(x,θ)(或f(x,θ))統一稱為理論分布。即
(1)

(2)
式中:h稱為窗寬;K(x)稱為核函數;I(x)為示性函數,當條件x為真時,其值為1,否則為0?,F有的很多理論和實踐都說明了上述估計的合理性和優良性[3-6,18-19]。
既然式(1)和(2)的非參數型擬合分布是理論分布的良好估計,那么在式(1)和(2)的基礎上,不妨逆向思考:擬合分布Fn(x)(或fn,h(x))中形式上已經不含未知參數θ了,不妨轉換其角色,將其當作最終的“理論分布”,而把含未知參數的理論分布F(x,θ)(或f(x,θ))當作“擬合分布”,然后通過合適的手段尋找最佳的θ,使得“擬合分布”逼近“理論分布”,這就誕生了求解未知參數點估計的一個方法,不妨稱之為“非參數逆向思維法”。
此處合適的手段指構建度量“擬合分布”和“理論分布”偏差的損失函數,通過最優化手段確定某個θ的值使得損失函數值達到最小,得到未知參數θ的最優估計,即,

(3)

(4)
上述約束優化模型的求解若需要一個合理的初值,不妨取參數θ為參數空間Θ的中間值。
從形式上看,式(3)適用于連續型分布,若直接套用給離散型分布是不行的。實際上對于離散型分布,取其密度估計為頻率即可,即

(5)
另外,若取函數K(x)=I(x),且窗寬h=1,即可將式(5)統一到式(1)中。
式(3)和(4)給出了通過密度函數和分布函數求解未知參數點估計的基本模型,而要實現參數的區間估計與假設檢驗,還需借助非參數的Bootstrap方法。利用Bootstrap理論和方法,以已知樣本產生足夠多的自助樣本并利用式(3)或(4)得到相應的θ估計序列,再基于Bootstrap方法體系中的求解區間估計的樞軸量和非樞軸量法[5],可計算出未知參數θ的置信度為1-α的區間估計,同理可計算未知參數θ相應的假設檢驗問題的檢驗p值。
考慮式(3)或(4)中的q值,若取q=2,則是基于最小二乘思想;若取q=1,則基于最小一乘思想,代入樣本觀測值,可得到θ的點估計。結合R軟件給出如下過程:
步驟1根據樣本選擇合適的窗寬h,由式(1)結合R軟件中的density函數得到核密度估計fn,h(x)或者由式(2)結合R軟件中的ecdf函數是得到經驗分布函數Fn(x)。
步驟2選定q值,構建式(3)或(4)的最優化目標函數g(θ)。

若希望得到未知參數θ更穩定可靠的估計,可以引入Bootstrap方法,以Bootstrap法估計的平均值作為最終的點估計值。

近似正態法:由
可得:

(6)
近似t分布法:由
可得:

(7)

(8)

近似正態法:

(9)
此處Φ(x)為標準正態分布的分布函數。同理可以實現近似t分布法的檢驗p值,不再贅述。
分位數法:

(10)
參照雙側檢驗的作法,同理可以得到右側檢驗或左側檢驗的檢驗p值。
例1設樣本X1,X2,…,Xn及觀測值x1,x2,…,xn來自柯西分布總體X~C(μ,γ),其概率密度函數為
求參數μ、γ的估計值。
柯西分布的各階矩均不存在,故無法使用矩法得到點估計,但可以使用極大似然估計法。而本文方法屬于非參數法,故不受此限制。下面給出本文方法(簡稱逆向法)和極大似然估計法(簡稱極大法)的隨機模擬比較結果。
從表1結果可知,式(3)和(4)定義的逆向法可用來求解分布參數的點估計。由于核密度估計受到窗寬的影響較大,故式(3)的逆向法與極大似然估計法的偏差較式(4)的逆向法大,而式(4)的逆向法與極大似然法估計的結果非常接近。再給出μ=0,γ=1時,100次模擬的各種基本統計指標,以便比較這3種估計,具體見表2。

表1 逆向方法與極大似然估計法對于參數μ和γ估計的模擬對比(q=2)
注:模擬用的隨機樣本采用函數set.seed(12)進行固定以便分析結果可驗證和可重現,具體程序見附錄1

表2 2個參數(μ,γ)100次隨機模擬估計值的統計指標(q=2)
注:具體程序見附錄2
從表2及圖1、2輸出可知,式(3)逆向法的各項統計指標弱于式(4)逆向法,且式(3)逆向法與極大似然估計有較大差異。但式(4)逆向法與極大似然估計法只有微小差異,由檢驗p值可知,兩者不存在統計意義上的顯著差異。

對表3給出的計算結果分析可知,式(4)逆向法結合Bootstrap方法得到的區間估計與檢驗p值與常規的方法相比沒有本質區別,這歸功于經驗分布函數的穩健性。而式(3)逆向法得到的區間估計與檢驗p值與常規方法相比大部分沒本質區別,但在指數分布上差異明顯,主要原因在于指數分布的密度最高值在邊界處達到,而核密度估計的劣勢正好在邊界處。但在分布密度基本對稱時,式(3)和(4)沒有本質差異。

圖1 參數μ的100次模擬結果的箱線圖

圖2 參數γ的100次模擬結果的箱線圖

表3 分布參數的雙側置信區間和檢驗p值(q=2)
注:各分布均采用函數set.seed(12)指定隨機數表,分別提取100個隨機數;Bootstrap方法的自助樣本依隨機數表1~2 000生成,見附錄3
總體而言,在實際應用中,式(4)逆向法優于式(3)逆向法。一方面式(4)逆向法估計精度有保障而且估計更穩??;另一方面式(4)逆向法的計算效率高于式(3)逆向法。
本文將非參數核密度估計和經驗分布函數這兩個實際上的擬合分布當作“理論分布”,而將密度函數和分布函數這兩個實際上的理論分布當作“擬合分布”,讓“擬合分布”逼近“理論分布”為基本思想逆向地構建了兩者之間逼近的損失函數,通過優化模型得到未知參數的估計,并給出了這兩套方法的R腳本程序。通過構建兩個例子和大量的隨機模擬過程,一方面將這兩套方法與極大似然估計方法進行比較;另一方面給出了解決統計推斷中的區間估計和假設檢驗兩類問題的基本過程,實現了統計方法教學所需要的完整過程,有利于學生充分思考、研究并掌握該方法,隨機模擬實驗設計與對應的R程序也便于實際教學及演示。隨機模擬結果表明基于分布函數的逆向法的普適性、精確性和穩健性優于基于密度函數的逆向法,基于分布函數的逆向法與極大似然估計方法沒有本質差別。在區間估計和假設檢驗方面,基于分布函數的逆向法表現與常規方法無顯著差異;除了在邊界處出現密度極端值外,其他情況下,基于密度函數的逆向法與常規方法也沒有太多差別。但在實際應用中,建議優先使用基于分布函數的逆向法。至于如何提高基于密度函數的逆向法的普適性,則需要在核密度估計的天生缺陷問題上進行改進,有一定難度。至于提高基于密度的逆向法的精度,需要樣本量、最優窗寬和核函數的綜合考量,有待進一步研究。