Two-Step versus Simultaneous Estimation of Survey-Non-Sampling 

Error and True Value Components of Small Area Sample Estimators 

a b P.A.V.B. Swamy , T.S. Zimmerman and J.S. Mehta c 

a,b Bureau of Labor Statistics, Room 4985, 2 Massachusetts Avenue, N.E., Washington, DC 20212, E-mail: paravastu_s@bls.gov 

and Zimmerman_tamara@bls.gov. 

c Department of Mathematics, Temple University, Philadelphia, PA 19122. E-mail: 

mehta1007@comcast.net. 

Abstract 

The precision of a small-domain sample estimator can be improved by estimating simultaneously its true value and sampling 

and non-sampling error components. In principle, this simultaneous estimation is superior to any two-step estimation of the 

true value and sampling error components, ignoring the non-sampling error component. In this paper, a time series model for 

state employment or unemployment is used to demonstrate the limitations of a two-step method. A cross-sectional model for 

state employment or unemployment is used to explain the advantages of simultaneously estimating a true value and the sums 

of sampling and non-sampling errors in two or more sample estimators of the true value. 

Key Words. Indirect model-dependent estimate; domain indirect; time indirect; domain and time indirect; coefficient driver. 

1. Introduction 

Accurate estimates of employment and unemployment at various levels of geographic detail are needed to formulate good 

regional policies, such as the determination of eligibility for and/or the allocation of Federal resources, with a good 

understanding of local economic conditions. Such estimates cannot be produced from survey data collected in each small area 

because sample sizes within those areas are often either zero or too small to provide reliable estimates. If an area-specific 

sample is available but is not large enough to yield “direct estimates” of adequate accuracy, then their accuracy can be 

improved, using “indirect” model-dependent estimators that “borrow strength” from auxiliary data collected in this and other 

small areas and/or at more than one time period. On a definition suggested by a passage in Rao (2003, p. 2) indirect 

estimators might be classified as “domain indirect”, “time indirect”, or “domain and time indirect” depending on whether 

they borrow strength cross-sectionally, over time, or both. In this paper, we develop domain indirect estimators that facilitate 

simultaneous estimation of the true value and sampling and non-sampling error components of sample estimators. In 

principle, simultaneous estimation is superior to two-step estimation of these components. To illustrate the proposed 

techniques, we consider the problem of improving the accuracy of current population survey (CPS) estimates of employment 

and unemployment for 51 “small areas” consisting of the 50 United States plus the District of Columbia. CPS state estimates 

are widely used and improving their accuracy is valuable. 

The remainder of this paper is divided into six sections. Section 2 describes the time series models used to estimate 

employment or unemployment for the states, major metropolitan areas, and corresponding balance of states in the U.S. 

Section 3 shows the limitations of a two-step method used to estimate these time series models. The consequences of 

incorrectly neglecting non-sampling errors are given in Section 4. Section 5 uses cross-sectional estimates based on two or 

more sample estimators of state employment (or unemployment) to estimate simultaneously the true values of employment 

(or unemployment) and the sums of sampling and non-sampling errors contained in those estimators. Seasonal adjustment of 

the true values of state employment (or unemployment) and estimation of the autocovariances of sampling errors are 

discussed in Section 6. Section 7 provides an example. Section 8 concludes. 

2. Time Series Models for the True Value Component of a Sample Estimator 

2.1 Sampling Model 

Acknowledgment. Comments by Christopher Bollinger are highly appreciated. 

Disclaimer. Any opinions expressed in this paper are those of the authors and do not constitute policy of the Bureau of Labor Statistics. 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ˆCPS Let Yit denote the current population survey (CPS) composite estimator of the true population value, denoted by Yit , of 

either unemployment or employment. 1

 Let i and j index states, major metropolitan areas, or balance of states and let t index 

ˆCPS months. The decomposition of Yit into its unobserved components is 

ˆCPS CPS CPS Y = Y + e + u , t = 1, ..., τ , i = 1, ..., m (1) it it it it i 

CPS CPS where eit and uit represent sampling (or survey) and non-sampling errors (such as non-response or measurement errors), 

CPS ˆCPS respectively. A value, denoted by yit , of Yit obtained from a specific CPS sample is an estimate of Yit . 

The annual averages of the components of a type of non-response rates for 1993-1996 and 2003 CPS national estimates are 

given in U.S. Department of Labor (2006, p. 16-4). These averages show that we cannot assume that uit

CPS = 0. The methods 

in Cochran (1977, Chapter 13) might have been already tried to reduce the magnitude of uit

CPS . Under the current 2000 CPS 

design, the number of assigned households in monthly CPS samples for different states ranges from 700 to 5,344. These 

ˆCPS 2

sample sizes are not large enough to yield Yit of adequate precision. It is necessary to use model-based approaches to 

ˆCPS improve the efficiency of Yit . 

The design-based approach to CPS inference treats the Yit as fixed quantities. Unlike this approach, the model-based 

approach to CPS sampling inference treats the Yit as a random sample from a “superpopulation” and assigns to them a

probability distribution implied by a model for Yit .

3

 Since the “true” model for Yit is unknown, we can use one of its 

approximations in place of Yit used in (1) but then we should be prepared to carefully investigate whether such a use 

ˆCPS 4

improves the precision of Yit . Let Ep 

and Vp  denote the expectation and variance operators with respect to the 

probability distribution induced by the 2000 CPS sampling design conditional on Yit , respectively. Let Em and Vm

 denote 

the expectation and variance operators with respect to a probability distribution assigned to the Yit by an assumed model for 

Yit , respectively (see Little, 2004, p. 547). 

2.2 Survey Error (SurE) Model 

CPS CPS Suppose that uit = 0 for all i and t. Then a model for eit implied by the CPS design is 

σ

CPS e

CPS

,it * 

e = e (2) it it

σ *

e i, 

CPS CPS CPS CPS 2 where for i, j = 1, ..., m and t = 1, ..., τ , E e( Y ) Cov e , ) | Y Y ] = 0 if i ≠ Ve | Y ) = = i p it | it = 0, p 

