王智光，閆在在，張婷婷

（1.內(nèi)蒙古醫(yī)科大學(xué)附屬醫(yī)院 神經(jīng)內(nèi)科，呼和浩特 010050；2.內(nèi)蒙古工業(yè)大學(xué) 理學(xué)院，呼和浩特 010051）

0 引言

不放回不等概率抽樣是抽樣調(diào)查中的一種重要抽樣形式，在實(shí)際中被廣泛應(yīng)用。Horvitz和Thompson(1952)[1]給出了不放回不等概率抽樣下總體總值的估計(jì)量，估計(jì)量的方差及其方差估計(jì)；如何構(gòu)造和實(shí)施πPS抽樣設(shè)計(jì)是一個(gè)長(zhǎng)期公開(kāi)問(wèn)題。Brewer和Hanif(1983)[2]總結(jié)了50種不放回不等概率抽樣方法；Zou和Feng(1995)[3]也給出了一種新的不等概率抽樣方法。而在實(shí)際應(yīng)用中，n=2時(shí)的πPS應(yīng)用最為廣泛，研究也最充分。Brewer(1963)[4],Durbin(1967)[5]提出了樣本單元數(shù)為n=2時(shí)的嚴(yán)格不放回πPS抽樣方法;Rao等(1962)[6]提出了隨機(jī)劃分總體為n個(gè)子總體，在每個(gè)子總體中實(shí)施大小為1的不等概率抽樣設(shè)計(jì)及其理論。在上述研究成果基礎(chǔ)上，本文探索這些方法的改進(jìn)或者結(jié)合使用以獲取新的具有潛在應(yīng)用性且精度更高的抽樣設(shè)計(jì)方法，提出一種新的近似πPS抽樣設(shè)計(jì)并建立該設(shè)計(jì)相應(yīng)的理論。對(duì)于一些難以得到解析表達(dá)式的結(jié)果，可以利用統(tǒng)計(jì)軟件R作數(shù)值模擬計(jì)算相應(yīng)的量，從而實(shí)現(xiàn)提出方法和已有一些經(jīng)典方法精度比較的目的。具體地，Rao等[6]提出的方法：將總體中的單元隨機(jī)地劃分成n組，每組的單元數(shù)記為N1，N2，…，Nn（預(yù)先確定的），在每組中按與單元大小Zi成比例的概率抽取一個(gè)單元入樣，得到固定容量為n的樣本。簡(jiǎn)言之，此方法是將總體劃分，然后在每個(gè)子總體中實(shí)施樣本量為1的不等概率抽樣。受此方法啟示，眾所周知，存在一些經(jīng)典的樣本量為2的嚴(yán)格的不放回不等概率抽樣設(shè)計(jì)，本文擬利用此思想并結(jié)合Brewer(1963)[4]或者Durbin(1967)[5]提出的n=2時(shí)的嚴(yán)格不放回πPS抽樣方法構(gòu)建新的抽樣設(shè)計(jì)。

1 經(jīng)典不等概率抽樣方法介紹

Rao等[6]提出一個(gè)簡(jiǎn)單而適用的方法?？傮wU={1，2,…，N}，Y是研究變量。Zt是總體U中抽取一個(gè)單元取到第t個(gè)單元的概率。將總體中的單元隨機(jī)地分成n，每組的單元數(shù)記為N2，…，Nn，在第g組中，每個(gè)單元對(duì)應(yīng)的Z值重新記為，對(duì)應(yīng)的Y值重新記為Yi(

g)，i=1，2，…，Ng,Z(g)是第g組中Z值的總和。在每組中按與Zj(g),j=1，2，…，Ng,Z(g)成比例抽取一個(gè)單元，最后合成樣本容量為n的樣本s。記第g組中抽到的樣本單元觀測(cè)值yg，相應(yīng)的Z值記為zg?？傮w總值Y的Rao-Hartley-Cochran估計(jì)量定義為：

Rao-Hartley-Cochran估計(jì)量是總體總值Y的無(wú)偏估計(jì),具有方差：

Hajek（1964）[7]設(shè)計(jì)了一種不放回的近似嚴(yán)格πPS抽樣方法,即泊松抽樣，設(shè)計(jì)如下：對(duì)每個(gè)總體單元賦予一個(gè)入樣概率πi,使得πi/Zi=ν,其中ν是一個(gè)常數(shù)。以πi為成功概率，作一次Bernoulli試驗(yàn)，若試驗(yàn)成功，則相應(yīng)的單元入樣，共做N次試驗(yàn)，實(shí)際樣本容量是一個(gè)隨機(jī)變量?？傮w總值Y的無(wú)偏估計(jì)量定義為：

