1

TESTING BIAS OF MEAN ESTIMATES

DUE TO NOT MISSING AT RANDOM NONRESPONSE

WITH APPLICATION TO THE CONSUMER EXPENDUTURE SURVEY

Michael Sverchkov

Bureau of Labor Statistics, Washington DC, USA

1. INTRODUCTION

There exists almost no survey without nonresponse. In practice most survey methods that 

adjust for nonresponse assume either explicitly or implicitly that the missing data are 

‘missing at random’ (MAR). That is, they assume response does not depend on the variable 

of interest (the outcome variable) given some auxiliary information known for the whole 

population. However, in many practical situations, this assumption is not valid, since the 

probability of responding often depends on the outcome value, even after conditioning on 

available covariate information. In such cases, the use of methods that assume nonresponse 

is MAR can lead to large bias of parameter estimators and distort subsequent inference.

The case where the missing data are not MAR (NMAR) can be treated by postulating a 

parametric model for the distribution of the outcomes before nonresponse and a model for 

the response mechanism. These two models define a parametric model for the observed 

outcomes, so that the parameters of these models can be estimated from the observed data. 

Once the parameters are estimated, the first model can be used for inference. See, for 

example, Pfeffermann and Sverchkov (2009) for details, with overview of related

literature.

Modeling the distribution of the outcomes before nonresponse is difficult since only the 

observed data are available. Sverchkov (2008) proposes an alternative approach, which

allows the parameters of the response model to be estimated without postulating a 

parametric model for the distribution of the outcomes before nonresponse. To account for 

the nonresponse, Sverchkov (2008) assumes a response model and estimates the response 

probabilities by applying the missing information principle (MIP), which consists of 

defining the likelihood as if there was complete response, and then integrating out the 
2

unobserved outcomes from the likelihood, employing the relationship between the 

distributions of the observed and unobserved data. We describe the main steps of this 

approach in Sections 2, 3 and 4.

In this paper, Section 5, we show how this approach can be applied to testing whether 

probability weighted estimators, corrected for MAR nonresponse, are still biased due to 

ignoring that nonresponse is NMAR actually. In Section 6 We illustrate this approach on 

real data (Consumer Expenditure Survey) example.

 2. NOTATION AND MODELS

Let 

{ , ; 1,..., } i i y i N x =

represent the data in a finite population of 

N

units, where 

i

y

is the 

value of the outcome variable for unit 

i

and

,1 , ( ,..., ) i i i K x

= x x

is a vector of 

corresponding 

K

covariates. Suppose that the population outcome values follow model

(2.1):

 | ~ ( | ), 1,..., i i i i y f y i N x x = (2.1)

The target is to estimate the population mean 

1

1

N

i

i

Y N y −

=

=  , based on a sample 

s

of 

n

units with inclusion probabilities 

Pr( ) i  = i s . Denote by

i

I

the sample indicator; 

1 i

I =

if unit

i

is selected in the sample and 0 otherwise. Let 

1/ wi i= 

denote the sampling 

weights.

In practice, not every unit in the sample responds. Define the response indicator;

1 Ri=

if 

unit 

i s 

responds and

0 Ri=

otherwise. The sample of respondents is thus

{ : 1, 1} R i I R = = =

i i

and the sample of nonrespondents among the sampled units is 

{ : 1, 0} c R i I R = = =

i i

. 

Assumption 1. The response process is assumed to occur stochastically, independently 

between units, and 

1

0

n

i

i

R

=   . (This assumption is used in the following Eq. 3.1 in 

Section 3 for simplification.)
3

The sample of respondents defines therefore a second, self-selected stage of the sampling 

process with unknown response probabilities. (Särndal and Swensson, 1987).

Under the population model (2.1), the observed data follow ‘respondents’ model:

 ( | ) ( | , ) R i i i i f y f y i R x x =  (2.2) 

The model (2.2) is again general and all that we state at this stage is that under informative 

sampling and/or NMAR nonresponse, the models for the respondents and for the 

population differ;

