Analysis of a modified gravity model

Transpn Res. Vol. 4, pp. 51-62. Pergamon Press 1970. Printed in Great Britain

ANALYSIS

OF A MODIFIED

GRAVITY

MODEL

H. J. EDENS De Leuw, Cather & Company

of Canada Limited, Don Mills, Ontario,

Canada

(Received 3 April 1969; in revised form 10 November 1969) INTRODUCTION THE GRAVITY

model has become a standard transportation planning procedure for estimating interzonal trip interchanges. The process of trip distribution using the gravity model is described in the Bureau of Public Roads (B.P.R.) manual (1965). A brief description of the basic procedures discussed in the manual is repeated in this paper in order to show the reasons for developing the modified gravity model. The latter part of this paper describes the modified gravity model and some of the results obtained from a calibration using work-to-home trips in the evening peak period of travel in Quebec City. The formula describing the gravity model is as follows: Tij =

Pi A, Fii Kii n C Aj Fii Kii

(1)

j=l

In this equation: Tsj = number

of trips between

Pi = total trips generated A, = total trips attracted

zone of origin i and zone of destination

j

by origin zone i by destination

zone j

Fii = travel time factor describing the effect of spatial separation activities on trip interchange (spatial separation can be measured time, distance or cost)

of in

Kij = a zone i to zonej adjustment factor to allow adjustments not otherwise accounted for (in the remainder of this paper this factor will be assumed equal to 1) n = the total number

of zones.

Based on experience of previous studies, it has been found that travel time functions Fij are significantly different for different trip purposes. Consequently, the trip distribution and calibration process is normally done for separate trip purposes. As explained in the B.R.P. manual, trip purposes are separated into two main groups-those trips that are home based and those that non-home based. A home-based trip has one trip end (either origin or destination) at the place of residence. Depending upon the amount of detail desired, homebased trips are normally separated into several purposes such as home-based work, homebased shopping and so on. The remainder of this paper will concern itself with one trip purpose, home-based work (HBW); the discussion applies equally to other trip purposes. The HBW trip production of a zone (Pi) includes all trips for the purpose of home-to-work and work-to-home leaving 51

52

H.J. EDENS

and entering the place of residence trips for the purpose work-to-home employment in that zone.

in that zone; the attraction of a zone (Aj) includes all and home-to-work leaving and entering the place of

CALIBRATION The equation

shown

(i) CT=1 Tii = Pi (ii) ~~=,Tij

MODEL These constraints

are

(column

totals are equal to zone attraction).

frequency of the surveyed origin and destination matrix should be within stated limits to the trip length frequency of the gravity model

The first constraint follows :

It should

GRAVITY

(row totals are equal to zone production).

= Aj

(iii) Trip length comparable distribution.

OF THE

above has to satisfy several constraints.

is achieved

automatically

due to the nature

be noted that cjnE1 A, Fij is a constant

of equation

(1) shown as

for each row and can be written

as

1

C$ = -

(2)

; Aj Fij j=l

The second constraint is achieved by the B.P.R. computer program by means of an iteration procedure. In this iteration procedure, the attraction value Aj is changed based on comparison of achieved column totals to desired column totals. Two iterations are followed through step by step below: STEP

1:

(Iteration

Tii = 7

1)

C A, Fij j=l

STEP 2:

Add over the columns

of the matrix

5 Tij = Sjcl, i=l

S icr, = achieved STEP 3:

Compare

Tii

=

:

Ai(2)Fii

X 4~

j=l n

STEP 5:

C Tij = Si(2j i=l

summation

for column j after iteration

S,,i, to desired trip attraction

In this equation sTEp 4:

column

Fii

Aicl, is attraction (Iteration

2)

Aj. Change attraction

1

A, for iteration

Aj used for the first iteration.

2 to

53

Analysis of a modified gravity model

STEP 6:

Aj(a) = $-

Ajcz,

J(2)

STEPS:

Ti*=

pi AjW

~

Fij

(Iteration

3)

c Aj(3) Fij j=l

The same equations

for the Kth iteration

etc. are

(3) 1

5Aj(k)Fii

j=l

is a constant

for row i and can be represented

Ajuc)=

as Ci(k). Also

Aj s

Aj(k-l)

j(k-1)

Substituting

Aj(k--l)= A

A.

