龐智強(qiáng),王朝旭,2,牛璽娟,2

(1.蘭州財(cái)經(jīng)大學(xué) 統(tǒng)計(jì)學(xué)院,甘肅 蘭州 730020;2.青海師范大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,青海 西寧 810008)

一、引言

在抽樣估計(jì)領(lǐng)域,由于受到小樣本乃至無(wú)樣本的挑戰(zhàn),小域估計(jì)方法得到了學(xué)者們的一致青睞。相比于傳統(tǒng)的直接估計(jì)方法,借助于輔助變量的小域估計(jì)方法能夠更好地解決小樣本和無(wú)樣本問題。實(shí)踐中,小域估計(jì)方法在政府統(tǒng)計(jì)、人口統(tǒng)計(jì)、醫(yī)學(xué)統(tǒng)計(jì)、農(nóng)業(yè)統(tǒng)計(jì)、貧困率估計(jì)等領(lǐng)域有廣泛的應(yīng)用[1-7]。在理論層面,小域估計(jì)也受到高度重視,從Fay提出區(qū)域水平的Fay-Herriot(FH)模型[1]、Battese等提出單元水平的誤差嵌套(NER)模型[2]、Datta和Ghosh提出的小域估計(jì)的貝葉斯預(yù)測(cè)方法[3],到Rao和Molina關(guān)于小域估計(jì)基礎(chǔ)理論的系統(tǒng)綜述[4]、Morales等提出的幾類混合模型小區(qū)域估計(jì)理論[5],小域估計(jì)已經(jīng)形成了較為完整的理論體系。在抽樣調(diào)查的范疇,小域估計(jì)方法也是解決多層次目標(biāo)量估計(jì)的有效途徑之一,金勇進(jìn)和趙雪慧對(duì)此曾做過專門討論[8]。

在小域估計(jì)方法方面,基于模型的小域估計(jì)方法受到學(xué)者們更多的關(guān)注?；谀Ｐ偷姆椒ú粌H能夠估計(jì)無(wú)樣本區(qū)域中感興趣的變量值,還能更好地?cái)M合觀測(cè)數(shù)據(jù)的結(jié)構(gòu)使估計(jì)具有更小的偏差。作為小域估計(jì)的基礎(chǔ)模型之一,單元水平模型能夠處理小域中每個(gè)單元級(jí)別的目標(biāo)變量估計(jì),并且能通過單元數(shù)據(jù)計(jì)算相應(yīng)的區(qū)域水平值。受限于數(shù)據(jù)的收集、輔助變量的獲取以及模型計(jì)算等方面的限制,單元水平模型并不像區(qū)域水平模型那樣受學(xué)者的關(guān)注。若能夠獲取單元級(jí)別的觀測(cè)數(shù)據(jù)和輔助信息,則建立單元水平模型是小域估計(jì)中更好的選擇。如Battese等利用誤差嵌套(NER)模型結(jié)合抽樣數(shù)據(jù)和衛(wèi)星數(shù)據(jù)估計(jì)縣級(jí)水平農(nóng)作物面積[2]。Datta和Ghosh推廣了NER模型,利用分層貝葉斯方法討論了廣義線性混合模型下的估計(jì)量[3]。在小域模型參數(shù)以及目標(biāo)變量的估計(jì)方面,經(jīng)驗(yàn)最佳線性無(wú)偏估計(jì)量(EBLUP)成為應(yīng)用最為廣泛的方法,該方法能夠很好地解決混合線性模型的估計(jì)問題。當(dāng)觀測(cè)變量為二分變量或計(jì)數(shù)變量時(shí),經(jīng)驗(yàn)貝葉斯(EB)方法受到了更廣泛的應(yīng)用[4]。在傳統(tǒng)的單元水平模型中,區(qū)域隨機(jī)效應(yīng)和模型隨機(jī)效應(yīng)部分的誤差均假定服從正態(tài)分布,但在實(shí)際應(yīng)用中,觀測(cè)數(shù)據(jù)中離群值很常見,使用傳統(tǒng)的EBLUP方法會(huì)產(chǎn)生較大的估計(jì)偏差。本文重點(diǎn)關(guān)注隨機(jī)效應(yīng)分布非正態(tài)及模型觀測(cè)存在離群值的單元水平模型,并提出針對(duì)這種模型的穩(wěn)健估計(jì)方法。

