Performance of Generalized Regression Estimator and Raking

Estimator in the Presence of Nonresponse

Daifeng Han1, Richard L. Valliant2, Jill M. Montaquila3, Keith Rust4

1Westat, 1600 Research Blvd, Rockville, MD 20850

2

Joint Program in Survey Methodology, University of Maryland, 1218 LeFrak Hall, College Park, MD 20742

3 Westat, 1600 Research Blvd, Rockville, MD 20850

4Westat, 1600 Research Blvd, Rockville, MD 20850

Abstract

Calibration weighting is widely used to decrease variance, reduce nonresponse bias, and improve the face validity of

survey estimates. In the purely sampling context, Deville & Särndal (1992) demonstrate that many alternative forms

of calibration weighting are asymptotically equivalent, so the generalized regression (GREG) estimator can be used

to approximate some general calibration estimators with no closed-form solutions. It is unclear whether this

conclusion holds when nonresponse exists and single-step calibration weighting is used to reduce nonresponse bias

(i.e., calibration is applied to the basic sampling weights directly without a separate nonresponse adjustment step).

In practice, poststratification (as a special form of the GREG estimator) and raking (as an example of general

calibration estimators) are commonly used calibration approaches, but decisions between these estimators are often

made ad-hoc based on sample sizes and availability of external data. In this paper, we compare the performance of

these estimators by examining their biases, variances, and effective coverage of the confidence intervals. The

theoretical work and simulation study demonstrate the need to consider models for both the outcome variable and

the response pattern. The model supporting the typical application of raking has main effects only while

poststratification (and more general forms of GREG) can include interactions. A framework involving both designbased and model-based thinking is developed to simultaneously evaluate the impact of sampling, outcome variable

structure, and nonresponse mechanism.

Since survey practitioners often lack the knowledge of the outcome variables and nonresponse mechanism in realworld surveys, we also develop a diagnostic that helps gauge the potential consequence of choosing an inappropriate

calibration estimator. The results of this research will provide guidelines for choosing between the commonly used

calibration estimators.

Keywords: calibration, GREG, poststratification, nonresponse adjustment

1. Introduction

Calibration weighting has originally been developed as a method for reducing sampling errors while retaining

randomization consistency. Deville and Särndal (1992) introduce calibration estimators using the distance function

approach. Later work by Särndal (2007) points out that there are two different approaches to take account of 

auxiliary information in the estimation – “calibration approach” and “regression approach”. The two approaches

generate the same estimator, the generalized-regression (GREG) estimator, in the situation where the general linear

squares (GLS) distance function is used in the calibration approach and linear regression model is used in the

regression approach. For the purpose of comparison, we use the term “general calibration estimators” to refer to the

other estimators in the calibration estimator family covered by Deville and Särndal (1992), as opposed to the GREG

estimator. 

Although almost all surveys in practice are subject to frame deficiencies and nonresponse, the theories in Deville 

and Särndal (1992) are developed for the ideal situation where non-sampling errors do not exist. In the purely

sampling context, many alternative forms of calibration weighting are asymptotically identical. This leads to a

breakthrough in our understanding of some commonly used calibration estimators that do not have closed-form 

1
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

solutions, such as raking. As a result, the GREG estimator is often considered a good approximation of the general

calibration estimators. However, non-sampling errors such as nonresponse almost always exist in real-world surveys.

In the past decade, Särndal and Lundström (1999, 2005), Kott and Chang (2006, 2008, 2010), and Thibaudeau and

Slud (2009) have proposed different methods for using calibration to correct nonresponse bias through one-step

weighting, yet we still lack understanding of the empirical properties of the calibration estimators generated by these

methods. For example, it is unclear whether (or under what conditions) the GREG estimator and the general

calibration estimators are asymptotically equivalent when nonresponse is present in a survey and calibration

weighting is used to reduce potential nonresponse bias. In practice, poststratification and raking are both widely

used in the US and European surveys, and the corresponding estimators are the special cases of the GREG estimator

and the general calibration estimator, respectively. Quite often survey practitioners choose between the two

estimators based on the availability of the benchmark totals and the case counts in the survey requiring calibration

with the hope that poststratification and raking reduce potential nonresponse bias to a similar extent and thus result

in “approximately equivalent” estimators. Sometimes only the marginal totals are available, and the practitioner may

have no choice but to use raking, even if poststratification might have done a better job in reducing nonresponse bias.

No systematic research has been conducted on comparing the performance of the poststratification estimator and the

raking estimator when calibration is used to correct nonresponse bias.

The first contribution we attempt to make to the literature is to release the assumption of no non-sampling error and

evaluate the properties of the calibration estimators when calibration is used to correct nonresponse bias through a 

one-step weighting approach. In the absence of non-sampling errors, the purely design-based properties of the

calibration estimators were assessed by Deville and Särndal (1992). When nonresponse exists, however, the 

properties of a calibration estimator may depend on the underlying outcome variable model and response model.

Our research evaluates the impacts of sampling, population structure, and response mechanism simultaneously 

through analytical work and simulation studies. The theoretical and empirical results provide survey practitioners

with guidance of how to evaluate different calibration estimators and choose between them under various population

structure models and response mechanism models.

The rest of the article is organized as follows: Section 2 summarizes the literature on the properties of various

calibration estimators in the purely sampling context as well as the research on using calibration for nonresponse

adjustment through single-step weighting. Sections 3 and 4 attempt to fill in some gaps in the literature through

analytical work and a simulation study, respectively. Section 5 summarizes the findings briefly and discusses the

direction of our future research.

2. Literature Review

2.1 Calibration in Absence of Non-sampling Error

There are two different approaches for incorporating auxiliary information in the estimation, labeled “regression

approach” and “calibration approach”. Under the umbrella of calibration approach, Section 2.1.1 explains how to

use the distance function method to obtain a calibration estimator, followed by the description of an alternative

method called the function form method in Section 2.1.2. Section 2.1.3 presents the existing theories developed by

Deville and Särndal (1992) on the comparison of the GREG estimator and the general calibration estimator. It is

important to note that all the theory presented in this section was developed for the situations absent of nonsampling error.

2.1.1 Two Approaches to Incorporate Auxiliary Information in Estimation

There are two systematic ways to take account of auxiliary information in the estimation. In their original definition

of the calibration estimator, Deville and Särndal (1992) require “minimum distance” between the calibration weights

and the original sampling weights, subject to satisfying the calibration equation. In general, the term “calibration

approach” often refers to creating estimators by benchmarking the auxiliary information to external controls.

For a sample s drawn from a population U , let yk

be the value of the variable of interest, y, for the kth population

element, which is associated with an auxiliary vector value, xk = ( , xk1 ...xkp ,..., xkP )

T

. For the elements k s ∈ , we

2
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

observe ( , ykxk) . For simplicity, the population total of x, t

x = ∑xk

, which is often referred to as the benchmark

U

control vector, is assumed to be accurately known. 

The objective is to estimate the population total t

y = ∑yk

. Let dkbe the basic sampling design weight calculated as 

U

the inverse of the inclusion probability π . The Horvitz-Thompson estimator is tˆ = y /π = d y . The k yπ ∑ k k ∑ k k

s s

calibration estimator is defined as tˆ w y , with weights w as close as possible, in an average sense based on yw = ∑ k k k

s

a distance function, to the basic design weights dk while respecting the calibration equation

wk k = t

x ∑ x (2.1)

s

Under a chosen distance function Gwd ( , ) , this becomes an optimization problem. The goal is to find a set of k k k

weights {w } that minimizes

∈ 

Gwd ( , )subject to (2.1). This leads to the Lagrange function k k s ∈ ∑k s k k k

Ψ= ( , ) + λ t −

∈ 

Gw d w (2.2) T ∑ x k k ( x ∑ k k ) k s ∈ k s

which is minimized to find the optimal set of weights{w } ∈ 

.

k k s

The calibration weights can be expressed as

w = d F (x

T

λ) (2.3) k k k k

where λ = ( , λ1 ...,λ1

,...,λ P)

T

is the vector of Lagrange multipliers determined from (2.2). λ corresponds to a

realized sample, but for simplicity we sometimes use λ as the shorthand for λ . F (x

T

λ) is the inverse function of

s k k

g wd (,) = ∂Gwd ( , ) / ∂w , the first derivative of the distance function taken with respect to the calibration k k k k k k k

T weight. F (x λ) uniquely corresponds to Gwd ( , ) . It is assumed that F is non-negative and convex, and that k k k k k k

Fk(1) = 0 , implying that when wk = dkthe distance between the basic design weights and calibrated weights is zero.

Moreover, it is required that F

'

is continuous, monotonic, and that F

'

(1) = 0 and F

''(1) > 0 , which makes w = d a k k k k k

local minimum.

The Horvitz-Thompson estimator of t isˆt d x , so the calibration equation can be expressed as:

x xπ = ∑ k k

s

(

T

) − d x = t − ˆ ∑ d F x λ x t (2.4) kk k k ∑ k k x xπ 

s s

Define

T Φ λ = d F x λ −1 x (2.5)

s ( ) ∑ k { k ( k ) } k

s

Then (2.4) can be written as

Φs (λ) = t

x − t

ˆ

xπ 

(2.6)

The task of obtaining wk

boils down to solving (2.6) for λ . The calibration estimator of t

y

is

tˆ = w y = d F (x λ) y (2.7)

yw ∑ k k ∑ k k

T

k k

s s

Depending on the distance function Gwd ( , ) , an iterative process may be required to obtain a solution for λ.

k k k

There is no non-sampling error in Deville and Särndal’s setup, so the Horvitz-Thompson estimator t

ˆ

yπ 

using basic

sampling weights dk

is an unbiased estimator of the true population total. If the calibration weights wkare as close 

as possible, according to Gwd ( , ) , to the basic sampling weights d , then a realistic expectation is that the k k k k

calibration weights will maintain near unbiasedness.

3
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Although several distance functions are discussed in Deville and Särndal (1992), most theoretical research has

2

focused on the GLS distance function (w − d ) / d q , where 1 qk ∑ is the positive weight associated with the kth k k k k

s

term and is unrelated to dk

. Under this distance function, the calibration equation has a closed-form solution. We

x

T

λ 1

T

obtain F ( )= + q x λ , and the calibration estimator is the GREG estimator k k k k

T T t

ˆ = d (1+ q x λ) y = t

ˆ + (t − t

ˆ ) Bˆ

(2.8)

yreg ∑ k k k k yπ x xπ s

s

where

−1

( λ = T t − ˆt ) (2.9) s x xπ 

ˆ

−1 B = T d q x y (2.10)

s s ∑ kkkk

s

T = ∑d q x x

T

(2.11)

s k k k k

s

An alternative method for obtaining the calibration estimator is referred to as the “regression approach”. With the

regression approach, estimators are calculated by using an assisting model that closely represents the relationship

between the outcome variable and the auxiliary variables. The assisting model is also referred to as the calibration 

model or the working prediction model by Kott (2006) to distinguish it from other models such as those used to

address response propensity. The assisting model can have linear or nonlinear forms. When the assisting model is a 

linear regression model, the weight happens to be calibrated to the auxiliary controls and the estimator (which is the

GREG estimator) is expressible as a linearly weighted sum with calibrated weights as a by-product. One advantage

of the GREG estimator is that the calibrated weights are independent of any particular outcome variable y and can

therefore be applied to all the variables of interest in a survey.