sj(k_2i

in this equation

Cj(k)

’ is constant

one arrives at

for column j. c.j(k)’ _

Substituting

Aj(k-2)

the equation

Ai

Aj

Sj(k--l)

Sj(k-2)

Aj Aj ‘.’

sj(2)

for Ci(k) and C,,,, ’ in the gravity

Tij= ci(k)

Cj(kJ'

& Aj ‘ij

(4)

sj(l)

model equation

results in (5)

The row and column constants Ci and Cj’ are adjusted in subsequent iteration runs. Both constants will be seen to be indicative of the relative accessibility of the zones concerned. The third constraint selected in the case of the B.P.R. gravity model is the trip length frequency comparison. Both the model and surveyed trip length frequencies are plotted and compared visually and also the average trip length of model results is usually required to be within 3 per cent of survey results. The calibration process for the Fij values is based upon this third constraint. It is essentially a trial-and-error process whereby the curve for Fij is plotted graphically and adjusted by comparing the two frequency curves. The gravity model is rerun with adjusted Fii values until an acceptable comparison is achieved. CONSIDERATIONS LEADING TO MODIFICATION OF THE GRAVITY MODEL Ci and Cj’ are constants related to the accessibility of zones. Ci has been commonly used as an accessibility index for modal split analysis. As shown on Fig. 1 for Quebec City, the 5

H. J. EDENS

54 Cj’ values arrived pattern.

at after the first iteration

c;= - “i rn Tij i= I FIRST

ITERATION

AS

OF

RESULT

GRAVITY

of the gravity

model also show a geographical

OF

MODEL

FIG. 1. Geographical

distribution

of Cj’ values.

A way of interpreting the significance of Ci and Cj’ is to consider these values as modifiers of the Fijvalues. The propensity against travel is then (5)

55

Analysis of a modified gravity model

pair. Citk) and Cick)’ are variables The value Cick) Cjck) ’ Fij is different for each zone-to-zone describing locational features of production zone i and attraction zone j. Fii is a function of the separation of the origin and destination zone, it is not dependent on the location of both zones. In other words, two suburban zones 5 min apart and two central business zones 5 min apart use the same Fij value. Distribution of trips for a future time period, for instance a design year 20 years hence, has been based on the assumption that the functional form of Fij will remain constant over this time period. As the pattern of urban development in the city or region can change considerably in this time period, it is questionable if the form of Fii will remain constant (Cowan and Mladinov, 1966). It would be a more acceptable assumption to consider Fii of a constant form over this time period for zones of similar accessibility. This makes it desirable to have time functions more closely tied to the locational characteristics of the origin and destination zone than equation (5) allows. One way of achieving this would be to give suburban-suburban trips one Fij curve and central business district to suburban districts another Fij curve and to calibrate the curves separately. This process is cumbersome and mathematically inelegant, especially if more than two geographical classifications are used. The modified gravity model discussed in the next section describes an alternative way of using time functions related to the accessibility of the zones. DESCRIPTION

The modified

OF THE

MODIFIED

gravity model equation

GRAVITY

MODEL

is

(6) E is a balancing

constant

for the matrix

E=

and is given

by

C m.

CPi Ajfpifaj

i=l j=l

C=

~Pi= i=l

(7) ~Aj j=l

i

C is the grand total of all trips distributed; n = total number of zones; fpi is a function of some measure of spatial separation such as time tij related to production zone i; f~i is applied to row i; J;lj is a function of some measure of spatial separation such as time tii related to attraction zone j; fai is applied to column j; each zone has its own fpi and Jaicurve assigned to it. For a 400 x 400 trip matrix, 400 fpicurves and 400 f,j curves would have to be calibrated and coded for input. This is an unacceptable number and would make the modified model impractical. Also, data on a zonal basis are not statistically meaningful to be used for calibration purposes, and grouping by districts therefore is required. The key to solving this problem in such a way that the use of equation (6) is practical is the following assumption. In the development of the computer program to calibrate and iterate the modified gravity model, it is assumed that the fpi and fai curves of zones in locations of similar accessibility will have similar shapes and values. This assumption allows the use of a family of curves rather than separate curves for each zone. The use of a family of curves is further explained in the next section of this paper.

56

H. J. ITERATION

PROCESS

OF THE

EDENS

MODIFIED

GRAVITY

MODEL