離群觀測(cè)是任何調(diào)查中都不可避免的事實(shí)[9],因離群值普遍存在,且會(huì)對(duì)估計(jì)量產(chǎn)生較大影響,因此眾多研究專門探討減小離群觀測(cè)對(duì)調(diào)查估計(jì)影響的方法,常見的有三種。第一種,刪除離群觀測(cè)值,用剩余觀測(cè)值進(jìn)行估計(jì)和預(yù)測(cè)。很顯然,這一方法當(dāng)樣本量非常大時(shí),刪除少數(shù)的離群觀測(cè)可能是可行的。但當(dāng)樣本量本來(lái)就比較小時(shí),這種方法不僅會(huì)帶來(lái)信息的丟失,同時(shí)還會(huì)導(dǎo)致估計(jì)量偏離真實(shí)值。第二種,用非離群觀測(cè)值替代離群觀測(cè)值進(jìn)行估計(jì),并在非抽樣總體中采用穩(wěn)健投影的方式進(jìn)行穩(wěn)健預(yù)測(cè)。然而正如Chambers等所言,若觀測(cè)值是準(zhǔn)確獲取的,不能認(rèn)為每個(gè)觀測(cè)值是唯一的,且沒有理由說(shuō)明在未抽樣總體中不包含離群觀測(cè)值[9]。第三種,通過構(gòu)建對(duì)于離群觀測(cè)值不敏感的穩(wěn)健估計(jì)量,來(lái)減小離群值對(duì)估計(jì)結(jié)果影響,這類方法也是學(xué)者們關(guān)注的主要方法,本文也采用了這類穩(wěn)健估計(jì)方法。

有關(guān)小域穩(wěn)健估計(jì)的研究最早見于Ghosh和Lahiri的成果,他們利用貝葉斯方法討論了分層均值的穩(wěn)健估計(jì)[10]。后來(lái),分層貝葉斯方法[3]、長(zhǎng)尾分布假設(shè)、分位數(shù)回歸等方法均用于研究穩(wěn)健小域估計(jì)問題,但只是處理特殊情形下的穩(wěn)健估計(jì)問題,且缺乏對(duì)估計(jì)參數(shù)性質(zhì)的研究。Ghosh等提出了針對(duì)區(qū)域水平模型的穩(wěn)健貝葉斯預(yù)測(cè)量,并用于克服由離群觀測(cè)造成的過度收縮問題[11],但并沒有考慮離群值對(duì)模型系數(shù)的影響。Sinha和Rao利用Huber-Ψ函數(shù),通過對(duì)由離群值引起的項(xiàng)減小權(quán)重的方式,建立了單元水平模型的小域穩(wěn)健估計(jì)量[12]。目前,該方法普遍應(yīng)用于小域穩(wěn)健估計(jì)。在穩(wěn)健估計(jì)的其他研究方面,Smith等用穩(wěn)健投影和分位數(shù)的方法研究了商業(yè)調(diào)查中的小域穩(wěn)健估計(jì)[13],Bertarelli等通過偏差修正的方式研究了小域穩(wěn)健估計(jì)問題[14],Jiang和Rao對(duì)穩(wěn)健小域估計(jì)的研究歷史進(jìn)行系統(tǒng)總結(jié)后指出,在以后的研究中,小域估計(jì)方法可應(yīng)用于大數(shù)據(jù)、關(guān)聯(lián)數(shù)據(jù)等領(lǐng)域,同時(shí)其他的現(xiàn)代統(tǒng)計(jì)方法、學(xué)科理論也可以用于小域估計(jì)[15]。由此可見,小域估計(jì)的穩(wěn)健估計(jì)研究仍需要進(jìn)一步的完善,探索具有普適性的估計(jì)方法成為該類研究中不可或缺的內(nèi)容。本文重點(diǎn)關(guān)注更普遍的假設(shè)條件下穩(wěn)健小域估計(jì),并給出估計(jì)方法和步驟。