Our research adopts the perspectives of both approaches. The weights are primarily justified by their consistency 

with the benchmark controls (which is the calibration approach). Although the calibration approach does not refer

explicitly to any assisting models, we demonstrate that the performance of a calibration estimator in the presence of

nonresponse depends on the choice of auxiliary vector and/or function form used in the calibration process, and this

requires a modeling effort in some sense.

2.1.2 Distance Function Method versus Function Form Method

Under the umbrella of the calibration approach, two methods are discussed in the literature. Deville and Särndal 

(1992) initially require that the set of calibration weights {w } minimize some distance function 

∈ 

Gwd (,) k k s ∈ ∑k s k k k

subject to satisfying the calibration equation – This is the “distance function method” described in Section 1.1. An

alternative approach is the “function form method” or “instrumental vector method” (Estevao and Särndal 2006,

Kott 2006), which can also generate many alternative sets of weights calibrated to the same auxiliary information.

The function form method removes the limitation that the calibration weights minimize a distance function, and

requires only that {w } ∈ 

satisfy the calibration equation and be of the function form w = d F(z

T

λ) , where d is k k s k k k k

the design weight, and zkis a vector with values defined for the units in the sample and sharing the dimension of the

specified benchmark control vector xk. The vector zkis called the instrumental vector for the calibration, which can 

be a specified function of xkor of other background data about unit k (Särndal and Lundström 2005). The vector λ

is determined from the calibration equation. The function F( )⋅ plays a similar role as Gwd ( , ) does in the k k k

distance minimization method. For easy reference, we refer to F( )⋅ as “weight adjustment function” or “adjustment

function” in our research. One possible form of the weight adjustment function is w = d (1+ z

T

λ) , and the k k k

corresponding calibration estimator is:

tˆ

ycal = ∑dk

(1+ z

T

k

λ) yk

s

(2.12)

where

4
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

T −1

λ = (∑d x z ) (t − d x ) (2.13)

k k k x ∑ k k

s s

The GREG estimator t

ˆ defined in (2.8) is a special case of (2.12) obtained for z = q x .

yreg k k k

When nonresponse exists in a survey, we think that it is more appropriate to understand the calibration process using

the function form method rather than the distance function method. This is because in the presence of nonresponse,

the Horvitz-Thompson estimator for the total of an outcome variable y using the basic design weights becomes

tˆ

yπ = ∑ k k d y , where r represents the responding set. This estimator is biased when r s ≠ . If the calibration process

r

aims to correct the nonresponse bias, it is neither necessary nor appropriate to require the calibrated weights to be

“as close as possible” to the basic design weights based on a distance function.

More discussions about the weighting adjustment function F( )⋅ are included in Section 2.2.1. When applying the

function form method in practice, survey practitioners face questions such as how to choose the variables to be

included in the vector zk

, or what is the advantage of making the form of zkdifferent from the calibration vector

xk. These questions have not been clearly answered by the existing literature. Särndal (2007) gives an example

showing that “even ‘deliberately awkward choices’ for zk

give surprisingly good results”. However, the conclusion

of the near-unbiasedness of the calibration estimator in this situation seems to depend on the assumption of no nonsampling error, which may not hold in the presence of nonresponse.

2.1.3 Relationship between GREG Estimator and General Calibration Estimators

As described in Section 1.1, various calibration estimators can be derived with the aid of different distance measures

under the same set of constraints on the auxiliary variables. Alternative distance functions are compared in Deville,

Särndal, and Sautory (1993), Singh and Mohl (1996), and Stukel, Hidiroglou, and Särndal (1996). In a purely

sampling context, there are usually very small differences between the point estimates corresponding to the various

distance functions, and changes in the distance function often have minor effects on the variance of the calibration 

estimator even if the sample size is rather small. In particular, the GREG estimator and the other members of the 

calibration estimator family (referred to as the “general calibration estimators”) are compared in Deville and Särndal

(1992). They conclude that the GREG estimator is a first approximation to the general calibration estimators; all the

general calibration estimators are asymptotically equivalent to the GREG; and the variance estimator for the GREG

can be used for the general calibration estimators. Although the GREG estimator is a special case of the calibration

T

=+ T

estimator to denote the GREG estimator and ˆ family when the function form is F(x λ) 1 q x λ , we use t

ˆ t

k k k yreg yw

to denote the other calibration estimators (i.e., the general calibration estimators) for the purpose of comparison.

Deville and Särndal (1992) consider a sequence of finite populations and sampling designs indexed by n, where n is

the sample size (for a fixed-sized sampling design) or the expected sample size (for a random-sized sample design).

The finite population size, N, tends to infinity with n. Several assumptions are made about the auxiliary vector x: (i)

1 −1 −1/2

lim N )

−

t exists; (ii) N (

ˆt − t ) = O (n , where the subscript p means probability induced by the sample; and x xπ x p

1/2 −1

(iii) n N (

ˆt − t ) converges in distribution to the multinormal N( , 0 A) . Two additional assumptions are also xπ x

made for proving their Results 3-5: (iv) max x = M < ∞ , where max is over n as well as over k; and (v) k

max Fk

''(0) = M

'

< ∞ . Assumptions (i) through (iii) have two practical implications. First, the components of

ˆ ˆ t

xπ − tx

are considered small and quantities on the order of

2

are considered negligible. Second, ˆt

xπ − tx

t

xπ − tx

−1

follows an approximately normal distribution with covariance matrix n N2A (where A can be viewed as a matrix

that describes an asymptotic effect of the sampling design used for the survey), and this is to justify the use of the

normal approximation in confidence intervals based on t

ˆ

yw

. Assumption (iv) is usually satisfied in practice.

Assumption (v) is verified for all the calibration estimators given in Deville and Särndal (1992).

Deville and Särndal (1992) show five results. Result 1 states that the calibration equation (2.6) has a unique solution 

belonging to an open neighborhood of 0, with probability tending to 1 as n → ∞. Results 2 and 3 are about the

5
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

−1 −1 −1/2 magnitude of the Lagrange multiplier. They prove that λs s = T t( x − ˆt

xπ ) +O n p p ( ) = O (n ) , so λ tends to 0 in s

design probability as n → ∞. Result 4 indicates that the general calibration estimators are design-consistent, and the

difference between the general calibration estimators and the Horvitz-Thompson estimator is asymptotically zero.

−1 −1/2

That is, N (tˆ − tˆ ) = O (n ). Result 5 compares the general calibration estimators with the GREG estimator. yw yπ p

For any weight adjustment function F (⋅) obeying their assumptions, t

ˆ

yw

given by equation (2.7) is asymptotically k

−1 −1

equivalent to the GREG estimator given by equation (2.8), in the sense that N (tˆ − tˆ

yw yreg ) = Op (n ). Results 4

and 5 together show that as n → ∞, the difference between the general calibration estimators and the GREG

estimator approaches to zero faster than the difference between the general calibration estimators and the HorvitzThompson estimator. The asymptotic variance of t

ˆ

yw

is thus the same as that of the GREG estimator.

The results above have important practical implications because many general calibration estimators do not have a

closed-form solution. For example, although the raking ratio method has a long history, the variance of the raking

estimator is difficult to derive even approximately. Deville and Särndal (1992) resolve the problem by using the 

property that the general calibration estimators and the GREG estimator are asymptotically equivalent. Thus the

large-sample variance of the raking ratio estimator can be calculated using the same formula as that for the GREG

estimator, given in Särndal, Swensson, and Wretman (1992).

It is important to note that all the results in Deville and Särndal (1992) are derived under the assumptions i)

−1 −1/2

)

1/2 −1 N (

ˆt − t ) = O (n ; and ii) n N (

ˆ

xπ x p t − t ) converges in distribution to a multinormal distribution with xπ x

mean of 0. That is, they require the Horvitz-Thompson estimator with the basic design weights to be unbiased and

consistent, which is true only in the purely sampling context. When non-sampling errors exist, it is unclear whether

the GREG estimator is still asymptotically equivalent to other calibration estimators.

2.2 Calibration for Nonresponse Bias Reduction

Our research focuses on the non-sampling error caused by nonresponse. To keep the picture simple, we assume that

the sampling frame has perfect coverage and there is no measurement error in the survey. In practice, calibration is

widely used to correct nonresponse bias in government-sponsored studies such as the National Health and Nutrition

Examination Survey and Medical Expenditure Panel Survey. There are several variations in the literature on how to

adjust for nonresponse and calibrate the weights to benchmark controls. The conventional approach uses auxiliary

information in two steps (Kalton and Flores-Cervantes 2003). In step (i), a response model is formed based on the

patterns of correlation between the response probabilities and available auxiliary variables. The aim is to derive

good proxies of the unknown response probabilities, so as to limit the nonresponse bias as much as possible. In step

(ii), the goal is to select the auxiliary variables that best meet the dual purpose of reducing the sampling variance and

of giving added protection against nonresponse bias. The conventional approach is embodied in the estimator

t

ˆ

y = ∑ dk

(1/ p y ˆ

k

)

k

, where dkis the design weight calculated as the inverse of the selection probability π k, pˆ

k

is

r

the estimated response propensity, and r is the set of respondents. Survey practitioners usually act (for the purpose 

of variance estimation, for example) as if π k

pˆ

k was the true selection probability of element k . An unavoidable bias

results from replacing the unknown pk with pˆ

k

based on limited auxiliary information. An alternative approach is

to skip explicitly estimating the response propensity, but use calibration for nonresponse adjustment directly. This

approach has the potential to simplify the derivation of the variance estimation formulas, so we adopt this approach

in our work. A single-step weighting approach through calibration was first proposed by Särndal and Lundström 

(1999, 2005). The literature has been expanded by Kott (2005) and Chang and Kott (2008, 2010) in the past decade.

2.2.1 Särndal and Lundström Method

In the Särndal and Lundström method, all the auxiliary controls, from either the population (with a control vector of

dimension J

*

) or the sample (with a control vector of dimension J

o

), are included in the calibration equation, with the

x

*

⎛ ˆ k

dual purpose of reducing both sampling error and nonresponse bias. The auxiliary vector x = ⎜  ̃ has dimension k ⎜  ̃ x

o

⎝ ̄ k

6
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

⎛ x

*

ˆ 

* o ∑ J J + . The corresponding information input is t = 

U k

⎜  ̃ . We seek a weighting system w for k r ∈ the

o

x

⎜∑ d x  ̃ k , 

⎝ s

k k  ̄ 

respondent set, that satisfies the calibration equation ∑ w x = t . The calibrated weights are w d = v , where v kk x k kk k r

corresponds to the weighting adjustment function F( )⋅ described in Section 2.1.2 and can take different forms.

Although the distance function method is used in Lundström and Särndal (1999) for obtaining wk

, their later work

adopts the function form method, which seems more appropriate when nonresponse is present and calibration is used

to correct nonresponse bias. The calibration equation poses only weak constraints on the weights. Depending on the

form vk

takes, there exist many sets of calibrated weights for a given auxiliary vector xk. Särndal and Lundström

(2005) discuss two alternative schemes for defining the function form for vk

: (i) as a function of the auxiliary vector

xk; and (ii) as a function of any vector z (referred to as instrumental vector) specified for k r k ∈ and with the same

dimension as xk

. Under scheme (i), vkshould reflect the known individual characteristics of the element k r ∈ ,

summarized by the vector value . The calibration equation can be expressed as d F(x

T

x λ ) = t , where λ is a k ∑ k kr x r r