A computer program was written in order to be able to test the practicability of the modified gravity model equation. The following describes the processes used in this program and some results of its test application on survey data collected for Quebec City. The following constraints have to be satisfied: (i) Row totals are equal to input value Pi. (ii) Column totals are equal to input value Aj. (iii) A restraint based on a comparison of model and survey on trip length frequencies in this case referring to groupings of zones. Constraints (i) and (ii) are satisfied using an iteration procedure that selectsf& and j$ curves from a family of curves. The process is best explained by using an example. A family of curves, of the general form shown in Fig. 2, is given as input. Each curve is identified by a number and up to 11 curves can be used. A curve numbered 2.60 would be determined by interpolation of values of curves 2 and 3, as shown in Fig. 2. The iteration could be started by assigning one fPi curve number and one fai curve number to all zones, say, curve No. 2.00 from thefPi family of curves and curve No. 340 from thef& family of curves. Other input is the value of each Pi and Aj and a skim tree matrix containing time or cost between all zones i-j. An origin and destination matrix containing the product of Pi. Aj .fpi .f&is calculated and row totals and column totals are determined. E as defined in equation (7) is calculated and row and column totals are multiplied by E. The resulting row and column totals are compared to input values Pi and Aj. Curve selection is done as follows. If a row or column total resulting from the first iteration is higher than input value, then a lowerf,, or&,. curve is selected by the computer program for that zone. If the achieved row or column total is lower than the input value, a higher curve is selected. The second iteration is then started with selected new curve numbers assigned to each zone. It is our experience that after three or four iterations a satisfactory balance is achieved. The iteration subroutine of the program stores and prints thef& and fai curve numbers selected for each zone. Geographical distribution of the selected curve numbers shows that the number is a measure of accessibility of the zone (see Fig. 5). The iteration process determines therefore which zones have similar accessibility. Just as described for the gravity model in the Introduction of this paper, a third constraint is needed to select from the many possible origin-destination matrices a realistic one. This constraint (iii) leads to curve calibration. In the gravity model this third constraint is based upon comparison of city-wide trip length frequencies. A similar procedure is used for the modified gravity model, but the comparison is based upon trip length frequencies for districts of the region. These districts, as will be shown later (Fig. 5), are made up from zones of similar accessibility. The calibration is done separately for trips produced by and trips attracted to the zones. The modified gravity model program, as it is written, allows 11 curves for each set of the a& and faj curve families. It is not necessary to use all 11 curves; the analyst can select to use a coarser description and use, for instance, a family of 5 curves. The calibration process follows the rule that all zones using, for instance, curves 2.51-3.50 are considered to belong to district 3 (see Fig. 3). In general terms, all zones using curves (m- 1) + 0.51 to curve m + 0.50 belong to district m. The curve calibration is started by coding a family of curves for fpi and for faj. Curves calibrated for a previous project could be used. The general form is shown in Fig. 2. The model is iterated as described earlier, to satisfy constraints (i) and (ii). The result of this iteration is that curve numbers are assigned to each zone. In order to calibrate an fpi

Analysis of a modified gravity model

57

50

time tij

70 60 Note

faj

Curve of

50

2.60

interpolotion a+0.60b

is

determined 0s

by

means

time

50 tij

follows

40

FIG.

2.

General form of family off,+

and f,,j curves.

curve No. m, all zones using fPi curves of group m are considered to belong to district m. For the zones of this district m, the trip length frequency of the trips produced by the model is determined. The same is done for all trips produced by the same zones of district m according to survey results. The computer program is written to do this automatically.

H. J.

58

EDENS

Note : All zones

using

curve

2.51

ore port of district used colibrotion of curve 3

time FIG.

T 0 0’ b .e 5%

E

3.

to curve

3.50

for

t p

Definition of districts used for curve calibration. -

’

1

Trip all of

length

frequency

model

\

Surveyed for

\ M

for

trips produced by zones group m resulting from

I

all

zones

trip trips of

length produced

group

time (trip

frequency by

m

t length)

FIG. 4. Comparison of trip length frequency. For each time interval the ratio R, of the two trip-length frequencies is determined as shown in Fig. 4. A new curvef& is calculated by multiplying for each time interval the value of the input curvefPi by a factor equal to &lit). The same process is repeated for each next curve of the

Analysis of a modified gravity model

59

family. The same procedure is also applied for& curves in that case using trips attracted and districts determined by the faicurve numbers. Each input curve and suggested new curve values are plotted and listed using the printer of the computer. A new family of smooth curves is graphically prepared on a large scale on graph paper and becomes the input for the next round of calibrations. The process was tested on work-to-home trip data collected for the Quebec City Transportation Study. Figure 5 shows the geographical distribution of the curve numbers

j$

OVER 6.0

/ ZONES USED FOR C/)‘,BRATION OF CVRVE Nc 6

OVER

6.0

FIG. 5. Contour

map of faj curve numbers.

60

H. J.

EDENS

assigned to the zones for the calibrated model. The location of zones contributing to the calibration of curve 6 are shown shaded on this exhibit, and it can be seen that the curve calibration for this city is basically for a grouping of zones located in rings around the central business district. SPECIAL

FEATURES