在統(tǒng)計(jì)推斷中,當(dāng)觀測(cè)數(shù)據(jù)存在嚴(yán)重的離群值時(shí),基于密度的最小距離方法成為解決問題的一種有效方法。Basu等提出用最小化密度冪散度(MDPD)的方法來(lái)得到參數(shù)的穩(wěn)健估計(jì),其核心思想就是給離群觀測(cè)項(xiàng)賦予更小的權(quán)重,從而獲得更加穩(wěn)健的估計(jì)量[16]。在密度冪散度(DPD)類方法中,通過一個(gè)調(diào)整參數(shù)α(0<α<1)來(lái)調(diào)整參數(shù)的穩(wěn)健性和估計(jì)有效性之間的平衡,隨著α增大,模型估計(jì)的效率降低而穩(wěn)健性增加。一般該類方法需要通過選取調(diào)整參數(shù),損失一點(diǎn)有效性而達(dá)到穩(wěn)健性的目的。Jones等提出了基于γ散度的穩(wěn)健估計(jì)方法[17]。該方法與DPD方法同屬基于密度的最小距離方法,應(yīng)用于穩(wěn)健估計(jì)和回歸等方面[18]。考慮到該類方法在穩(wěn)健估計(jì)中的優(yōu)異性質(zhì),本文將γ散度引入小域估計(jì)中用于解決小域估計(jì)中的穩(wěn)健性難題。

在本文中,將γ散度應(yīng)用于單元水平的小域估計(jì)模型,考慮模型誤差有偏和存在離群值情形下,討論模型的隨機(jī)誤差和區(qū)域誤差被污染時(shí)的穩(wěn)健估計(jì)問題。首先,提出基于γ散度的單元水平模型的參數(shù)估計(jì)方法,給出參數(shù)估計(jì)方程和估計(jì)算法,討論了參數(shù)的漸近性質(zhì)。其次,給出了有限總體區(qū)域均值的估計(jì)量,并根據(jù)估計(jì)參數(shù)的漸近性質(zhì),給出了穩(wěn)健調(diào)節(jié)參數(shù)γ的選取算法。進(jìn)一步利用Bootstrap方法給出了穩(wěn)健估計(jì)量MSE的估計(jì)。最后,通過模擬數(shù)據(jù)和實(shí)際數(shù)據(jù)對(duì)本文提出的估計(jì)方法進(jìn)行了模擬驗(yàn)證,并與Sinha和Rao的穩(wěn)健估計(jì)方法進(jìn)行了對(duì)比[12],驗(yàn)證了本文穩(wěn)健估計(jì)方法的有效性。

二、模型及參數(shù)估計(jì)

(一)模型

當(dāng)抽樣得到的數(shù)據(jù)和輔助信息都是單元級(jí)別的數(shù)據(jù)時(shí),單元水平模型成為了小域估計(jì)的有力工具。誤差嵌套模型(Nested Error Regression model,NER)是最基本的單元水平模型,由Battese等提出[2]。其具體模型如下:

(1)

不失一般性,本文假設(shè)抽樣得到的數(shù)據(jù)滿足式(1)所示的模型,則由樣本數(shù)據(jù)構(gòu)成的模型可以寫成如下的矩陣形式:

(2)

(3)

利用Sherman-Morrison公式可求出Vi的逆矩陣,

(4)

(二)γ散度及參數(shù)估計(jì)

假設(shè)某一總體未知的真實(shí)密度函數(shù)為g(x),用參數(shù)模型fθ(x)=f(x;θ)估計(jì)g(x),其中θ∈Θ,Θ為θ的所有可能的取值構(gòu)成的p維參數(shù)空間。γ散度用于衡量?jī)蓚€(gè)概率密度之間的差異程度,γ散度的定義如下:

(5)

其中,γ>0是調(diào)整參數(shù),θ為待估的未知參數(shù)。根據(jù)Kawashima和Fujisawa關(guān)于γ散度的論述[18],當(dāng)γ→0時(shí),式(5)中的γ散度退化為Kullback-Leibler(KL)散度,且最優(yōu)的參數(shù)通過最小化式(5)得到。式(5)中的第一項(xiàng)與未知參數(shù)無(wú)關(guān),因此最小化γ散度的問題變?yōu)?

(6)

式(6)中,當(dāng)γ→0時(shí),最小化γ散度的方法就退化為極大似然估計(jì)(MLE)。一般用抽樣數(shù)據(jù)的條件經(jīng)驗(yàn)分布函數(shù)來(lái)替代式(6)中的密度函數(shù)g(x),則上述目標(biāo)函數(shù)可以寫成:

(7)

稱式(7)為γ似然函數(shù)。將式(3)的條件密度函數(shù)帶入式(7),便可得:

(8)

對(duì)式(8)的γ似然函數(shù)關(guān)于θ的每個(gè)分量求偏導(dǎo)數(shù),便可得參數(shù)估計(jì)方程:

(三)估計(jì)參數(shù)的性質(zhì)

為了便于標(biāo)記,定義:

(1)存在參數(shù)空間Θ中的開集Λ包含最佳擬合參數(shù)θ,使得θ∈Λ時(shí),J是正定矩陣。

若θg是θ的真實(shí)值,對(duì)上式在θ=θg處進(jìn)行泰勒展開,化簡(jiǎn)可得:

其中,Rm是泰勒展式中的余項(xiàng)。利用弱大數(shù)定律可知:

該定理說(shuō)明了估計(jì)參數(shù)的一致性和漸近分布,在漸近分布的計(jì)算中,J(i)、K(i)中均含有未知的真實(shí)密度函數(shù)g(yi|xi),因此直接進(jìn)行計(jì)算是難以做到的。在基本假設(shè)條件中,假定密度族fθ(yi|xi)中包含真實(shí)的密度函數(shù),因此在J(i)、K(i)、ξ(i)的定義中,用fθ(yi|xi)替換g(yi|xi)進(jìn)行估計(jì)。

從J(i)、K(i)、ξ(i)的定義中不難發(fā)現(xiàn),當(dāng)把矩陣J(i)認(rèn)為是γ的函數(shù)時(shí),矩陣K(i)可以表示為:K(i)=J(i)(2γ)-ξ(i)ξ(i)′,因此,矩陣K可以表示為:

通過上述協(xié)方差矩陣的計(jì)算可知,γ增大時(shí)參數(shù)方差增大,這說(shuō)明MDPD估計(jì)量的效率隨γ增大而減小。這一點(diǎn)進(jìn)一步驗(yàn)證了優(yōu)化參數(shù)γ用于控制MDPD估計(jì)量的效率和魯棒性之間的權(quán)衡,當(dāng)γ增加時(shí)穩(wěn)健性增強(qiáng)效率降低。然而,后面的模擬實(shí)驗(yàn)表明這種效率的損失并不大。

(四)調(diào)整參數(shù)的選取