( | ) ( | ) R i i i i f y f y x x  .

Remark 1. The respondents’ model refers to the observed data and hence can be estimated 

and tested by standard methods.

Let

( , ) r i i p y x = Pr( 1| , , ) R y i s i i i =  x

. If the probabilities 

( , ) r i i p y x

were known, the 

sample of respondents could be considered as a sample from the finite population with 

known sampling probabilities 

( , ) i i r i i   = p y x . In this case, the population mean 

Y

can 

be estimated, for example, by probability weighted estimator, 

1

1

/

ˆ

1/

n

i i

i

n

i

i

y

Y





=

=

=





. Also, if 

known, the response probabilities could be used for imputation of the missing data within 

the selected areas, by applying the relationship between the sample and samplecomplement distributions, (Sverchkov and Pfeffermann, 2004);

( | , )c

i i f y i R x  =

1

1

[ ( , ) 1] ( | , )

{[ ( , ) 1]| , }

r i i i i

r i i i

p y f y i R

E p y i R

−

−

− 

− 

x x

x x

 (2.3)

(here and in what follows 

1

a a 1/ −=

). 

Remark 2. Equation (2.3) is a key relationship for derivation of the main result of the next 

section, Equation (3.2). It allows the distribution of the unobserved outcomes to be written 

as a function of the distribution of the observed outcomes and response and/or selection 

model.

3. ESTIMATION OF RESPONSE PROBABILITIES
4

Unlike the sampling probabilities, the response probabilities are generally unknown. We 

assume therefore a parametric model, which is allowed to depend on the outcome and the 

covariate values;

Pr[ 1| , , ; ] R y i s i i i =  x γ ( , ; ) r i i= p y x γ

, where 

γ

is a vector of 

unknown coefficients. We assume that 

( , ; ) r i i p y x γ

is differentiable with respect to 

γ . In 

Section 4 we suggest a method for selection the response model.

Under these assumptions and Assumption 1, if the missing outcome values were observed, 

γ

could be estimated by solving the likelihood equations:

 

log ( , ; ) r i i

i R

p y





+





x γ

γ k

log[1 ( , ; )] 0

c

r k k

R

p y



 −

=





x γ

γ

. (3.1)

In practice, the missing data are unobserved and hence the likelihood equations (3.1) are 

not operational. However, one may apply in this case the missing information principle:

Missing Information Principle (MIP, Cepillini et al. 1955, Orchard and Woodbury, 

1972): Let 

{ , , ; , t 1,..., } O y n i R n =  = i t x

represent the known observed data used 

below. Since no observations are available for 

c

i R 

, solve instead,

k

log ( , ; ) log[1 ( , ; )]

c

r i i r k k

i R R

p y p y E O

 

        −   +             x γ x γ

γ γ

(2.3) log ( , ; ) by

r i i

i R

p y





=





x γ

γ

 

1

1

log[1 ( , ; )] [ ( , ; ) 1] | ,

0.

c {[ ( , ; ) 1]| , }

r k k

r k k k

k R r k k k

p y E p y k R

E p y k R

−

−



   −   − 

   + =

−  

x γ

x γ x

γ

x γ x

 (3.2) 

Notice that the expectations in the last expression are with respect to the model holding for 

the observed data for the respondents. This result is possible because of key relation (2.3), 

see Remark 2 See Sverchkov (2008) for detailed derivation of (3.2). 

Remark 3. When the response probabilities 

( , ; ) r i i p y x γ

depend on only 

i

x , they are 

referred to as propensity scores, and the missing data are missing at random. The estimating 

equations in (3.2) reduce in this case to the common log-likelihood equations,

 

log ( ; ) r i

i R

p





+





x γ

γ

log[1 ( ; )] 0,

c

r k

k R

p



 −

=





x γ

γ

 (3.3) 
5

where 

( ; ) Pr(R 1| ; ) r i i i p x γ x γ = = .

Imagine an unrealistic situation when we know response probabilities, 

*

Pr[ 1| , , ; ] R y i s i i i =  x γ

*

