武雅萱,劉曉宇
(1.中國人民大學統(tǒng)計學院,北京 100872;2.首都經(jīng)濟貿(mào)易大學統(tǒng)計學院,北京 100070)
抽樣調(diào)查中的小域估計問題指的是,在對有限總體的子總體進行估計時,所存在的樣本量無法滿足估計精度要求的問題[1]。由于組織時間、調(diào)查費用等因素的限制,擴大樣本量或補充樣本的解決方式收效較低,欲從根本上解決小域估計問題需要依賴對估計方法的改進[2]。傳統(tǒng)基于設計的推斷方式,例如以樣本目標變量觀測值與抽樣權(quán)數(shù)的乘積為基礎構(gòu)造的HT估計,需要較大樣本量來保證估計量的漸近性質(zhì),但小域估計中的域樣本量通常不足以獲得可靠的基于設計的估計量,基于設計的估計量方差可能會大到難以接受的水平,此時需要采用基于模型的方法進行小域估計。
基于模型的小域估計方法通過模型刻畫目標變量和輔助變量之間的關系,并將不同小域聯(lián)系起來,借助其他相似小域的樣本信息對本小域目標變量進行估計,其推斷效果僅依賴于模型假設,適用于數(shù)據(jù)缺失嚴重、樣本量極小的小域估計問題[3]。一般根據(jù)輔助信息的層次將模型分為域?qū)哟文P蚚4]和單元層次模型[5]。有限總體和小域之間存在層次結(jié)構(gòu),多層次模型既能對域間異質(zhì)性和域內(nèi)相關性予以考察,又能利用域?qū)哟蔚妮o助信息,還能利用單元層次的輔助信息,在小域估計中有著重要應用[6—8]。
測量誤差是測量值和真實值之間的差異,導致測量誤差產(chǎn)生的原因較多,如問題設置不當、受社會期許影響、受訪者受教育程度有限等。最常用的校正思路是,假定觀測值和真實值之間滿足某種函數(shù)關系,通過預測變量真實值來校正測量誤差[9]。目前,鮮有研究直接對目標變量的測量誤差進行校正,對基于模型小域估計方法的測量誤差的考察集中于輔助變量,當納入模型的輔助變量存在測量誤差時,會導致參數(shù)估計出現(xiàn)偏差,推斷效果有可能劣于基于設計的估計[10]。貝葉斯方法常被用于測量誤差校正[10—13]。
我國對貝葉斯小域估計模型的理論研究較少,主要集中于應用上[14]。本文在前人研究的基礎上進行拓展,采用多層次模型進行小域估計。第一層模型用于刻畫域間異質(zhì)性,對各小域的變量分布進行建模;第二層模型用于刻畫域間相關性,借助其他域的樣本單元實現(xiàn)對指定小域的估計,并對第二層模型的參數(shù)建立測量誤差模型。同時,考察抽樣誤差和非抽樣誤差的影響,直接對目標變量建立抽樣誤差模型和測量誤差模型。針對所提出的模型,本文提出具體的參數(shù)估計與誤差估計方法,通過模擬實驗驗證所提方法的具體效果,并將其應用于實際數(shù)據(jù)集。
記有限總體為U,總體規(guī)模為N,總體可劃分為m個小域(U1,U2,…,Um),Ui的規(guī)模為Ni。經(jīng)一定抽樣設計,從總體中抽取樣本,記樣本為s,樣本量為n,第i個小域中抽得的樣本量為ni,。以第i個小域中第j個單元的比率估計Rij為目標變量,Rij=Yij Mij,i=1,…,m;j=1,…,ni。其中,變量Mij是變量Yij相應的總量度量,Rij∈( 0,1) ,以各小域均值為最終估計目標,,記輔助變量為Xij。對目標變量建立抽樣誤差模型和測量誤差模型。
抽樣誤差模型為:
其中,uij為第i個小域中第j個單元的抽樣誤差;為抽樣誤差方差的估計值,可由抽樣設計得到。
測量誤差模型為:
以小域為單位建立第一層模型,考慮抽樣誤差與測量誤差,結(jié)合式(1)、式(2)可得第一層模型如下:
其中,p(·)為連接函數(shù),用于將數(shù)值轉(zhuǎn)化至0~1范圍內(nèi),以滿足Rij∈( 0,1) 的要求。在第一層模型中,各個小域的回歸系數(shù)βi及方差膨脹因子ψi不同,由此刻畫了域間異質(zhì)性。
考慮不同小域之間的相似性,對域?qū)哟蔚膮?shù)建立第二層模型:
其中待估計參數(shù)θ=(β,Σ,ψ,Φ )。由于無法獲取βi和ψi的真實值,只能從第一層模型中獲得βi和ψi的估計值和,故進一步考慮其測量誤差,假定:
其中,和為相應的方差估計,本文采用三明治公式,根據(jù)第一層模型估計結(jié)果和單元層次數(shù)據(jù)計算得到。三明治估計量的優(yōu)點在于,即使在擬合參數(shù)模型不成立或甚至沒有指定的情況下,它也能為參數(shù)估計提供一致的協(xié)方差矩陣估計值[15]。根據(jù)貝葉斯公式,由式(7)、式(8)可得βi和ψi后驗分布分別為:
對于第一層模型,先采用最小二乘法,計算各個小域的回歸系數(shù)βi及方差膨脹因子ψi的估計值和:
將估計結(jié)果作為初始值和,采用偽極大似然估計法進行迭代更新。給定第t次迭代值和,計算第t+1次迭代值和:
對于第二層模型,先采用三明治公式計算測量誤差模型中的方差和,化簡得:
再利用EM 算法求解待估參數(shù)θ=(β,Σ,ψ,Φ ),具體過程如下:
首先,求解(β,Σ)。取初始值重復如下E 步和M 步迭代至收斂,結(jié)果記作
(1)E步:給定第t次迭代值和,計算聯(lián)合分布的條件概率期望:
其中,f(βi;β),Σ為βi的概率密度函數(shù)。
(2)M步:最大化聯(lián)合分布的條件概率期望,記第t+1次迭代值為和,有:
解得:
其次,類似于求解(β,Σ)的過程,再次利用EM 算法求解待估計參數(shù)(ψ,Φ )。取初始值。重復如下E步和M步迭代至收斂,結(jié)果記作
(1)E步:給定第t次迭代結(jié)果和,計算聯(lián)合分布的條件概率期望:
其中,f(ψi;ψ,Φ )為ψi的概率密度函數(shù)。
(2)M步:最大化聯(lián)合分布的條件概率期望,記第t+1次迭代值為和,有:
解得:
最后,利用βi和ψi的貝葉斯后驗分布,獲得估計值分別為:
根據(jù)Prasad和Rao(1990)[16]對最優(yōu)無偏線性預測的討論,可得第i個小域中第j個單元目標變量Rij的最優(yōu)估計為:
為估計的誤差,定義,有:
為了檢驗多層次模型在小域估計中的應用效果,本文模擬真實應用場景,按照如下步驟生成樣本數(shù)據(jù)。
考慮兩種不同大小的域數(shù)量m=10、m=20,假設子總體規(guī)模Ni=2000,i=1,…,m。首先,對于第i個小域中第j個單元(i=1,…,m;j=1,…,Ni)生成變量Mij~50+Bin(100,0.5) ,輔助變量Xij1~Mij·Beta(2,5+0.1×i),Xij2~Mij·Beta(5,5+0.1×i)。然后,記Xij=(1,Xij1,Xij2)T,取β=(0,0.4,0.6)T,Σ=diag(σ2,σ2,σ2),σ2=0.01,生成服從多元正態(tài)分布的域回歸系數(shù)βi~N(β,Σ),由此可計算得到比率。假設方差膨脹系數(shù)服從正態(tài)分布ψi~N(0.5,0.16 ),生成測量誤差eij|pij~N(0,ψi pij(1-pij)),進而得到帶測量誤差的比率為Rij=pij+eij,假設抽樣誤差的方差為,則可以進一步生成各單元的抽樣誤差,進而得到觀測值。最后,根據(jù)輔助變量Xij2的中位數(shù)和上下四分位數(shù)將總體劃分為四層,采用與規(guī)模成比例的分層抽樣抽得樣本量n=100 的樣本。將模擬生成的觀測值,輔助變量Mij、Xij1和Xij2,以及抽樣誤差的方差作為已知信息,估計目標變量Rij。
重復實驗1000次,采用HT估計和本文提出的多層次模型進行估計,結(jié)果見表1。對于域個數(shù)較多的情況(m=20),表1中僅列示各小域估計結(jié)果的均值。其中,相對有效性表示模型估計誤差平方與HT估計誤差平方的比值,比值在0 到1 之間說明模型估計比HT 估計有效性更高,比值越小,表明模型校正效果越好。
從表1可知,整體而言,HT估計的相對偏差大于模型估計的相對偏差,模型估計僅有一個小域的相對偏差略大于HT估計,且HT估計的MSE均不小于模型估計的MSE,說明模型估計的結(jié)果更穩(wěn)定。相對有效性均小于1,進一步說明了模型估計的優(yōu)勢。
本文以恩格爾系數(shù)的估計為例,基于2018 年中國綜合社會調(diào)查(CGSS)數(shù)據(jù)進行實證分析。恩格爾系數(shù)可由年家庭食品消費支出除以年家庭總支出計算得到?;凇笆》荨薄熬游瘯虼逦瘯薄俺錾攴荨薄笆芙逃潭取薄?017年家庭食品支出”“2017年家庭總支出”“樣本權(quán)重”構(gòu)造完整數(shù)據(jù)集,以不同省份、直轄市、自治區(qū)為小域,得到有效樣本3579個。
將樣本按地理分布劃分為華北、東北、華東、中南、西南、西北5 個地區(qū)。各省份發(fā)展情況存在差異,但同一地區(qū)的省份較為接近,因此,假設多層次模型中域間回歸系數(shù)存在差異,同一地區(qū)內(nèi)的域回歸系數(shù)獨立同分布。以城鄉(xiāng)、年齡、受教育程度為輔助變量,采用本文所提多層次模型進行估計,并與傳統(tǒng)基于設計的HT估計進行比較,各域有效樣本量及估計結(jié)果見表2。
表2 各域有效樣本量及恩格爾系數(shù)估計結(jié)果
國家統(tǒng)計局公布的2017 年全國居民恩格爾系數(shù)為0.293,城鎮(zhèn)居民恩格爾系數(shù)為0.286,農(nóng)村居民恩格爾系數(shù)為0.312。從表3中可以發(fā)現(xiàn),HT估計整體偏高,特別是北京的恩格爾系數(shù)估計受異常值影響,估計結(jié)果為2.710,誤差遠超可承受范圍。模型估計可以較好地解決抽樣誤差和測量誤差較大的問題,對于大部分域而言,模型估計結(jié)果更接近國家統(tǒng)計局公布的數(shù)據(jù)。
模型估計不僅可以給出域級別估計,還可以充分利用小域之間的關聯(lián),自下而上得到地區(qū)、國家層面的估計。對于華北地區(qū),受到北京數(shù)據(jù)異常值的影響,HT估計結(jié)果為1.514,遠超全國平均水平0.293,誤差依然處于不可忽略的程度;而多層次模型由于利用了地區(qū)內(nèi)省份的相似性,因此極大地降低了估計誤差,模型估計結(jié)果為0.214,相較于全國平均水平0.293略低,比HT估計更為合理。類似地,對于東北、華東、中南等地區(qū),多層次模型估計結(jié)果相比HT估計均有略微下降,起到了誤差校正的作用。對于恩格爾系數(shù)低于全國平均水平的西南和西北地區(qū)而言,HT估計結(jié)果與模型估計結(jié)果無較大差異。由于樣本代表性較差,HT估計用于估計全國恩格爾系數(shù)不再可靠,估計結(jié)果為0.536,與國家統(tǒng)計局公布的數(shù)據(jù)差距較大,但本文所提模型估計結(jié)果為0.294,與國家統(tǒng)計局公布的數(shù)據(jù)幾乎一致,更加體現(xiàn)了本文所提模型的優(yōu)勢。
利用一套樣本進行子總體估計,往往存在劃分后的各小域樣本量較小且分布不均勻的問題,尤其是當所關注的子總體在抽樣設計階段未考慮到時,需要依賴較少甚至為0 的樣本量推斷子總體特征?;谀P偷耐茢嘈Ч軜颖玖恐萍s較小,相比基于設計的推斷更加適用于小域估計。本文方法采用多層次模型進行小域估計,盡可能挖掘樣本信息,第一層模型用于刻畫域間異質(zhì)性,第二層模型用于刻畫域間相關性,借助其他域樣本估計指定小域特征,并在此基礎上對抽樣誤差和非抽樣誤差中的測量誤差進行校正。模擬結(jié)果驗證了本文方法的效果:本文方法在樣本自身質(zhì)量較差時仍可作出較為合理、可靠的估計,且在自下而上進行多層次估計時仍能保證估計效果。
雖然基于模型能得到較好的估計結(jié)果,但不能因此而忽視數(shù)據(jù)采集階段的抽樣設計和質(zhì)量把控,高質(zhì)量的樣本對保證小域估計的精度至關重要。我國大型抽樣調(diào)查一般以省為單位進行樣本抽取和指標估計,調(diào)查樣本僅能實現(xiàn)國家層面和省級單位的估計,難以滿足市級、縣級等細分單位的低層次估計需要,此類多層次估計需求可以通過小域估計實現(xiàn),本文所提多層次模型為此提供了實踐方案。事實上,本文模型不僅可以估計比率,而且可以直接對其他分布的目標變量建模,或通過比率間接估計目標變量,兩種思路的實際效果如何需要另行討論??傮w而言,我國基于模型小域估計的理論研究與相關實踐還不多,本文作了一些嘗試,小域之間空間相關性的建模、域?qū)哟文P汀卧獙哟文P偷倪x擇等都有待進一步研究。