具有方差：

雖然泊松抽樣設(shè)計(jì)實(shí)施簡(jiǎn)單，但存在一大缺點(diǎn)，即樣本量n是隨機(jī)的。為了克服泊松抽樣設(shè)計(jì)上的缺點(diǎn)，Hajek(1981)[8]討論了一種固定樣本量n的泊松抽樣方法即條件泊松抽樣，它的具體實(shí)施方法如下:以pi，i∈U(滿足作為一組工作概率連續(xù)進(jìn)行泊松抽樣，直到出現(xiàn)容量正好等于預(yù)定的固定樣本量n的樣本，則抽樣結(jié)束,否則,繼續(xù)上述抽樣。條件泊松抽樣是基于泊松抽樣得到的一種不等概率抽樣設(shè)計(jì)，此設(shè)計(jì)是樣本量n固定，嚴(yán)格不放回的，但包含概率πi與單元大小不是嚴(yán)格成比例的近似πPS抽樣設(shè)計(jì)。之后，基于泊松抽樣、條件泊松抽樣的不放回不等概率抽樣設(shè)計(jì)有大量的研究，如Grafstrom(2009)[9]提出的Repeated Poisson抽樣和Laitila等（2011）[10]提出的一種二相πPS抽樣設(shè)計(jì)。條件泊松抽樣的一階包含概率計(jì)算有相應(yīng)的遞推公式，即：

二階包含概率的遞推公式為：

2 抽樣方法的提出

令總體U={1，2，…，N},y是研究變量。Zt是總體U中抽取一個(gè)單元取到第t個(gè)單元的概率,例如按與單元大小Xt成正比抽樣預(yù)定的樣本容量n,提出的抽樣方法由以下兩步組成:

(1)將總體中的單元隨機(jī)地分成m=n/2（若i=1，2，…，Ng為偶數(shù)）組或者m=(n+1)/2（若n為奇數(shù)）組,每組的單元數(shù)記為N1，N2，…，Nm,其中N1+N2+ … +Nt=N。

（2）如果n為偶數(shù),在每組中按照Brewer方法或者Durbin方法抽取兩個(gè)單元,即可得到容量為n的樣本s；如果n為奇數(shù),同樣在前面m-1個(gè)組中按照Brewer方法或者Durbin方法抽取兩個(gè)單元,在最后一個(gè)組中按單元大小成比例抽取一個(gè)單元,即可得到容量為n的樣本s。

對(duì)于上面提出的抽樣設(shè)計(jì),構(gòu)造總體總值估計(jì)量,并計(jì)算其均值、方差及其方差估計(jì),給出解析表達(dá)式。如解析表達(dá)式不易得到或過(guò)于復(fù)雜,可以通過(guò)數(shù)值模擬來(lái)計(jì)算相應(yīng)的量,進(jìn)而比較本文提出的方法和Rao-Hartley-Cochran方法的精度。

2.1 n為偶數(shù)的情形

假定將總體隨機(jī)劃分為m(=n/2)組：記為G1，G2，…，Gm。在第g組Gg中,每個(gè)單元對(duì)應(yīng)的Z值在此重新記為,對(duì)應(yīng)的Y值重新記為將每組中的Z值歸一化記為

利用Brewer方法,在第g組Gg中兩個(gè)樣本單元的抽取方法是：

提出的抽樣方案是容易實(shí)施的,但是抽樣機(jī)制相對(duì)復(fù)雜,可以想到估計(jì)量（6）的理論分析是困難的。下面從理論和數(shù)值上評(píng)價(jià)提出的估計(jì)（6）。本文記E1和V1分別表示隨機(jī)分組的數(shù)學(xué)期望和方差；E2和V2分別表示在固定分組條件下抽樣設(shè)計(jì)的數(shù)學(xué)期望和方差。

可以證明估計(jì)量是總體總值Y的無(wú)偏估計(jì)。

定理1：在提出抽樣設(shè)計(jì)下,E()=Y。

由式(1)和式(3)可以得到：

上述證明得到估計(jì)量是總體總值Y的無(wú)偏估計(jì)。下面給出估計(jì)量的方差表達(dá)式。

定理2：在提出抽樣設(shè)計(jì)下：

證明：

根據(jù)Horvitz-Thompson估計(jì)量的性質(zhì)，有：

進(jìn)一步化簡(jiǎn)因此：

2.2 n為奇數(shù)的情形

將總體隨機(jī)劃分為m(=(n+1)/2)個(gè)組：記為G1，G2，…，Gm。在第g組Gg中,每個(gè)單元對(duì)應(yīng)的Z值在此重新記為Zi(g),對(duì)應(yīng)的Y值重新記為Yi(g),i=1，2，…，將每組中的Z值歸一化記為

在每組Gg(1≤g≤m-1)中,按Brewer方法,抽取兩個(gè)樣本單元；在第m組Gm中按與，j=1，2，…，Nm成正比抽取一個(gè)單元,最后合成容量為n的樣本s。

按Brewer方法,構(gòu)造每組Gg(1≤g≤m-1)的總體總值的估計(jì)：

按Rao-Hartley-Cochran方法,構(gòu)造最后一組Gm的總體總值的估計(jì)：

進(jìn)而構(gòu)造總體總值Y的估計(jì)：

利用Rao-Hartley-Cochran的結(jié)果和n為偶數(shù)時(shí)的結(jié)果直接可得：

與條件泊松抽樣、Rao-Hartley-Cochran方法比較，提出方法精度的改進(jìn)，方差減少量的大小將通過(guò)數(shù)值模擬展示。

3 方差估計(jì)

條件泊松抽樣下,對(duì)于總體總值Y的估計(jì)量Y?CP,根據(jù)Horvitz-Thompson估計(jì)量可以得到方差的一個(gè)無(wú)偏估計(jì)：

其中πi,πj,πij在第上文已給出相應(yīng)的遞推公式。

Rao-Hartley-Cochran針對(duì)估計(jì)量,提出了方差的一個(gè)無(wú)偏估計(jì)：

其中Z(g)是第g組Ng個(gè)單元Zi值的總和,yg是第g組抽到單元的觀測(cè)值,相應(yīng)的Z為zg。

對(duì)于本文提出的估計(jì),提出其方差的一個(gè)無(wú)偏估計(jì)：

其中和是第g組抽到兩個(gè)單元的觀測(cè)值,相應(yīng)的Zg*值為和

4 數(shù)值模擬

為了說(shuō)明提出方法的優(yōu)良性，本文基于不同超總體模型下的總體，利用Monte-Carlo模擬計(jì)算提出方法方差的第二部分或進(jìn)而計(jì)算提出方法的方差和相對(duì)標(biāo)準(zhǔn)差。同樣可以利用Monte-Carlo模擬計(jì)算條件泊松抽樣方法的方差和相對(duì)標(biāo)準(zhǔn)差。

Rao-Hartley-Cochran方法的方差和相對(duì)標(biāo)準(zhǔn)差可以直接計(jì)算。調(diào)查變量根據(jù)八個(gè)不同的模型產(chǎn)生,每一個(gè)模型利用單變量回歸函數(shù)產(chǎn)生E(yk|x)=fk(x)，k=1，…，8。本文考慮下列回歸函數(shù):

Linear:y1=1+2(x-0.5)+?

Quadratic:y2=1+2(x-0.5)2+?

Bump:y3=1+2(x-0.5)+exp(-200(x-0.5)2)+?

Jump:y4=1+2(x-0.5)I(x≤0.65)+0.65I(x≥0.65)+?

CdF:y5=Φ((0.5-2x)/0.02)+?，其中Φ是標(biāo)準(zhǔn)正態(tài)分布函數(shù)

Exponential:y6=exp(-8x)+?

Cycle1:y7=2+sin(2πx)+?

Cycle4:y8=2+sin(8πx)+?,其中x∈(0，1),參見(jiàn)文獻(xiàn)[12]。

本文考慮x為一個(gè)有偏分布,實(shí)施模擬從Beta(66/49，165/49)分布獨(dú)立同分布地產(chǎn)生。擾動(dòng)變量?～N(0，,為了保持調(diào)查變量y與輔助變量x的主要回歸關(guān)系，正態(tài)擾動(dòng)變量的方差的大小的選擇為調(diào)查變量y的方差的1/5。對(duì)每一個(gè)模擬，按照提出的抽樣設(shè)計(jì)和RHC設(shè)計(jì)以及條件泊松抽樣設(shè)計(jì)，10000個(gè)樣本針對(duì)不同的總體容量和樣本容量產(chǎn)生。估計(jì)量和他們的方差估計(jì)被計(jì)算。通過(guò)下面的數(shù)量計(jì)算方法評(píng)價(jià)提出方法的性能。

的 Monte-Carlo 模擬為隨機(jī)劃分的一個(gè)實(shí)現(xiàn)。因此估計(jì)量

的方差的Monte-Carlo模擬為：

估計(jì)量的相對(duì)標(biāo)準(zhǔn)差為：

本文提出的方法與RHC方法、條件泊松抽樣方法的效率比較情況如表1所示。

表1 RHC方法、條件泊松抽樣的效率比較(N=60,n=6)

本文根據(jù)方差估計(jì)的相對(duì)方差和置信區(qū)間覆蓋百分比對(duì)本文提出的方差估計(jì)和Rao-Hartley-Cochran提出的方差估計(jì)以及條件泊松抽樣設(shè)計(jì)相應(yīng)的方差估計(jì)進(jìn)行評(píng)價(jià)。

一個(gè)方差估計(jì)$v$的相對(duì)方差：

置信區(qū)間覆蓋百分比為：

本文提出的方法與RHC方法、條件泊松抽樣方法的方差估計(jì)比較情況如表2所示。

表2 提出方法、RHC方法、條件泊松抽樣方法的方差估計(jì)比較(N=60,n=6)

5 結(jié)論

本文應(yīng)用隨機(jī)劃分和Brewer方法,成功地設(shè)計(jì)了一種不等概率抽樣設(shè)計(jì)，與經(jīng)典的Rao-Hartley-Cochran的方法相比，調(diào)查精度能夠達(dá)到RHC方法，具有實(shí)施方便簡(jiǎn)單的優(yōu)點(diǎn)，主要解決了構(gòu)造高精度的方差估計(jì)的難點(diǎn)問(wèn)題。并且將本文提出的方法和條件泊松抽樣方法進(jìn)行了詳細(xì)地比較，結(jié)果發(fā)現(xiàn)本文所提出的方法在精度上明顯優(yōu)于已有的條件泊松抽樣方法，而且實(shí)施也較條件泊松抽樣簡(jiǎn)單很多，方差估計(jì)也較條件泊松抽樣簡(jiǎn)單，所以本文所提出的方法是一種很有實(shí)際應(yīng)用價(jià)值的不等概率抽樣方法。

參考文獻(xiàn)：

[1]Horvitz D G,Thompson D J.A Generalization of Sampling Without Re?placement From a Finite Universe[J].Journal of the American Statisti?cal Association,1952,47（260）.

[2]Brewer K,Hanif M.Sampling With Unequal Inclusion Probabilities[M].New York：Springer-Verlag，1983.

[3]Zou G H,Feng S Y.A New Unequal Probability Sampling Design[C].Contributed Papers of 50th ISI,1995.

[4]Brewer K R W.A Model of Systematic Sampling With Unequal Proba?bilities～[J].Austral.J.Statist.,1963,（5）.

[5]Durbin J.Design of Multistage Surveys for the Estimation of Sampling Errors[J].Applied Statist,1967,（16）.

[6]Rao J N K,Harteley H O,Cochran W G.On a Simple Procedure of Un?equal Probability Sampling Without Replacement[J].Jour.Roy.Stat.Soc,1962,24（2）.

[7]Hajek J.Asymptotic Theory of Rejective Sampling With Varying Prob?abilities From a Finite Population[J].Annals of Mathematical Statis?tics,1964,（35）.

[8]Hajek J.Sampling From a Finite Population[M].New York:Marcel Dekker,1981.

[9]Grafstrom A.Repeated Poisson Sampling[J].Statistics and Probability Letters,2009,（79）.

[10]Laitila T,Olofsson J.A Two-Phase Sampling Scheme and πPS De?signs[J].Journal of Statistical Planning and Inference,2011,（141）.

[11]李苗苗,閆在在.條件泊松抽樣下二階包含概率的遞歸計(jì)算[J].應(yīng)用數(shù)學(xué)學(xué)報(bào),2014,37（1）.

[12]Montanari G E,Ranalli M G.Nonparametric Model Calibration Esti?mation in Survey Sampling[J].Journal of the American Statistical As?sociation,2005,100（472）.