從上述結(jié)果中不難發(fā)現(xiàn),利用DPD方法迭代得到的參數(shù)θ的估計(jì)表達(dá)式中,含有未知的調(diào)整參數(shù)γ。而γ的選取決定著穩(wěn)健性和有效性的平衡,當(dāng)γ的取值越接近于1,則估計(jì)得到的參數(shù)具有更強(qiáng)的穩(wěn)健性,反之,穩(wěn)健性越弱,有效性越強(qiáng),因此選取適當(dāng)?shù)恼{(diào)整參數(shù)γ是穩(wěn)健估計(jì)中很關(guān)鍵的因素。目前關(guān)于調(diào)整參數(shù)選取的方法主要有兩種,一種是根據(jù)有效性和穩(wěn)健性的占比關(guān)系,研究者根據(jù)自己的需求確定兩者之間的占比,在此基礎(chǔ)上,選取最優(yōu)的調(diào)整參數(shù)。另一種方法是基于數(shù)據(jù)驅(qū)動(dòng)的參數(shù)選擇方法,Warwick和Jones通過最小化估計(jì)參數(shù)的MSE得到最優(yōu)的調(diào)整參數(shù)[19]。但是,該方法依賴于初始值的選取,不同的參數(shù)初始值可能會(huì)產(chǎn)生不同的調(diào)整參數(shù)。Basak等在Warwick和Jones的研究基礎(chǔ)上,提出了不依賴于初始值的調(diào)整參數(shù)選擇方法[20]。本文采用Basak等提到的參數(shù)選取方法,選擇最優(yōu)的調(diào)整參數(shù),用于構(gòu)建穩(wěn)健小區(qū)域估計(jì)量。

對(duì)于未知參數(shù)θ的真實(shí)值θ*最優(yōu)的調(diào)整參數(shù)γ通過最小化MDPD估計(jì)量的MSE來(lái)獲得,即:

(9)

算法 IWJ算法

重復(fù):1.利用WJ算法最小化正文中的式(9),并在γ的取值區(qū)間Iγ內(nèi)更新γ:

2.固定γ(i+1),將其帶入最小化γ散度的估計(jì)方程并更新得到θ*(i+1);

3.重復(fù)1、2,直至|γ(i+1)-γ(i)|<ε或|θ*(i+1)-θ*(i)|<ε*,其中ε、ε*為估計(jì)精度。

輸出:γ(i+1)。

三、穩(wěn)健小域估計(jì)量

(一)穩(wěn)健經(jīng)驗(yàn)貝葉斯預(yù)測(cè)量

結(jié)合上述單元水平模型的假設(shè),在抽樣數(shù)據(jù)yis的條件下,yir的條件分布為:

yir|yis～N(μir|s,Vir|s)

(10)

其條件均值和條件協(xié)方差矩陣分別為:

(11)

(二)區(qū)域均值的EBP

(三)區(qū)域EBP的MSE

估計(jì)單元水平模型的均方預(yù)測(cè)誤差(MSPE)通常是比較困難的,原因有二:其一,單元水平模型中誤差項(xiàng)和未抽中樣本單元所服從的真實(shí)分布并不知道,從而無(wú)法獲取密度函數(shù)進(jìn)行MSPE的計(jì)算;其二,有時(shí)即便知道未抽樣單元的分布,由于單元水平模型涉及的是單元層級(jí)的數(shù)據(jù),在計(jì)算其期望的過程中會(huì)出現(xiàn)多重積分,這對(duì)MSPE的計(jì)算造成了極大的挑戰(zhàn)。

在本文中,利用Sinha和Rao中提到的參數(shù)bootstrap方法來(lái)估計(jì)MSPE[12]。具體估計(jì)步驟如下:

四、數(shù)值模擬與實(shí)例分析

(一)數(shù)值分析

在模擬實(shí)驗(yàn)中,重點(diǎn)比較傳統(tǒng)的極大似然方法(ML)、Sinha和Rao提到的穩(wěn)健估計(jì)方法(RML)[12]以及本文提出的取不同調(diào)整參數(shù)(γ=0.1,0.2,0.3)時(shí)的穩(wěn)健最小化γ散度的方法(RMG)。在RMG方法中,通過兩種方式選取了調(diào)整參數(shù)γ,第一類為固定選取調(diào)整參數(shù)γ為0.1、0.2、0.3;第二類為使用IWJ算法選取的調(diào)整參數(shù)。將選取的調(diào)整參數(shù)帶入?yún)?shù)估計(jì)方程,則可得模型參數(shù)的數(shù)值估計(jì),其估計(jì)的平均結(jié)果見表1。首先比較四種污染情形下的模型參數(shù)的估計(jì)。表1給出了幾種穩(wěn)健估計(jì)方法在不同污染情形下的參數(shù)估計(jì)的偏差和均方誤差,其中每種方法對(duì)應(yīng)的第一列表示估計(jì)的偏差,第二列表示相應(yīng)的MSE。

