基于Brewer設計的不放回不等概率抽樣方法

王智光，閆在在，張婷婷

（1.內蒙古醫科大學附屬醫院 神經內科，呼和浩特 010050；2.內蒙古工業大學 理學院，呼和浩特 010051）

0 引言

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

1 經典不等概率抽樣方法介紹

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

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

Rao-Hartley-Cochran估計量是總體總值Y的無偏估計,具有方差：

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

具有方差：

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

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

2 抽樣方法的提出

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

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

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

對于上面提出的抽樣設計,構造總體總值估計量,并計算其均值、方差及其方差估計,給出解析表達式。如解析表達式不易得到或過于復雜,可以通過數值模擬來計算相應的量,進而比較本文提出的方法和Rao-Hartley-Cochran方法的精度。

2.1 n為偶數的情形

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

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

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

可以證明估計量是總體總值Y的無偏估計。

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

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

上述證明得到估計量是總體總值Y的無偏估計。下面給出估計量的方差表達式。

定理2：在提出抽樣設計下：

證明：

根據Horvitz-Thompson估計量的性質，有：

進一步化簡因此：

2.2 n為奇數的情形

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

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

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

按Rao-Hartley-Cochran方法,構造最后一組Gm的總體總值的估計：

進而構造總體總值Y的估計：

利用Rao-Hartley-Cochran的結果和n為偶數時的結果直接可得：

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

3 方差估計

條件泊松抽樣下,對于總體總值Y的估計量Y?CP,根據Horvitz-Thompson估計量可以得到方差的一個無偏估計：

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

Rao-Hartley-Cochran針對估計量,提出了方差的一個無偏估計：

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

對于本文提出的估計,提出其方差的一個無偏估計：

其中和是第g組抽到兩個單元的觀測值,相應的Zg*值為和

4 數值模擬

為了說明提出方法的優良性，本文基于不同超總體模型下的總體，利用Monte-Carlo模擬計算提出方法方差的第二部分或進而計算提出方法的方差和相對標準差。同樣可以利用Monte-Carlo模擬計算條件泊松抽樣方法的方差和相對標準差。

Rao-Hartley-Cochran方法的方差和相對標準差可以直接計算。調查變量根據八個不同的模型產生,每一個模型利用單變量回歸函數產生E(yk|x)=fk(x)，k=1，…，8。本文考慮下列回歸函數:

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)+?，其中Φ是標準正態分布函數

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

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

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

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

的 Monte-Carlo 模擬為隨機劃分的一個實現。因此估計量

的方差的Monte-Carlo模擬為：

估計量的相對標準差為：

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

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

本文根據方差估計的相對方差和置信區間覆蓋百分比對本文提出的方差估計和Rao-Hartley-Cochran提出的方差估計以及條件泊松抽樣設計相應的方差估計進行評價。

一個方差估計$v$的相對方差：

置信區間覆蓋百分比為：

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

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

5 結論

本文應用隨機劃分和Brewer方法,成功地設計了一種不等概率抽樣設計，與經典的Rao-Hartley-Cochran的方法相比，調查精度能夠達到RHC方法，具有實施方便簡單的優點，主要解決了構造高精度的方差估計的難點問題。并且將本文提出的方法和條件泊松抽樣方法進行了詳細地比較，結果發現本文所提出的方法在精度上明顯優于已有的條件泊松抽樣方法，而且實施也較條件泊松抽樣簡單很多，方差估計也較條件泊松抽樣簡單，所以本文所提出的方法是一種很有實際應用價值的不等概率抽樣方法。

參考文獻：

[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].應用數學學報,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）.