( , ; ) = p y r i i x γ , where 

*

γ

is a true value of

γ

, but we still 

want to estimate 

γ

by solving (3.2). Then the equations (3.2) can be solved by maximizing 

the log-likelihood,

( ) log ( , ; ) [1 log ( , ; )] log ( , ; )

c

r i i r i i r i i

i R i R i R

l E p y p y O p y

  

 

= + − =    γ x γ x γ x γ    

1 *

1 *

{[ ( , ; ) 1]log[1 ( , ; )]| , }

c {[ ( , ; ) 1]| , }

r k k r k k k

k R re r k k k

E p y p y k R

E p y k R

−

−



− − 

+

− 



x γ x γ x

x γ x

. (3.4)

Although such scenario is unrealistic, it suggests the following iteration algorithm of 

solving (3.2): starting with some initial estimate of parameter 

γ

, say 

0

γ

, for example 

0

γ

is a solution of (3.3), maximize in the (q+1) iteration the expression,

( 1) log ( , ; ) q

r i i

i R

p y

+



 x γ

1 ( ) ( 1)

1 ( )

{[ ( , ; ) 1]log[1 ( , ; )]| , }

c {[ ( , ; ) 1]| , }

q q

r k k r k k k

q

k R r k k k

E p y p y k R

E p y k R

− +

−



− − 

+

− 



x γ x γ x

x γ x

 (3.5) 

with respect to 

( 1) q+

γ . The maximization can be carried out, for example, by SAS Proc 

NLIN.

Remark 4. A fundamental question regarding the solution of the MIP equations is the 

existence of a unique global solution or more generally, the identifiability of the response

model. Riddles et al. (2016) propose a similar approach to deal with NMAR nonresponse 

in the general context of survey sampling inference and establish the following 

fundamental condition for the response model identifiability: the covariates 

x

can be 

decomposed as 

1 2 x x x = ( , )

with 

2

dim( ) 1 x 

, such that

Pr( 1| , ) R y i i i= x

= = Pr( 1| , ) R y i i i 1

x

. In other words, the covariates in 

2

x

that appear in the outcome model 

do not affect the response probabilities, given the outcome and the other covariates.

Riddles et al. (2016) prove asymptotic normality of the estimate 

γ

ˆ

under general regularity 

conditions.
6

4. SELECTION OF A RESPONSE MODEL

When following the approach proposed in Section 3, the likelihood (3.4) suggests at least 

two procedures for the selection of the response model under NMAR nonresponse. One 

procedure involves comparing different models based on information criteria such as the 

Akaike information criterion, AIC= 2 ( ) 2dim( ) − + l γ γ , or Schwarz information 

criterion, 

BIC= 2 ( ) dim( )log( ), i

i s

l n n n



− + = γ γ 

; a second procedure involves testing

a saturated versus a nested model based on the likelihood ratio test. In Section 5 we 

illustrate via a simulation study how the likelihood (3.4) can be used for the application of

these selection procedures.

 5. TESTING BIAS OF MEAN ESTIMATES DUE TO NOT MISSING AT 

RANDOM NONRESPONSE

Let 

ˆ

MAR

wi

be a good estimate of 

1 1

( ) Pr( 1| , ) i i

r i i i

w w

p R i s

=

x x = 

, in other words 

ˆ

MAR

wi

is a good estimate of final inclusion weight if (or “assuming”) the nonresponse is 

MAR. Then the obvious estimate for the mean will be 

1

1

ˆ

ˆ

ˆ

n

MAR

i i

MAR i

n

MAR

i

i

w y

Y

w

=

=

=





and if the 

response is really MAR then this estimate is consistent.

How one can check in real situation that the latter estimate is not biased? Formally one can 

estimate response probabilities 

Pr[ 1| , , ; ] R y i s i i i =  x γ ( , ; ) r i i= p y x γ

following the 

same approach as in Section 3 and then compare model 

( , ; ) r i i p y x γ

with 

( ; ) r i p x γ

by 

some information criterion as in Section 4 or check whether MAR and NMAR estimates, 
7