表1 四種污染情形下模型參數(shù)估計(jì)的偏差和MSE

從表1可以得出結(jié)論:當(dāng)隨機(jī)誤差未被污染時(shí),ML方法在參數(shù)估計(jì)中表現(xiàn)最好,但是RML以及具有較小調(diào)整參數(shù)的RMG方法和ML估計(jì)結(jié)果非常相近,說(shuō)明在這種情形下,具有較小調(diào)整參數(shù)的RMG方法、RML與ML方法是幾乎一樣有效的,其偏差和MSE相差都不多。具有較大調(diào)整參數(shù)的RMG方法表現(xiàn)不好,恰好表明調(diào)整參數(shù)的選取至關(guān)重要,而最優(yōu)調(diào)整參數(shù)可根據(jù)我們提供的算法得到。

在區(qū)域效應(yīng)和模型誤差均被污染時(shí),ML方法對(duì)方差部分的估計(jì)受離群值觀測(cè)影響嚴(yán)重,以致于產(chǎn)生很大的偏差和MSE,RML方法減小了離群值帶來(lái)的影響,但表現(xiàn)并非最好,RMG方法顯著好于RML方法,只需要選取合適的調(diào)整參數(shù)即可。

表2 四種污染情形下均值估計(jì)的偏差及其MSE

圖1 模型系數(shù)的MSE隨污染比例變化圖

圖2 估計(jì)方差的MSE隨污染比例變化圖

圖3 模型系數(shù)的MSE隨污染方差變化圖

圖4 估計(jì)方差的MSE隨污染方差變化圖

接下來(lái)考慮有限總體中小域均值的估計(jì)。比較在有限總體中用上述穩(wěn)健估計(jì)方法得到的區(qū)域平均值的估計(jì)結(jié)果。設(shè)區(qū)域個(gè)數(shù)m=40,考慮區(qū)域總體數(shù)量Ni分別為40、80、200的情形下總體均值估計(jì)的表現(xiàn)。從第i個(gè)區(qū)域的Ni個(gè)單元中選取ni=4個(gè)單元作為隨機(jī)觀測(cè)樣本。對(duì)每次模擬的數(shù)據(jù)集,可以利用上述模擬中提到的穩(wěn)健方法獲取每個(gè)區(qū)域上觀測(cè)變量的均值,最后比較500次模擬之后區(qū)域均值的平均估計(jì)效果。

從表3的模擬結(jié)果上看,總的來(lái)說(shuō)RMG方法在大部分情形下具有更小的MSE,在某些情形下雖不及RML方法的估計(jì)效果,但相差不大。另外,當(dāng)區(qū)域尺度變大時(shí)如Ni=200,RMG方法得到的估計(jì)量具有更小的MSE,估計(jì)效果好于傳統(tǒng)的估計(jì)方法。

表3 不同污染情形下有限總體均值估計(jì)的偏差及其MSE

(二)實(shí)例分析

下面采用Battese等給出的數(shù)據(jù)(1)該數(shù)據(jù)由Battese等(1988)給出,主要結(jié)合了農(nóng)作物種植面積的抽樣數(shù)據(jù)和衛(wèi)星觀測(cè)數(shù)據(jù),用于估計(jì)愛德華州12個(gè)縣的玉米和大豆的種植面積。來(lái)驗(yàn)證本文提出的穩(wěn)健估計(jì)方法[2]。該數(shù)據(jù)可以從R包sae中獲得,其中包含了來(lái)自于12個(gè)縣區(qū)的玉米和大豆的面積樣本數(shù)據(jù)37個(gè),以及每個(gè)區(qū)域上的玉米和大豆的像素值。用農(nóng)場(chǎng)采訪收集的數(shù)據(jù)作為因變量,衛(wèi)星數(shù)據(jù)作為輔助變量,建立單元水平模型。

