PREDICTING BIRTHS IN THE CURRENT EMPLOYMENT STATISTICS SURVEY

Steve Woodruff, Bureau of Labor Statistics

2 Massachusetts Ave. N.E., Suite 4985, Washington D.C. 20212

KEY WORDS: Birth and Death Process, Maximum

Likelihood Estimation

1) Introduction

The Current Employment Statistics survey is a

monthly survey of about 400,000 business

establishments. It is used to estimate total national

non-farm employment along with several other less

important economic statistics. The universe for this

survey (and sample) is constantly changing due to

establishments going out of business and others

starting up. The effect of these “births” and

“deaths” remains an important source of nonsampling error.

This problem of births and deaths was studied

by, Lent and Grezsiak , (1985). The operational

conclusions of their work was that births and deaths

tend to cancel and that the effect of these on

estimation was minimal. The sample automatically

picks up deaths, they report zero or just stop

reporting although, in this case, nonrespondents

must still be separated from deaths. Births can go

unnoticed for extended periods of time (up to a

year) between sample replenishments. With

implementation of random sampling and possibly

rotation sampling, the effect of births and deaths

can be measured and these measurements used to

improve the accuracy of the estimates.

The strata in this survey are industry/size-class

groups. Employment in these strata may be

shrinking or growing at different rates and in

different directions. By observing births and deaths

in the immediate past for each stratum, it is possible

to estimate the “infinitesimal” birth and death rates

for a Markov process that

describes this birth and death behavior. This

Markov model is used to predict the number of

births for the current estimation period.

Modeling the number of establishments in the

sampling frame, known in the Bureau as the

Universe Data Base (UDB), as the realization of a

birth and death process is a useful means of

predicting population size when actual numbers of

births is unknown.

ARIMA models were also fitted to the data on

these Birth/Death processes. Autoregressive

models within ARIMA seemed to provide the best

fit of the various Time Series techniques. These

models were used to forecast the number of

reporting units in an estimation cell and these

forecasts compared to the birth/death forecasts.

Any views expressed here are those of the

author and not the Bureau of Labor Statistics.

2) A Markov Birth and Death Process

For a sampling stratum, let X(t) be the total

number of establishments in that stratum for month

t. X(t) can be thought of as a stochastic process as

each month, business establishments migrate to

other strata or go out of business and others start

up or migrate into the stratum. Karlin (1973, Ch

7), outlines several “Birth and Death Processes”.

One possible candidate and a method of finding it’s

moments is outlined next.

Let X(t) be a stationary Markov Process:

P(X(t+h)εA| X(t)εB) = P(X(h)εA| X(0)εB) ∀ t,

and let P h P X h j P i

i, j

( ) = ( ( ) = | (0) = ) . Let the

transition probabilities satisfy the following

conditions where l m i iand are positive numbers:

P

i i

h

i

h o h

as h i

P

i i

h

i

h o h as h i

P

i i

h

i i

h o h as h i

P

i j

ij

,

( ) ( )

,

,

( ) ( ) ,

,

( ) ( ) ( ) ,

,

( )

+

= + =

↓ ≥

−

= + ↓ ≥

= − + + ↓ ≥

=

1

0 0

1

0 1

1 0 0

0

l

m

l m

d

P(X(h) = i +1|X(0) = i)

where δ is the kronecker-δ (*)

Let M(t) = E[X(t)], be the expected value of

X(t) under the stochastic structure just described

and V(t) be the variance of X(t), (both assuming

X(0)=M(0) is known). For a stochastic process
called Linear Growth with Immigration, the

infinitesimal birth rate when the population contains

n members is: l l n= n + a and the infinitesimal

death rate is: m m n= n where a>0, λ>0, and μ>0.

M(t) and V(t) are derived by solving the “forward

Kolmogorov differential equations” which are

themselves derived from (*), the prime denotes

derivative with respect to t; these equations are:

Pit

0

'

( ) = − l + m

0 0 1 1

P

i

t P

i

( ) (t)

P

ij tj

P

i j

t

j j

P

ij t

j

P

i j

t

'

( )

,

( ) ( ) ( )

,

( )

=

− −

− + +

+ +

l l m

m

1 1

1 1

 

for j ≥ 1

When X(0) = i then be definition of expected value,

M(t) is:

M t jPij j( ) = (t)

=

∞

∑

1

.

Differentiate both sides of this equation and

substitute the forward Kolmogorov differential

equations to get the first order linear differential

equation in M(t):

M’(t)=a +(λ -μ)M(t), M(0)=X(0)=i .

The solution of this is:

M(t)= 

a

e M e

t t

l m

l m l m

−

− +

− −

{ } ( )

( ) ( )

1 0

for λ ≠ μ and

M(t) = at + M(0) for λ =μ

Using the same technique, E(X(t)2) can be