1

1

ˆ

ˆ

ˆ

n

MAR

i i

MAR i

n

MAR

i

i

w y

Y

w

=

=

=





and 

1

1

ˆ

ˆ

ˆ

n

i i

NMAR i

n

i

i

w y

Y

w

=

=

=





, 

1 ˆ

( , , )ˆ

i i

r i i

w w

p y

=

x γ

, are significantly 

different (use 

ˆ

var( ) ˆ

MAR Y

for the test, the latter usually available).

As one can see, estimation procedure (iterative maximization of the likelihood (3.5)) is 

complicated optimization problem if dimension of vector of covariates 

i

x

is big, which is 

usually the case in real applications. Moreover, often correction to the response mechanism 

is made through different calibration procedures and so that response model, 

Pr[ 1| , , ; ] R y i s i i i =  x γ ( , ; ) r i i= p y x γ , can be hard to identify. Assuming that NMAR 

estimate is consistent, 

ˆ ˆ 1

( )

j

N

MAR MAR j

y

Bias Y Y

N

=

= −

 1 1

1

ˆ

( ) ˆ

ˆ

j

N n

MAR MAR

i i

j i

n

i

i

y Y w

N

w



= =

=

−

=  −

 



,

where

1

1

ˆ

ˆ

n

MAR

i i

MAR i

i i n

MAR

i

i

w y

y

w



=

=

= −





, and therefore one can test whether 

ˆ

MAR Y

is biased or not 

by testing if 

Pr( 1| ) MAR Ri i

= 

depends on 

MAR

i



or not which can be done by the approach 

suggested in Sections 3 and 4. Therefore dimension of optimization procedure can be 

reduced dramatically. Note that although 

Pr( 1| ) MAR Ri i= 

does not involve any 

covariates, we still need a covariate 

2i

x

such that 

MAR

i



correlates with 

2i x

but 

Pr( 1| , ) Pr( 1| ) 2

MAR MAR R R i i i i i

= = =   x

, otherwise the solution of (3.5) can be not 

unique, see Remark 4.

 6. APPLICATION TO CONSUMER EXPENDUTURE SURVEY

We consider estimation of the following expenditures, Total, Food, Housing, Health, based

on 2019, 2020 and 2021 samples. We assume truncated logit model for response 

probability
8

0 1

0 1

exp( + ) 0.0001

1+exp( + ) Pr[ 1| , ; ] ( ; ) 1.0001

MAR

i

MAR

MAR MAR i R i s p i i r i

  

  

 

 

  +

  =  = = γ γ

where 

MAR

i



is defined as in Section 5 by a particular expenditure and MAR final sample 

weight, FINLWT21, used currently by Consumer Expenditure Survey. We use “median 

income” variable as instrumental variable

2i x

(see Remark 4), correlation coefficients 

between 

2i x and Total, Food, Housing, Health are equal to 0.31, 0.26, 0.35, 0.08 

respectively for year 2021, and simple linear model explain relation between the 

expenditures and 

2i x

reasonably although with low R-square statistics, 0.15, 0.07, 0.16, 

0.20. Parameter estimates for response model, 

Pr[ 1| , ; ] ( ; ) MAR MAR R i s p i i r i =  =   γ γ , 

are as follows (all significant with p-value<0.0001, p-values are produced by SAS Proc 

NLIN used in maximization of the likelihood (3.4)-(3.5), this suggests that the response 

model is NMAR not MAR):

Total: 

0



= -0.513, 

1



= 0.0000085.

Food: 

0



= -0.510, 

1



= 0.0000086.

Housing: 

0



= -0.496, 

1



= 0.0000303.

Health: 

0



= -0.619, 

1



= 0.0007000.

(Similar results for these relationships for samples of year 2019 and 2020).

AIC and BIC criteria for all 4 expenditures and all years considered in this research prefer 

response model that include 

MAR

i



as a covariate (NMAR nonresponse) versus response 

model that assume MAR response, results summarized in the table below. 
9

AIC

 MAR NMAR

BIC