[( it ejt it , jt j, p ( it it σ e

CPS 

,it 

2 2 2 σ b eCPS + D CPS σ CPS with σ b eCPS = between-Primary Sampling Unit (PSU) variance component of the total design 

, ,it e ,it 1e ,it , ,it 

CPS 2

variance of e arising from selecting a subset of all non-self-representing (NSR) PSUs in state i, D CPS CPS σ = within- it e ,it 1e ,it 

PSU variance component of the total design variance of eit

CPS arising from sampling of housing units within selected PSUs of 

state i, D 

it 

 = the ratio of the within-PSU variance component of the total design variance of eit

CPS 

CPS , assuming the CPS 

e , 

sample design, to the total design variance, σ

2 

CPS = (Nit / nit )Yit (1 − (Yit / Nit )) , of an unbiased estimator of Yit , assuming 

1e ,it 

simple random sampling (design effect), Nit = the civilian non-institutionalized population 16 years of age and older for the 

*

state, (Nit / nit ) = the state sampling interval (see U.S. Department of Labor, 2006, p. 3-6) and eit follows a mixed 2 nd order 

1

 Here employment is defined on a place-of-residence basis, since the CPS is a household survey. 

2 We strictly adhere to Cochran’s (1977, p. 16) distinction between “precision” and “accuracy”. 

3

 Rao (2003, p. 77) calls this model a linking model. 

4

 The “true” model for Y is an unknown function, Y = f (x ,..., x , x ,..., x ) , of all of the determinants of Y with the it it it 1it Kit K +1, it Lit ,it it

correct but unknown functional form and without any omitted determinants or mismeasured variables. Any of its estimable 

approximations is misspecified if it has an incorrect functional form and suffers from omitted-variable and measurementerror biases (see Freedman and Navidi, 1986). 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

* 2

autoregressive and 17th order moving average, or ARMA(2, 17), model with m

( ) = σ * V e , as shown by Tiller (1992, 2005). it e i,

This ARMA(2, 17) model is approximated by a 15th order autoregression or AR(15) model. 

2.3 Linking Model 

Durbin and Koopman (2001) work with Harvey’s (1989) structural time series models. Following them, it is assumed that the 

finite population is generated according to the superpopulation model 

Y = T + S + I (3) it Y it Y it Y it it ,

it ,it ,

where T is the trend-cycle, S the seasonal, and I the irregular component, which are treated as the unobserved Y it Y it Y it it ,

it ,it ,

components of Yit .

5

 All these components are assumed to be independent of each other for all i and t. This assumption may 

not be true if the decomposition of Yit in equation (3) is not unique. To account for this non-uniqueness, Havenner and 

Swamy (1981) assume that T , S , and I are correlated with each other for all i and t.

Y it Y it Y it it ,

it ,it ,

It is assumed that TY it follows a random walk model with one period lagged value of a random walk drift. Durbin and 

it ,

Koopman (2001, p. 39) call this the local linear trend model. The seasonal component, SY it , is expressed in a trigonometric 

it ,

form to make it follow Durbin and Koopman’s (2001, (3.6), p. 40) quasi-random walk model. The irregular component, 

I

Y it , is assumed to be serially independent and normally distributed with mean zero and constant variance. 6 

it ,

It is known that employment has a strong tendency to move cyclically, downward in general business slowdowns and upward 

in expansions. These cycles are asymmetric because employment decreases at a faster rate than it increases. The behavior of 

unemployment over time is the opposite of employment’s behavior and is called asymmetric counter-cyclical behavior. 

Montgomery, Zarnowitz, Tsay and Tiao’s (1998, p. 487) results imply that the time series models assumed above for the 

components of Yit in (3) may be misspecified because they may not exhibit the asymmetric cycles of employment and 

unemployment. More specifically, their results are: (i) A first order autoregressive model in first differences, denoted by 

ARIMA(1, 1, 0), that fitted the U.S. quarterly unemployment rate series for 1948-1993 quite well was not able to accurately 

represent the asymmetric cycles of unemployment during this period; (ii) An ARIMA (1, 1, 0)(4, 0, 4) model with a 

multiplicative seasonal factor, denoted by ARIMA(4, 0, 4), under-predicted the U.S. unemployment rate during the rapid 

increase of 1982 and exhibited forecasts that fluctuated a great deal more during stable periods of unemployment. 7 

2.4 Covariate Model 

Consider an extension of model (3) to allow one or more of the unobserved components of Yit to be related to corresponding 

components in another series called a covariate. The unemployment insurance claims (U.I. Claims) from the Federal-State 

Unemployment Insurance System are used as a covariate series if Yit represents unemployment, and nonagricultural payroll 

employment estimates from the Current Employment Statistics (CES) survey are used as a covariate series if Yit represents 

employment. 

3. A Two-Step Method of Estimating the Parameters of Survey Error and Linking Models 

We call the model consisting of equation (1) without uit , equation (2), and equation (3) with covariates “the combined 

model”. For each i, the unknown parameters of model (2) are: the design variances, σ

2

, t = 1, ..., τ i , and 17 parameters 

e

CPS 

,it 

of AR(15) model. For each i, the unknown parameters of model (3) are the error variances of the models for T , S , and Y it Y it it ,

it ,

I

Y it , which are four in number (see Tiller, 2005, p. 7). When the yit

CPS for t = 1, ..., τ i , are augmented with the available

it ,

5 Tiller (2005) also allows for temporary or permanent shifts in the level of the series, { Yit }, t = −∞ , +∞ . 

6 Someimes it is assumed that I

Y it follows a lower-order autoregressive process (see Tiller, 1992, p. 151). 

it ,

7 It can be shown that Yit contains a random walk component if it contains a stationary component after being first differenced d times. It is 

precisely this condition that is violated when the nonlinear models that accurately represent the asymmetric cycles of employment and 

unemployment are considered. 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

observations on the covariates, the number of the unknown parameters of (3) with covariates increases to 10 for each i (see 

Tiller, 2005, p. 10). All these parameters of the combined model are not identifiable on the basis of the available CPS 

estimates, yit

CPS 

, t = 1, ..., τ i , i = 1, ..., m, and observations on the covariates alone. Consequently, for each i, a two-step 

procedure is used: First step. For each i, the available CPS estimates ( yit

CPS 

, t = 1, ..., τ i ) are used to estimate the parameters 

of (2) independently of (3) with covariates; Second step. Given the CPS estimates ( yit

CPS 

, t = 1, ..., τ i ), the parameters of (3) 

with covariates are estimated, holding the parameters of (2) and AR(15) fixed at their estimated values. 

3.1 Survey Error Variance Estimates 

CPS 2 Let i and t be fixed so that the components of Y in (3) are fixed. Let be absent. Then the design variance, σ CPS , in (2) it uit e ,it 

can be estimated using the generalized variance function (GVF) (see Wolter (1985, p. 203), Lent (1991, 1994), and U.S. 

Department of Labor (2006, p. 3-6)), 

CPS 2

( VY ˆ| Y ) = −(b / N )Y + bY (4) p it it it it it it it 

2 2 2 where = ( / n )D CPS , = 1/(1 , = σ CPS / σ CPS , and D CPS , , , σ b e

kN CPS b , and it it it it e 

kit − pit ) pit it , e , 

Nit nit 

,it be , , e it it , , it 

2 ˆCPS σ (= VY | ) are as defined below (2). 

p 

( it Yit )

e

CPS 

,it 

A method of estimating the bit is summarized in the following steps: 

1. All survey statistics are divided into several groups with model (4) fitted independently in each group. If bit differs 

among the statistics in a group, then they are not consistently estimable because the number of unknown parameters 

increases with the size of the group (see Lehmann and Casella, 1998, p. 481). To avoid this difficult incidental 

parameter problem, care is taken to group together all survey statistics that follow a common model such as (4) with 

fixed i and t. This may involve grouping together the statistics with similar design effects; the same characteristics for 

selected demographic or geographic subgroups. Examples of such groupings are given in Wolter (1985, p. 209). 

2. Let g index the statistics in a group formed in Step 1. The design variances, VY CPS ( ) | Y g )

p (

ˆ 

it g it ( ) , are computed for 

ˆ ˆCPS 2 ˆ ˆCPS ˆ ˆ ˆ ˆ

several statistics of this group using −( / )( ( b N Y g )) ( + bY g) , where b D = k (N / n ) , k =1/(1 − pˆ ) ,

it it it it it it it it it CPS it it e ,it 

2 2 2 2 2 2

pˆ = σˆ CPS / σˆ CPS , σˆ CPS = σˆ

b eCPS + Dˆ 

CPS σˆ CPS , σˆ

b eCPS is obtained by multiplying the 2000 decennial census it be it ,it , 

ˆCPS () 

, ,it e ,it e , , ,it e 1e it , ,it 

estimate of between-PSU variance by the square of the ratio of the annual average of Yit g for the current year to the 

2000 census estimate of Yg , Dˆ is computed from the national within-PSU design effects by adjusting for it ( )

e

CPS

,it 

differences in a certain noninterview rate, see the memorandum (2007) from Khandaker Mansur, and σˆ

2 

= σ

2 

CPS CPS 1e ,it 1e ,it 

ˆCPS () ˆ ˆCPS (

ˆCPS when Yg VY g g g g .

it ( ) = Y g . Let ( ( ) | Y ( ) ) denote the estimator of VY ( ) | Y ( ) )

it p it it p it it 

ˆCPS () ˆ ˆCPS Yit p ( it g it ( ) )from Step 2 gives the estimate of it 3. Fitting model (4) to the data ( g , VY ( ) | Y g ) b for the group. The 

model fitting technique is an iterative weighted least squares procedure, where the weight is the inverse of the square of 

ˆ ˆCPS ˆCPS 2 b

#

it the predicted value of p ( it ( )| Yit ( )) /( it VY g g Y g ( )) (see Wolter, 1985, p. 207, (5.4.2)). Let denote an iterative 

weighted least squares estimator of bit . 

ˆCPS 4. An estimate of the design variance of a survey statistic, say Y ( ) it g , for which the successive difference replication 

CPS method of Fay and Train (1995) is not applied is now obtained by evaluating model (4) at the point ( Y ( ) ˆ

it g ; b

#

it ). It is 

called a GVF estimate. 

The estimator, b

#

it , obtained in Step 3 may not have any desirable statistical properties because of the effects of the errors, 

ˆCPS Y ( ) it g - Yg , and can be very imprecise if all statistics within a group do not behave according to model (4). These it ( )

problems may get resolved if model (4) is replaced by the less restrictive model, 

ˆ ˆCPS ˆCPS ˆCPS 2 

p ( it ( )| Yit ( )) = bitg { it () − (1/ Nit )( it VY g g Y g Y g ( )) } , g = 1, 2, ..., G (5) 

with 

ˆ ˆ bitg = δ 01kitg (Nit / nit )De

CPS 

,itg 

+ξ 0itg (6) 
 

 

 

 

 

 

 

 

 

 

 

 

where Dˆ 

e , 

and k

ˆ 

itg are the estimators of De 

and kitg given in Step 2, and ξ 0itg CPS CPS is a random variable with mean zero 

itg ,itg 

and constant variance. Chang, Swamy, Hallahan and Tavlas’ (2000) method can be used to estimate model (5) under 

ˆCPS assumption (6). Estimator (5) is imprecise, since Y ( ) it g is an imprecise estimator of Yg .

it ( )

3.2 Survey Error Autocorrelation Estimates 

In Zimmerman and Robison (1995), the variances and autocovariances of the eit 

* in (2) are estimated from the separate panel 

(rotation group) estimates of Yit . A panel is defined as the set of sampling units joining and leaving the sample at the same 

time. Each panel is a representative sample of the population. The CPS sample consists of 8 such panels in every month. To 

estimate model variance and autocovariances of the CPS composite estimator from the estimates of model variances and 

autocovariances of the panel estimators, Zimmerman and Robison (1995) and Zimmerman (2007) consider the following 

model: 

ˆCPS CPS Y = μ θ+ + β + e (7) ijt i ij it ijt 

where Yˆ

ijt

CPS = CPS estimator of Yit from panel j, μi = overall mean, θij = month-in-sample effect, βit = time effect, and

CPS ˆCPS eijt = sampling error of Yijt . 

CPS CPS Let eijt denote eit in (2) for panel j. Model (7) is misspecified if the time series models for the components of Yit in 

CPS CPS ˆ ˆ (3) are correctly specified. If this is so, model (7) incorrectly estimates eijt as Yˆ 

ijt − μˆ

i −θij − βit instead of as 

ˆCPS ˆ ˆ ˆ CPS Yijt −Tyijt ,ijt − Syijt ,ijt − I

yijt ,ijt . The incorrectly estimated eijt cannot give accurate estimates of the parameters of model (2) 

other than σ

2 even if model (5) gives accurate estimates of σ2 

CPS CPS . Thus, the First step--of estimating the survey error 

e ,it e ,it 

model parameters from design based information independently of the time series models for the components of Yit in (3)--

can lead to very inaccurate estimates of the survey error model parameters. If the condition, “holding the parameters of the 

survey error model fixed at their estimated values”, means the condition, “the parameters of the survey error model are set 

equal to the estimates obtained in the First step”, then the inaccurate estimates obtained in the First step can lead to inaccurate 

estimates of the parameters of model (3) with covariates in the Second step. There will be no consistency between the

estimates obtained in the First step and those obtained in the Second step. Furthermore, the random walk models assumed for 

some of the components of Yit lead to a predictor of Yit that is unconditionally inadmissible relative to quadratic loss 

functions because the predictor does not possess finite unconditional mean. They do not provide the predictors of Yit with 

good conditional and unconditional properties. Brown (1990, p. 491) shows the importance of working with such predictors. 

3.3 State-Space Form of the Combined Model 

Let { xit } be a covariate series. Let yxit = ( yit

CPS 

, xit )′ . Let the combined model be written in state-space form, as in Durbin 

and Koopman (2001, p. 38). Hold the parameters of the survey error model fixed at the estimates obtained in the First step. 

For calculating in the Second step loglikelihood function and the maximization of it with respect to the parameters of model 

(3) with covariates, the joint density of the sample observations, yxi1

, yxi2, ..., yxiτ

, implied by the combined model is 

i 

written as 

τ i 

p y( ,..., y ) = ∏ py y ( | y ,..., ) (8) xi1 xi τ xit xit −1 xi 1 i t =1 

py (

1 where it is assumed that ( | y ) = p y ) with finite mean (see Durbin and Koopman, 2001, p. 30). The truth of this xi1 xi 0 xi

assumption is usually unknown. For this reason, it is convenient to assume that Y#

Xit is distributed with density 

( | y ,..., y 1) for all t;

8

 unfortunately, this assumption contradicts the assumption that ( | y ) = p y i py py ( ) with xit xit −1 xi xi1 xi0 x 1 

finite mean. This is because the former assumption says that under the random walk models assumed for some of the 

components of Y , Y# d Xit oes not possess finite unconditional mean for all t and the latter assumption says that Y#

Xit possesses it

8 Here we distinguish a random variable from its value by a tilde. For example, yxit is the value taken by the random variable 

CPS Y#

Xit = (Yˆ

, X )′ .

it it 
 

 

 

 

 

 

 

 

 

 

 

 

finite unconditional mean for t = 1. Even if it is assumed that Y# is Xit distributed with density py ( | y ,..., y )

xit xit −1 xi 1 for all t > 1,

py (

1

the assumption that ( | y ) = p y )

xi

 with finite mean cannot be true for all those data sets for which the value t = 1 xi1 xi 0 

occurs at different points on the time axis. 

One may use l (> 0) non-stationary elements in the state vector which determines the number of observations required to 

form priors of these elements. Without this or any other initialization, it is not possible to apply the Kalman filter to statespace form of the combined model in (1)-(3) (see Durbin and Koopman (2001, pp. 27-30, 99-104) and Maddala and Kim 

(1998, pp. 475-477)). The initial value yxi0

 is usually unknown and any assumption about it may be questioned. Incorrect 

assumptions about the yxi0 can lead to the estimates of Yit , t = 1, ..., τ i , with incorrect initial values and hence with incorrect 

time profiles. 

4. Consequences of Incorrectly Neglecting Non-Sampling Errors 

CPS CPS In the previous section, we have considered equation (1) with uit suppressed even though uit is not equal to zero with 

probability 1. Adding to this equation a linking model of the general linear mixed (GLM) model’s type for Yit , which 

implies that Yit possesses finite second moment for all i and t, the best linear unbiased predictor (BLUP) of Yit can be 

found. 9

 This BLUP can be expressed as a weighted average of Yˆ

it

CPS and the regression-synthetic estimator of Yit implied by 

ˆCPS the linking model. Under certain regularity conditions given in Rao (2003, p. 117), the BLUP of Yit coincides with Yit that 

ˆCPS ˆCPS is not affected by the misspecifications in the linking model as the design variance of Yit goes to 0. That is, when Yit is 

design-consistent, the BLUP of Yit can also be design-consistent. This result is the basis of Little’s (2004, p. 551) statement 

that one way of limiting the effects of model misspecification is to restrict attention to models that yield design-consistent 

ˆCPS estimators (see also Rao, 2003, p. 148). This observation is of no use to us when Yit is subject to non-sampling errors, in 

ˆCPS which case Yit is not design consistent (see Little, 2004, p. 549). 

In Rao’s (2003) terminology, the BLUP of Yit becomes an empirical BLUP (EBLUP) when the variance parameters involved 

in the BLUP are replaced by their respective sample estimators. A second-order approximation to the mean square error 

(MSE) of EBLUP involves three terms, as Rao’s (2003, p. 104, (6.2.31)) very elegant derivation shows. 10 This derivation 

further shows that an unbiased estimator of the MSE of EBLUP to a desired order of approximation may involve four terms 

(see Rao, 2003, p. 105, (6.2.37)). The value of the sum of these four terms obtained from a specific sample can exceed the 

ˆCPS ˆCPS value of the design variance of Yit , showing that the EBLUP of Yit can be less efficient than Yit which is itself imprecise 

because of the smallness of the sample on which it is based. The MSE of the BLUP of Yit contains two terms and only one of 

ˆCPS these terms is smaller than the design variance of Yit , as shown by Rao (2003, p. 117, (7.1.7)). These results raise the 

question: in what sense does the linking model for Yit “borrow strength” from its explanatory variables in making an estimate 

of Yit ? Under certain regularity conditions, two of the three terms of a second-order approximation to the MSE of an EBLUP 

ˆCPS of Yit go to zero and the remaining term remains below the design variance of Yit as the number of observations on the 

explanatory variables of the linking model goes to infinity (see Rao, 2003, p. 117). Thus, in the limit the MSE of an EBLUP 

of Yit can involve only one term that is smaller than the design variance of Yˆ

it

CPS 

.

11 Any linking model for Yit can borrow 

strength from its explanatory variables in the limit as the number of observations on its explanatory variables goes to infinity 

ˆCPS if Yit is design-consistent and some regularity conditions are satisfied. This result holds even when the linking model for 

9 Swamy, Zimmerman, Mehta and Robison (2005) extend Rao’s (2003, pp. 96-98) derivation of the BLUP from a GLM model for a 

sample estimator to take account of omitted-variable, measurement-error and incorrect functional-form (or simply specification) biases in 

the GLM model. 

10 This derivation has been extended to take account of specification biases in linking models of the GLM type in Swamy, Yaghi, Mehta 

and Chang (2006).

11 The gain in efficiency associated with the use of Rao’s (2003, pp. 116-117) BLUP may get reduced if the linking model that provides the 

BLUP suffers from specification biases (see Swamy et al., 2005). 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ˆCPS Yit is misspecified, provided Yit does not contain non-sampling errors. It is of limited use for us because large sample 

theories of estimation are irrelevant to small area estimation and the sources of error other than random sampling variation 

given in Cochran (1977, p. 359) are present in the CPS and CES survey. There is usually no justification to ignore the 

possible misspecifications in the time series models for the components of Yit in equation (3) with covariates. 

5. Simultaneous Estimation of Sums of Survey and Non-Sampling Errors and the Corresponding True Values 

ˆCPS ˆ A

Let us now reconsider Yit . Suppose that an additional estimator, denoted by Yit , of Yit is available. This additional 

estimator will be given by U.I. Claims data if Yit represents unemployment and by data from the CES survey or by 

extrapolated data from the Quarterly Census of Employment and Wages (QCEW) program if Yit represents employment. We 

ˆCPS ˆ A ˆCPS CPS assume that each of Yit and Yit is subject to both sampling and non-sampling errors and write Yit = Yit + ε it with 

CPS CPS CPS ˆ A A A A A CPS A CPS A

ε it = eit + uit and Yit = Yit + ε it with ε it = eit + uit , where the vectors, ( eit , eit )′ and ( uit , uit )′ , denote the 

ˆCPS A ˆ A

sampling and non-sampling errors of ( Yit , Yˆ 

it

A )′ , respectively.12 Let yit denote a value of Yit obtained from a specific 

sample. 

We now fix CPS A CPS At and let i vary to focus on cross-sectional variations in ε it and ε it . These types of variations in eit and eit 

are different from their random sampling variations that are present because only parts of the population have been measured 

13 CPS A

using some sampling designs. Modeling assumptions are needed to analyze the non-sampling errors, uit and uit , and 

CPS A ˆCPS ˆ A

cross-sectional variations in the sampling errors, eit and eit . The problem of choosing between Yit and Yit can only be 

solved if we can derive from them a single estimate that is closer to Yit than either estimate. We show in this section that such 

A CPS a single estimate can be found even when | yit - Yit | > | yit - Yit |. 

5.1 Random Coefficient Regression Model 

A A ˆ A ˆ A A ˆ Writing A A α1it 1( it / Yit ) , a function of Yit and ε it , we transform the sampling model, Yit = Yit + = − it ε ε , into the linking 

A ˆ A A ˆ A model, Yit = α1itYit . Replacing Yit in the sampling model in (1) by α1itYit gives the model: 

ˆCPS CPS A A Y =α +α Yˆ

, i = 1, ..., m (9) it 0it 1it it 

CPS CPS where α = ε ,

0it it 

CPS ca ca ca ca α = π + π BP + π HP + ζ (10) 0it 00 01 it 02 it 0it 

A ca ca BPit ca HPit ca α = π +π ( ) +π ( ) +ζ (11) 1it 10 11 TP 12 TP 1it it it 

where BP = the black population 16 years of age and older for the area, HP = the Hispanic population 16 years of age and 

older for the area, TP = the civilian non-institutionalized population 16 years of age and older for the area, all the π s are 

ˆCPS ˆ A (

BP fixed, the superscript, ca, of the π s and ζ s is shorthand for “regression of Y it on Y it ”, and the variables, BP, HP, , ) TP 

ˆ A ˆ 12 A * A A * Note that if Yit is an estimator of employment, then its correct decomposition is Yit = Yit + ( Yit - Yit ) + eit + uit , where Yit is a 

ˆCPS “place-of-work” employment that is different from the place-of-residence employment, Yit , the true-value component of Yit . To adopt 

* A ˆ A A

the “place-of-residence” concept, the difference, Yit - Yit , is added to uit so that Yit = Yit + ε it . Mathiowetz and Duncan (1988) adopt the 

“place-of-work” concept to study response errors in retrospective reports of unemployment. 

13 CPS CPS Here we ignore the autocorrelations of eit implied by the CPS sampling design because with fixed t and varying i, eit exhibits only 

cross-sectional variation. Our purpose in this section is to model this variation in the presence of non-sampling errors. An advantage of 

having cross-sectional data alone on two or more estimators of Yit is that they can provide information about the sum of the sampling and 

ˆCPS ˆCPS non-sampling errors of Yit , whereas time-series data alone on Yit muddle the two errors, with no prospect of estimating even their 

sum. The former result follows from the analysis of this section and the latter result follows from the analysis of Section 3. 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

HP

and ( ) , are called “the coefficient drivers”.14 Additional coefficient drivers may be included in (10) and (11) if they are TP 

thought to be appropriate. 

We estimate equations (9)-(11) under the following assumptions: For i = 1, ..., m and fixed t, 

ca ca (A1) the ( , ) ζ ζ ′ are conditionally and independently distributed with mean vector 0 and constant covariance matrix 0it 1it 

σ ζ

2 Δζ

 given the coefficient drivers. 

CPS A ˆ A CPS A (A2) P(α ,α | Y BP , , HP ,TP = P(α ,α | BP , HP ,TP ) where P( )⋅ is a joint probability distribution function 0it 1it it it it it )

0it 1it it it it 

of the variables in (9)-(11). 

The costs of assumptions (A1) and (A2) are considerably less than those of assumptions underlying the combined model in 

Section 3, as the following discussion shows: 

(i) Design-based inference is strictly inapplicable to situations where uit

CPS ≠ 0 (see Little, 2004, p. 540). Hence model (2) 

that is needed for design-based inference is not employed, since such situations are considered in this section. 

(ii) What we have done in (9) is that we retained the sampling model (1) and avoided all the misspecifications in, and the 

identification problems with, models (2) and (3) by replacing the linking model (3) by its alternative, Yit = α1 

A

itYˆ 

it

A . Not 

much is known about the true model for Yit in footnote 4. Therefore, the correctness of any specified model for Yit may 

be questioned. However, the misspecifications in the time series models for the components of Yit in (3) seemed to us to 

cry out for alternative linking models. A careful comparison of the results produced by model (3) and its alternative, 

A ˆ A

Y = α Y , might expose the weaknesses of both models. it 1it it 

(iii) The conditional independence assumption (A2) is weaker than the unconditional independence assumption, 

CPS A ˆ A CPS A ˆA A A A P(α ,α | Y ) = P( , α α ) , which is false because of the nonzero correlations between ( , , Y e u )′ and α . We 0it 1it it 0it 1it it it it 1it 

include the coefficient drivers in (10) and (11) to avoid this false assumption. 

(iv) Since (9) is a relation between two estimators of employment or unemployment, the range of its intercept is much wider 

than that of its slope. To account for this difference in ranges, its intercept is made to depend on BPit and HPit , which 

BP it it HP

take values over wider intervals than TP and TP , on which its slope is made to depend. 

it it 

ˆCPS (v) Simultaneous estimation of the intercept and slope of (9) permits simultaneous estimation of all the components of Yit

CPS CPS A ˆ A

in (1), since α = ε and Y = α1itYit . In principle, this simultaneous estimation is superior to the two-step 0it it it

CPS CPS estimation of Yit and eit , ignoring uit , as in Section 3. The simultaneous estimation based on (9)-(11) eliminates the 

inconsistencies inherent in the two-step estimation of Section 3. The Chebychev inequality shows that if the variance of 

ca CPS ζ is small, then the distribution of α is tightly concentrated about the right-hand side of (10) with the error 0it 0it 

suppressed (see Lehmann, 1999, p. 52). In this case, the right-hand side of (10) with the error suppressed gives a good 

ˆCPS CPS CPS determination of the sum of sampling and non-sampling errors in Yit . Subtracting a good estimator of eit + uit 

ˆCPS ca from Yit gives a good estimator of Yit .

15 Thus, under a restriction on the variance of ζ 0it , equations (9)-(11) permit 

ˆCPS accurate simultaneous estimation of all the components of Yit . These results reveal the advantages of using equations 

(9)-(11) instead of the combined model in Section 3. 

(vi) Swamy, Mehta and Chang (2006) show that superior estimates of Yit are unlikely to be obtained when the dependent and 

explanatory variables of (9) are interchanged. Therefore, we do not make such an interchange in this paper. 

A ˆ A

For i = 1, ..., m and fixed t, the Y it are assumed to be a random sample from the “superpopulation” defined by Y it = α1it it Y

ˆCPS CPS A ˆ A ˆ A A and are assigned a model conditional probability distribution of (Y - α ) = α Y given Y = y and the coefficient it 0it 1it it it it 

14 Swamy, Mehta and Chang (2006) extend model (9) to include more than two estimators of Y and also to sub-state areas it

for which CPS data are either too sparse or unavailable. 

15 The importance of having a good estimate of employment can be seen from Rowthorn and Glyn ((2006) and Magnac and Visser (1999). 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

CPS CPS A ˆ A A ˆ A

drivers. Assumptions (10) and (11) imply that ε (= α ) and ε = Y - Y = (1 - α )Y follow a bivariate conditional it 0it it it it 1it it

ˆ A A

frequency distributions given Yit = yit and the coefficient drivers, as i varies for fixed t. Let the model means of this 

CPS ˆ A A CPS A ˆ A A BPit HPit A CPS A

distribution be denoted by E m it (ε | Y it = y it it it , BP , HP ) = μit and E m it it it (ε | Y = y , , ) = μ . Then μ and μ TP TP it it it 

ˆCPS ˆ A CPS A ) 

it it 

represent the biases in Y and Y , respectively. The connection between ( μ , μ ′ and the terms in equations (10) and it it it it 

(11) is as follows: 

π

ca ˆCPS CPS ˆCPS = Constant bias in Y that affects all i or the constant term of the bias, μ , in Y ;

00 it it it

ca ca CPS ˆCPS 

π BP + π HP = Variable component of the bias, μ , in Y ;

01 it 02 it it it

ca CPS ζ = Fluctuating component of the error, ε ;

0it it 

ca ca BPit ca HPit A A ˆ A

[(1 −π ) −π ( ) −π ( )]y = Variable component of the bias, μ , in Y ;

10 11 TP 12 TP it it it it it 

ca A A

−ζ 1it yit = Fluctuating component of the error, ε it . 

A ˆ A A Thus, α presents no difficult problems even though it is a non-linear function of Y and ε if assumption (A2) is true. The 1it it it 

ˆCPS above decomposition of the sampling and non-sampling errors of Yit is related to Cochran’s (1977, p. 378) decomposition 

of an error of measurement on a unit based on his model of measurement error. 

16 It should be noted that with fixed t, the time 

ˆCPS ˆ A CPS A

series shocks of Y and Y are fixed and get subsumed into μ and μ , respectively. it it it it 

Inserting (10) and (11) into model (9) gives 

ˆCPS ca ca ca ca Y = (π + π BP + π HP + ζ ) it 00 01 it 02 it 0it 

ca ca BPit ca HP it ca ˆ A + {π +π () +π ( ) +ζ }Y , i = 1, ..., m. (12) 10 11 TP 12 TP 1it it it it 

Note that in this model, the interactions between some of the coefficient drivers and the explanatory variable of model (9) 

ca ca ˆ A

appear as additional explanatory variables and the disturbances, ζ +ζ Y , are heteroscedastic. This means that the 0it 1it it

formulation in (12) is much richer than model (9) with fixed intercept and fixed slope which is misspecified. 17 In (12), crossˆCPS ˆ A 

π

ca πca sectional variation in the pair, ( Y , Y ), is modeled by getting the national intercept, , and slope, ,--which are it it 00 10 

ˆCPS ˆ A

common to all areas--modified by individual area components. In other words, the m pairs, ( Yit , Yit ), are modeled 

π

ca πca simultaneously so that model (12) for each area consists of the common national intercept, , and slope, , and their 00 10 

ˆCPS ˆ A

deviations from the national level. Model (12) recognizes major aspects of the sampling designs yielding Yit and Yit

because the coefficient drivers are selected to explain high proportions of spatial and temporal variations in Y and the sums it

ˆCPS ˆ A

of sampling and non-sampling errors of Yit and Yit . On a definition suggested by a sentence in Hwang and Dempster 

(1999, p. 298) attempting to recognize major aspects of sample design, temporal variation, and spatial variation amounts to 

getting the science right. Model (12) gets the science right in this sense. The misspecifications in (12) with these features can 

be less serious than those in (3). 

ˆ A A Under assumptions (10) and (11), model (12) implies that for i = 1, ..., m and fixed t, given Yit = yit and the coefficient 

ˆCPS CPS A A drivers, the Yit are independently and conditionally distributed with model mean equal to μit - μit + yit and model 

ca ca ˆ A ˆ A A variance equal to V {(ζ +ζ Y ) | Y = y , BP , HP ,TP }. The basic fitting algorithm for (12) is an Iteratively Re-Scaled m 0it 1it it it it it it it 

Generalized Least Squares (IRSGLS) method of Chang, et al. (2000) where in every iteration, the weighting of the sums of 

squares and cross products of observations on the dependent and explanatory variables in (12) by the elements of the inverse 

16 Estimation of this model requires replicated data which we do not have. 

17 The treatment of the coefficients of model (9) as constants is a misspecification because it ignores (i) variations in the sampling and nonˆ A ˆCPS A

sampling errors of Y , (ii) variations in the sampling and non-sampling errors of Y , and (iii) the correlations between α and the it it 1it 

ˆ A A

elements of ( Yit , eit , uit

A )′ . 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

ca ca A of the covariance matrix of the heteroscedastic disturbances, ζ 0it + ζ 1it yit , will be performed. These model-based weights as 

ˆCPS ˆ A well as the design-based weights embedded in Yit and Yit affect the estimates of the coefficients and the error terms of 

ˆCPS ˆ A

(12). This is how the sampling mechanisms generating Yit and Yit are modeled in (12). A model-based approach that 

ignores the sampling mechanism is not valid unless the sampling distribution does not depend on the survey outcomes (see 

Little, 2004, p. 548). 

5.2 Estimation of Model (12) 

The number of the unknown parameters of model (12) is 10: six π s, three distinct elements of Δζ

, and one σ ζ

2 . Given the 

CPS A

cross-sectional data, yit , yit , BPit , HPit , TPit , for i = 1, ..., m and fixed t, an IRSGLS method can be used to obtain good 

approximations to the minimum variance linear unbiased estimators (MVLUEs) of the coefficients and the BLUPs of the 

ˆ

ca ˆ

ca ˆca errors of model (12). These approximations are denoted by ( π , k = 0, 1, h = 0, 1, 2, ζ , k = 0, 1). For k = 0, 1, ζ is an kh kit kit 

ca EBLUP of ζ kit . Using these approximations in place of their population counterparts used in (10) and (11) gives the 

CPS A ˆ CPS 2

estimators of α and α , denoted by α and αˆ A , respectively. The corresponding estimator of σ Δ is denoted by 0it 1it 0it 1it ζ ζ

2 σ ˆ ˆ

ζΔζ

. Appropriate formulas for computing the standard errors of these estimates have been worked out in the Appendix. 

Kariya and Kurata (2004, pp. 42 and 73) prove that the IRSGLS estimators of the coefficients and the EBLUPs of the errors 

of model (12) possess finite second-order moments under very general conditions. To these conditions the condition that m – 

6 ≥ 6 should be added. With this additional condition, two or more degrees of freedom remain unutilized after the estimation 

of the coefficients of (12), σ ζ

2 , and the elements of Δζ

. Sufficient conditions for the consistency and asymptotic normality of 

the IRSGLS estimators of the coefficients of model (12) are given in Cavanagh and Rothenberg (1995). We need 

assumptions (A1) and (A2) because a necessary condition for the consistency of the IRSGLS estimators of the coefficients of 

ca ca ˆ A A

(12) is that the model conditional expectations of ζ 0it and ζ 1it , given Yit = yit and the coefficient drivers are zero. 

5.3 Bias- and Error-Corrected Version of the CPS Estimator 

CPS CPS ˆCPS CPS Since α = ε represents the sum of sampling and non-sampling errors in Y = Y + ε , satisfying equation (10), it 0it it it it it 

follows that the formula 

BECCPS CPS CPS CPS ca ca ca ca Yˆ = Yˆ −αˆ = Yˆ −πˆ −πˆ BP −πˆ HP −ζ

ˆ (13) it it 0it it 00 01 it 02 it 0it 

BECCPS gives “bias- and error-corrected CPS (BECCPS) estimator” of Yit . A value of this estimateor is denoted by yit . 

ˆ

BECCPS The standard error (SE) of Yit is the square root of an approximately unbiased estimator of the model MSE, 

BECCPS 2

( EY ˆ −Y ) . Its derivation is given in the Appendix. Cochran’s (1977, p. 14) computations show that the effect of the m it it 

ˆCPS CPS ˆCPS bias in Yit due to uit on the probability of the error, Yit - Yit , of more than (or less than) 1.96 (or -1.96) times the 

ˆCPS standard deviation of Yit is appreciable if the absolute value of the bias is greater than one tenth of the standard deviation 

ˆCPS ˆ

BECCPS of Yit . In these cases, the bias correction in (13) is desirable if the absolute value of the bias remaining in Yit after the 

ˆ

BECCPS correction is less than one tenth of the standard deviation of Yit . Such desirable bias corrections may frequently occur 

ˆ

BECCPS with Yit . 

Applications. Use the successive difference replication method described in U.S. Department of Labor (2006, Chapter 14) to 

CPS ˆCPS obtain a new estimate of the within-PSU variance contribution to the total design variance, ( ˆ VY | Y ) , of Y in (2). p it it it

BECCPS Apply an IRSGLS method to equations (5) and (6) after inserting this new estimate and yit every place an estimate of the 

within-PSU variance and a CPS estimate of Yit are used in these equations, respectively. This application can give improved 

estimates of design effect and bit used in equation (4). 

5.4 Improved Additional Estimator 
 

 

 

 

 

 

 

 

 

 

 

ˆ A A A ˆ A A ˆ A Using α1it and yit in place of α1it and Yit used in Yit = α1itYit , respectively, gives 

IA A A ca ca BPit ca HPit ca A

Yˆ = αˆ y = (πˆ +πˆ +πˆ +ζ

ˆ

) y (14) it 1it it 10 11 TP 12 TP 1it it it it 

where Yˆ

it

IA denotes an improved additional (IA) estimator of Yit . 

ˆ

IA (

ˆ

IA 2

The SE of Y is the square root of an approximately unbiased estimator of the model MSE, EY −Y ) . Its derivation is it m it it 

given in the Appendix. 

5.5 Comparison of BECCPS and IA Estimators 

ˆ

BECCPS ˆIA A proof of the result, Yit = Yit with probability (w.p.) 1, is given in the Appendix. The following prior information 

ˆCPS ˆ A ˆ

BECCPS helps us assess the relative accuracies of Y , Y , Y , and YˆIA 

.

it it it it

Prior information. (i) The overall CPS sample size is sufficient to produce national-level monthly employment or 

unemployment estimators that satisfy prespecified precision requirements. (ii) However, relatively small CPS state sample 

sizes do not permit the production of reliable monthly employment and unemployment estimates for the states. (iii) The CPS 

and CES survey being the household and establishment surveys provide information about “place-of-residence” and “placeˆ A

of-work” employment, respectively (see U.S. Department of Labor, 1997, p. 45). Consequently, Yit may contain a larger 

ˆCPS magnitude of non-sampling error than Yit when it is viewed as an estimator of “place-of-residence” employment, 

particularly for metropolitan areas like Washington, DC where a good many suburbanites work (see footnote 12). (iv) 

Extrapolated QCEW data generated from a time series model may contain large errors because the best nonlinear models that 

accurately represent asymmetric cycles of state employment are unknown. (v) Both at the national and state levels U.I. 

Claims are very inaccurate estimates of unemployment. 

This prior information implies that at the national level, CPS estimates are more accurate than the corresponding additional 

ˆCPS estimates; that may or may not be the case for state level data. For this reason, Yit being a state level estimator may not be 

ˆ

BECCPS ˆIA ˆ A ˆCPS * more precise than Yit or Yit and Yit may not be more precise than Yit . Let yit be an estimate of Yit based on either 

ˆ

BECCPS ˆIA

Yit or Yit . Then 

* CPS Restriction I. y shall lie within an estimated confidence interval of .

it yit 

Let m be the number of states or areas which geographically exhaust the entire U.S. Then 

Restriction II. the sum ∑

m 

yit 

* shall equal the CPS estimate of the national employment or unemployment, ∑

m Yit .

i=1 i=1 

* CPS The idea of Restriction I is to limit the deviation of y from . This limiting reduces the loss in efficiency due to the it yit 

misspecifications in model (12) for some states without losing the benefits of correctly specified model (12) for other states 

*

(see Jiang and Lahiri, 2006, p. 36). Prior information (i) stated above justifies Restriction II. The y ’s do not satisfy it 

Restrictions I and II if the coefficient drivers in (10) and (11) are inappropriate. Swamy, Mehta and Chang (2006) found that 

ˆ A when Yit based on the American Community Survey was used in (9), the coefficient drivers in (10) and (11) were adequate 

ˆ

BECCPS in the sense that they produced the values of Yit that satisfied Restriction I for all i = 1, 2, ... m. It is a good practice to 

experiment with all possible combinations of coefficient drivers on which data are available. After everything that can be 

done in this way has been done, if some of the resulting yit 

* ’s do not satisfy Restriction I, then this restriction may be 

* * CPS imposed externally. That is, use yit if yit lies within an estimated confidence interval of yit and use an estimated lower or 

CPS *

upper confidence limit of yit whichever is closer to yit otherwise (see Rao, 2003, p. 118). Hopefully, major changes in 

*

these restricted y ’s are not needed to satisfy Restriction II. it 

6. Time Series of Cross-Sectional Estimates 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

BECCPS CPS Estimating model (12) separately for each month denoted by t = 1, ..., τ (= min( τ 1

, ..., τ m )) gives ( yit , a0it )′ , i = 1, 

BECCPS CPS ˆ

BECCPS ˆ ..., 

CPS m, t = 1, ..., τ , where yit and a0it , are the values of Yit and α0it , respectively, obtained from cross-sectional 

BECCPS data in (12) for each t = 1, ..., τ . To produce seasonally adjusted estimates of the Yit , the estimates, yit , t = 1, ..., τ , 

will be adjusted externally by X-12 ARIMA for each i = 1, ..., m. Let σˆ

2 

CPS be an estimate of σ

2 

CPS . Such estimates are 

e ,it e ,it 

computed by the Bureau of Labor Statistics (BLS). An AR (15) model may be fitted to CPS /σˆ CPS a0it , t = 1, ..., τ , for each i = e ,it 

1, ..., m, to estimate the variances and autocovariances of eit

CPS implied by the CPS sampling design. 

7. Example 

ˆ A ˆCES In this section, Yit denotes total employment, Yit = Yit , the CES survey estimator of Yit , t denotes January 2006, i indexes 

the 50 United States and the District of Columbia, and m = 51. Following the IRSGLS method and using the data on the 

variables in (12) with these values of i, t, and m, we obtain the estimates of the coefficients of (12) given below: 

CPS ˆca

y = (0.00000043158 + 0.021745 BP − 0.008542 HP +ζ )

it it it 0it 

(16.99) (0.42495) ( 0.18712) −

+(1.0403 + 0.016003 

BPit 

− 0.005132 

HP it +ζ

ˆ

ca ) yCES (15) 1it it (71.028) (0.53134) TP it it 

( 0.18362) − TP

where yit

CES is the CES survey estimate of Yit , all the other variables are as defined in (12), and the figures in parentheses 

below the coefficient estimates are the t-ratios. 

π

ca Equation (15) shows that the estimate, 0.00000043158, of the common national intercept, 00 , is very close to zero and the 

π

ca 51 CPS CES estimate, 1.0403, of the common national slope, 10 , differs slightly from (∑ yit / ∑

51 

yit ) = 1.0664. Both these 

i=1 i=1 

estimates are significant. The 51 estimates of the errors of (10) and (11) and the coefficients of (9) lie within the ranges: -

ˆ

ca ˆca ˆ CPS ˆ CES 417140 ≤ ζ ≤ 207120, -0.019683 ≤ ζ ≤ 0.0097733, -412470 ≤ α ≤ 207870, and 1.03 ≤ α ≤ 1.1238. The 0it 1it 0it 1it 

ca ca estimates of the variances of ζ and ζ are: 9.185E+08 and 0.00002045, respectively. These estimates show that the 0it 1it 

coefficient drivers, BP and HP , included in (10) explain only a very small proportion of the variation in α

CPS 

, i = 1, ..., it it 0it 

51. We need to include additional coefficient drivers in (10) to reduce the variance of ζ 0 

ca

it to a small number. The capital 

equipment and the ratio of capital to labor for each area in period t can serve as additional coefficient drivers. Factors that 

adjust the “place-of-work” nonfarm employment estimates, yit

CES , to a place-of-residence basis, as in the CPS, can also serve 

as additional coefficient drivers. Unfortunately, data on these variables are not available. 

The CPS estimate, yit

CPS , for each of the 51 areas is given in the column, labeled “CPS”, of Table 1. The values of the 

BECCPS estimator in (13) for the 51 areas implied by the estimates in (15) are given in the column, labeled “BECCPS”, of 

ˆCPS CPS Table 1. They are the estimates of the true value component, Yit , of Yit in (1). The estimate of the error component, eit + 

CPS ˆCPS 

uit , of Yit in (1) for each of the 51 areas is given in the column, labeled “S&NSE”, of Table 1. Note that the estimates in 

ˆCPS the columns, labeled “BECCPS” and “S&NSE”, are the results of simultaneously estimating the components of Yit using 

model (12). The coefficient drivers included in equations (10) and (11) are adequate for 36 areas in the sense that the 

unrestricted BECCPS estimates for these areas statisfy Restriction I. We had to impose Restriction I externally on the 

BECCPS estimates for the remaining 15 areas. 

CPS ˆCPS CPS The two-step estimates of the components, Yit and eit , of Yit for each state when uit is ignored, are given in the 

columns, labeled, “Signal” and “SurE”, of Table 1, respectively. We have now three different estimates of Yit . They are 

denoted by “CPS”, “BECCPS”, and “Signal”, respectively. The CPS estimates are not satisfactory because their accuracy is 

inadequate, as we have already pointed out in Section 2.1. The Signal estimates are also not satisfactory because they are 

implied by the incorrect estimates of the parameters of models (2) and (3), as we have already shown in Section 3. Sufficient 

conditions for the MSE of the BECCPS estimator in (13) to be smaller than the design variance (4) of the CPS estimator, 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ˆCPS Yit , in (1) are given in the Appendix. These conditions can be satisfied in large samples even when an approximately 

unbiased estimate of the MSE of (13) is larger than a sample estimate of the design variance in (4), as we have explained in 

Section 4. If these conditions are not satisfied by model (12) in small samples, then they can be satisfied when we expand the 

set of coefficient drivers in (10). 

ˆCPS The design variance of Yit given in (4) involves unknown quantities and hence is unknown. Its estimate was obtained using

CPS 

y and the BLS estimates of and D CPS in place of , , and D CPS used in (4), respectively. It is given for each of it kit e ,it 

Yit kit e ,it 

the 51 areas in the column, labeled “varcps” of Table 1. An alternative estimate of (4) was obtained using the 

BECCPS CPS estimate, y ( ) it , of Yit and the IRSGLS estimates of δ 01 and ξ 0itg in place of Yˆ

it g , δ 01 , and ξ 0itg used in (5), 

respectively. It is given for each of the 51 areas in the column, labeled “varbeccps” of Table 1. The “varbeccps” estimates can 

ˆCPS be more accurate than the “varcps” estimates because the BECCPS estimator in (13) corrects for bias and error in Yit . 

Each value given in the column, labeled “VarSignal”, of Table 1 is an estimate of the model variance of the two-step 

estimator of Yit for a state. It is an underestimate of the true model variance because the parameters of model (3) are 

estimated, holding the parameters of model (2) fixed at their incorrectly estimated values. 

8. Conclusions 

Simultaneous estimates of the true value and sampling and non-sampling error components of small area sample estimators 

presented in this paper are based on weaker assumptions than their two-step estimates. Misspecifications in linking models 

can result in misestimated components of sample estimators and design variances. Specification errors in the linking models 

used in simultaneous estimation can be negligible compared to those in the linking models used in two-step estimation. 

Appendix 

A. Derivation of the MSE of the BECCPS Estimator 

A BP HP ) it it Π

ca Let xit be the 2-vector, (1, Yˆ 

it ′ , zit be the 5-vector, (1, BPit , HPit , , )′ , be the (2×5) matrix having TP TP 

ca ca ca ca ca ca ca ca ca 

it it 

(π , π , π , 0, 0 ) and (π , 0, 0, π ,π ) as its first and second rows, respectively, and ζ be the 2-vector, (ζ ,ζ )′ ,

00 01 02 10 11 12 it 0it 1it 

where the transpose of a matrix is denoted by a prime. Using these definitions the m equations in (12) are written

ˆCPS π

caLong ca Y = X + D ζ (A1) t zt xt t 

ˆCPS ˆCPS ˆCPS where Yt

 is the m-vector, ( Y1t , ..., Ymt )′ , Xzt is the (m×10) matrix having a Kronecker product between zit ′ and xit ′ , 

18 π

caLong Π

ca denoted by (zit ′ ⊗ xit ′ ) , as its ith row, is the 10-vector given by a column stack of , Dxt is the (m×2m) matrix, 

ca ca ca π

caLong diag[ x′ , ..., x′ ], ζ is the 2m-vector, (ζ ′ , ..., ζ ′)′ , and there are zero restrictions on the elements of . These 1t mt t 1t mt 

π

caLong restrictions can be stated as R = 0, where R is the 4 ×10 matrix of full row rank having ones as its (1, 4)-th, (2, 6)-th, 

(3, 7)-th, and (4, 9)-th elements and zeros elsewhere, and 0 is the 4-vector of zeros. Now a (6×10) matrix C of full row rank 

ca ca can be found such that R C ′ = 0. Under assumptions (10) and (11), ED ( ζ | X = 0 and VD ζ t ) ( | X ) = m xt t z m xt t zt 

Dxt (Im ⊗σ ζ

2 Δζ

)Dxt ′ = Σω

, where I denotes an identity matrix and a subscript is included to indicate its order, and ω is the 

4-vector having σ ζ

2 and the distinct elements of Δζ

 as its elements. Swamy and Tinsley (1980) explain how we go from 

equation (12) to equation (A1). 

Identification in the sense of Lehmann and Casella (1998, p. 24). The coefficient vector, π

caLong , is identifiable if Xzt has 

full column rank. The error vector, ζ t

ca , is unidentifiable because Dxt does not have full column rank. This result implies that

ζ t

ca is not consistently estimable (see Lehmann and Casella, 1998, p. 57). The coefficient drivers are used in (10) and (11) to 

18 The definition of a Kronecker product we use is given by Greene (2003, p. 824), among others. 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

reduce the unidentifiable portions of the coefficients of (9). However, Dxt ζ t

ca is identifiable and its predictor can be used to 

obtain a consistent estimator of ω . 

For known ω , applying the derivation in Greene (2003, p. 100) with appropriate modifications to model (A1) gives the 

π

caLong πcaLong MVLUE of that satisfies the restriction R = 0. This estimator is 

caLong −1 −1 −1 ˆCPS 

ˆ

R ω = C C( Ψω C ) CX 

zt ′ Σω Yt π () ′ ′ (A2) 

ˆ caLong −1 − where the subscript of π R is shorthand for “restricted” and Ψω

 = ( X zt ′ Σω X zt ) 1 . It follows from C.R. Rao (1973, p.77, 

−1 −1 −1 1 caLong Problem 33) that ′( Ψ C ) C = Ψω

 - Ψ R RΨ R′) Ψ

− CC ω

′ ′

ω

( ω RΨω

. Post-multiplying both sides of this equation by ω π

−1 −1 −1 caLong π

caLong caLong gives C C′ ′ ( Ψω C ) CΨω π = , since Rπ = 0. This result is needed to prove that estimator (A2) is unbiased 

with the model covariance matrix 

caLong −1 −1 V (π ω ˆ ( ) | X ) = C ′(CΨω C ′) C . (A3) m R zt 

For known ω , the BLUP of ζ t

ca is 

ˆ

ca 2 −1 ˆCPS ζ (ω ) = (I ⊗σ Δ )D′ Σ M Y (A4) Rt m ζ ζ xt ω ω t

−1 −1 −1 where = I - XC ′(C Ψ C ′) CX ′ Σ . The matrix, Mω M , is idempotent (though not symmetric) with the property that ω m zt ω zt ω

caLong ca ca Mω X C′ = 0. It can be shown that E ( ζ Co π ζ ( )) | V ζ ( ) ) = ˆ

ca ( ω)| X ) = 0, v [( ˆ( ), ω ˆ ω X ] = 0, and (

zt 

ˆ ω | X m Rt zt m R Rt zt m Rt zt 

2 −1 2 (I ⊗σ Δ )D′ Σ M D (I ⊗σ Δ ) because M Σ M ′ = . m ζ ζ xt ω ω xt m ζ ζ ω ω ω MωΣω

ˆ

ca () ˆCPS caLong () ˆ

BECCPS Pre-multiplying both sides of equation (A4) by D gives D ζ ω = Y - X πˆ ω which proves that Y in xt xt Rt t zt R it

ˆ

IA

(13) is equal to Yit in (14) with probability 1 for all i = 1, ..., m and fixed t. For known ω case, the BECCPS estimator of 

ˆ

BECCPS () ˆCPS CPS ˆCPS caLong ˆ

ca Y in (13) becomes Y = −α ( ) = − z′ ′) ˆ ω − j D ( ) ω , where is the m-vector it it ω Yit ˆ

0it ω Yit ( it ⊗l

1 π R () i ′

l ζ Rt j

i 

having 1 as its ith element and zeros elsewhere, l

1 is the 2-vector, (1, 0)′ , and Dl1 

 is the m×2m matrix, (Im ⊗l

1 ′) . 

1 

ˆ

BECCPS The MSE of Y ( ) it ω is 

BECCPS 2 ˆCPS CPS ˆCPS CPS 2 E [{ ˆ Y ( ) ω −Y } | X ] = E [{ Y −αˆ ( ) ω −Y +α } | X ] m it it zt m it 0it it 0it zt 

= g1( ) + ω g 2 ( ) ω (A5) 

where 

2 2 −1 2

g1( ) ω = j Di l m ′ [(I ⊗σ ζ ζ Δ ) − (I m ⊗σ ζ ζ Δ )Dxt ′ Σω D xt m ( I ⊗σ ζ ζ Δ )]D j l i ′ + 1 1 

2 −1 −1 −1 −1 2 [ j D′ (I ⊗σ Δ )D ′ ′ ′ Σ X − 2( z ⊗l )] C′( C Ψ C ′ ′ ) CX Σ D (I ⊗σ Δ ) D j ′ (A6) i l m ζ ζ xt ω zt it 1 ω zt ω xt m ζ ζ l i 1 1 

and 

−1 −1

g ω ′ ′ ′ ′ Cz 2( ) = (zit ⊗l

1

)C (CΨω C ) ( it ⊗l

1 ) . (A7) 

caLong As in Rao (2003, p. 99, (6.2.11)), the second term in (A5) arises as a direct consequence of using πˆ ( ) R ω instead of 

π

caLong 2

 in (A4). In the cases where the coefficient drivers in (10) and (11) reduce the magnitudes of the elements of σ ζΔζ

to 

small values, the first term of g1

( ) ω can be much smaller than Rao’s (2003, pp. 99 and 117) g δ 1 ( ) which, in turn, is smaller 

CPS CPS than the design variance of eit when uit is absent. Also, g1

( ) ω is smaller than its first term if its second term is negative. 

Now we can elaborate on our discussion in Section 4. Rao (2003, p. 117) proved that the first term in (A5) is smaller than the 

CPS CPS design variance of eit if (i) the non-sampling error, uit , is zero with probability 1 for all i and t, (ii) the sampling errors, 

eit

CPS 

, i = 1, ..., m, are independently distributed with known design variances, (iii) for i = 1, ..., m and fixed t, the Yit follow 

a GLM model of Rao’s (2003, p. 116) type with no constraints on its coefficients, (iv) the errors of the GLM model (or the 

vi ’s in Rao’s notation) are identically and independently distributed with known model variance, (v) the GLM model error is 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

CPS ˆ

BECCPS independent of eit for all i and t, and (vi) the estimator, Yit , is replaced by the BLUP of Yit given by the GLM model. 

Even when these conditions hold, the sum of the two terms in (A5) may not be smaller than the design variance of eit

CPS 

unless m is sufficiently large and the regularity conditions given in Rao (2003, p. 117) are satisfied. In any case, ignoring the 

non-sampling error, uit

CPS , can result in an inconsistent and inefficient predictor of Yit . 

We now assume that the error term of model (A1) is normally distributed. 

ca CPS 2 − CPS E (ζ | M Y ′

ˆ ) = (I ⊗σ Δ )D M ′ (M Σ M ′ ) M Yˆ (A8) m t ω t m ζ ζ xt ω ω ω ω ω t

where (M Σ M ′ )

−

 is a generalized inverse of M Σ M ′ and this generalized inverse is defined as in C.R. Rao (1973, p. ω ω ω ω ω ω

24). Swamy and Mehta (1975, p. 596) proved that the right-hand side of equation (A8) is equal to the BLUP in (A4). Thus, 

ca ca ˆCPS ca when ζ is normal, its BLUP is equal to its best unbiased predictor (BUP), since E (ζ | M Y ) is the BUP of ζ .

t m t ω t t 

Let A be a × ( − 6) matrix of full column rank such that m m A zt ′X C′ = 0. Then 

ca CPS 2 −1 CPS E ( ′

ˆ D A A ′ A ′ˆ ζ | A Y ) = (I ⊗σ Δ ) ′ ( Σ ) A Y . (A9) m t t m ζ ζ xt ω t 

−1 −1 −1 −1 −1 −1

It follows from C.R. Rao (1973, p. 77, Problem 33) that A A( ′Σω A) A′ + Σω X ztC ′(CX zt ′ Σω X ztC ′) CX zt ′ Σω = Σω

. 

Inserting this identity into (A9) shows that (A9) is equal to the BLUP of ζ t

ca . This derivation extends Jiang’s (1997) proof to 

the cases where the coefficients of model (A1) are subject to linear restrictions. 

ˆCPS 2 2 We now turn to the case where ω ˆ

 is unknown. We use A′Y to estimate σ Δ so that our estimator, denoted by σˆ

ζΔ ,

t ζ ζ ζ

2 ˆCPS of σ Δ is a function of A′Y . Let ωˆ be the 4-vector having σˆ

ζ

2

 and the distinct elements of Δˆ

ζ

 as its elements. The ζ ζ t 

ˆ

BECCPS ˆBECCPS () ˆCPS CPS () ˆCPS caLong BECCPS estimator, Y (

it , of Yit in (13) can be written as Yit ωˆ = Yit - αˆ

0it ωˆ = Yit - (zit ′ ⊗l

1 ′)πˆ

R ωˆ) -

ˆ ˆ

BECCPS ( ) ca j D′ ζ () ωˆ . Sufficient conditions for E Y[ ωˆ −Y ] = 0 are given in Kariya and Kurata (2004, pp. 42 and 73). The i l Rt m it it 

ˆ

BECCPS () 

1 

MSE of Yit ωˆ is 

BECCPS 2 BECCPS 2 E Y([ ˆ ( ) ([ ωˆ −Y ] | X ) = EY ( ) Y ] | X ) m it it zt m it it zt 

ˆ ω −

ˆ

BECCPS ˆBECCPS 2 + EY ([ ( ) ωˆ −Y ( )] | ω X ) m it it zt 

ˆ

BECCPS ˆ BECCPS ˆBECCPS + 2 EY ([ ( ) ω −Y ][ Y ( ) ωˆ −Y ( )] ω | X ) . (A10) m it it it it zt 

In (A5), we have already evaluated the first term on the right-hand side of this equation. To show that the third term on the 

ˆ

BECCPS () ˆ

BECCPS () CPS () CPS right-hand side of this equation vanishes, we first note that Y ωˆ - Y ω = -αˆ ωˆ +αˆ ( ) ω = it it 0it 0it 

caLong ( ) ˆ ca () caLong ˆ ca () −1 −1 −1

-(z ′ ⊗l ′)πˆ ωˆ - j D′ ζ ωˆ + (z′ ⊗l′)πˆ (ω) + j D′ ζ ω = (z ′ ⊗l ′){ C ′(CΨ C ′) CX ′ Σ - it 1 R i l Rt 1

it 1 R i l Rt 1 it 1 ω zt ω

−1 −1 −1 ˆCPS −1 −1 − caLong 2 −1 CPS caLong C C( Ψ C ′) CX Σ }[Y - XC (C Ψ C ′ CΨ1

πˆ ( ) ] + jD {(I σ Δ )D′ [ ′ ′ ′ ) ω ′ ⊗ Σ Y ˆ − X πˆ ( )] ω - ωˆ zt ωˆ t zt ω ω R i l m 1 ζ ζ xt ω t zt R

2 −1 CPS caLong ˆCPS ˆCPS −1 −1 − caLong ⊗ ˆ Δ (

ˆ

)D′ Σ [ ˆ − X ˆ (ω)]} is a function of A′Y because ωˆ , - XC (C Ψ C ′ Ψ

1

(I σ Y π ˆ Y ′ ) C πˆ ω) , m ζ ζ xt ωˆ t zt R t t zt ω ω R 

ˆCPS caLong () ˆCPS caLong () ˆCPS BECCPS Y - X πˆ ω , and Y - X πˆ ωˆ are all functions of A′Y . Furthermore, the equation, Y ( ) ˆ ω −Y = t zt R t zt R t it it 

caLong caLong ca ca

− π ω [

ˆ (z′ ⊗l ′)[ ˆ ( ) −π ] − j D′ ζ () ω −ζ ] , is such that the first term on its right-hand side is independent of it 1 R i l Rt t 1 

ˆCPS A′Y because of the condition that A′X C′ = 0, and the last term on its right-hand side can be shown to be equal to t zt 

ca CPS ca jD ′ [E ( | ζ AY ′

ˆ

) −ζ ] using the result in (A9). Hence, the third term on the right-hand side of equation (A10) is i l m t t t 1 

BECCPS BECCPS caLong caLong CPS equal to 2 E Yˆ ( ) −Yˆ ω E [ ( z ⊗l ){ ˆ ( ) −π }| AY ˆ

m it it m it ([ ωˆ ( )] − ′ ′1 π ω R t zt ′ , X ]) +

ˆ BECCPS ˆ BECCPS ca ˆCPS ca ˆCPS 2 E ([ ( ) Y ωˆ −Y ( )][ ω − jD ′ (E {E ( | ζ AY ′ ) −ζ } | AY ′ , X )]) which vanishes. This proof extends the m it it i l m m t t t t zt 1 

proofs of Swamy and Mehta (1969), Jiang (2001), and Rao (2003, p. 114) to the cases where the coefficients of model (A1) 

are subject to linear restrictions. 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ˆ

BECCPS Because of the second term on the right-hand side of equation (A10), the MSE of Y ( ) it ωˆ is always larger than that of 

BECCPS () BECCPS () BECCPS Y ( ) ˆ ω in the normal case. The method of approximating the MSE of Yˆ ωˆ by the MSE of Yˆ ω could, it it it

therefore, lead to serious underestimation. 

Unfortunately, the exact evaluation of the second term on the right-hand side of equation (A10) is generally not possible 

except in some special cases, as Rao (2003, p. 103) has pointed out. It is therefore necessary to find an approximation to this 

term. We begin the derivation of such an approximation with the assumptions (see, e.g., Lehmann and Casella, 1998, p. 430) 

CPS () CPS that permit an expansion of αˆ ωˆ about αˆ ( ) ω with bounded coefficients. Using a Taylor approximation, we obtain 0it 0it 

CPS CPS αˆ ˆ () ω α− ˆ ( ) ω ≈ d ( )( ωωω ′ ˆ − ) (A11) 0it 0it 

where d( ) ω = ∂α

CPS 

ˆ ( )/ ˆ − are of lower order relative 0it ω ∂ω and it is assumed that the terms involving higher powers of ω ω

2 −1 CPS () caLong caLong to d ( )( ω ′ ω ω ˆ − ) . Let b1 ′ = jD i l m ′ (I ⊗σ ζ ζ Δ )D xt ′ Σω

. Then αˆ 0it it ω = [(z l ′ ⊗ − 1 1 ′ ′ ) bX zt R ][ ( ) πˆ ω π− ] + 1 

caLong CPS caLong (zit ′ ⊗l

1 ′)π + b Y 1 ′(

ˆ

t − X zt π ) . Under normality, 

ˆCPS caLong d( ) ω ≈ ( ∂ ∂ b′ / ω)( Y − X π ) = d *

( ) ω , (A12) 1 t zt 

caLong () caLong since the terms involving the derivatives of πˆ

R ω - π with respect to ω are of lower order. Therefore, 

2 * 2 * * 

m

[ ( ) ( ′ ˆ− )] ≈ E d m ω ′ ω ω− )] ≈ tr[ E d Ed ω ω ω [ ( )( ˆ 

m

( ( ) ( ) ) ( )] ω d ω ′ V ωˆ

∂ Σ ∂ ′ ( ) ] g3 = tr[( ∂b′ / ω) ( b / ∂ω)′V ωˆ = ( ) ω (A13) 1 ω 1 

where V ( ) ωˆ is the asymptotic covariance matrix of ωˆ , and the neglected terms are of lower order. It now follows from 

(A11)-(A13) that 

ˆ

BECCPS ˆ BECCPS 2 CPS CPS 2 EY ([ ( ) ωˆ −Y ( )] ω | X ) = E ([αˆ ( ) ωˆ −αˆ ( )] | ω ω X ) ≈ g ( ) . (A14) m it it zt m 0it 0it zt 3 

BECCPS Inserting (A5) and (A14) into (A10) gives a second-order approximation to the MSE of Y ( ) ˆ 

it ωˆ as 

BECCPS 2

([ ˆ EY ( ) ωˆ −Y ] | X ) ≈ g ( ) ω + g ( ) ω + g ( ) ω (A15) m it it zt 1 2 3 

where the terms, g2( ) g3 ω and ( ) ω , arise as a direct consequence of using the estimators rather than the true values of

caLong BECCPS () BECCPS π and ω , respectively, in Y ( ) ˆ ωˆ and are of lower order than g ( ) ω . The estimator Yˆ ωˆ cannot be rejected it 1 it

ˆCPS CPS in favor of Yit even when the design variance of eit is smaller than the MSE in A(15). The reason is that this design 

ˆCPS CPS variance understates the MSE of Yit by ignoring uit . 

B. Estimation of the MSE of the BECCPS Estimator 

It follows from Rao (2003, p. 104) that Eg ωˆ ≈ g ω and Eg ( ) 2

() ( ) ωˆ ≈ g ω to the desired order of approximation, 3( ) m 2 m 3 

but g1( ) g1( ) g2( ) g3 ωˆ is usually a biased estimator of ω with a bias that is generally of the same order as ω and ( ) ω . To 

evaluate this bias, we make all the assumptions that permit a Taylor expansion of g 1 ( ) ωˆ about g ω with bounded 1( )

coefficients (see Lehmann and Casella (1998, p. 430)). Under these assumptions, 

1 2

g ( ) ωˆ = g ( ) ω + (ωˆ − ∇ ω) ′ g (ω) + (ωˆ − ∇ ω)′ g ( ω)( ωˆ −ω) (A16) 1 1 1 1

2 

where g1

() is the vector of first-order derivatives of 1 ω with respect to ω and 2 ∇ ω g ( ) ∇ g 1

( ) ω is the matrix of second-order 

derivatives of g1( ) ω with respect to ω . The estimator ωˆ is generally a biased estimator of ω and hence the model 

expectation of the second term on the right-hand side of equation (A16) is generally nonzero. Consequently, 

1 2 Eg ( ) ωˆ ≈ g () ω + E (ωˆ −ω) ′∇g () ω + tr[ ∇ g ( ) ( )] ω ω V ˆ . (A17) m 1 1 m 1 1 2 

If Σω

 has a linear structure, then (A17) reduces to 

E g () () ωˆ ≈ g ω + E (ωˆ − ∇ ω)′ g () ω − g ( ) ω . (A18) m 1 1 m 1 3 

ˆ

BECCPS This result shows that an estimator of the MSE of Y ( ) it ωˆ to the desired order of approximation is given by 
 

 

 

 

 

 

 

 

 

 

 

 

g ( ) ωˆ ˆ − estimate of [ E ( ω ω− ) ′∇g ( )] ω + g () ˆ ˆ + 2g ( ) 2 ω ω . (A19) 1 m 1 3 

BECCPS The model expectation of (A19) is approximately equal to the MSE of Y ( ) ˆ 

it ωˆ . The second term in (A19) can be ignored 

if it is of lower order than - g3( ) ω .

C. Derivation of the MSE of the IA Estimator 

caLong A

Yˆ

IA () ′ ˆ j D ˆ

ca When ω is known, the IA estimator of Y () ,

it in (14) can be written as it ω = [(zit ⊗l

2 ′)π R (ω) + i ′

l ζ Rt ω ] yit 2 

where l

2 is the 2-vector, (0, 1)′ , ji is the m-vector having 1 as its ith element and zeros elsewhere, and Dl2 

 is the ( m × 2m ) 

matrix, (I ⊗l′) . The MSE of Yˆ

IA ( ) ω is m 2 it

IA 2 A A 2 A 2

[{ E Y ˆ ( ) ω −Y } | X ] = E [{α ω ˆ ( ) −α } ( y ) | X ] m it it zt m 1it 1it it zt 

= f

1

( ) f

2 ω + ( ) ω (A20) 

where 

2 2 −1 2 A 2

f ( ) ω = jD ′ {(I ⊗ Δ− σ ) (I ⊗ Δ σ )D ′ Σ D (I ⊗ Δ σ )} D j ′ ( y ) +{ j D′

1 i l m ζ ζ m ζ ζ xt ω xt m ζ ζ l i it i l 2 2 2 

2 −1 −1 −1 −1 2 A

×(I ′ C′( ′ ⊗ Δ D j ( m ⊗σ ζΔζ

)Dxt ′ Σω X 

zt − 2( zit ⊗l2 ′)} CΨω C′) CX zt Σω Dxt (Im σ ζ ζ

) l ′

i yit )

2 (A21)

2 

and 

−1 −1 A 2 ω = ′ ′ ′

ω

′

it 2

f .

2

( ) (zit ⊗l

2

)C (CΨ C ) Cz ( ⊗l )( yit ) (A22) 

IA A second-order approximation to the MSE of Y ( ) ˆ 

it ωˆ is 

ˆ

IA 2 EY ([ ( ) ωˆ −Y ] | X ) ≈ f

1

( ) ω + f

2

( ) ω + f

3 ( ) ω (A23) m it it zt 

2 2 −1 where f ( ) ω = tr[ (∂b ′ / ) ( ∂ Σ∂b′ / ∂ω) ′V ( ωˆ )]( y with = jD ( σ ) ′

A ) b′ ′ I ⊗ Δ Σ .

3 2 ω ω 2 it 2 i l m 2 ζ ζ Dxt ω

D. Estimation of the MSE of the IA Estimator 

ˆ

IA An estimator of the MSE of Y ( ) it ωˆ to the desired order of approximation is given by 

1

ˆ − estimate of [ Em ( ˆ − ∇f

1

( )] ω + f

2

( ) + 2 f

3

f ( ) ω ω ω)′ ωˆ ˆ ( ) ω . (A24) 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

References 

Brown, L. D. (1990). An Ancillarity Paradox Which Appears in Multiple Linear Regression (with discussion). The Annals of

 Statistics, 18, 471-538. 

Cavanagh, C. L. and Rothenberg, T.J. (1995). Generalized Least Squares with Nonnormal Errors, Advances in Econometrics 

 and Quantitative Economics, G.S. Maddala, P.C.B. Phillips and T.N. Srinivasan (eds). Oxford, UK: Blackwell. 

Chang, I., Swamy, P.A.V.B., Hallahan, C. and Tavlas, G.S. (2000). A Computational Approach to Finding Causal Economic

 Laws. Computational Economics, 16, 105-136. 

Cochran, W.G. (1977). Sampling Techniques, 3rd edition. New York: John Wiley & Sons. 

Durbin, J. and Koopman, S.J. (2001). Time Series Analysis by State Space Methods. Oxford: Oxford University Press. 

Fay, R. and Train, G. (1995). Aspects of Survey and Model-Based Postcensal Estimation of Income and Poverty

 Characteristics for States and Counties. Proceedings of the Section on Government Statistics, American Statistical 

Association. 

Freedman, D.A. and Navidi, W.C. (1986). Regression Models for Adjusting the 1980 Census (with discussion). Statistical

 Science, 1, 3-39. 

Greene, W.H. (2003). Econometric Analysis, 5th edition. Upper Saddle River, New Jersey: Prentice Hall. 

Harvey, A.C. (1989). Forecasting, Structural Time Series Models and the Kalman Filter. Cambridge: Cambridge University 

 Press. 

Havenner, A. and Swamy, P.A.V.B. (1981). A Random Coefficient Approach to Seasonal Adjustment of Economic Time

 Series. Journal of Econometrics, 15, 177-209. 

Hwang, J.S. and Dempster, A.P. (1999). A Stochastic System for Modeling Labor Force Series of Small Areas. Statistica 

 Sinica, 9, 297-324. 

Jiang, J. (1997). A Derivation of BLUP – Best Linear Unbiased Predictor. Statistics & Probability Letters, 32, 321-324. 

Jiang, J. (2001). Private Communication to J.N.K. Rao. 

Jiang, J. and Lahiri, P. (2006). Mixed Model Prediction and Small Area Estimation (with discussion). Sociedad de Estadistica 

 e Investigacion Operativa Test, 15, 1-96. 

Kariya, T. and Kurata, H. (2004). Generalized Least Squares. Hoboken, New Jersey: John Wiley & Sons, Ltd. 

Lehmann, E.L. (1999). Elements of Large-Sample Theory, New York: Springer. 

Lehmann, E.L. and Casella, G. (1998). Theory of Point Estimation. 2nd edition. New York: Springer. 

Lent, J. (1991). Variance Estimation for Current Population Survey Small Area Labor Force Estimates. Proceedings of the 

 Section on Survey Research Methods, American Statistical Association. 

Lent, J. (1994). Memorandum for Michael Welch, Bureau of Labor Statistics, Washington, DC. 

Little, R.J. (2004). To Model or Not To Model? Competing Modes of Inference for Finite Population Sampling. Journal of

 the American Statistical Association, 99, 546-556. 

Maddala, G.S. and Kim, I. (1998). Unit Roots, Cointegration, and Structural Change. Cambridge: Cambridge University

 Press. 
 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Magnac, T. and Visser, M. (1999). Transition Models with Measurement Errors. The Review of Economics and Statistics, 

 81, 466-474. 

Mathiowetz, N.A. and Duncan, G.J. (1988). Out of Work, Out of Mind: Response Errors in Retrospective Reports of 

Unemployment. Journal of Business & Economic Statistics, 6, 221-229. 

Memorandom from Khandaker Mansur. “Documentation of the CPS History of Design Effects, Sampling Intervals, and Type 

A Rates. Procedures and Sources Used to Obtain These Estimates (Var2000-3),” dated February 22, 2007. 

Montgomery, A.L., Zarnowitz, V., Tsay, R.S. and Tiao, G.C. (1998). Forecasting the U.S. Unemployment Rate. Journal of 

 the American Statistical Association, 93, 478-493. 

Pfeffermann, P. and Tiller, R. (2006). Small-Area Estimation with State-Space Models Subject to Benchmark Constraints. 

 Journal of the American Statistical Association, 101, 1387-1397. 

Rao, C.R. (1973). Linear Statistical Inferences and Its Applications, 2nd edition. New York: John Wiley & Sons. 

Rao, J.N.K. (2003). Small Area Models. Hoboken, New Jersey: John Wiley & Sons, Inc. 

Rowthorn, R. and Glyn, A.J. (2006). Convergence and Stability in U.S. Emplyment Rates. Contributions to Macroeconomics, 

6, 1-42. 

Swamy, P.A.V.B. and Mehta, J.S. (1969). On Theil’s Mixed Regression Estimator. Journal of the American Statistical

 Association, 64, 273-276. 

Swamy, P.A.V.B. and Mehta, J.S. (1975). Bayesian and Non-Bayesian Analysis of Switching Regressions and of Random

 Coefficient Regression Models. Journal of the American Statistical Association, 70, 593-602. 

Swamy, P.A.V.B. and Tinsley, P.A. (1980). Linear Prediction and Estimation Methods for Regression Models with 

 Stationary Stochastic Coefficients. Journal of Econometrics, 12, 103-142. 

Swamy, P.A.V.B., Zimmerman, T.S., Mehta, J.S. and Robison, E.L. (2005). New Improved Small Area Models. Federal

 Committee on Statistical Methodology Research Conference (Conference Papers). 

Swamy, P.A.V.B., Yaghi, W., Mehta, J.S. and Chang, I. (2007). Empirical Best Linear Unbiased Prediction in Misspecified

 and Improved Panel Data Models with an Application to Gasoline Demand. Computational Statistics & Data Analysis, 51, 

3381-3392. 

Swamy, P.A.V.B., Mehta, J.S. and Chang, I. (2006). An Efficient Method of Estimating the True Value of a Population 

 Characteristic from its Discrepant Estimates. Unpublished manuscript. 

Tiller, R.B. (1992). Time Series Modeling of Sample Survey Data from the U.S. Current Population Survey. Journal of 

Official Statistics, 8, 149-166. 

Tiller, R.B. (2005). State Model-Based Labor Force Estimation. Unpublished manuscript. 

Wolter, K.M. (1985). Introduction to Variance Estimation. New York: Springer. 

Zimmerman, T.S. and Robison, E.L. (1995). Estimation of Autocorrelations for Current Population Survey Labor Force 

 Characteristics. Proceedings of the Section on Survey Research Methods, American Statistical Association, 414-419. 

Zimmerman, T.S. (2005). Reliability of CPS Estimators Under Current (2000) Design. Memorandum. 

Zimmerman, T.S. (2007). Sampling Error Correlation Status Report. Memorandum. 
 

 

 

 

U.S. Department of Labor, Bureau of Labor Statistics (1997). BLS Handbook of Methods, Bulletin 2490. Washington, DC: 

Government Printing Office. 

U.S. Department of Labor, Bureau of Labor Statistics (2006). Design and Methodology: Current Population Survey, 

 Technical Paper 66. Washington, DC: Government Printing Office. 
 

Table 1. Two-Step versus Simultaneous Estimates of True Value and Error Components

of Sample Estimates of State Employment for January 2006 

CPS Sample Estimates Simultaneous Estimates Two-Step Estimates

St CPS varcps BECCPS S&NSE varbeccps Signal VarSignal SurE 

AL 2038673 3679861578 2022994 15679 3679900196 2064406 1677393936 -25733 

AK 308055 95009042 304142 3913 95009433 312698 34951744 -4643 

AZ 2719607 4940016080 2710402 9206 4940037094 2742883 1872639076 -23276 

AR 1335773 931402110 1275956 59817 931377568 1279302 386712225 56471 

CA 16568973 19713973309 16565410 3563 19713785649 16815860 19757394721 -246887 

CO 2417786 2594841766 2339775 78011 2594846962 2428008 1521468036 -10222 

CT 1728484 1184483549 1722675 5809 1184456506 1709671 577200625 18813 

DE 430941 71573417 441423 -10482 71573611 416387 27405225 14554 

DC 274617 37453913 286613 -11996 37453577 292316 25030009 -17699 

FL 8524044 15431076249 8506356 17688 15431008706 8415622 7395140025 108422 

GA 4418295 5567633630 4272046 146249 5567646053 4399394 2121523600 18901 

HI 637847 136983872 627152 10695 136983835 626053 78021889 11794 

ID 723132 314448720 688376 34756 314459782 701278 118766404 21854 

IL 6102978 8049548107 6106455 -3477 8049551901 6074746 4563948249 28232 

IN 3116148 4327930599 3050431 65717 4327891224 3051273 2294601604 64875 

IA 1588221 811523654 1532386 55835 811511184 1553633 377408329 34588 

KS 1400258 822429885 1375147 25111 822416077 1386284 295083684 13974 

KY 1811409 2988915564 1881199 -69790 2988889434 1856101 896403600 -44692 

LA 1862668 2733058634 1858502 4166 2733026297 1869696 1643572681 -7028 

ME 666650 186489439 639884 26766 186490357 660612 88717561 6038 

MD 2844012 2088972084 2754430 89582 2088966050 2817955 929579121 26057 

MA 3117108 5040669842 3256264 -139156 5040668540 3182095 1567684836 -64987 

MI 4664360 7475499516 4501730 162630 7475551167 4663248 3227944225 1112 

MN 2714723 2432139574 2810682 -95959 2432132934 2765644 1109889225 -50921 

MS 1205433 1491094340 1154127 51306 1491091754 1195684 725709721 9749 

MO 2858183 4227831478 2827705 30478 4227860017 2831368 2008832400 26815 

MT 470794 169403823 445284 25510 169404383 458981 45131524 11813 

NE 912968 256986524 944388 -31420 256987727 925837 120253156 -12869 

NV 1172640 694617481 1224297 -51657 694624454 1176211 238702500 -3571 

NH 712864 131084819 690423 22441 131084294 703637 54730404 9227 

NJ 4186935 4993905767 4199985 -13050 4993927587 4226153 2227462416 -39218 

NM 890814 616281651 854395 36419 616295929 877161 110250000 13653 

NY 9051645 16731126410 9024711 26934 16731098425 8935803 10327437376 115842 

NC 4125019 6143432242 4111840 13179 6143465977 4117926 2625537600 7093 

ND 350767 70109708 355857 -5090 70109997 332271 32821441 18496 

OH 5464068 10148532036 5561621 -97553 10148624788 5447573 2539353664 16495 

OK 1647282 1897191571 1585666 61616 1897214340 1628011 651985156 19271 

OR 1747033 1377585050 1732349 14684 1377595302 1737507 462121009 9526 

PA 5905700 6871837932 5902957 2743 6871921628 5879965 3680484889 25735 

RI 543188 89395889 524657 18531 89395489 540791 54405376 2397 

SC 1932052 2196475459 1931836 216 2196491622 1922926 548074921 9126 

SD 403813 70320546 400649 3164 70320508 401614 26081449 2199 

TN 2799013 4799435890 2840960 -41947 4799471524 2727933 1528731801 71080 

TX 10675536 22253332600 10753712 -78176 22253467249 10647979 7967704644 27557 

UT 1213992 810295800 1216001 -2010 810306676 1214375 356643225 -383 

VT 353711 40295008 341269 12442 40295305 340964 18198756 12747 

VA 3737584 5901222652 3816349 -78765 5901175918 3793403 1311453796 -55819 

WA 3143887 2853273738 3039191 104696 2853224390 3110833 1535307489 33054 

WV 777019 700320633 769822 7197 700338137 747780 268304400 29239 

WI 2914222 3758736287 2915293 -1072 3758743232 2845595 1782106225 68627 

WY 270567 42070126 273908 -3341 42069787 265398 15856324 5169 

Tot 141481487 140969680 511806 141118844 362643 

-511806 -362643 

CPS ˆ ˆCPS | BECCPS Note: CPS = yit , varcps = p 

(

it Yit ) y VY in Section 3.1, BECCPS = in Section 5.3, S&NSE = it

CPS 2 Estimate of ε it in Section 5.2, varbeccps = σˆ CPS in Section 7, Signal = Two-step estimate of e ,it 

Yit , VarSignal = estimated model variance of Signal, SurE = Two-step estimate of eit

CPS .