derived and this used to find V(t). The result

is:

V t M a e

a M e a

t

t

( )

( )

{( ( )( ) )

[( ) ( )( )] }

( )

( )

=

−

− +

− + + − +

−

−

1

0

0

2

2 2 2

2 2

l m

l m l

l m l m m

l m

l m

for λ ≠ μ and

V(t) = aμt + M( )μ t + at for λ = μ

2

2 0

M(t) gives the expected value of X(t) given

M(0), λ , μ, and a. The mean, M(t), will be used to

predict X(t) given X(t-m) is known, m>0. In this

case the prediction for X(t) is M(m) given

M(0)=X(t-m).

X(t) = X(t-m) + (b(t)-d(t)) where b(t) is the

number of births between months t-m and t and d(t)

is the number of deaths in that same interval. Then

b(t) = X(t) - X(t-m) + d(t). Since X(t-m) and d(t)

are known, predicting X(t) is equivalent to

predicting b(t).

If in addition to the above assumptions about

{X(t)}, the {X(t)} are normally distributed (neither

independent or identical ), the likelihood function

for this Markov process, {X(t)}, is:

L

v

Exp

X i m

v i i

k

i

i

= −

−

=

∏

1

2 2 1

2

π

(

( ( ) )

) where

m and v are i i M(1) and V(1) evaluated with

M(0) = X(i-1) for i=1,2,.....,k. That is, mi=

a

e X i e

λ μ

λ μ λ μ

−

− + −

− −

{ } ( )

( ) ( )

1 1 and

v X i a e

a X i e a

i

=

−

− − +

− + + − − +

−

−

1

1

1

2

2 2 2

2 2

( )

{( ( )( ) )

[( ) ( )( )] }

( )

( )

l m

l m l

l m l m m

l m

l m

for i=1,2,3,......,k and λ ≠ μ .

For λ=μ,

m a X i and

v a X i a

i

i

= + −

= + − +

( )

( )

1

m 2 1 m

The terms in L are the conditional densities of

X(i) given the outcome of X(i-1) for i=1,2,.....,k.

L(X(1),X(2),X(3),........,X(k)) is the joint density of

(X(1),X(2),X(3),........,X(k)) for this Markov

process. Finding the values of λ , μ, and a which

maximize this likelihood function is done by

computer rather than algebraic manipulation.

Using historical data for 31 months (Apr 1990 to

Apr 1993) in each of 10 industries in Iowa, the

values of λ , μ, and a that maximize this likelihood

were found.

This initial testing indicated that the relatively

simple expressions for first and second moments

when μ = λ may be most appropriate, with the

common value of λ and μ of about: .02.

When λ = μ, the problem reduces to estimating

“a”.

3) ARIMA Modeling and Other Alternatives

In addition to the MLE for “a” described

above, three other estimators are tested: a

regression estimator, ar

, a simple slope estimator,

a n, and an ARIMA based estimator, aAR. Let a m
denote the maximum likelihood estimator of “a”

(maximizes L).

ar= 

iX i M i

i

i

n

i

n

i

n

( ) − ( )

= =

=

∑ ∑

∑

0

1 1

2

1

 is the value of a

that minimizes: (X(i) ai M( ))

i

n

− −

=

∑ 0

2

1

.

a n= 

X n X

n

( ) − ( )

−

1

1

 is the average month-tomonth change from month 1 to month n.

In spite of good estimates for the mean monthto-month change, the variation around this mean of

the realized change is large enough to make

predictions of the change for any particular month

very problematic.

The correlations between the three estimates of

“a” are all about .9 This implies that little is to be

gained by using a composite estimator (with

optimal weights only a variance reduction of about

10% and actual weights would necessarily be only

approximately optimal). .

The Bureau may be sampling for births and

estimating birth employment periodically. Let yb

denote this proposed the sample estimate of birth

employment. Let nb nb d nd

= + − and let

xb nbAEt k= − be the estimate of birth

employment projected from the historical data.

These two estimators are stochastically independent

and there respective variances are available

(estimable). V(y ) b

 is estimated from the

sampling distribution.

V(x ) AE V(n ) b t k b d= − −

2

, and

( ) ( ) [

]

V n V ka k a k AE k

ak

b − d t − k

= = + +

2 2

m 2 m

k=3 in case births are estimated for a quarter at a

time and updated from the last quarter UDB

employment figures.

The composite estimator for birth employment is:

C

y V x x V y

V x V y

b

b b b b

b b

=

+

+

( ) ( )

( ) ( )

. and the

variance of this composite estimator is:

V C

V x V y V y V x

V x V y

b

b b b b

b b

( )

[ ( )] ( ) [ ( )] ( )

[ ( ) ( )]

=

+

+

2 2

2

and this is estimated with

( )

( ) ( )

( ) ( )