該模型可以作為模型(1)中kij=1,xij=(1,xij1,xij2)′的特殊情形。其中,yij表示第i個(gè)縣第j個(gè)區(qū)域的玉米(大豆)的畝數(shù),xij1和xij2分別表示第i個(gè)縣第j個(gè)區(qū)域的玉米、大豆的像素值。

Battese等將Hardin縣的值識(shí)別為離群值,并且在預(yù)測(cè)玉米和大豆的面積時(shí)只是將這一觀測(cè)值簡(jiǎn)單地刪除。Sinha和Rao則利用穩(wěn)健估計(jì)方法對(duì)這一數(shù)據(jù)進(jìn)行了分析,并給出了存在離群值時(shí)相應(yīng)的預(yù)測(cè)值[12]。在這里,運(yùn)用本文提出的方法對(duì)該數(shù)據(jù)進(jìn)行建模分析,對(duì)每個(gè)區(qū)域上的玉米種植面積進(jìn)行估計(jì)和預(yù)測(cè)。在用提出的方法進(jìn)行預(yù)測(cè)估計(jì)之前,首先選取適當(dāng)?shù)膮?shù)γ,根據(jù)IWJ算法得出適用于該數(shù)據(jù)的最優(yōu)γ=0.01。由于本數(shù)據(jù)中僅存在一個(gè)離群觀測(cè),因此在提出的穩(wěn)健估計(jì)方法中,選取γ=0.01、0.05兩種情形進(jìn)行估計(jì)。在表4中給出了利用ML方法、RML方法以及本文的RMG方法估計(jì)得到回歸系數(shù)和隨機(jī)誤差的方差。同時(shí)表中括號(hào)內(nèi)給出了根據(jù)漸近分布得到的每個(gè)參數(shù)的標(biāo)準(zhǔn)誤。結(jié)合表中數(shù)據(jù)比較而言,在調(diào)整參數(shù)γ=0.01時(shí),用RMG方法估計(jì)得到的系數(shù)介于ML方法和RML方法得到的估計(jì)值之間。當(dāng)調(diào)整參數(shù)γ增加到0.05時(shí),模型估計(jì)的系數(shù)有了顯著的變化。通過比較表4中展示的參數(shù)的標(biāo)準(zhǔn)誤,可見用本文的方法估計(jì)得到的參數(shù)具有更小的標(biāo)準(zhǔn)誤。

表4 模型參數(shù)的估計(jì)及其標(biāo)準(zhǔn)誤

為了體現(xiàn)估計(jì)的效果,比較了幾種估計(jì)方法對(duì)每個(gè)區(qū)域上玉米種植面積平均畝數(shù)的預(yù)測(cè)。由于每個(gè)區(qū)域上玉米種植面積的真實(shí)值并不知道,因此重點(diǎn)分析估計(jì)值的MSPE。采用上述估計(jì)方法得到的EBLUP值在表5中,其中括號(hào)內(nèi)展示了用500個(gè)bootstrap樣本得到的估計(jì)量的MSPE的估計(jì)值。

表5 區(qū)域玉米種植面積的預(yù)測(cè)及其Bootstrap MSPE

首先從估計(jì)結(jié)果上來(lái)看,RMG對(duì)沒有離群值區(qū)域的估計(jì)更接近ML的結(jié)果,并且對(duì)區(qū)域Hardin的預(yù)測(cè)有了一定的改進(jìn)。其次,通過比較bootstrap MSPE,提出方法RMG1和RMG2比其他兩種估計(jì)方法得到的MSPE值更小,且最優(yōu)調(diào)整參數(shù)γ=0.01時(shí)的RMG1表現(xiàn)更好?？梢?RMG方法是有效的。

五、結(jié)論與研究展望