The calibration and development of a trip distribution model should not be restrained by an inflexible computer program. The program and model equations used should be as versatile as possible in order to allow the analyst to account for features specific to the problem at hand. In order to express special distribution characteristics of a zone or groups of zones, five separate sets of fp4 and .faicurves can be applied in the above-described program. This allows, for instance, the input and calibration in one operation of one set of curves for zones representing external cordon stations and one set for internal zones. In addition, up to three other groups of zones could be given a separate set of curves. As an example, this would be useful to analyse trips to and from a large airport separately in a “home-based other” model. The quality of a model depends partly on the manner in which the separation of zones is measured. Time, distance and cost are units most frequently used. It is normal practice to add terminal penalties at central business district zones to account for parking costs and congestion delays. Penalties are also often used on river crossings to express psychological or sociological separation not expressed by the speeds coded in the network. The determination of such penalty times or costs necessitates reruns of the gravity model on a trialand-error basis. As a tool to determine such penalties quicker, the following procedures are carried out in the modified gravity model and incorporated as options in the computer program. Two 10 x 10 district trip matrices of survey and model results respectively are calculated by compression of the zonal matrices. Values in each of the 100 cells in the two matrices are compared. For the cell with the largest difference in number of trips, a time (or cost) penalty is calculated that would make survey and model results equal. This correction is distributed over the remainder of the district-district trip matrix of the model to retain desired row and column totals. The comparison is then repeated for the next cell with the largest difference. Significant time penalties can be plotted in order to determine a geographical pattern and the information is then used to improve the model further. The facility to use K factors as used in the gravity model is also included in this computer program for the modified model. COMPARISON

For models not the The

OF GRAVITY

MODEL

AND

MODIFIED

GRAVITY

MODEL

the test case of Quebec City, an investigation was carried out to determine if both might result in the same distribution. The following comparison shows that this is case. gravity model equation is [repeating equation (5) for convenience]

2 = Ci i

The modified

gravity

model equation

Ci’ Fii

3

is [repeating

T.. --c = Efpifaj Pi Aj

equation

(6) for convenience]

61

Analysis of a modified gravity model

Equation

(5) could be rewritten

as ~, = Ci F~~jCj’ Faij 2 3

(8)

In this equation Fpij = Fiid Fag = Fiie, AS Ci and Cj’ are different

d+ e = 1

for each zone, equations

(6) and (8) are identical

Ci Fpij = f,i JE,

that is

Ci Fp<*_ Ci Fiid f~i = - JE dE

Cj Faij = faj JE,

that is

CjFaij _ CjFij” faj = - JE JE

if

The conclusion is that the gravity model is a special case of the more general modified model. The modified gravity model program would simulate the results of a gravity model if one used as input a family of curves for fp and f, of identical shape. Figure 6a shows such a family off, curves. The two curves f,,and f,, shown in Fig. 6a are the Fij” curve multiplied by two different values of Cj/,/E. A constant ratio C,, as shown on this Figure applies for each value of time t. 40

0.8

30

0.6

20

0.4

IO

0.2

0 0

IO

20

30

40

50

TIME

FIG. 6a. Family off,

curves simulating

the normal gravity model.

This ratio was also determined for the curves calibrated for the modified gravity model in the Quebec test case and here it was found that this ratio is not constant but shows a significant difference between short and long trip lengths as shown in Fig. 6b. The f(t) characteristics [f(t) curve shape] are therefore considerably different in our test case for zones with different accessibility. Consequently, the modified model describes more accurately the existing trip distribution than the gravity model does.

62

H.J.

EDENS

The shape of the single Fij curve of the gravity model can be considered to be the result of a city-wide weighted average of fl,i and f& curves. This consideration could point to the reason for the different Fii curve forms found for different cities for the same trip purposes. Also the Fij curve should change form, if the pattern of development of the city significantly

0

IO

20

30

$0

50

TIME

FIG. 6b. Family off,

curves calibrated

by the modified gravity model.

changes over time. Estimating the new Fij curve, however, is difficult and it is assumed for that reason that the curve remains the same. On the other hand, in the modified model the trip distribution is determined by curves calibrated for zones with similar locational characteristics. Application to future conditions of the modified model assumes that the fP and f, curves for zones of a given accessibility remain the same. The iteration process determines the new (future) accessibility of each zone in relation to the size and location of the other zones and assigns the curve number accordingly. It therefore takes the new pattern of development fully into account. This feature is especially important for large regional studies and even more so if landuse configurations considerably different from today’s are to be tested realistically in these regions. REFERENCES BUREAU

OF PUBLIC ROADS (1965). Calibrating

a Gravity

ment Printing Office, Washington. COWAN G. R. and MLADINOV J. K. (1966). Discussion Res. Rec. 141, 43.

Model for any Size

Urban Area.

U.S. Govern-

of paper: Factors in work trip lengths.

Highw.