vector to be determined through the calibration equation. A simple function form is recommended that depends 

T T T −1

linearly on x : F(x λ ) = +1 x λ , where λ = ( d x x ) (t − d x ) . Although nonlinear forms can also be k k r k r r ∑ k k k x ∑ k k

r r

T T considered such as F(x λ ) = exp(x λ ) , Särndal and Lundström (2005) suggest that the linear form will suffice k r k r

due to its considerable computational advantage and fits the routine production environment. However, little

theoretical or empirical justification is provided to support this statement. 

An alternative scheme is to define the weighting adjustment function using the instrumental vector zk

, specified for

k r xk. The vector zkxk ∈ and with the same dimension as can be a specified function of or any background data

about k . Only linear function form based on zk

is considered by Särndal and Lundström (2005). The calibrated

T T −1 weights are w d = (1+ z λ ) , where λ = ( d x z ) (t − d x ) . Besides alerting the reader to this generality k k kr r ∑ k k k x ∑ k k

r r

of the calibration approach, Särndal and Lundström (2005) give little information about how to choose zk

except to

suggest that zk = xk

is the “standard choice”.

Särndal and Lundström (2005) claim that their approach meets the double objective of reducing sampling error and

nonresponse error in the presence of powerful auxiliary information. A variance estimator is also developed to take

into account the increased variance caused by nonresponse, as described in Section 2.2.

We can see that scheme (ii) is the generalization of scheme (i) in Särndal and Lundström (2005). When zk = xk

, the

two schemes give identical estimators. Furthermore, when r s xk = 

*

k = (indicating full response) and x (meaning

that the auxiliary vector contains information only from external benchmarks and not from the sampling frame), the

calibration estimator w d (1+

T

= x λ ) and the GREG estimator defined in equation (2.8) are identical. k k k r

2.2.2 Kott and Chang Method

Recent development by Kott (2006) and Chang and Kott (2008, 2010) emphasize that the set of variables modeling

the response mechanism (referred to as “model variables”) can be different from the benchmark variables in the

calibration equation. The vector for the benchmark controls in the calibration equation is still xk

, with known

population totals t

x

. Unit nonresponse is viewed as an additional phase of Poisson sampling. Using the quasirandomization perspective, each element k in the original sample is assumed to have a response probability pk

( )⋅ , 

which is a function of the response model covariate vector zk

. The statistician can specify the function form for 

pk( )⋅ and the unknown parameters in the function can be estimated. Kott (2006) explains why it may be desirable to

use a vector zk

that is different from xk– Sometimes the variables the response mechanism depends on are known

7
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

only for respondents, not for the whole sample. For example, in an agriculture survey, the benchmark variables can 

be previous-census frame variables known for every farm in the population while the model variables are currentperiod variables known only for survey respondents. Kott (2006) still requires that the dimensions of zk

and xk

coincide. Chang and Kott (2008) expand the method so that it allows the number of benchmark variables (i.e., the

dimension of xk

) to exceed the number of response model covariates (i.e., the dimension of zk).

Based on Kott and Chang, the response propensity for each responding unit k can be specified as p(z

T

kβ) , an

unknown but estimable linear combination of the response model covariate vector zk

. The input weight for the 

calibration equation is calculated as the product of basic design weight dk

and 1 (z

T

p β , and then the vector β can k

)

dk

be estimated from the data using the calibration equation t

x = ∑k r T ∈ x . This equation is sufficient to

p(z β)

k

k

determine βˆ

if the dimension of xkequals the dimension of zk(Kott 2006). On the other hand, when the dimension

of xkexceeds the dimension of zk, the calibration equation can be modified into a nonlinear regression-type

dk

equation t = ∑k r ∈ x + ε , where zk

and xkdenote the vectors for response model covariates and benchmark x T k

p(z β)

k

variables respectively, t

x

is the calibration target values consisting the known population totals, and ε is the error

term between the calibrated estimates and the population controls of the auxiliary variable (Chang and Kott 2008).

Although in theory, the response propensity p( )⋅ can take different forms, the discussions in Chang and Kott (2008)

are restricted to the situation where the response propensity is the linear function of the response model covariates.

2.3 Summary of Gaps in the Literature

In Sections 2.1 and 2.2, we provided a review of the calibration literature. In this section, we note gaps in the

literature that remain to be addressed.

First, there is little research on evaluating the asymptotic properties of different calibration estimators when

nonresponse is present and calibration is used for correcting nonresponse bias. For example, the raking ratio 

estimator and poststratification estimator are special cases of the general calibration estimator and GREG estimator

respectively and widely used in practice. Based on Deville and Särndal (1992), these two estimators are

asymptotically equivalent in absence of non-sampling errors. However, non-sampling errors such as nonresponse

error almost always exist in surveys. It is important to re-examine conclusions in Deville and Särndal (1992) in the

context of using calibration for nonresponse adjustment. 

Second, if the conclusions in Deville and Särndal (1992) do not hold when calibration is used for nonresponse

adjustment, then the existing literature provides neither a good framework for comparing the performances of

different calibration estimators, nor practical guidance for choosing the appropriate auxiliary vectors and/or function

forms for calibration weighting. To answer these questions, we need to go beyond the purely design-based approach

used in Deville and Särndal (1992) and examine the underlying models for population structure (i.e., what variables

are correlated with the key outcome variable) and response mechanism (i.e., what variables are corrected with

response). Survey practice calls for guidelines for how to select variables to be included in the auxiliary vectors xk

and response model covariate vector zk when calibration is used for nonresponse adjustment, but there is very little

research in this area.

3. Comparison of the Design-Based Properties of GREG Estimator and General Calibration Estimator in the

Presence of Nonresponse

This section attempts to fill in one gap in the literature by comparing the asymptotic properties of the general

calibration estimator and GREG estimator when calibration is used for nonresponse adjustment through a single-step

weighting approach. We use the Särndal and Lundström calibration method described in Section 2.2.1 and focus on

the situation where the instrumental vector zk

coincides with the calibration vector xk. In the presence of

nonresponse, the Horvitz-Thompson estimator of the total for the auxiliary vector using the basic design weights is a

8
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

function of the respondent set and can therefore be “far” from the benchmark control total. This violates one of the

key assumptions in Deville and Särndal (1992), so it is unclear whether their conclusions still hold.

The setup and analytical work in this section largely follow the approach taken by Deville and Särndal (1992),

which is purely design-based. We use the terms “new assumption” and “new result” to differentiate our assumptions

and findings from those in Deville and Särndal (1992). Our theoretical results are applicable to a family of general

calibration estimators. At the end of this section, we point out the limitations of the purely design-based approach

and emphasize the importance of examining the underlying models for the outcome variable and response rate when

comparing different calibration estimators. 

3.1 Scope and Assumptions

First, we assume the analytic survey (i.e., the survey requiring calibration) and benchmark survey come from the

same population U of size N . Although the benchmark control totals are often estimated and subject to sampling

and non-sampling errors in practice, we assume that the total for the auxiliary vector x is known and equal to the true 

population total t

x = ∑xk

. 

U

Second, we assume that the analytic survey has no coverage or measurement error, but may suffer from nonresponse

error that can bias the estimates of parameters such as population totals. In the presence of nonresponse, the survey

has a respondent set r of size nr

. We assume that no separate nonresponse adjustment is conducted prior to

calibration, so that in the absence of calibration, the population estimates are calculated using only the basic design

weights dk

, i.e., using the Horvitz-Thompson estimators of the population totals of the auxiliary vector and outcome

variable ˆt d x and tˆ d y respectively. The proofs in this section require no explicit specifications r = ∑ k k r = ∑ k k

r r

xπ xπ 

for the sample design of the analytic survey, but it is reasonable to assume that the units in the survey are selected

with a method that results in unbiased estimates of the totals for various variables in the absence of nonresponse. 

Finally, although survey nonresponse is caused by a random mechanism, for the simplicity of theoretical derivation,

we assume that each population member has fixed response propensity of either 1 or 0. That is, a population member

either always agrees or always refuses to participate in the analytic survey. In the presence of nonresponse, the

design-based expectation of the Horvitz-Thompson estimator reflects the characteristics of the “responding

population” U of size N . We define E (t

ˆ ) = t and E t( ) ˆ = t , where E means design-based expectation,

r r π r r xπ xπ π r π

yπ 

ryπ 

and t and t are the totals of the auxiliary variables and outcome variable for the respondent population U ,

r r r x y

respectively. It is reasonable to assume that the size of the responding population, N , is large so that N N = O(1) .

r r

The theoretical derivation in this section requires the following assumptions. We refer to these as “new assumptions”

in contrast of those in Deville and Särndal (1992).

−1

t

−1 −1 New assumption (i): lim N exist, but in general, lim N t ≠ N t .

r r r r x x x

− −1 −1/2 New assumption (ii): N

1

(

ˆt − t ) → 0 in design probability. N ˆt − t = O (n ).

r r r r ( r r ) p r xπ x xπ x

1 2 −1 New assumption (iii). n N (

ˆt − t ) converges in distribution to the multinormal N( , 0 A) , where A can be r r r r xπ x

viewed as a matrix that describes an asymptotic effect of the sampling design used for the analytic survey.

Recall that one of the key assumptions in Deville and Särndal (1992) is that in the purely sampling context, the

Horvitz-Thompson estimators of the population totals of the auxiliary vector approach to the true values of the

−1 ˆ

−1/2

population as the sample size increases. That is, N (t )

π − t