MAR NMAR

ˆ

ˆ ( )

ˆ

ˆ ( )

MAR

MAR

Bias Y

std Y

Total 2019

Total 2020

Total 2021

68372.9 66794.6

70538.0 69337.4 

71907.6 71082.4 

68380.9 66810.6

70545.9 69353.3 

71915.6 71098.2 

1.00

0.49

1.67

Food 2019

Food 2020

Food 2021

68372.9 66872.5 

70538.0 69337.4 

71907.6 71159.3 

68380.9 66888.4 

70545.9 69353.3 

71915.6 71175.1

0.71

2.50

0.65

Housing 2019

Housing 2020

Housing 2021

68372.9 66858.2 

70538.0 69366.4 

71907.6 71049.2

68380.9 66874.2 

70545.9 69382.0 

71915.6 71065.1

0.01

2.86

3.15

Health 2019

Health 2020

Health 2021

68372.9 66872.1 

70538.0 69337.4

71907.6 71142.1 

68380.9 66888.2 

70545.9 69353.3 

71915.6 71157.9

-1.11

-0.66

0.28

The last column of the table presents an estimate of normalized bias of MAR estimator, 

ˆ

ˆ ( )

ˆ

ˆ ( )

MAR

MAR

Bias Y

std Y

, where 

1

1

ˆ

ˆ ˆ ˆ

ˆ ( ) ( )

ˆ

n

MAR

i i

MAR NMAR MAR i

n

i

i

w

Bias Y Y Y

w



=

=

= − − = −





, and 

ˆ

ˆ ( ) MAR std Y

is 

a BRR standard error estimate of 

ˆ

MAR Y

. Note that the later statistics can be used only as 

crude approximation of normalize bias 

ˆ

( )

ˆ

[ ( )]

MAR

MAR

Bias Y

std Bias Y

, since it does not include errors 

due to 

ˆ

NMAR Y

estimation.
10

The above three statistics (naïve approximations) suggest that the biases of the estimate 

due to ignoring NMAR nonresponse are not (or slightly, see Housing) significant based on 

the standard errors of the final estimates. Note that standard error of 

ˆ

ˆ ( ) MAR Bias Y

is a 

standard error of the difference of MAR and NMAR estimators, therefore the standard error 

of 

ˆ

MAR Y

can be used only for crude comparison of two estimators.

Conclusion: Although current Consumer Expenditure Survey nonresponse 

adjustment of sample weights does not remove nonresponse bias (MAR estimates for 

the mean total expenditures are still biased after the adjustment), the bias is not or 

only slightly significant based on the standard error of the final estimates.

REFERENCES

Cepillini, R., Siniscialco, M., and Smith, C.A.B. (1955). The estimation of gene 

frequencies in a random mating population. Annals of Human Genetics, 20, 97-115.

Orchard, T., and Woodbury, M.A. (1972). A missing information principle: theory and 

application. Proceedings of the 6th Berkeley Symposium on Mathematical Statistics and 

Probability, 1, 697-715.

Pfeffermann, D., and Sverchkov, M. (2009). Inference under Informative Sampling. In: 

Handbook of Statistics 29B; Sample Surveys: Inference and Analysis. Eds. D. Pfeffermann 

and C.R. Rao. North Holland. pp. 455-487.

Riddles, K.M., Kim, J.K. and Im, J. (2016). A propensity-score adjustment method for 

nonignorable nonresponse. Journal of Survey Statistics and Methodology, 4, 215-245.

Särndal, C.E. and Swensson B. (1987). A general view of estimation for two phases of 

selection with applications to two-phase sampling and nonresponse. International 

Statistical Review, 55, 279-294.
11

Sverchkov, M. (2008). A new approach to estimation of response probabilities when 

missing data are not missing at random. Joint Statistical Meetings, Proceedings of the 

Section on Survey Research Methods, 867-874.

Sverchkov, M., and Pfeffermann, D. (2004). Prediction of Finite Population Totals Based 

on the Sample Distribution. Survey Methodology, 30, 79-92.