V C

V x V y

V x V y

b

b b

b b

=

+

 ≤







( )

( )

V x

V y

b

b

. The

variance of the composite estimator is minimized

when V(x ) V(y ) b b= and is half their common

value.

 The Time Series programs that were used to

fit ARIMA models, test them, and predict the

industry sizes (number of reporting units) from a

fitted model for this study was done using S-PLUS

software, Venables and Ripley (1994), on a PC. As

was noted above, when birth and death processes

were used to model this problem, the noise in these

stochastic processes overwhelmed the signal. The

same was true when ARIMA models were fitted to

the same data.

AR(1) to AR(3) processes seemed to best

describe the behavior of the data as measured by

the Akaike Information Criteria. The

Autoregression coefficients were estimated via

maximum likelihood assuming a Gaussian

distribution. The ARIMA estimates of “a” denoted

aAR (or multiples of “a” as tabulated in the table

above) appear uncorrelated with the other

estimates. This fact suggests that a composite of the

ARIMA estimates and one of the other three may

give an improved estimate of “a”.

4) An Empirical Study

Historical data for Iowa from SICs 16, 23, 24,

25, 28, 30, 53, 63, 73, and 75 for 37 months

(March 91 to March 94) was used to test the three

estimators 6a (half year ahead projections of SIC

growth in number of reporting units). The MLE

was generally best (possibly because it tends to be

more conservative, smaller) and is based on a more

specific stochastic description of birth and death

processes than the other predictions. For these six

month ahead projections, data for the first 31

months was used to estimate change between month

31 and month 37. These projections are given below

in Table 1.
Table 1a. Root MSEs (across SICs) of Six month

ahead changes (smallest in bold). These come from

Table 1, root of average of ([change row] -

[predictor row])2.

MSE →

Est ↓

All SICs W/O SIC 16

aAR 8.13 7.75

a m

11.31 3.37

ar

10.43 4.12

a n

10.37 4.34

SIC 16 is highly seasonal and the ARIMA

predictor did a much better job of picking up this

seasonality than the other methods as shown in

Table 1a. When SIC 16 was removed (second

column of Table 1a, the MLE was the best method.

These estimated changes in number of

reporting units must be translated into change in

employment. Since births and deaths generally

predominate in the smaller size classes, it would

give biased estimates of employment change to

simply multiply average reporting unit employment

in an SIC by the estimated change in the number of

reporting units. The historical data can also be

used to determine month-to-month change in

employment as well as month-to-month change in

the number of reporting units. The ratio of the one

sequence to the other then estimates the change in

employment per unit change in reporting units.

Some weighted average of these ratios can

be multiplied by the projected change in number of

reporting units to estimate projected change in

employment.

5) Conclusions

It may be futile to try to predict changes in

number of reporting units in an estimation cell

based on past numbers of births and deaths. A

random sampling and estimation plan is a necessary

adjunct to the indirect techniques studied in this

paper. One current plan is to sample for births on a

“just in time” basis. Look for births between frame

refinements and include them in the estimation

process in a more timely fashion than is now

occurring. At present, the effect of births is only

accounted for during annual benchmark revisions.

We propose getting birth employment by sampling

new Unemployment Insurance accounts and

screening these to remove pre existing

establishments while asking new establishments for

their employment.

The indirect techniques studied here may be

useful in a composite estimator where the other

component is a direct sample based estimator for

birth employment. The MLE of “a” that is derived

from the Birth-Death process seems to be better

(except in highly seasonal industries) than the other

alternatives considered: the regression estimates,

naive estimates, and the Time Series estimates. The

MLE can be improved if seasonality is first

removed from the historical series before “a” is

estimated, then this estimate of “a” is used to

predict the 6-month ahead nonseasonal component

of the time series, and then the seasonal component

added back to get the final prediction.

References

Grzesiak T. J.and Lent J. (1988) Estimating

Business Birth Employment in the Current

Employment Statistics Program, Proceeding of the

American Statistical Association, Section on

Survey Research Methods.

Karlin Samuel (1973). A First Course in Stochastic

Processes, Academic Press.

Venables W. N. and Ripley B. D. (1994). Modern

Applied Statistics with S-Plus, Springer-Verlag.

Table 1. Six month ahead changes for four predictors and the observed change.

SIC→

Est↓

16 23 24 25 28 30 53 63 73 75

a n

6.1 1.1 3.1 2.2 -2.4 3.1 -3.6 0 3.6 2.7

ar

6.6 2.1 4.1 3.5 -1.3 3.2 -.7 1.1 4.3 3.5

a m

10.3 .4 4.5 1.7 .1 2.9 0 0 4.1 4.5

aAR -13 0 -1 -1 4 -3 3 -1.5 -6.5 -11

change -24 -6 7 1 3 -3 0 2 3 7