x ) = O (n . Based on our new assumption (ii), the x p

Horvitz-Thompson estimators from the respondent set approach only to t = E (t

ˆ ) . We know that t ≠ t in the r π r r x x xπ x

presence of nonresponse. This has important implications in the theoretical derivation below.

9
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

3.2 Analytical Results 

In this section we re-examine the results in Deville and Särndal (1992) in the context of using calibration for

nonresponse adjustment through single-step weighting. The input weights for the calibration equation are the basic 

design weights d . In this setup, the Horvitz-Thompson estimator t

ˆ using the basic sampling weights d is k r k yπ 

biased due to nonresponse (unless the nonresponse is generated by a missing completely at random mechanism), so

calibration is used to reduce such bias to the extent possible. It is appropriate to conduct calibration using the 

function form method instead of the distance function method, and the calibration weights may not be as “close” to

the basic design weights as required in Deville and Särndal (1992). We suspect that whether the calibration equation

has a solution may depend on the overall response rate as well as how the response rates differ by subgroups formed

by the variables used as benchmark controls in the calibration. We show that the vector for the Lagrange multiplier

determined from the calibration equation, λr

, consists of a term that is driven by the difference between the HorvitzThompson estimator of the auxiliary vector (using the basic design weights) for the respondent population total

(denoted by ˆt

r

) and the benchmark control total (denoted by t

x

). Unless nonresponse is negligible, this term does

xπ 

not decrease as the survey sample size increases, so λr may tend to a non-zero constant vector in design probability.

Our analytical work results in the formulas for: (1) the difference between a general calibration estimator and

Horvitz-Thompson estimator; and (2) the difference between a general calibration estimator and the GREG

estimator, in the presence of nonresponse. We prove that when nonresponse exists and calibration is used to reduce

nonresponse bias through single-step weighting, the general calibration estimators and the GREG estimator are not

asymptotically equivalent in general situations. 

In the presence of nonresponse, the calibration equation is ∑ w x = t and the calibration estimator is tˆ w y .

kk x yw = ∑ k k

r r

Equations (2.5) and (2.6) in Section 2 should be modified into

Φ ( ) λ = d F{ (x λ) −1} x

T

(3.1)

r ∑ k k k k

r

Φ (λ) = t − ˆt = (t − t ) + (t − t

ˆ ) (3.2) r x r x r r r xπ x x xπ 

Note that in this section, we sometimes use λ as a short-handed form of λrfor convenience.

We know that t

r ≠ tx

in the presence of nonresponse, so the right-hand side of (3.2) contains a non-zero term that

x

does not exist in equation (2.6) of Section 2. This non-zero term plays an important role in the discussions below.

We have five new results in parallel to the ones in Deville and Särndal (1992).

New Result 1. As nr → ∞, whether equation (3.2) has a solution may depend on the difference between t

r

and t

x x

as well as the function form Fk (⋅) used in the calibration.

For this result, we give intuitive explanations instead of strict proof. In the presence of nonresponse, the calibration

equation can be written as:

−1 −1 T −1 −1 −1 N Φ (λ ) = N d F ( ) − N d x = N (t − t ) + N (t − ˆ x λ x t ) (3.3)

r r r ∑ kk k k r ∑ k k r x r r r r

r r

x x xπ 

For the right-hand side of (3.3), the second term is similar to that in Deville and Särndal (1992). 

−1 −1/2 N (t − t

ˆ

r r r ) = Op r (n ) and is asymptotically 0. However, when nonresponse exists, t

x r ≠ t and the first term x xπ x

N

−1

(t − t ) = O(1) . Due to this additional term, the right-hand side of (3.3) does not tend to 0, but becomes a nonr x rx

zero constant vector as nr

increases.

A more intuitive way to understand this result is that in Deville and Särndal (1992), only “small” adjustments need

to be made to the basic design weights to obtain the calibration weights, and that is essentially why the calibration

10
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

equation almost always has a solution for large samples. When nonresponse exists, the Horvitz-Thompson estimator

∑d x may be “far” from the benchmark controls t and therefore “large” adjustments on the basic design weights k k x

r

may be required to satisfy the calibration constraints. In this situation, whether the calibration equation has a

solution may depend on the difference between t and t as well as the function form F (⋅) used in the calibration. r x k x

An empirical example is that for the same calibration constraints and respondent set r , poststratification always has

a solution but raking does not always converge. 

New Result 2. Let λ be the solution to equation (3.3) if one exists. If t − t ≠ 0 , then λ = O ( )1 in general

r x r r p x

situations. This means that λr

tends to a non-zero vector in design probability. 

−1 −1 −1 −1

Proof: Define z = N (t − t ) and z = N (t − ˆt ) , so λ = (N Φ ) (z + z ) if a solution to (3.3) exists. 1 r x r 2 r r r r r 1 2 x x xπ 

−1 −1

1 −1 −1

Since N

− Φ 0 = 0 , we have λ − = 0 N Φ z + z − N Φ 0 . Following the notations in Deville and r ( ) r ( r ) ( 1 2 ) ( r ) ( ) 

Särndal (1992), we can obtain:

≤ z + z K(1−β )

−1

≤ z K(1−β )

−1

+ z K(1−β )

−1

λ (3.4)

r 1 2 1 2

1

where K is defined in Section 1.1 of Appendix in Deville and Särndal (1992) and 0 << β .

2

−1 2 −1 2

= (

1 2 p r

) , inequality (3.4) implies that λr = O(1) + p( Since r

z = O(1) and z O n O n ) . The second term tends to 0

as nrincreases. However, the first term is a non-zero constant vector in general situations, and does not decrease as

nrincreases.

−1 T

λr = Trt

rx

− ˆt

rxπ

+Op(1) New Result

r ∑ k k k k

3. In general situations, ( ) , where T = d q x x .

r

T

Proof: We use Fk k (x λ) to denote the adjustment function for a general calibration estimator. For the GREG

estimator, the adjustment function takes the form 1+ q x

T

λ . The difference between the two adjustment functions is k k

expressed as:

T T T

θ (x λ) = ( F (x λ) − 1+ q x λ) (3.5) kk kk k k

From (3.1), (3.2) and (3.5), we obtain

T T (t − t ) + (t − ˆt ) = d x q x λ +θ (x λ ) (3.6)

x r r r ∑ k k { kk r k kr } 

r

x x xπ 

Multiplying both sides of (3.6) by Tr

−1

and rearranging the terms, we obtain

−1 −1 −1 T λ − T (t − t ) − T (t − ˆt ) = −T d x θ (x λ ) (3.7)

r r x r r r r r ∑ k kk k r

r

x x xπ 

An important assumption in Deville and Särndal (1992) is that Fk

''(0) is uniformly bounded, which is equivalent to

2 T T T −1/2 θ ( ) x λ = maxθ ( ) ( ) x λ = O ( x λ ) . Note that this assumption requires the condition that λ = O n( ) ,

k r k kr k r r p

T which does not necessarily hold when t ≠ t . But given that max | x λ |< ∞ , when nonresponse in the analytic x r k r x

'' T

survey is not extremely severe, we can still assume that for anyε > 0 , there exists K such that, for all k, | x λ |<ε k

T '' T will imply that θ (x λ) ≤ K (x λ)

2

.

k k k

11
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

For λrsufficiently small (which happens when nonresponse in the analytic survey is not extremely severe),

−1

3  ̧

ˆ

−1 −1 '' ⎧ −1 −1

λ −T t( − t ) ≤ (NA Tr) K N d x

2

 ̋ + T t( − t ) (3.8) r r r r ⎨ A ∑ k k

λr r x r ⎩ r

x xπ x

 ̨ 

−1 −1

(Nr Tr) = Op(1) and Nr

−1∑d We k

know that x

3

= Op(1) . Based on the New Result 2, 2

= Opλ ( ) 1 , so the r

r

k

first term of the right-hand side of (3.8) is Op

(1) . The second term of the right-hand side of (3.8) is also Op(1) . So

−1 −1

λ = T t − t

ˆ +O (1) T t − t we have ˆ

r r r r p ( ) . Although r r r ( ) tends to 0 as n → ∞ , the magnitude of λ is O (1) in

x xπ x xπ A r p

general situations. Unless t = t , λ does not tend to 0 as n → ∞ .

x rx

r A

New Result 4. The difference between the general calibration estimator and the Horvitz-Thompson estimator can

be expressed in two ways.

In terms of totals:

) ) T T T

ˆ tˆ = ( − t ) +B t − t

ˆ ) + ( − ˆ t − B t ( Y Y D) θ (3.9) r r rx r rr r r r r r yw yπ x x xπ 

where

−1

−1 T ⎛ T ˆ

Bˆ = T X DQY = ⎜ d q x x d q x y r r r r r r ∑ kkkk  ̃ ∑ k k k k

p×1 ⎝ r  ̄ r

T T T = X DQX = d q x x r r r r r ∑ k k k k

p p × r

⎛ x ... x ˆ

11 1p

⎜  ̃

X M O M r = ⎜  ̃ 

n p ⎜  ̃ r × 

x L x

⎝ n 1 n p r r  ̄ 

T T T T ( ,...x ,...,xn = x )

1 k r

⎛ ˆ q1

⎜ 0  ̃ 

Qr = ⎜ O  ̃ 

r × n nr ⎜  ̃ ⎜0 q  ̃

⎝ nr  ̄ 

⎛ d 0ˆ ⎜ 1

 ̃ 

Dr = ⎜ O  ̃ 

n n r r × ⎜  ̃ ⎜ 0 d  ̃

⎝ nr  ̄ 

( ,...y ,..., y )

T Y = y r 1 k nr

Yˆ = X Bˆ

r r r

( ,...θ ,...,θ 

r

)

T

θr = θ 1 k n

In terms of means:

)

−1 −1 T −1 T −1 T N tˆ − N tˆ = B

) 

μ 1− (μ μ ) p −μˆ (1− p) + N B t − t

ˆ + N (Y − Y D ˆ ) θ (3.10) r r r r x ( r x ) r r ( r r ) r r r r yw yπ x yπ x xπ 

where μ is the mean of the auxiliary vector for the respondent population, μ is the true population mean, μˆ is r x r x yπ 

the Horvitz-Thompson estimator of the mean for the outcome variable estimated from the respondent set, and p is

the response rate of the analytic survey. In the special situation where μr = μx

(indicating ignorable nonresponse),

x

the difference between the two estimators becomes:

) )

−1 −1 T −1 T −1 T N tˆ − N tˆ = (B μ −μˆ )(1− p) + N B t − ˆt + N ( − ˆ Y Y Dθ ) (3.11) r r r r x r r ( r r ) r r r r yw yπ yπ x xπ 

12
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Proof: If the calibration equation has a solution λ , then from (3.5) the difference between the general calibration

estimator and the Horvitz-Thompson estimator can be written as:

ˆ ˆ dy q T T

t − t = { x λ +θ (x λ )} (3.12)

r r ∑ k k k k r k kr

r

yw yπ 

From (3.7), 

−1 −1 −1 T λ = T t − t + T t − ˆt − T d x θ x λ (3.13)

r r ( x r ) r ( r r ) r ∑ k k k ( k r ) 

r

x x xπ 

Replacing λrin (3.12) by the right-hand side of (3.13), we obtain

tˆ

r − tˆr yw yπ 

T T

= ∑ dy q{ x λ +θ (x λ )} kk kk r k kr

r

⎧ ⎧  ̧  ̧ T −1 −1 −1 T T

= d y q x T ˆ ∑ (t − t ) + T (t − t ) − T d x θ (x λ ) +θ (x λ )

k k ⎨ k k ⎨ r x r r r r r ∑ kkk k r  ̋ k k r  ̋ 

r ⎩ ⎩ r  ̨  ̨

x x xπ 

T −1 T −1 T −1 T T

= dyq x T ˆ ∑ (t − t ) + d yq x T (t − t ) − d yq x T d θ (x λ x) + d y θ (x λ )

k kkk r x r ∑ k kkk r r r ∑ k kk k r ∑ k k kr k ∑ k kk k r

r r r r r

x x xπ 

T T T

T −1 T −1 T −1 T T

= ˆ (X DQY ) T (t − t ) + (X DQY ) T (t − t ) − (X DQY ) T X D θ + Y Dθ r r rr r x r r r rr r r r r r rr r r rr r rr x x xπ 

T T )

T −1 T −1 T T T

= ˆ (X DQY r r rr ) Tr x r (t − t ) + (X DQY r r rr ) Tr r (t − t

r ) − BXD r r rr θ + Y Dθ r rr x x xπ 

) ) T T T T

= B t r ( x r − t ) + B t r ( r − t

ˆ

r ) − Y D ˆ

r rr θ + Y Dθ r rr x x xπ 

= B t − t + B t − t

ˆ + (Y − Y ) D θ

) T

( ) ) T( ) ˆ

T

r x r rr r r r r r x x xπ 

(3.14)

Then the difference between two means is

−1 −1 N tˆ − N tˆ

r r r yw yπ 

−1 T T −1 −1

= N ∑ dy q{ x λ +θ (x λ )}+ (N − N t)ˆ

kk kk r k kr r ryπ 

)

r

)

−1 T −1 T ˆ−1 ˆT

= N B t( − t ) + N B t( − t ) + N (Y − Y ) D θ + (t

ˆ

r

)( N Nr N −1)

r x r r r r r r r r r x x xπ yπ 

) )

T −1 T ˆ

−1 ˆT

