Kapteyn, A.J.

Wansbeek, T.

English

Research Papers in Economics

r~.,;~
CBM
R ~
7626
1984
154
I T1j'-?S~'.~??rTLNBUREAU
Bestemmingl B;SLI.~:~~'í?'~:-~n
Ii:`.'i~IC:~.:`-~i:~
~C~ HOG~S;:H~Oi,
TiLBURG
subfaculteit der econometrie
RESEARCH MEMORANDUM
i~uuiiiMiiiiiiniiinNiiNUNuiiin~iuiiui ~ .~,.oo.:,.
TILBURG UNIVERSITY
DEPARTMENT OF ECONOMICS
Postbus 90153 - 5000 LE Tilburg
NetherlandsFEW
154
Statistically and Computationally Efficient
Estimation of the Gravity Model
Tom Wansbeek~)
and
Arie Kapteyn~~)
~) Groningen University, Econometric Institute, P.O. Box 800,
9700 AV Groningen, The Netherlands.
~~) Tilburg University, Department of Econometrics, P.O. Box 90153,
5000 LE Tilburg, The Netherlands.Abstract
Least squares estimation of the gravity model in log-linear form
is, under certain assumptione, statistically efficient. Computational
efficiency ( essential when the number of nodes distinguished is large)
is obtain2d by estimation via rhe an-~alled covariance transformatíon.
This transformation can only be applied when flows are observed between
all nodes. It therefore breaks down in the generic case where there are
no flows observed from node i to itself. In the paper, an alternative
transformation i s presented that allows for efficient computation of the
estimates in the generic case.3
1. Introduction
In this paper we present a computationally efficient calibration
methodl) for the widely used gravity model. In log-linear form the model
can be written as follows:
(1.1)
yij : ai } Yi }~ dii } Ea, ~ i~j ' 1,...,R,
--~ ~J
where yij is the logarithm of the traffic flow from node i to node j; a
i
and Yj are generation and attraction factors respectively, specific to
nodes i and j; dij is some measure of the distance between i and j, E
ij
is a random disturbance term, which we assume to be normally distributed
with mean Zero and constant variance, a2, for all i and j. The parame-
ters ai, Yj and S are unknown and have to be estimated.
Depending upon the particular context to whích the model is ap-
plied, a different operationalisation of dij will be used. A fairly
general specification would be to replace b dij by
Ek~lskxijk' where
the sk are parameters and xijk are different measures of the distance
between 1 and j. This transforms model ( 1.1) into
K
(1.2) Yij ~ ai f Yj f kEl Skxijk f eij , i,j ~ 1,...,R
If we would have observations on yij and xijk for all 1, j and k, esti-
mation of (1.2) by means of ordínary least squares would be statistical-
ly efficient, i.e. satisfy the Gauss-Markov theorem. Although the number
of parameters may be large in practice, namely 2R t K, use of the well-
known covariance transformation reduces the computational burden of
calibration to manageable proportions. For later reference, we give
these computationally efficíent formulas for the estimates
and Sk in a Lemma.
Define:
of ai, Yj
1) We will use the terms "calibration" and "estimation" interchangeably,
but the latter more frequently.4
(1.3) x - 1 E x , x - 1 E x - 1
i.k - R j ijk .jk - R i ijk' x..k - R2 i~ xijk '
k
and y, , y „ y analogously. Next define
(1.~) xijk - xijk -
xi.k - Y.jk } x1c.. ' k` 1,,...K ,
Yij - Yij - Yi~ - Y~j f Y~~ ~
Furthermore, let bk be the OLS-estima[e of gk, k~ 1,...,K, in the
transformed model
(1.5) Yij
K
~ E sk xijk f eij ,
k~l
where eij - eij - ei~ - e.j f e~., in obvious notation.
Lemma. Least squares estimates of Sk, ai and yj i n model ( 1.2) are bk,
ai and cj, respectively, where
(1.6) ai - (yi.-y..) -
E(xi.k-x..k)bk k
(1.7) cj -(Y.j-y..) -
E(x.jk x..k)bk } y.. - E
x..kbk '
k k
The estimatea of Sk are unique, those of ai and y j are unique up to an
additive constant.
Proof: e.g. Scheffé (1959), Cesario (1975), Judge et al. (1980, Ch. 8).
The non-uniqueness of the estimates of ai and yj is already
clear from model (1.2) where we can add a constant to all ai and sub-
tract the same constant from all yj without affecting the yij. This lack
of uniqueness can be used to obtain a slightly more elegant expression
for cj. Notice that the ai sum to zero. We can introduce a constant
term SO in (1.2) and impose the restriction that both the ai and the cjadd up to zero. In that case the last two terms on the right hand side
of (1.7) drop out, giving (1.6) and (1.7) a similar structure. The esti-
mate of SO is then
(1.8) b0 -
Y.. - E
x..kbk ' k
The obvious advantage of the least squares solution here is its
computational simplícity. Rather than having to invertl) a matrix of
order (2RtK) x(2RfK), when applying OLS to (1.2) dírectly, we can carry
out a simple transformation of the data, first estimate the sk which
requires inversion of a KxKrmatrix only and then compute the estimates
of the ai and Yj straightforwardly.
A typical characteristic of traffic flow data ie that there will
be no observations for i- j, i.e. one does not observe a flow from node
i to itself. As a result, the lemma cannot be applied. Although one
could still try to apply least squares to (1.2) directly, this will turn
out to be non-feasible for many realistic problems where R may large (we
return to this in section 3). Hence, we require results similar to those
presented in the lemma, for the model where observations are missing
whenever i- j. These results are the main objective of thís paper and
are collected in the theorem below.
Define:
1 ~ 1 1 (1.9)
xi.k - R-1 ~ xijk' x.jk - R-i i xíjk' x..k - R(R-1) i~ xijk '
k a 1,...,K,
and
Yi.' y.j' Y..' E i.' E.j' e.. analogously. To appreciate the defini-
tions, notice that the total number of
stead of R2 and that for given k and i
are R-1 observations xijk. Next define
observations is now R(R-1) in-
(and for given k and j), there
(1.10) x - x -(R-1 2 X -(R-1 2 - R-1 ~
ijk - ijk R(R-2) i.k R(R-2 x.jk R(R-2) xj.k -
1) Since this matrix has at most rank 2R f K-1, a generalized inverse
has to he taken, of course.6
- R(R-2) x.ik } RR2 x..k . k ~ 1,...,K
and-yij and eij analogously. Note that for large R, xijk will be close
to xijk defíned by (1.4). Furthermore, let bk be the OLS-estimate of
Bk, k~ 1,...,K, in the transformed model
K
(1.11) yij ~ E Sk xijk f Ei kal j ~
then we
typical
tity á2
var(Eij)
define St as ft ~ o2(X'X)-1, where X ís the
element xijk and denote a typical element
is defined as the "usual" best quadratic
in this model, a nd
(1.12) rij - Yij - k bkxijk
R(R-1) x K matrix with
of R by wkR . The quan-
unbiased estimator of
1 1 1
, r - - E r 1.13)
ri. - R-1 j rij R-1 i ij' r.. - R(R-1) i~ rij '
In addition we define:
2
(1.14)
qik - R(R12) xi.k } R(R-2) x.ik - R12 x..k
R-1 (R-1 2 R-1 (1.15) sjk - R(R-2) xj.k } R(R-~ x.jk - R-2 x..k '
Theorem. Let in model (1.2) observations be míssinq if and only if i~ j.
Least squares estimates of gk, ai and yj are bk, ai and cj, where
(1.16) " R-1 2
ai - R(R-2) ri.
R-1 (1.17) cj -
R(g-2) rj.
R-1
} R(R-2) r.i - R-2 r..
(R-1)2 R-1
t R(R-2) r.j - R-2 r..
The estimates of Sk are uníque, those of ai and yj are unique up to an
additive constant. The variance covariance matrix of the estimates bk
is n, The best quadratic unbiased estimator of o2 is
17
(1.18) 02 - á2 R(R-1)-K
R(R-1)-(2RfK-1) '
The coefficient of determination is
( 1.19) R2 - 1 - {R(R-1)-(2RfK-1)}a2
E E yi j-y. )Z---
ij '
The variances of ai and cj are
(1.20) var(a )- 02 2R2-3R-1 } E E q q w, í L 1. i R(R-1)(R-2) k R ik iR kR '"'R
(1.21) var(c ) ~ a2 (~ f E E s s w , j R(R-2) k~ jk jR kR j ~ 1,...,R .
Section 2 will be devoted to a proof of this theorem and section
3 to a discussion of its practical importance.
2. Proof of the Theorem
Rewrite model ( 1.2) in matrix form as follows,
(2.1) y- Zd t XB t e,
where
(2.2) Y - (Y12.Y13....,y1R,y21,y23,...,y2R,...,yR1,...,YR~R-1)~
(2.3) - ,e ,...,e ,e ,e ,... e ~ E E12 13 1R 21 23 ' 2R'"''ERl'"''ER,R-1)
X is a R(R-1) x K matrix of which the k-th column xk is
- 2.4) xk -
(x12k'x13k' " ''x1Rk'x21k'x23k'" ''X2Rk'" ''xRlk' " ''xR,R-l,k)r
(2.5) Z - (Z1,Z2) , R(R-1) x 2R,8
with
R-1)
1
(2.6) Z1 -
1
1
1
~~( R-1) Z2 -
0 1
1 0
1
1
KR-1)
1
(2.7)
R R
S - (61~62....,Sk)' , d - (a',Y')' ~
1
1
1
with
(2.8) a - ( a1,a2,...,aR)'~ Y - (Y1~YZ~...,yR)'
In a slightly more compact form, (2.1) can be wrítten as
(2.9) Y - Wp t E,
1 0
(R-1)
(R-1)
(R-1)
with
(2.10) W - (Z.X). P - (B'~d')'
The matrix Z defined by (2.5) has rank 2R-1 and hence W will
generally be of rank 2R-1fK. Consequently, there is no unique leaet
squares solution for the (2RtK)-vector p. It is well-known however,
(cf., e.g., Searle, 1971, Ch. 5), that all least squares solutions for
p are generated by9
(2.11) P - (W'W)- W'Y~
where (W'W)- is a g-inverse of W'W. By varying over the set of all pos-
sible g-inverses of W'W we generate all possible solutions p. Further-
more,
(2.12) var(p) ~ a2(W'W)- W'W(W'W)-
and the best quadratic unbiased estimator of o2 is
(2.13) 02 ~ y'[I-W(W'W)-W']y ~
R(R-1)-rank(W)
The coefficient of determination is
(2.14) R2 - 1 -
~~[I-W(W~W)-W,~~Y
y'y-R(R-1)yL .
The proof of the theorem consists of an elaboration of these
formulas, with respect to the model at hand. We will start with ( 2.11).
I Z' Z Z' X~
(2.15) W'W -
X' Z X' X
so that a g-inverse of W'W is (cf. Rohde, 1965)-
(z.16) (w'w)- 3
(z'z)-f(z'z)-z'x(x'Px)-lx'z(z'z)- -(z'z)-z'x(x'px)
I-(x'Px)-lx'z(z'z)- (X,Px)-1
where (Z'Z)- is a g-inverse of Z'Z and p- I-Z(Z'Z)-Z; P is idempotent
and invariant under the choice of (Z'Z)-. Using (2.16) to elaborate
(2.11) yields
„ d (Z'Z)-Z'{y-(X'PX)-1X'Py}
(2.17) p - w t
b (X'PX)-1X'PY10
Notice that P and hence b is invariant under the choice of a g-inverse
of Z'Z, i.e. b is unique, as claimed in the theorem. The non-uniqueness
of d can be investigated formally by employing the theory of estimable
functions. For [he sake of brevity, this part of the proof is omitted
here.
The next step is to elaborate the matríx P, to obtain a formula
for b that is computationally efficient. First, notice that
(2.18) Z'Z ~
(R-1)IR JRIR
~
JRIR (R-1)IR
where IR is the identity matrix of order R and JR is a square RxR-matrix
of ones. It can be verified dírectly that the following matrix is a g-
inverse of Z'Z:
(2.19) (Z'Z)- a
0 0
R-1
} R(R-2)
(R~ ER R11 ER
1
R-1 ER ER
with ER - IR - R JR. Next, it is a matter of straightforward manipula-
tion to show that
(2.20) Z(Z'Z)-Z' - R(R12){(R-1)Z1Zi f Z1Z2 f Z2Zi t(R-1)Z2Z2 - RR1 JN},
where JN i s an NxN-matrix of ones, N- R(R-1). For P we thus obtain
(2.21) P - I-Z(Z'Z)-Z' -
- I - 1(z z'tz z') t 1 (z -z )(z -z )' f 1 J N R-2 1 1 2 2 R(R-2) 1 2 1 2 (R-1)(R-2) N
To see how this works out for b, consider as an example the expression
for Py. Define nl - Ziy, n2 - ZZy, n3 - iNy, with 1N an N-vector of
ones. This yields
(2.22) py 3 y - R12 (ZlnltZ2n2) f R(R12) (zl-z2)(nl-n2)
f(R-1)(R-2) tN
Thus the (i,j)-element of the vector Py is11
(2.23) (PY) ~ Y - 1(n fn ) f 1 ~ ij ij R-2 li 2j R(R-2)(nli-nlj-n2i~2j) }( R-1)(R-2)
or
(2.24) (PY)ij - Yij ,
as defined on the previous section (cf. (1.10) and below). Premultipli-
cation of X by P transforms the columns of X analogously, so that xijk
is replaced by xijk, cf. (1.10). As a result we recognize the formula
for b as the least squares estimate of g in model (1.11), as claimed in
the theorem.
n w Now consider d-(a',c')', where a-(a1,a2,...,aR)',
c~(c1,c2,...,cR)'. The expression for d in (2.17) can be rewritten as
(2.25) d - ( z'Z)-Z'r ,
where r- y- Xb. Defíne
(2.26) P1 - Zir, P2 - ZZr. P- tNr - tÁ P1 ~ t
Using (2.19), we can write
(2.27) á R-1
R(R-2)
R P2
R11 (R11 P1 f P2 - Á tR)
1 - ~
R-1 pl } p2 R-1 tR
J
R(R-2)
(R-1)
P1 } P2 - R-1 tR
pl t(R-1) p2 - p.tR
J O
Inserting the expressíons for pl, p2 and p in (2.27) directly yields
(1.16) and (1.17) of the Theorem.
Now consider the variances of the estimators. It is easy to
verify that (Z'Z)- given by (2.19) is a reflexive g-inverse of Z'Z and12
tha[ (W'W)- given by (2.16) is a reflexive g-inverse of W'W. Hence the
variance covariance matrix of p is a2(W'W)-. Consequently a2(X'PX)-1 is
the variance covariance matrix of b and thís is nothing else than S2
defined below (1.11), as claimed in the Theorem.
From (2.16), we have that the variance covariance matrix of d is
equal to
a
(2.28) var(d) - var
c
- a2{(z'z)-f(z'z)-z'x(x'Px)-lx'z(z'z)-}
To prove (1.20) and (1.21) in the Theorem, we will show that these ex-
pressions correspond to diagonal elements of [he matrix on the right
hand side of (2.28). First consider (Z'Z)-. From (2.19) it is clear that
its diagonal elements in the upper left (RxR)-block are equal to
2
(2R -3R-1)~R(R-1)(R-2). The diagonal elements in the lower right (RxR)-
block are equal to (R-1)2~R2(R-2). Next, consider (Z'Z)-Z'X(X'PX)-1.
X'Z(Z'Z)-. To evaluate this expression, we note that the k-th column of
(Z'Z)-Z'X has the same structure as d, given in (2.25). By analogy with
(2.27) and (1.16) and (1.17) it is clear that
(2.29) (Z'Z)-Z'R - I S~ ,
where the (RxK)-matrices Q en S have typícal elements qik and sj~, de-
fined in (1.14) and (1.15). Denoting the (k,R)-element of a2(X'PX) by
wkR, as in section 1, we obtain as typical diagonal elements of
(Z'Z)-Z'X(X'PX)-1X'Z(Z'Z)-,
EkXRqikqiRwkR' 1~ 1,...,R, 1n the upper
left (RxR)-block and EkE~s~ks~RwkR' 3- 1,...,R, in the lower right
(RxR)-block. Combining this with the results for (Z'Z)- yields (1.20) and
(1.21) of the Theorem.
To arrive at the best quadratic estimate for o2 we note that in
(2.13) y'[I-W(W'W)-W']y is equal to y'P(I-PX(X'PX)-1X'P)Py, which is
simply the residual sum of squares for model (1.11). However, when ap-
plying least squares to (1.11), we would divide the residual sum of
squares by R(R-1)-K to arrive at the best quadratic estimator of o2 in
that model, whereas in (2.13) we divide by R(R-1)-(2RtK-1), which ex-
plains the "degrees of freedom correctiona" in (1.18) of the Theorem.13
Finally, (2.14) immediately implies (1.19).
3. Discussion
As with the standard model, we can use the indeterminateness of
the estimates of ai and y~ to introduce an intercept gD, and to restrict
ai and c~ to sum to zero. It is easily seen from (1.17) that the c~
already add up to zero, so that we only have to adapt the r..-term in
(1.16); (1.16) is hence replaced by
2
(1.16)' ai - R(R12) ri. } R(R-2) r.i - R-2 r..
and the estimate of gD is then
(3.1) b~ - r.. ,
which is analogous to (1.8).
A comparison of the formulas in the Theorem with the correspond-
ing ones in the Lemma, shows that for R i m they become pairwise iden-
tical. This is what one might have expected, because for increasing R the
missing observations make up a smaller proportion of the total number of
observations.
The formulas in the Theorem look a bit more complicated than
those in the Lemma, but their computational complexity is the same. Given
the estimate for ~, b, the computation of all other quantities
(ai,ci, etc.) requires at most 0(N) tíme (N - R(R-1)). The computation
of b itself requires simple manipulations of the data to arrive at model
(1.11). These manipulations and the computation of b also require 0(N)
time. In contrast, if least squares would be applied to (1.2) directly,
one would have to invert a(2RtK) X(2RfK)~atrix, which requires
0((2R-~K)3) time. Thus, application of the Theorem reduces the computing
time required approximately by a factor of R. In practice, where R may be
in the hundreds, this may very well be the difference between feasibílity
and non-feasibility of the estimation task.14
References
Cesario, F.J. (1975), "Least Squares Estimation of Trip Distribution
Parameters", Transportation Research, 9, pp. 13-18.
Judge, G.G., W.E. Griffiths, R.C. Hill, T.C. Lee (1980), The Theory and
Practice of Econometrics, Wiley, New York.
Rohde, C.A. ( 1966), "Some Results on Generalized Inverse ", SIAM Review,
8, pp. 201-205.
Scheffé, H. (1959), The Analysis of Variance, Wiley, New York.
Searle, S.R. (1971), Línear Models, Wiley, New York.i
IN 1983 REEDS VERSCHENEN
126 H.H. Tigelaar
Identification of noisy iinear systems with muïtiple arma inputs.
127 J.P.C. Kleijnen
Statistical Analysis of Steady-State Simulations: Survey of Recent
Progress.
128 A.J. de Zeeuw
Two notes on Nash and Information.
129 H.L. Theuns en A.M.L. Passier-Grootjans
Toeristische ontwikkeling - voorwaarden en systematiek; een selec-
tief literatuïroverzicht.
130 J. Plasmans en V. Somers
A Maximum Likelihood Estimation Method of a Three Market Disequili-
brium Model.
131 R. van Montfort, R. Schippers, R. Heuts
Johnson SU-transformations for parameter estimation in arma-models
when data are non-gaussian.
132 J. Glombowaki en M. Kriiger
On the R81e of Distribution in Different Theories of Cyclical
Growth.
133 J.W.A. Vingerhoets en H.J.A. Coppens
Internationale Grondstoffenovereenkomsten.
Effecten, kosten en oligopolisten.
134 W.J. Oomens
The economic interpretation of the advertising effect of Lydia
Pinkham.
135 J.P.C. Kleijnen
Regression analysis: assumptions, alternatives, applications.
136 J.P.C. Kleijnen
On the interpretation of variables.
137 G. van der Laan en A.J.J. Talman
Simplicial approximation of solutions to the nonlinear complemen-
tarity problem with lower and upper bounds.ii
IN 1984 REEDS VERSCHENEN
138 C.J. Cuypers, J.P.C. Kleijnen en J.W.M. van Rooyen
Teeting the Mean of an Asymetric Population:
Four Procedures Evaluated
139 T. Wanebeek en A. Kapteyn
Estimation in a linear model with serially correlated errors when
observations are missing
140 A. Kapteyn, S. van de Geer, H. van de Stadt, T. Wansbeek
Interdependent preferences: an econometric analysis
141 W.J.H. van Groenendaal
Discrete and continuous univariate modelling
142 J.P.C. Kleijnen, P. Cremers, F. van Belle
The power of weighted and ordinary least squares with estimated
unequal variances in experimental design
143 J.P.C. Kleíjnen
Superefficient estimation of power functions ín simulation
experiments
144 P.A. Bekker, D.S.G. Pollock
Identifícation of linear stochastic models with covariance
restrictions.
145 Max D. Merbis, Aart J. de Zeeuw
From structural form to state-apace form
146 T.M. Doup and A.J.J. Talman
A new variable dimension simplicial algorithm to find equilibria on
the product space of unit simplices.
147 G. van der Laan, A.J.J. Talman and L. Van der Heyden
Variable dimension algorithms for unproper labellings.
148 G..1.C.Th, van Schijndel
Dynamic firm behaviour and financial leverage clienteles
149 M. Plattel, J. Peil
The ethico-political and theoretical reconstruction of contemporary
economic doctrines
150 F.J.A.M. Hoes, C.W. Vroom
Japanese Business Policy: The Cash Flow Triangle
an exercise in sociological demyatification
151 T.M. Doup, G. van der Laan and A.J.J. Talman
The (2n}1-2)-ray algorithm: a new aimplicial algorithm to compute
economic equilibriaiii
IN 1984 REEDS VERSCHENEN (vervolg)
152 A.L. Hempeniue, P.G.H. Mulder
Total Mortality Analysis of the Rotterdam Sample of the Kaunas-
Rotterdam Intervention Study (KRIS)
153 A. Kapteyn, P. Kooreman
A disaggregated analysis of the allocation of time within the
household.uiiiui'uui~uhMiiiiiii~~uiNïWi

A new variable dimension simplicial algorithm to find equilibria on the product space of unit simplices.

de Zeeuw From structural form to state-apace form

der Heyden Variable dimension algorithms for unproper labellings.

Groenendaal Discrete and continuous univariate modelling

Heuts Johnson SU-transformations for parameter estimation in arma-models when data are non-gaussian.

Kapteyn Estimation in a linear model with serially correlated errors when observations are missing

Kleijnen en J.W.M. van Rooyen Teeting the Mean of an Asymetric Population: Four Procedures Evaluated

Kooreman A disaggregated analysis of the allocation of time within the household.uiiiui'uui~uhMiiiiiii~~uiNïWi

Kriiger On the R81e of Distribution in Different Theories of Cyclical Growth.

Least Squares Estimation of Trip Distribution Parameters&quot;,

Mulder Total Mortality Analysis of the Rotterdam Sample of the KaunasRotterdam Intervention Study (KRIS)

On the interpretation of variables.

Peil The ethico-political and theoretical reconstruction of contemporary economic doctrines

Pollock Identifícation of linear stochastic models with covariance restrictions.

Regression analysis: assumptions, alternatives, applications.

Schijndel Dynamic firm behaviour and financial leverage clienteles

Somers A Maximum Likelihood Estimation Method of a Three Market Disequilibrium Model.

Statistical Analysis of Steady-State Simulations: Survey of Recent Progress.

Superefficient estimation of power functions ín simulation experiments

Talman Simplicial approximation of solutions to the nonlinear complementarity problem with lower and upper bounds.ii IN

The (2n}1-2)-ray algorithm: a new aimplicial algorithm to compute economic equilibriaiii IN

The Analysis of Variance,

The economic interpretation of the advertising effect of Lydia Pinkham.

Theuns en A.M.L. Passier-Grootjans Toeristische ontwikkeling - voorwaarden en systematiek; een selectief literatuïroverzicht.

Tigelaar Identification of noisy iinear systems with muïtiple arma inputs.

van Belle The power of weighted and ordinary least squares with estimated unequal variances in experimental design

Vingerhoets en H.J.A. Coppens Internationale Grondstoffenovereenkomsten. Effecten, kosten en oligopolisten.

Vroom Japanese Business Policy: The Cash Flow Triangle an exercise in sociological demyatification

Wansbeek Interdependent preferences: an econometric analysis

Zeeuw Two notes on Nash and Information.

Statistically and computationally efficient estimation of the gravity model

Statistically and computationally efficient estimation of the gravity model

Abstract

Similar works

Full text

Available Versions

Research Papers in Economics