- Email: [email protected]
Least-squares estimation of trip distribution parameters: A note
Tmnrpn Rex, Vol. Il. w. 4ZW31.
Pqamon
Press
1!?77. Printed h Gmal Britain
LEAST-SQUARES ESTIMATION OF TRIP DISTRIBUTION PARAMETERS: A NOTE T. J. WANSBEEK Economic Institute, Leyden University, Leiden, The Netherlands (Received 12 December 1975;in revised
form
1977)
10 January
Abstrset_-Inthe context of the linearized trip distribution model, simple formulae are derived for the least-squares estimators of the parameters and their covariance matrix, using a generalized inverse to solve the normal equations. An extension to the case of a non-linear distance function is given. A method for solving the zero observation problem is proposed.
In this note, estimation of the parameters of a standard trip distribution model by means of the ordinary leastsquares (OLS) method is analyzed. Simple formulae for the various quantities of interest are derived. A standard trip distribution model can be written: x; = u,P@/(@,
i, j = 1, . . . , n,
log-distances from region i: di,=$gdij,
i=l,...,
n
d; = (d,,, . . . , 4.)
(1)
and d2 is the n-vector of arithmetic means d.i of distances to region j, j = 1,. . . , n. where xi and 4, denote the observed number of trips and Clearly, D’D is singular, and as a result, the OLS the distance from region i to region j respectively; a:, /Z?? normal equations and y are the (2n + 1) trip distribution parameters. Taking logarithms, the model can be rewritten: D’S
x,,=cQ+/~~+&,
i,j=l,...,
n;
(2)
the variables appearing in (2) are the logarithms of the corresponding superscripted variables in (1). The parameters ai, Bi and y can be estimated by OLS. To this end, we add a linear disturbance ejj to the right-hand side of (2). We impose the usual assumptions of zero expectation, constant variance and independence on the e+ The model (2) can be rewritten in matrix form: x=fi+o,
= 1yx
(6)
do not yield a unique solution. To overcome this problem, we may proceed along either of two lines. One possibility is to impose some restriction on the parameters (e.g. Cesario, 1973). The alternative is to proceed with a generalized inverse of 1YD.t In the present context, the latter approach is the simpler one, We construct a g-inverse of D’D by repeated application of a result due to Rohde (1965). Define
(3)
M= I+
(7)
N = M@M;
(8
where x’ = (x,1,. . . , Xl”,. . . ,x.1.. * * 9xlul) d’ = (d,,, . . . , d,., . . . , d,,, . . .v da.) E’= (En,. . . , Cl,, . . . ,419. . . 9E”.) n’=(a’,p’,y)“(a I,... , amB,....,#LY) D = [I@+@Ijd];
M is the familiar average operator (The& 1971, p. 12). N may be interpreted as a “double” average operator; for example:
I denotes the identity matrix of order R, and L is the
n-vector of unit elements. The symbol @ represents the Kronecker product of two matrices. In the sequel, we derive a simple formula for the OLS estimator li, using the typical structure of D as given in (4). We observe
D’D=[3
2;
;]
(5)
d’Nx = d’(M @ M)x = d’x - nd;x,
- nd$Ux,
= d’x - nd;Mx, - nd;x, = zip = xx
ii
Cd, - di.)(xii
- x,j)
(d,
- xi.),
- d.i)(xii
as can easily be verified. Next, define: A, = [ZjMjO]
where d, is the n-vector of the arithmetic means di. of
(10)
A2 = [d;jd;Mil]
tA generalized inverse (or g-inverse) of a matrix H is any matrix G satisfying HGH = H (see e.g. Searle, 1971). 429
(9)
1
A = ,A;A,
(11) 1
+ mA;A2.
(12)
430
T. J.
WANSBEEK
The following properties of A are of interest: D’DAD’D = D’D p AD’D=[;
8’ of u2 is given by: (13)
1 d2 = -x’(I n(n -2)
;LL~ 01 ;
;I
= &{x’Nx
(14)
AD’DA = A
(15)
DAD’ = I @ I - N t (d’Nd)-‘Ndd’N.
(16)
A proof of (13H16) is omitted. They are readily verified by straightforward multiplication, in which process the result d’Nd = d’d - ndld,- nd$Ud,
- DAD’)x - (d’Nd)-‘(x’Nd)*}
in view of (16); R2 = (x’x)-‘{x’x -x’Nx t (d’Nd)-‘(x’Nd)‘}. An unbiased
predictor
(23)
.f of x is given by
(17)
(see (9)) is frequently applied. Our results are now at hand. From (13), it follows that A is a g-inverse of D’D. We may use A to obtain a solution of the normal eqns (6). After some algebraic manipulations, we have:
(22)
(24 or .
xi/= xl. +
X.j
-
X..-
%C+ d.j- d.. - di,}, i,j= 1,. . . ,n.
(25)
1
These results lead to the conclusion that the approach presented in this note yields simple analytical exXI -md, pressions for all quantities of interest. As an illustration, the method will be applied to the 6 = AD’x = M(x2-$d2) hypothetical trip distribution example given by Cesario d’Nx (1975, p. 16). The data are contained in the first two : d’Nd columns of Table 1. After estimating the trip distribution parameters on basis of these data, predictions can be or, alternatively, made of the numbers of trips between the various origins and destinations. These predictions, obtained both by 6, = xi. - +di., i=l,..., n Cesario’s method and by the method presented here, are /$ = x.j - x., - +(d.j - d..), j = 1,. . . , n given in columns 3 and 4 of Table 1 respectively. As expected, the predictions obtained by both methods alf=zx (d:j- 4 .I( xi1 - x./l most coincide. We note in passing that, due to the ij transformation back to the multiplicative model, our where x.. and d.. denote the arithmetic mean of all x,, and predictions are no longer unbiased; however, if the varidti respectively. Note that, from (2), (I and ,3 can be ance of the disturbances is known, a corrective method determined up to an additive constant only. This implies can be applied to restore unbiasedness (Teekens and that Q and /3 are not estimable. The normal equations Koerts, 1972). solution li derived above is an unbiased estimator of (see Theoretically, it is important to stress the fact that, Searle, 1971, p. 169): under our stochastic specification, our method yields, a.o., an unbiased and efficient estimator for the parameter of greatest interest, the distance elasticity ‘y. at/_% EG=AD’Dr= /l-/t%, (20) Surprisingly, we obtained an estimate 9 = -0.025 as Y [ d'NX
(‘*I
1
in view of (14), where
fi=fl’Jn. Clearly, the normalization rule implied by our choice of a g-inverse appears to be: @L = 0. Comparable estimators were proposed by Cesario (1973, p. 243, formula 24). A daerence between his approach and ours is that, apart from a different normalization rule, he does not employ the tools of matrix algebra. The following properties of the OLS estimator can be derived immediately (see Searle, 1971,pp. 169 ff): var(7i) = AD’DAd = AU’,
(21)
from (IS), while A is given in (12). An unbiased estimator
Table 1. Hypothetical trip distribution data and predictions of the number of trips 1
Relation l-l l-2 l-3 2-l 2-2 2-3 3-1 3-2 3-3
Distance 5 1 3 10 2 1 1 7 5
3 2 Number of Cestio’s predictions IliPS 4 2 1 5 6 1 10 3 4
3.325 2.779 1.075 5.476 4.577 2.539 9.964 4.069 2.101
4 Our predictions 3.697 2.150 1.006 5.660 3.291 1.611 9.558 5.088 2.468
431
Least-squares estimationof trip distributionparameters:a note
compared to -0.202 in Cesario (1975). In view of the standard error of our estimate (0.334) however, the interregional distance does not contribute siicantly to the explanation of the trafhc flows in this example. Another advantage of our approach is its suitability for generalization. A frequently occurriug alternative to the linearized standard trip distribution model (2) is: xu=a,+&tf(4,y),
i,j=l,...,n,
(26)
where d,, now denotes a k-vector of impedance factors between i and j (e.g. distance, travel cost and road quality) and y is now an I-vector of parameters to be estimated. Usually, k = I, but this is not necessary. When determining the OLS estimators along the lines given above, we find that 9 is the solution of the following equations system of order 1: (x - fyNjT7 = 0,
i.e. a linear impedance function,
w=~toz.
Q* = w’w = (2 + 8, -
(28)
(29)
(34)
where p = ~‘(1 - DAD’)z.
The proof is straightforward and is omitted. In view of (9) and (17), (34) can be rewritten:
6=x,-{, (30)
(33)
rij= -jz’DALYi &=ADylit;z)
and when, moreover, C = d, (29) reduces to (18). The OL!I estimators aEand /$ in (26) are given by:
s^= M(X2 - f)3,
D'l;)'(f t oz - D’li).
This yields a system of equations in li and P. A solution of this system is:
(27) reduces to:
7 = (C’nrC)_‘C’IVX
(32)
The least-squares estimator 3 of o can be found by putting to zero the derivatives with respect to 0 and li of the sum of squares Q’:
(n)
where fV is the derivative of f with respect to 7, evaluated iu 7 = $, and f = f(d,, 9). In the special case: f(4, Y) = c7,
servation”-problem arises in all large-scale practical research work. To conclude this note we discuss the following heuristic canonical-correlation type solution: to all zero observations, a small positive constant W” is added. This constant is to be estimated along with the other parameters. Let o = In (04, and let z be au n2vector whose (i,]j)th element vanishes when the trafhc flow from i to j is nonzero and equals unity otherwise. Let x’ be the moditied n*-vector of logarithms of trathc flows, with zeros inserted when the original observation is zero. Define the n*-vector w:
1 @=
d’Nz.d’Nf - d’Nd.z’Nx’ (CNzp - d’Nd.z’Nz *
(35)
where
(31)
Again, these formulae are straightforward generalizations of (18). The advantage of our decomposition ap preach to the estimation of the parameters of (26) is that the size of the (possibly non-linear) equations system to be solved is reduced from (2n + ! - 1) to I only. The main drawback of the methods introduced so far lies in the requirement that the tratIic flows between any two regions should be positive. This familiar “zero ob-
After inserting G in the relevant places of the logobservations vector, the estimation method presented in this note can be applied.
Cesario F. J. (1973) A generalized trip distributionmodel. L Regional Sci. 12.233-248. Cesario F. J. (1975)Least-squares estimationof trip distribution parameters. Tmnspn Res-9. 13-18. Rohde C. A. (1965)Generalized inverses of partitionedmatrices. J. Sot. Ind. App. Math. 13, 1033-1035. Searle S. R. (1971)Linear Models. Wiley, New York. Teekens R and Koerts J. (1972) Some statisticalimplications of the log transformationof multiplicativemodels.Econometrica 40.793-819. Theil H. (1971) Principles of Econometrics.Wiley,New York.
Pqamon
Press
1!?77. Printed h Gmal Britain
LEAST-SQUARES ESTIMATION OF TRIP DISTRIBUTION PARAMETERS: A NOTE T. J. WANSBEEK Economic Institute, Leyden University, Leiden, The Netherlands (Received 12 December 1975;in revised
form
1977)
10 January
Abstrset_-Inthe context of the linearized trip distribution model, simple formulae are derived for the least-squares estimators of the parameters and their covariance matrix, using a generalized inverse to solve the normal equations. An extension to the case of a non-linear distance function is given. A method for solving the zero observation problem is proposed.
In this note, estimation of the parameters of a standard trip distribution model by means of the ordinary leastsquares (OLS) method is analyzed. Simple formulae for the various quantities of interest are derived. A standard trip distribution model can be written: x; = u,P@/(@,
i, j = 1, . . . , n,
log-distances from region i: di,=$gdij,
i=l,...,
n
d; = (d,,, . . . , 4.)
(1)
and d2 is the n-vector of arithmetic means d.i of distances to region j, j = 1,. . . , n. where xi and 4, denote the observed number of trips and Clearly, D’D is singular, and as a result, the OLS the distance from region i to region j respectively; a:, /Z?? normal equations and y are the (2n + 1) trip distribution parameters. Taking logarithms, the model can be rewritten: D’S
x,,=cQ+/~~+&,
i,j=l,...,
n;
(2)
the variables appearing in (2) are the logarithms of the corresponding superscripted variables in (1). The parameters ai, Bi and y can be estimated by OLS. To this end, we add a linear disturbance ejj to the right-hand side of (2). We impose the usual assumptions of zero expectation, constant variance and independence on the e+ The model (2) can be rewritten in matrix form: x=fi+o,
= 1yx
(6)
do not yield a unique solution. To overcome this problem, we may proceed along either of two lines. One possibility is to impose some restriction on the parameters (e.g. Cesario, 1973). The alternative is to proceed with a generalized inverse of 1YD.t In the present context, the latter approach is the simpler one, We construct a g-inverse of D’D by repeated application of a result due to Rohde (1965). Define
(3)
M= I+
(7)
N = M@M;
(8
where x’ = (x,1,. . . , Xl”,. . . ,x.1.. * * 9xlul) d’ = (d,,, . . . , d,., . . . , d,,, . . .v da.) E’= (En,. . . , Cl,, . . . ,419. . . 9E”.) n’=(a’,p’,y)“(a I,... , amB,....,#LY) D = [I@+@Ijd];
M is the familiar average operator (The& 1971, p. 12). N may be interpreted as a “double” average operator; for example:
I denotes the identity matrix of order R, and L is the
n-vector of unit elements. The symbol @ represents the Kronecker product of two matrices. In the sequel, we derive a simple formula for the OLS estimator li, using the typical structure of D as given in (4). We observe
D’D=[3
2;
;]
(5)
d’Nx = d’(M @ M)x = d’x - nd;x,
- nd$Ux,
= d’x - nd;Mx, - nd;x, = zip = xx
ii
Cd, - di.)(xii
- x,j)
(d,
- xi.),
- d.i)(xii
as can easily be verified. Next, define: A, = [ZjMjO]
where d, is the n-vector of the arithmetic means di. of
(10)
A2 = [d;jd;Mil]
tA generalized inverse (or g-inverse) of a matrix H is any matrix G satisfying HGH = H (see e.g. Searle, 1971). 429
(9)
1
A = ,A;A,
(11) 1
+ mA;A2.
(12)
430
T. J.
WANSBEEK
The following properties of A are of interest: D’DAD’D = D’D p AD’D=[;
8’ of u2 is given by: (13)
1 d2 = -x’(I n(n -2)
;LL~ 01 ;
;I
= &{x’Nx
(14)
AD’DA = A
(15)
DAD’ = I @ I - N t (d’Nd)-‘Ndd’N.
(16)
A proof of (13H16) is omitted. They are readily verified by straightforward multiplication, in which process the result d’Nd = d’d - ndld,- nd$Ud,
- DAD’)x - (d’Nd)-‘(x’Nd)*}
in view of (16); R2 = (x’x)-‘{x’x -x’Nx t (d’Nd)-‘(x’Nd)‘}. An unbiased
predictor
(23)
.f of x is given by
(17)
(see (9)) is frequently applied. Our results are now at hand. From (13), it follows that A is a g-inverse of D’D. We may use A to obtain a solution of the normal eqns (6). After some algebraic manipulations, we have:
(22)
(24 or .
xi/= xl. +
X.j
-
X..-
%C+ d.j- d.. - di,}, i,j= 1,. . . ,n.
(25)
1
These results lead to the conclusion that the approach presented in this note yields simple analytical exXI -md, pressions for all quantities of interest. As an illustration, the method will be applied to the 6 = AD’x = M(x2-$d2) hypothetical trip distribution example given by Cesario d’Nx (1975, p. 16). The data are contained in the first two : d’Nd columns of Table 1. After estimating the trip distribution parameters on basis of these data, predictions can be or, alternatively, made of the numbers of trips between the various origins and destinations. These predictions, obtained both by 6, = xi. - +di., i=l,..., n Cesario’s method and by the method presented here, are /$ = x.j - x., - +(d.j - d..), j = 1,. . . , n given in columns 3 and 4 of Table 1 respectively. As expected, the predictions obtained by both methods alf=zx (d:j- 4 .I( xi1 - x./l most coincide. We note in passing that, due to the ij transformation back to the multiplicative model, our where x.. and d.. denote the arithmetic mean of all x,, and predictions are no longer unbiased; however, if the varidti respectively. Note that, from (2), (I and ,3 can be ance of the disturbances is known, a corrective method determined up to an additive constant only. This implies can be applied to restore unbiasedness (Teekens and that Q and /3 are not estimable. The normal equations Koerts, 1972). solution li derived above is an unbiased estimator of (see Theoretically, it is important to stress the fact that, Searle, 1971, p. 169): under our stochastic specification, our method yields, a.o., an unbiased and efficient estimator for the parameter of greatest interest, the distance elasticity ‘y. at/_% EG=AD’Dr= /l-/t%, (20) Surprisingly, we obtained an estimate 9 = -0.025 as Y [ d'NX
(‘*I
1
in view of (14), where
fi=fl’Jn. Clearly, the normalization rule implied by our choice of a g-inverse appears to be: @L = 0. Comparable estimators were proposed by Cesario (1973, p. 243, formula 24). A daerence between his approach and ours is that, apart from a different normalization rule, he does not employ the tools of matrix algebra. The following properties of the OLS estimator can be derived immediately (see Searle, 1971,pp. 169 ff): var(7i) = AD’DAd = AU’,
(21)
from (IS), while A is given in (12). An unbiased estimator
Table 1. Hypothetical trip distribution data and predictions of the number of trips 1
Relation l-l l-2 l-3 2-l 2-2 2-3 3-1 3-2 3-3
Distance 5 1 3 10 2 1 1 7 5
3 2 Number of Cestio’s predictions IliPS 4 2 1 5 6 1 10 3 4
3.325 2.779 1.075 5.476 4.577 2.539 9.964 4.069 2.101
4 Our predictions 3.697 2.150 1.006 5.660 3.291 1.611 9.558 5.088 2.468
431
Least-squares estimationof trip distributionparameters:a note
compared to -0.202 in Cesario (1975). In view of the standard error of our estimate (0.334) however, the interregional distance does not contribute siicantly to the explanation of the trafhc flows in this example. Another advantage of our approach is its suitability for generalization. A frequently occurriug alternative to the linearized standard trip distribution model (2) is: xu=a,+&tf(4,y),
i,j=l,...,n,
(26)
where d,, now denotes a k-vector of impedance factors between i and j (e.g. distance, travel cost and road quality) and y is now an I-vector of parameters to be estimated. Usually, k = I, but this is not necessary. When determining the OLS estimators along the lines given above, we find that 9 is the solution of the following equations system of order 1: (x - fyNjT7 = 0,
i.e. a linear impedance function,
w=~toz.
Q* = w’w = (2 + 8, -
(28)
(29)
(34)
where p = ~‘(1 - DAD’)z.
The proof is straightforward and is omitted. In view of (9) and (17), (34) can be rewritten:
6=x,-{, (30)
(33)
rij= -jz’DALYi &=ADylit;z)
and when, moreover, C = d, (29) reduces to (18). The OL!I estimators aEand /$ in (26) are given by:
s^= M(X2 - f)3,
D'l;)'(f t oz - D’li).
This yields a system of equations in li and P. A solution of this system is:
(27) reduces to:
7 = (C’nrC)_‘C’IVX
(32)
The least-squares estimator 3 of o can be found by putting to zero the derivatives with respect to 0 and li of the sum of squares Q’:
(n)
where fV is the derivative of f with respect to 7, evaluated iu 7 = $, and f = f(d,, 9). In the special case: f(4, Y) = c7,
servation”-problem arises in all large-scale practical research work. To conclude this note we discuss the following heuristic canonical-correlation type solution: to all zero observations, a small positive constant W” is added. This constant is to be estimated along with the other parameters. Let o = In (04, and let z be au n2vector whose (i,]j)th element vanishes when the trafhc flow from i to j is nonzero and equals unity otherwise. Let x’ be the moditied n*-vector of logarithms of trathc flows, with zeros inserted when the original observation is zero. Define the n*-vector w:
1 @=
d’Nz.d’Nf - d’Nd.z’Nx’ (CNzp - d’Nd.z’Nz *
(35)
where
(31)
Again, these formulae are straightforward generalizations of (18). The advantage of our decomposition ap preach to the estimation of the parameters of (26) is that the size of the (possibly non-linear) equations system to be solved is reduced from (2n + ! - 1) to I only. The main drawback of the methods introduced so far lies in the requirement that the tratIic flows between any two regions should be positive. This familiar “zero ob-
After inserting G in the relevant places of the logobservations vector, the estimation method presented in this note can be applied.
Cesario F. J. (1973) A generalized trip distributionmodel. L Regional Sci. 12.233-248. Cesario F. J. (1975)Least-squares estimationof trip distribution parameters. Tmnspn Res-9. 13-18. Rohde C. A. (1965)Generalized inverses of partitionedmatrices. J. Sot. Ind. App. Math. 13, 1033-1035. Searle S. R. (1971)Linear Models. Wiley, New York. Teekens R and Koerts J. (1972) Some statisticalimplications of the log transformationof multiplicativemodels.Econometrica 40.793-819. Theil H. (1971) Principles of Econometrics.Wiley,New York.