= B t / N − (t / N )(N N)) + (t

ˆ

r r

)( r ( x r r r N Nr N −1) + N B t r ( r r − t ) + N (Yr r − Y ) Drθr

) T

x yπ x xπ 

)

−1 T −1 T

= B μ − μ (N N)) + μˆ (N N −1) + N B t − ˆt + N (Y − Yˆ

r ( x r r r r r ( r r ) r r ) D θr r yπ 

) 

x xπ 

) T

x

−1 T −1 T

= B (μ − ( ˆ ( ˆ + N (Y − ˆ

r x μ μ x r μx

) p) − μr(1− p) + N B t rr r − t ) r r rr Y ) D θ

yπ x xπ 

)

−1 T −1 T ) T

x

= Br x μ ( 1− (μ μ r x) p) − μˆ

r r (1− p) + N B t( r r − t

ˆ ) + N (Yr r − Yˆ

) Dr r θ

x yπ x xπ 

(3.15)

For the right-hand side of (3.15), the first two terms do not cancel out except in some special situations such as

μ = μ (indicating ignorable nonresponse) and B

) T

μ = μˆ (meaning that the assisting linear regression model has r x r x r x yπ 

−1/2 perfect predicting power). The third term is O n( . Based on the New Result 3, we know θ = O (1) , so the p r ) r p

fourth term does not necessarily diminish as nrincreases. Instead, its magnitude seems to depend on the variation of 

the outcome variable, the predicting power of the regression model underlying the GREG estimator, and the form of

the weight adjustment function used in calibration. In general, the difference between general calibration estimator

and Horvitz-Thompson estimator does not decrease as nr

increases.

13
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

New Result 5. The difference between the general calibration estimator and the GREG estimator is

−1 −1 T N (tˆ − tˆ ) = N Y D θ = O (1). r r r r r p yw yreg

Proof: From (3.5) and (3.14), the general calibration estimator can be expressed as:

T T tˆ = d y + d q x λ y + d θ (x λ ) y r ∑ k k ∑ k k k r k ∑ k k k r k

r r r

yw

T ⎧ −1 −1 −1 T  ̧ T

= d y + d q x T t − t + T t − ˆ ∑ t − T d x θ x λ y + d θ x λ y k k ∑ k k k ⎨ r ( x r ) r ( r r ) r ∑ k k k ( k r ) ̋ k ∑ k k ( k r ) k

r r ⎩ r  ̨ r

x x x

T −1 T −1 T −1 T T

= d y + d q x T t − t y + d q x T t − ˆ ∑ t y − d q x T d x θ x λ y + d θ x λ y k k ∑ k k k r ( x r ) k ∑ k k k r ( r r ) k ∑ k k k r ∑ k k k ( k r ) k ∑ k k ( k r ) k

r r r r r r

x x x

) ) T T T T

= ˆ ˆ ˆt + B t − t + B t − t − Yr D θ + Y D θ r r ( x r ) r ( r r ) r r r r r yπ x x x

(3.16)

But the first four terms of the right-hand side of (3.16) is the GREG estimator

) ) T T ˆ

T

tˆ = tˆ + B t − t + B t − t

ˆ − Y D θ (3.17) r r r ( x r ) r ( r r ) r r r yreg yπ x x x

−1 −1 T So N (tˆ − tˆ ) = N Y D θ = O (1). r r r r r p yw yreg

The term θ captures the difference between the weight adjustment function for any general calibration estimator and

the weight adjustment function for the GREG estimator. When calibration is used for nonresponse adjustment,

θ = Op

(1) in general situations and does not tend to zero as the sample size nr

increases. As a result, the GREG

estimator and the general calibration estimator are not asymptotically equivalent. 

The results in this section are purely design-based and provide some initial insight on the difference between the

general calibration estimator and GREG estimator when calibration is used for nonresponse adjustment. To gauge

the magnitude of the difference, we need to go beyond the design-based approach and examine the underlying

models for the population structure and response mechanism. For example, a set of variables may be correlated with

the outcome variable of interest. Another set of variables may be correlated with response propensity. The question

is how to incorporate these covariates in the calibration process to reduce potential nonresponse bias without

increasing variance significantly.

4. Comparison of Three Commonly Used Calibration Estimators for Nonresponse Adjustment through

Simulation Study

In this section, we focus on three commonly used calibration estimators in the situation where the auxiliary

information is in the form of counts in a frequency table in two or more dimensions. We examine raking (as an

example of the general calibration estimators), poststratification (as a special form of the GREG estimator that

accounts for the interaction effects of the auxiliary variables), and the GREG estimator that accounts for only the

main effects of the auxiliary variables. In practice, the choice between these estimators is often based on the

distribution of respondents in the analytic survey and the availability of external data. This section develops a

systematic approach for evaluating the performance of these estimators through a simulation study. We compare the

unconditional and conditional empirical biases, empirical variances, and the coverage rates of 95 percent confidence 

intervals of these estimators. The findings demonstrate the importance of accounting for the outcome variable model

and response model when choosing the appropriate calibration estimator. A framework involving both design-based 

and model-based thinking is developed to simultaneously evaluate the impact of sampling, outcome variable

structure, and nonresponse mechanism.

Since survey practitioners often lack the knowledge of the outcome variables and nonresponse mechanism in a realworld survey, we also develop a diagnostic method that helps gauge the potential consequence of failure to

incorporate significant covariates in the calibration process. The results of this section provide survey practitioners 

with guidelines for choosing between these commonly used calibration estimators.

14
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4.1 Outcome Variable Model, Response Model, and Covariates in Weighting

The results in Section 3 indicate that the GREG estimator and general calibration estimators are not necessarily

asymptotically equivalent when calibration is used for nonresponse adjustment. To gauge how much a general

calibration estimator may diverge from a GREG estimator in a particular setting, it is necessary to go beyond the

design-based approach and examine the underlying models for the outcome variable and response mechanism. For

example, a set of covariates X1 may determine the outcome variable of interest, while a set of covariates X2 may

drive the response propensity. The relationship between X1and X2can fall into one of the three situations: i) X1

and X2are exactly the same; ii) X1and X2are different but have overlapping components; and iii) X1and X2

are completely different with no overlapping components. In practice, we often face the situations that some the

covariates corresponding to the response propensity are not correlated with the outcome (situation ii), so the

question is whether and how to incorporate these covariates in the calibration process to reduce nonresponse bias

without increasing variance significantly. 

Based on Little and Vartivarian (2006), “the most important feature of variables for inclusion in weighing

adjustments is that they are predictive of survey outcome; prediction of propensity to respond is secondary, although

useful.” As shown in Table 1 (which is reproduced from Table 1 in their paper), Little and Vartivarian (2006) assess

four scenarios based on the association of the auxiliary variables with response and outcome, and draw the following

conclusions:

L&V (i): Substantial bias reduction requires adjustment cell variables that are related both to nonresponse and to the 

outcome of interest.

L&V (ii): If the adjustment cell variables are unrelated to nonresponse, then weighting tends to have no impact on

bias, but reduces variance to the extent that the adjustment cell variables are good predictors of the outcome.

L&V (iii): If adjustment cell variables are good predictors of nonresponse but unrelated to the outcome variable, then

weighting increases variance without any reduction in bias.

L&V (iv): If the adjustment cell variables are related to neither outcome nor nonresponse, then weighting affects

neither bias nor variance.

Table 1. Summary of Little and Vartivarian (2006) Conclusions

Scenario Association with Outcome Association with Response Bias Variance

L&V (i) High High ↓ ↓

L&V (ii) High Low -- ↓

L&V (iii) Low High -- ↑

L&V (iv) Low Low -- --

SOURCE: Little and Vartivarian (2006), Table 1.

In the single-step weighting approach, calibration is applied to the basic sampling weights directly without a

separate nonresponse adjustment step, so Little and Vartivarian (2006) offer a useful framework for choosing

auxiliary variables and corresponding calibration estimator (e.g., poststratification versus raking). However, the

messages in Little and Vartivarian (2006) are not quite clear to the reader sometimes. For example, on the one hand, 

they assert that “[a] covariate for a weighting adjustment must have two characteristics to reduce nonresponse bias –

it needs to be related to the probability of response, and it needs to be related to the survey outcome.” On the other

hand, they state that “the most important feature of variables for inclusion in weighting adjustment is that they are

predictive of survey outcome; prediction of propensity to respond is a secondary, though useful, goal.” The former

statement seems to suggest that the outcome variable model and response model should play equally important roles

in determining the appropriate covariates for nonresponse adjustment, while the latter seems to indicate that the

outcome variable model should be the dominant factor. Moreover, the wordings in their text and their Table 1 are

not quite consistent. The text seems to address extreme conditions where the variables are either “related” or “not

related” to the outcome and/or response, while Table 1 shows “high” and “low” correlations, which are the middleground conditions that we are more likely to see in reality. Finally, Little and Vartivarian (2006) address only main

15
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

effects and do not provide any guidance about how to handle the interaction effects. Since the interaction terms of

the main effect variables are not completely new variables, the conclusions in Little and Vartivarian (2006) cannot

be applied directly to the comparison between poststratification, raking, and the GREG with only the main effect 

terms. In our research, we attempt to address the issues related to the interaction terms and refine the conclusions in

Little & Vartivarian (2006) through a simulation study.

4.2 Poststratification, Raking, and the GREG without Interaction Effects

In the simulation study, we focus on three widely used calibration estimators: (1) poststratification estimator as a

special case of the GREG estimator, where both main and interaction effects of the categorical auxiliary variables

are taken account of; (2) raking ratio estimator as an example of the general calibration estimators; and (3) the

GREG estimator when only the main effects of the auxiliary variables are accounted for. For simplicity, we refer to

the GREG estimator accounting for only the main effects as “GREG_Main”.

The poststratification estimator is generated under the group-mean assisting model. The auxiliary information

consists of known cell counts in a frequency table in any number of dimensions. For simplicity, we consider a twoway table with r rows and c columns, and thus r c × mutually exclusive cells. The auxiliary vector x is composed k

of rc −1 entries of 0 and a single entry of 1 indicating the cell to which k belongs. The population cell Uij contains

r c

Nij elements, i=1, ..., r; j=1, ..., c. So N = ∑∑Nij . Sometimes we do not have all the cell counts Nij , but only

i=1 j=1

marginal counts for the benchmark controls. One way to utilize the auxiliary information is to calibrate on known 

marginal counts, referred to as generalized raking (Deville and Särndal, 1993). Deming and Stephan (1943) suggests

an iterative proportional fitting procedure that adjusts one marginal at a time until convergence is achieved. An

alternative approach to take advantage of the marginal counts is to use the GREG estimator by accounting for only

the main effects of the auxiliary variables (i.e., GREG_Main). We conduct pairwise comparison of these three

estimators through simulation.

Raking versus GREG_Main: These two calibration estimators share the same set of auxiliary variables, and their

difference lies in the form of distance function G( )⋅ and corresponding adjustment function F( )⋅ . In the pure

sampling context as discussed in Deville and Särndal (1992), these two estimators are asymptotically equivalent. 

That is, conditioning on the same set of auxiliary variables, the particular form of the distance function has

negligible impact on the asymptotic property of the calibration estimator if non-sampling error does not exist.

However, the conclusion in Deville and Särndal (1992) does not always hold when nonresponse exists and

calibration is used to reduce nonresponse bias. The theoretical results in Section 3 suggest that the difference