本文提出了一種針對(duì)存在離群觀測(cè)值的單元水平模型小域穩(wěn)健估計(jì)方法。通過引入MDPD方法,給出了解決具有離群觀測(cè)和非正態(tài)分布誤差的穩(wěn)健估計(jì)方法。首先,給出了單元水平模型參數(shù)的估計(jì)方程和漸近性質(zhì)。其次,結(jié)合參數(shù)的漸近分布,給出了最優(yōu)調(diào)整參數(shù)的選擇程序。再次,給出了有限總體中單元和區(qū)域均值的EBLUP值和估計(jì)量的MSE。最后,通過模擬數(shù)據(jù)和實(shí)際數(shù)據(jù)驗(yàn)證了本文提出方法的優(yōu)越表現(xiàn)。在模擬部分,模擬了分布被污染時(shí)的穩(wěn)健估計(jì),討論了三種污染情形下幾類穩(wěn)健估計(jì)方法的效果,還特別討論了污染比例變化以及污染分布的方差變化時(shí)幾類估計(jì)方法的MSE的變化情況。同時(shí)模擬結(jié)果表明,本文提出的方法能更好地解決這種離群情況。實(shí)際數(shù)據(jù)中,用一個(gè)小區(qū)域估計(jì)中很經(jīng)典的數(shù)據(jù)說(shuō)明本文提出的方法不僅十分有效,還能夠很好地處理離群觀測(cè)這一特殊情況。

進(jìn)一步的驗(yàn)證表明,本文的方法針對(duì)隨機(jī)效應(yīng)服從其他有偏分布時(shí)也是有效的。當(dāng)隨機(jī)誤差的分布被污染,且污染概率大于0.3時(shí),本文的方法表現(xiàn)一般,比Sinha和Rao研究中的穩(wěn)健估計(jì)方法差[12],但是在這種情形下,幾類方法得到的MSE均很大,穩(wěn)健估計(jì)結(jié)果都不太有價(jià)值。下一步,應(yīng)嘗試將本文方法進(jìn)一步推廣應(yīng)用于指數(shù)分布情形的小區(qū)域估計(jì)問題。在本文的模擬實(shí)驗(yàn)部分,僅展示了模型誤差項(xiàng)分布來(lái)自于混合正態(tài)分布時(shí)的估計(jì)結(jié)果。實(shí)際上,當(dāng)模型誤差項(xiàng)分布屬于其他情形的有偏估計(jì)時(shí),本文提出的估計(jì)方法都是適用的。例如當(dāng)模型誤差來(lái)自于t分布、對(duì)數(shù)正態(tài)分布、伽馬分布等情形時(shí),本文提出的穩(wěn)健估計(jì)方法具有更小的偏差和均方誤差,是一致有效的。這一點(diǎn)進(jìn)一步說(shuō)明該方法的普適性,能夠在更一般的模型假設(shè)條件或存在離群觀測(cè)時(shí)使用,并得到較為理想的估計(jì)結(jié)果。

在本文的研究基礎(chǔ)上,還可以做如下工作:有關(guān)小域穩(wěn)健估計(jì)的應(yīng)用方面,Tang等提出了一種全局—局部收縮限樣的方法來(lái)刻畫小域估計(jì)中的隨機(jī)效應(yīng),在貝葉斯分析的框架下給出了參數(shù)和目標(biāo)變量的后驗(yàn)估計(jì)[21]。與本文提出的穩(wěn)健估計(jì)不同,該方法主要適用于區(qū)域數(shù)量較大的區(qū)域水平模型。Kurisu等利用了和本文同樣的估計(jì)方法,給出了區(qū)域水平模型的置信區(qū)間的估計(jì)[22]。類似地,可以在本文的基礎(chǔ)上給出單元水平模型的置信區(qū)間估計(jì)。在本文的穩(wěn)健估計(jì)方法中,調(diào)整參數(shù)的選取采用了IWJ算法,而Sugasawa和Yonekura給出了最小化密度冪類方法的調(diào)整參數(shù)選擇的另外一種方法[23]。因此,在本文的研究框架下,同樣可以采用Sugasawa和Yonekura提出的方法來(lái)選擇調(diào)整參數(shù),從而達(dá)到穩(wěn)健小域估計(jì)的目的。