between raking and GREG_Main could be as large as Op

(1) . The question is in what situation the two estimators

tend to give very similar results and in what situations they tend to diverge significantly.

GREG_Main versus poststratification: These two estimators both belong to the GREG estimator family, although

with poststratification there is a unique solution to the calibration equations regardless of the distance function

involved. GREG_Main accounts for only the main effects of the auxiliary variables while poststratification accounts

for the interaction effects as well. This comparison demonstrates the impact of outcome variable model and response

model in the choice of calibration estimator. We relate our simulation results to Little and Vartivarian (2006) and

provide more refined guidelines for choosing auxiliary variables in nonresponse adjustment weighting.

Poststratification versus raking: These are probably the two most commonly used calibration estimators in US

government surveys. From the practical perspective, the key difference between poststratification and raking seems

to be that the former fits a fully saturated model with both main and interaction effects of the auxiliary variables,

while the latter fits a model including only the main effects. On the other hand, Deville, Särndal, and Sautory (1993)

refer to poststratification as complete poststratification and raking ratio as incomplete poststratification. Does this

imply that the raking estimator accounts for the interaction effects to some extent? We attempt to investigate to what

extent raking can get closer to poststratification compared to what GREG_Main does.

16
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4.3 Scope and Conceptual Framework for Simulation Study

The simulation study aims to evaluate the empirical properties of the poststratification estimator, the raking ratio

estimator, and the GREG_Main estimator for finite population totals and means when calibration is used for

nonresponse adjustment in a one-step weighting approach. We measure the magnitude of their differences in terms

of empirical bias, variance, mean square error (MSE), and coverage rate of the 95 percent confidence interval, under

different model assumptions for the outcome variable and nonresponse mechanism. The research is conducted in the

following scope.

First, we evaluate estimates for population totals and means for a single outcome variable. In the presence of

nonresponse, calibration is used to reduce the bias, variance, and mean squared error (MSE) of the estimate for this

single outcome variable.

Second, although Section 2 points out that it is possible to use a covariate vector zk

for the calibration adjustment

function F( )⋅ that is different from the auxiliary vector xkin the calibration equation, our evaluation focuses on the

situation where zk = xk, which is the “standard choice” in Särndal and Lundström (2005).

Third, the outcome variable model and response model contain the same main effect covariates. We also assume that

there are only two main effect covariates and they are both categorical variables. Sometimes a substantively and 

statistically significant interaction term between the two main effect variables is also included in either or both

models. Depending on the particular calibration estimator, the interaction term may or may not be included in the 

calibration equation.

Fourth, for the response mechanism, we assume missing at random (MAR). This means that the probability of

response does not depend on the outcome variable once we control for the known covariates. The classes or cells

defined by the covariates are response homogeneity groups.

Finally, the results focus on overall estimates in the context of simple random sampling (SRS). Although practical

surveys almost always involve complex sample designs, the SRS assumption allows us to focus on the impact of

population structure and response mechanism on the performance of a calibration estimator. The findings about

how to choose auxiliary variables and calibration estimators apply in general to complex designs, although the

technical details become more complicated.

In the simulation study, the models for the outcome variable and response propensity are varied by either including

or excluding the interaction term of the main effect covariates. We then compare the properties of the

poststratification estimator, raking ratio estimator, and GREG_Main estimator in one-step calibration weighting. The

alternative models for the outcome variable Y and response propensity R are specified as below. Depending on

whether the interaction term is included, we refer to the models as “Y_Main” and “R_Main” (meaning only main

effects are in the model) versus “Y_Interaction” and “R_Interaction” (meaning that the interaction effects are

included in the model as well).

Y_Main: Yijk =μY +α Yi +β 

Yj+ε Yijk , i = 1, 2; j = 1, 2; k = 1, ... , Nij (4.1)

Y_Interaction: Yijk =μ Y +α Yi +β 

Yj +γ Yij +ε Yijk , i = 1, 2; j = 1, 2; k = 1, ... , Nij (4.2)

R_Main: Rijk =μ R +α Ri +β Rj +ε Rijk , i = 1, 2; j = 1, 2; k = 1, ... , Nij (4.3)

R_Interaction: Rijk =μ R +α Ri +β Rj +γ Rij +ε Rijk , i = 1, 2; j = 1, 2; k = 1, ... , Nij (4.4)

~ (0, 2, and N σ 

2 where N )

ij is the population size in cell ijfor the survey, ε 

Yijk N σ Y

) ε Rijk ~ (0, R.

We are aware that in theory, it would be more appropriate to use logistic models for the response propensity, but

choose to use linear models in order to manipulate the simulation parameters easily.

17
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Table 2 summarizes the four simulation scenarios and the corresponding models governing the outcome variable and

response propensity. EM means expectation in terms of the underlying superpopulation model that generates the

finite population, and ER means expectation in terms of the response model. A hypothetical example can be that we

are interested in a single outcome variable income. Both the outcome variable y and the response indicator r can

be explained by two dichotomous variables, education (high and low) and age (young versus old), and possibly an 

interaction effect term. In our models, αmeasures the main effect of education, β measures the main effect of age,

and γ measures the possible interaction effect between education and age. In poststratification, the weights are

adjusted by four cells defined by education and age. In raking, the weighting adjustment is conducted iteratively by

using age and education as marginal controls until convergence is achieved. In GREG_Main, the calibration

estimator is a function of the regression coefficient as the result of modeling the outcome variable income by only

the main effects of education and age.

Table 2. Simulation Scenarios

Scenario Outcome Variable Model Response Propensity Model

Y_Main & R_Main: neither outcome model nor

response model includes interaction effect term E Y M

( ) ijk = μ 

Y +α Yi + β Yj E r R

( ) ijk = μ R +α Ri + β Rj

Y_Main & R_Interaction: only response model

includes interaction term

E Y M( ) ijk = μ 

Y +α Yi + β Yj E r( ) = μ +α + β +γ R ijk R Ri Rj Rij

Y_Interaction & R_Main: only outcome model

includes interaction effect term E Y( ) = μ +α + β +γ M ijk Y Yi Yj Yij E r R

( ) ijk = μ R +α Ri + β Rj

Y_Interaction & R_Interaction: both outcome

model and response model include interaction

effect term

E Y( ) = μ +α + β +γ M ijk Y Yi Yj Yij E r( ) = μ +α + β +γ R ijk R Ri Rj Rij

Through the simulation study, we examine the performance of poststratification, raking, and GREG_Main under 

different outcome variable model and response model combinations. We evaluate the consistency between our

results and those in Little & Vartivarian (2006) and refine their conclusions. At the same time, we attempt to expand

Deville and Särndal (1992) and shed light on the empirical difference between the GREG estimators (i.e.,

GREG_Main and poststratification) and the raking estimator (as an example of the general calibration estimator) in

the presence of nonresponse.

4.4 Simulation Steps and Parameters

The implementation of the simulation is based on the framework shown in Table 2 and involves the following steps:

1. Generate an artificial finite population of size N that contains 4 subpopulations defined by the categories of the

two auxiliary variables (2x2). The subpopulation size, Nij , is generated using Poisson distribution and

approximately equal across the four subpopulations.

2. Generate two sets of values for the outcome variable, one set corresponding to the outcome model Y_Main

specified in (4.1) and the other set corresponding to the outcome model Y_Interaction specified in (4.2).

3. Select a simple random sample of size n from the finite population. 

4. From the SRS sample, draw subsamples of respondents using the response models R_Main and R_Interaction

as specified in (4.3) and (4.4). This is achieved through Poisson sampling and the final sample size (of

respondents) is nr

.

5. For each of the four respondent samples corresponding to the different outcome variable and response model

combinations shown in Table 2, conduct calibration using poststratification, raking, and GREG_Main, 

18
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

respectively. Obtain the estimates for the outcome variable associated with the three calibration estimators, and

then compare the empirical results using the evaluation criteria specified in Section 4.5. 

Several factors may affect the properties of and differences between these three calibration estimators, including: (1)

the number of simulation samples; (2) the substantive and statistical significance of the interaction effect in the

outcome variable model; (3) the substantive and statistical significance of the interaction effect in the response

model; and (4) the overall sample size for the respondent sample and the distribution across the four subpopulations.

During the initial investigation, we choose the simulation parameters in a way to minimize the factors that could

cloud our comparison.

First, the number of simulation samples is very large ( s = 1, 2, ..., 10000 ). That is, for each outcome variable model

and response model combination, we draw 10,000 SRS samples from the finite population and then 10,000

corresponding respondent subsamples. The large number of repetitions ensures that any observed differences in the

three calibration estimators are not due to the random errors in the simulation samples. This is particularly important

when we partition the samples into groups and evaluate the conditional properties of the calibration estimators.

Second, for all the outcome variable and response propensity models (regardless of whether the interaction effect is

included or not), the random error terms are set to be very small, so the explanatory power of the overall model is

strong. At the same time, the interaction terms in the Y_Interaction model and R_Interaction model are set to be

substantively and statistically significant. Such setup ensures that any significant impact of the outcome variable

model and response model on the properties of the calibration estimators can be detected. Under these criteria,

several sets of outcome model parameters and response model parameters are used. The results associated with these

different parameters tend to lead to the same conclusions, so we choose to present the results based on only one set

of parameters, as shown below.

The parameters corresponding to the outcome variable models (4.1) and (4.2) are

μY = 1000

α = ( , α α ) = (−200, 300) Y Y 1 Y 2

βY = (β Y 1

, β Y 2) = (− 100, 150)

γY = (γ 

Y11, γ Y12,γ Y 21,γ Y 22 ) = (100, 300, 700, 1200)

ε 

Yijk ~ (0,900) N

The parameters corresponding to the response propensity models (4.3) and (4.4) are

μ R = 0.05

α = (α , α ) = (0.05, 0.1) R R1 R2

β = (β , β ) = (0.2, 0.4) R R1 R2

γR = (γ R11, γ R12,γ R21,γ R22 ) = (0.05, 0.15, 0.2, 0.4)

ε Rij ~ (0,0.0025) k N

Tables 3 and 4 show the expected cell means corresponding to the models with and without interaction terms, for the

outcome variable and the response propensity respectively.

Table 3. Expected Cell Means for Outcome Variable

Outcome Variable Model E y( ) M 11k E y( ) M 12k E y( ) M 21k E y( ) M 22k

Y_Main 700 950 1,200 1,450

Y_Interaction 800 1,250 1,900 2,650

Table 4. Expected Cell Means for Response Propensity

Response Propensity Model E r( ) R 11k E r( ) R 12k E r( ) R 21k E r( ) R 22k

R_Main 0.30 0.50 0.35 0.55

R_Interaction 0.35 0.65 0.55 0.95

19
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Finally, the respondent sample size is determined by the SRS sample size n as well as response propensity. Under a

particular set of parameters for the response propensity model, we vary the respondent sample size by changing the

SRS sample size n . Comparing the results corresponding to large ( n = 8,000 ), medium ( n = 2,000 ), and small

( n = 200 ) respondent sample sizes allows us to evaluate asymptotic properties of these three calibration estimators.

The simulation is conducted in R (Lumley, 2005; R Development Core Team, 2005) because of its efficiency in

handling matrix calculations and extensive capacity for analyzing survey data.

4.5 Evaluation Criteria

We first examine the empirical properties of the three calibration estimators using repeated sampling approach (i.e.,

averaging across the 10,000 simulation samples). Then we compare the results conditioning on the types of samples

defined by a proposed distance measure. In real-world survey practice, only one sample can be fielded and all the

estimates are based on that particular sample. The conditional properties of these calibration estimators help shed

light on the importance of choosing the appropriate estimator based on the particular sample a survey practitioner

may obtain.

The empirical results for the three calibration estimators are compared under the four outcome variable and response

model combinations using several measures across the simulation samples. The measures are described below in

terms of totals. We also evaluate the properties of the means using a similar set of measures.

1. Relative bias RelBias(t

ˆ

yw ) = (1/ S)∑

S

(t

ˆ

yw − ty

) t

y s s=1 s

where sindicates a particular sample, S is the total number of samples included. t

y

is the true population

total, and t

ˆ

yw is the estimate from sample s using one of the three calibration estimators.

s

S

2. Relative standard error RelSE(tˆ

yw

) = var(tˆ

yw

) t

y = (1/ S)∑ (tˆyw − Ep

(tˆ

yw

))

2

t

y s s s=1 s s

where E t(

ˆ ) = (1/ S)

S

t

ˆ , the expected value of tˆ over repeated sampling.

p yws ∑ s=1ywsyws

3. “Relative” square root of MSE

S 2 S

RelRMSE(t

ˆ ) = MSE(t

ˆ ) t = (1/ S)∑ (tˆ − E (tˆ )) + ((1/ S)∑ (tˆ − t ))

2

t

yw yw y yw p yw yw y y s s s=1 s s s=1 s

4. Coverage rate of the 95 percent confidence interval

(1/ S)∑s

S

=1

I ( zˆ ≤ z ), whereα = 0.05 and zˆ = (tˆ − t ) var(tˆ )

j 1−α 2 j yws s y yw

where α = 0.05 and zˆ

j = (t

ˆ

yw − ty

) var(t

ˆ

yw

) .

s s

5. Bias ratio: calculated as the ratio of Bias(t

ˆ ) and square root of var(tˆ ) :

ywsyws

S S

BiasRatio(tˆ

yw

) = (1/ S)∑ (tˆ

yw − ty

) (1/ S)∑ (tˆ

yw − Ep

(tˆ

yw

))

2

s s=1 s s=1 s s

4.6 Simulation Results

We evaluate the performance of poststratification, GREG_Main and raking first over repeated sampling, and then

conditioning on sample types defined by a distance measure. The results for the total and those for the mean

demonstrate almost exactly the same pattern, so we focus on discussing the properties of the estimators for the total.

20
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4.6.1 Properties over Repeated Sampling

The empirical bias and variance used for calculating the measures in Table 5 were estimated by averaging across the

10,000 simulation samples. For each of the two outcome variable models (Y_Main and Y_Interaction), we first

obtain the results corresponding to 100 percent response as the baseline, and then evaluate the properties of the three

calibration estimators under the R_Main and R_Interaction models.

Impact of Outcome Variable Model and Response Propensity Model

Table 5 shows that the relative biases associated with all the three calibration estimators are very small, indicating

that calibration helps reduce nonresponse bias significantly so long as the main effects of the key auxiliary variables

are accounted for. At the same time, the outcome variable model seems to be the driving factor that determines the

performance of the calibration estimators, not the response propensity model. When the outcome variable model

does not contain the interaction effect of the auxiliary variables (in the Y_Main scenarios), poststratification (which 

accounts for the interaction effect) does not seem to reduce nonresponse bias or variance further than raking or

GREG_Main, and this is true regardless of the response model (R_Main or R_Interaction). For example, in the

“ n = 8,000 , Y_Main & R_Interaction” scenario, the three calibration estimators have almost the same relative bias

(approximately 1.8E-5), relative standard error (approximately 4.0E-4), relative square root of MSE (approximately

0.00040), and coverage rate of 95% confidence interval (93 percent).

On the other hand, when the outcome variable model contains the interaction effect term, poststratification performs

better in terms of both relative bias and relative standard error regardless of the response propensity model. In the

“ n = 8,000 , Y_Interaction & R_Interaction” scenario, the relative biases for raking and GREG_Main are 65 times

and 300 times as large as that for poststratification, and the relative standard errors for raking and GREG_Main are

approximately 3 times as large as that for poststratification. 

We think the key to understanding this pattern of results is the following:If an auxiliary variable is correlated only to

nonresponse but uncorrelated to the outcome variable, the auxiliary variable does not cause any nonresponse bias. In

the Y_main scenarios, although the interaction term affects nonresponse, such nonresponse does not introduce any 

nonresponse bias in addition to the nonresponse bias that is already caused by the main effects (because the

interaction effect itself is not correlated with the outcome variable). Including the interaction term in calibration does

not help because no bias is caused by the interaction term in the first place. This is why raking and GREG_Main

performed almost as well as poststratification in the “Y_main” scenarios. The interaction effect in the response

model does seem to play a role in the performance of the three calibration estimators, conditioning on the fact that

the interaction effect is corrected to the outcome variable. Although GREG_Main performs poorly in the

“Y_Interaction & R_Main” scenarios, its worst performance is observed in the “Y_Interaction & R_Interaction”

scenarios when the interaction term is related to both the outcome variable and the response propensity.

Coverage Rate of 95% Confidence Interval and Bias Ratio

Despite the differences among the three calibration estimators in the Y_Interaction scenarios, the relative biases are

very small even for the calibration estimator that fails to account for the interaction effect that is correlated to the

outcome variable. However, the small relative bias can be misleading because the coverage rate of 95% confidence

interval can be poor for GREG_Main and raking even with very small relative bias. A very interesting pattern in

Table 5 is that as the SRS sample size increases, the confidence interval coverage rates actually become worse for

raking and GREG_Main. For example, in the “Y_Interaction & R_Interaction” scenario, the coverage rate of 95%

confidence interval for GREG_Main is only 55% with the sample size of 2,000, and it drops further to 3% when the

sample size increases to 8,000!

How can it be possible that the coverage rate of 95% confidence interval becomes almost unacceptable even when

the relative bias is very low? Why does the increase in sample size hurt the confidence interval coverage rate? The

answer lies in the asymptotic property of the bias ratio. We can re-write the t-statistic into the summation of two

terms:

21
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Table 5. Properties of Poststratification, GREG_Main, and Raking over Repeated Sampling

Relative Bias

( RelBias ˆt ) in 10-5

s

yw

Relative Standard Error

( RelSE ˆt ) in 10-4

s

yw

Relative Square Root of

MSE

( RelRMSE ˆt ) in 10-4

s

yw

Coverage of 95%

Confidence Interval

No calibration

PoststratificationGREG_MainRakingNo calibration

PoststratificationGREG_MainRakingPoststratificationGREG_MainRakingPoststratificationGREG_MainRaking

SRS sample n = 8,000

Y_Main &100% response -0.9 -0.1 -0.1 -0.1 26.1 2.8 2.8 2.8 2.2 2.2 2.2 95% 95% 95%

Y_Main & R_Main 4100.0 1.1 1.2 1.2 42.0 4.8 4.8 4.8 3.8 3.8 3.8 93% 93% 93%

Y_Main & R_Interaction 7390.0 1.8 1.9 1.8 33.6 4.0 4.0 4.0 3.2 3.2 3.2 93% 93% 93%

Y_Interaction & 100% response 0.8 -0.1 -0.3 -0.3 39.4 1.8 4.9 4.9 1.5 3.9 3.9 95% 95% 95%

Y_Interaction & R_Main 8660.0 1.1 -135.0 -67.9 65.6 3.1 9.0 8.3 2.5 14.0 8.7 94% 63% 88%

Y_Interaction & R_Interaction 13900.0 1.0 -304.0 -64.5 52.8 2.6 8.2 7.0 2.1 30.4 7.8 94% 3% 90%

SRS sample n = 2,000

Y_Main &100% response -0.8 -0.7 -0.7 -0.7 56.5 6.1 6.1 6.1 4.9 4.9 4.9 95% 95% 95%

Y_Main & R_Main 4100.0 0.4 0.5 0.4 88.3 9.8 9.8 9.8 7.8 7.8 7.8 95% 95% 95%

Y_Main & R_Interaction 7370.0 2.1 2.3 2.2 70.0 8.5 8.4 8.5 6.7 6.7 6.7 94% 94% 94%

Y_Interaction & 100% response 1.5 -0.2 -0.6 -0.6 87.0 3.9 10.8 10.8 3.2 8.6 8.6 95% 95% 95%

Y_Interaction & R_Main 8670.0 2.2 -136.0 -68.0 136.0 6.4 18.6 17.3 5.1 18.6 14.8 95% 88% 95%

Y_Interaction & R_Interaction 13900.0 1.4 -304.0 -63.9 110.0 5.4 17.0 14.7 4.3 30.9 12.8 95% 55% 96%

SRS sample n = 200

Y_Main &100% response -30.2 2.8 2.7 2.7 185.0 20.0 20.0 20.0 15.9 15.9 15.9 94% 94% 94%

Y_Main & R_Main 4130.0 -2.1 -1.9 -2.1 280.0 31.9 31.7 31.7 25.5 25.3 25.3 94% 94% 94%

Y_Main & R_Interaction 7400.0 3.1 3.7 3.2 230.0 27.2 26.9 27.0 21.6 21.4 21.4 94% 94% 94%

Y_Interaction & 100% response 3.6 0.0 1.9 2.0 284.0 12.9 35.2 34.8 10.3 28.1 27.8 95% 95% 95%

Y_Interaction & R_Main 8710.0 3.4 -149.0 -77.0 433.0 21.0 62.2 56.4 16.8 50.7 45.2 94% 94% 96%

Y_Interaction & R_Interaction 13900.0 -3.5 -314.0 -72.7 356.0 17.8 56.3 48.4 14.2 50.8 38.7 94% 92% 98%

22
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

tˆ − t tˆ − E E (tˆ ) EE (tˆ ) − t

yw y yw p M yw p M yw y s s s s

t statistic - = = + (4.5)

var(tˆ ) var(tˆ ) var(tˆ )

yw yw s s yws

The first term on the right-hand side of (4.5) is asymptotically N(0, 1). The second term is the standardized bias or

bias ratio. As the sample size increases, the denominator of the second term decreases. However, if the calibration 

estimator is model-biased as in the situation of GREG_Main and raking, the numerator in the second term of does

not decrease with increase sample size, but rather stays constant. As a result, the increase in sample size leads to an

increase in the bias ratio and therefore makes the coverage rate of 95% confidence interval worse. An important

lesson here is that increasing sample sizes does not help improve the performance of a calibration estimator that is 

model-biased.

Why Does Raking Perform Better Than GREG_Main?

Another interesting pattern in Table 5 is that in the “Y_Interaction” scenarios, the raking estimator has much smaller

relative bias and better confidence interval coverage than GREG_Main. Neither approach takes account of

interaction effects during the calibration process, so the question is why raking performs better than GREG_Main.

An important feature of raking is that the algorithm forces the weights to conform to the marginal totals without

perturbing the associations in the unadjusted table (Haberman, 1979), so the method preserves the interaction effect

that already exists in the data before calibration. That is, raking retains the cross-product ratios or odds ratios of the

cell totals in the observed data (Brick, Montaquila, and Roth 2003). In our situation, we conduct raking on the

respondent sample using iterative proportional fitting, and this procedure does not change the existing association

between the two main effect variables. As shown in Table 6, the odds ratios in the raking column are the same as

those in the respondent sample column. In contrast, GREG_Main fits a regression model using only the main effect

variables, and it destroys the existing correlation between the two main effect variables in the respondent sample.

Raking and GREG_Main both account for only main effects, but are associated with different distance functions.

This comparison shows that the form of distance function matters.

Table 6. Impact of Poststratification, GREG_Main, and Raking on Odds Ratios of Cell Totals

Odds Ratio of Cell Totals

SRS sample size = 8,000 Population Respondent

sample Poststratification GREG_Main Raking

R_Main scenarios,

regardless of Y model 0.99 0.93 0.99 0.88 0.93

R_Interaction scenarios,

regardless of Y model 0.99 0.93 0.99 0.76 0.93

6.2 Properties Conditioning on Sample

The relative biases of raking and GREG_Main seem acceptable over repeated sampling. However, a survey

practitioner can obtain one and only one sample in the real world, so it is important to understand how a calibration

estimator may perform for a particular sample. We define a measure that helps survey practitioners gauge the

potential consequence for choose the inappropriate calibration estimator for a particular sample. The discussions

below focus on the comparison between raking and poststratification. A similar comparison can be conducted

between GREG_Main and poststratification.

Assume that the model for the outcome variable Y contains an interaction effect term, as specified in equation (4.2).

In the SRS setting, the calibration estimator for the total can be expressed as:

2 2 nij

t

ˆ

yw = ∑∑wij∑ yijk (4.6)

i=1 j=1 k=1

where wij is the calibrated weight for a unit in cell ij .

23
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

The correct calibration method should be poststratification. If raking is used, then the model expectation of the

estimator is:

E (t

ˆ ) M raking

2 2 nij

= ∑∑ w E ( y )

ij∑ M ijk

i=1 j=1 k=1

2 2 nij

= ∑∑wij∑(μY +αYi + β Yj +γ Yij)

i=1 j=1 k=1

2 2 2 2 2 2 2 2

μ w n +∑α w n w n +∑∑ γ = Y ∑∑ ij ij Yi∑ ij ij +∑β Yj∑ ij ij Yijw nij ij

i=1 j=1 i=1 j=1 j=1 i=1 i=1 j=1

2 2 2 2

= μY Nˆ +∑αYiNˆ

i⋅ +∑β YjNˆ⋅ j +∑∑γ YijNˆij

i=1 j=1 i=1 j=1

(4.7)

The model bias of the raking estimator t

ˆ

raking is:

EM(tˆ

raking − t)

2 2 2 2 (4.8)

= μ ˆ ˆ ˆ ˆ

Y

(N N) +∑αYi(Ni⋅ − Ni⋅) +∑β Yj(N⋅ j − N⋅ j) +∑∑ γ Yij(Nij − Nij − )

i=1 j=1 i=1 j=1

The raking process forces the estimated row totals and column totals to be equal to the control totals, so the first

three terms of the right-hand side of equation (4.8) should be approximately zero, leaving the fourth to be the driving

term.

In practice, we do not know the values for μY

, α Yi , β Yj , and γ Yij . However, if we have the cross-classification and

corresponding cell totals for the population, we can compute the difference between Nˆ and N for each of the four ij ij

cells defined by the two auxiliary variables. Then we can define a distance measure as the square root of the

summation of the differences across all the cells:

2 2

2 Draking = ∑∑(Nˆ

ij − Nij) (4.9)

i=1 j=1

This distance measure is computable for each particular sample and can help predict the model bias of the raking

estimator for a particular sample when poststratification is supposed to be the appropriate estimator. We use the

“ n = 8,000 , Y_Interaction & R_Interaction” scenario to demonstrate how this distance measure can be used – First,

we compute the distance measure for each of the 10,000 simulation samples. Then we partition the 10,000 samples

into 20 groups based on the magnitude of the distance measure, with approximately 500 samples in each group.

Finally, we calculate the average relative biases and average coverage rate of 95% percent confidence interval for

each of the 20 groups. The results are presented in Table 7.

Table 7 shows that different samples defined by the magnitude of the distance measure behave differently. As the

distance measure increases, the relative bias increases and the coverage rate of 95% confidence interval decreases.

The coverage rate of 95% confidence interval becomes drops to below 88 percent for a sample with a distance

measure above the 80th percentile. If a survey practitioner happens to obtain a sample with distance measure above

the 95th percentile, then coverage rate of 95% confidence interval is only 10 percent. So the impact of choosing the

wrong estimator can be detrimental for some “unlucky” samples.

In practice, only one sample is fielded for a survey, so we can calculate a single value for the distance measure

instead of multiple values associated with many possible samples. Then the question is how to interpret and make

24
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

use of this single value. The next step in our research is to refine the form of the distance measure and derive the

statistical distribution for it, so survey practitioners can determine the statistical significance of the value associated

with a particular sample based on the distribution of the distance measure.

Table 7. Relative Bias and Confidence Interval Coverage Rate for Raking Estimator Conditioning on Sample, in the

“ n = 8,000 , Y_Interaction & R_Interaction” Scenario

Range of Draking

Relative Bias in 10-5 Coverage Rate of 95%

Confidence Interval

Poststratification

GREG_

Main Raking

Poststratification

GREG_

Main Raking

0

th

– 5

th percentile: (0.04,29] 1.5 -228.2 1.0 93% 1% 100%

5

th

– 10th percentile: (29,59] 0.6 -231.7 -1.4 94% 1% 100%

10th– 5

th percentile: (59,88] 2.1 -231.6 -1.3 92% 2% 100%

15th– 20th percentile: (88,117] 1.9 -237.7 -6.4 95% 1% 100%

20th– 25th percentile: (117,148] 0.2 -245.9 -14.5 94% 1% 100%

25th– 30th percentile: (148,180] 2.6 -247.7 -16.6 92% 3% 100%

30th– 35th percentile: (180,210] 3.5 -259.6 -23.9 93% 3% 100%

35th– 40th percentile: (210,240] 0.0 -269.9 -35.3 93% 4% 100%

40th– 45th percentile: (240,270] -0.4 -272.1 -36.8 94% 5% 100%

45th– 50th percentile: (270,301] 1.1 -287.1 -49.7 93% 6% 100%

50th– 55th percentile: (301,332] 1.0 -296.2 -58.3 96% 4% 100%

55th– 60th percentile: (332,366] 0.1 -307.2 -66.6 93% 6% 99%

60th– 65th percentile: (366,401] 1.1 -318.2 -77.4 93% 3% 99%

65th– 70th percentile: (401,439] 0.6 -325.8 -84.7 93% 5% 98%

70th– 75th percentile: (439,480] 0.1 -340.7 -98.2 94% 3% 97%

75th– 80th percentile: (480,528] -0.3 -354.2 -109.9 95% 2% 93%

80th– 85th percentile: (528,582] 1.1 -367.5 -119.8 95% 2% 88%

85th– 90th percentile: (582,650] 1.4 -385.3 -136.0 95% 1% 70%

90th– 95th percentile: (650,750] 0.7 -410.4 -157.4 92% 0% 42%

95th– 100th percentile:(750,1450] 0.2 -457.9 -197.7 94% 0% 10%

5. Summary

Our research proves that in general, the GREG estimator and general calibration estimators are not asymptotically

equivalent in the presence of nonresponse. Through a simulation study, we demonstrate the importance of

accounting for the outcome variable model and response model when choosing the appropriate estimator. 

Furthermore, we point out that the outcome model should be the dominant factor in determining what covariates

should be included in the calibration process.

One of the interesting findings is that small relative bias associated with the inappropriate calibration estimator

could result in very poor confidence interval coverage. Increasing sample size only makes the situation worse

because the bias tends to remain constant while the variance decreases with increasing sample sizes.

In the real-world survey practice where only a single sample is obtained, a distance measure could help gauge the

potential consequence of choosing inappropriate estimator. We plan to conduct further research on this and provide

more details in our future work.

References

Brick J.M, Montaquila, J., and Roth, S (2003). Identifying problems with raking estimators. Proceedings of the

Section on Survey Research Methods, American Statistical Association, 710-717.

25
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Chang, T., and Kott, P. (2008). Using Calibration Weighting to Adjust for Nonresponse under a Plausible Model.

Biometrika, 95, 555-571.

Deming, W.E. (1943). Statistical Adjustment of Data. New York: John Wiley & Sons, Inc.

Dever, J. (2008). Sampling Weight Calibration with Estimated Control Totals. (Dissertation)

Dever, J., and Valliant, R. (2010). A Comparison of Variance Estimators for Poststratification to Estimated Control

Totals. Survey Methodology, 36, 45-56.

Deville, J-C., and Särndal, C-E. (1992). Calibration Estimators in Survey Sampling. Journal of the American

Statistical Association, 87(418), 376-382.

Deville, J-C., Särndal, C-E., and Sautory, O. (1993). Generalized Raking Procedures in Survey Sampling. Journal

of the American Statistical Association, 88(423), 1013-1020.

Estevao, V. M., and Särndal, C-E. (2000). A Functional Form Approach to Calibration. Journal of Official Statistics, 

16, 379-399.

Estevao, V.M., and Särndal, C.E. (2006). Survey Estimates by Calibration on Complex Auxiliary Information.

International Statistical Review, 74, 127-147.

Fuller, W. A. (2000). Two-phase Sampling. Proceedings of the Section on Survey Research Methods, American

Statistical Association, 23-30.

Haberman, S. (1979). Analysis of Qualitative Data. Academic Press, Inc. New York.

Kalton G., and Flores-Cervantes, I. (2003). Weighting Methods. Journal of Official Statistics, 19, 81-97.

Kott, P.S. (1994). A note on Handling Nonresponse in Surveys. Journal of the American Statistical Association, 89,

693-696.

Kott, P. S. (2006). Using Calibration Weighting to Adjust for Nonresponse and Coverage Errors. Survey

Methodology, 32, 133-142.

Kott, P. S., and Chang, T. (2010). Using Calibration Weighting to Adjust for Nonignorable Unit Nonresponse. 

Journal of the American Statistical Association, 105(491), 1265-1275.

Lumley, T. (2010). Complex Surveys: A Guide to Analysis Using R. John Wiley & Sons: New Jersey.

Lundström, S., and Särndal, C-E. (1999). Calibration as a Standard Method for Treatment of Nonresponse. Journal

of Official Statistics, 15, 305-327.

R Development Core Team (2005). R: A Language and Environment for Statistical Computing. R Development

Core Team, Vienna, Austria. URL http://www.R-project.org

Särndal, C-E. and Lundström, S. (2005). Estimation in Surveys with Nonresponse. John Wiley & Sons: Chichester.

Särndal, C-E., Swensson, B., and Wretman, J. (1992). Model Assisted Survey Sampling. New York: Springer.

Särndal, C-E. (2007). The Calibration Approach in Survey Theory and Practice. Survey Methodology, 33, 99-119.

Singh, A.C., and Mohl, C.A. (1996). Understanding Calibration Estimators in Survey Sampling. Survey

Methodology, 22, 107-115.

26
 

 

 

Stukel, D.M., Hidiroglou, M.A. and Särndal, C-E. (1996). Variance Estimation for Calibration Estimators: A

comparison of Jackknifing versus Taylor Linearization. Survey Methodology, 22, 117-125.

Thibaudeau, Y. and Slud, E. (2009). Simultaneous Calibration and Nonresponse Adjustment. Proceedings of the

Survey Research Methods Section of the American Statistical Association, 2263-2272.

27