category: literaturenote citekey: roqueinvestigatingrelationshipbetweenrun2014 title: Investigating the relationship between run-off-the-road crash frequency and traffic flow through different functional forms authors: "Roque, Carlos; Cardoso, João Lourenço" year: 2014 date: 2014-00-00 2014 doi: 10.1016/j.aap.2013.10.034 publication: "Accident Analysis & Prevention" url: "https://doi.org/10.1016/j.aap.2013.10.034" zotero_key: MDPWCMPP zotero_storage: H87N43SS collections: transpordi_planeerimine folder: 001_artiklid/Road Safety Indicators firstAuthor: "Roque, Carlos"

status: converted

Contents lists available at ScienceDirect

Accident Analysis and Prevention

journal homepage: www.elsevier.com/locate/aap

Investigating the relationship between run-off-the-road crash frequency and traffic flow through different functional forms

Carlos Roque∗, João Lourenc¸ o Cardoso

Laboratório Nacional de Engenharia Civil, Departamento de Transportes, Núcleo de Planeamento, Tráfego e Seguranc¸ a, Av do Brasil 101, 1700-066 Lisboa, Portugal

a r t i c l e i n f o

Article history: Received 5 July 2013 Received in revised form 24 September 2013 Accepted 31 October 2013

Keywords: Crash prediction models Run-off-the-road crash Functional form Confidence intervals

a b s t r a c t

Crash prediction models play a major role in highway safety analysis. These models can be used for various purposes, such as predicting the number of road crashes or establishing relationships between these crashes and different covariates. However, the appropriate choice for the functional form of these models is generally not discussed in research literature on road safety. In case of run-off-the-road crashes, empirical evidence and logical considerations lead to conclusion that the relationship between expected frequency and traffic flow is not monotonously increasing.

1. Introduction

Crash prediction models play a major role in highway safety analysis. These models can be used for various purposes, such as predicting the number of road crashes or establishing relationships between these crashes and different covariates. However, the appropriate choice for the functional form of these models, relating crash frequency and traffic flow, is generally not present from research literature on road safety.

As suggested by Lord et al. (2005a), it may be preferable to begin to develop models that consider the fundamental crash process rather than making efforts for the most-fitted model.

In case of run-off-the-road (ROR) crashes, empirical evidence and logical considerations lead to conclusion that the relationship between expected frequency and traffic flow is not a linear one. For low traffic flows one may expect the number of ROR crashes per unit of time to be proportional to traffic flow. But, as traffic flows increase it becomes more and more difficult notto hit another car. Hence, for ROR crashes, proportionality cannot be expected to hold at high flows. In traffic congestion, ROR crashes are, in fact, impossible, except for very low skid resistance conditions, such as ice and snow covered pavements. This is a reflex of the fact that drivers behave differently in sparse, heavy or congested traffic and

Two main approaches have been followed to modelthe relationship between roadway and roadside characteristics and ROR crash risk (Pardillo-Mayora et al., 2010).

One method, usually referred to as encroachment-based, uses a set of conditional probabilities ofthe sequence of events thatlead to a ROR crash following the encroachment of an errant vehicle on the roadside (Mak, 1995; Mak et al., 1998; Ray et al., 2012). The main obstacle in the development of this type of models is the shortage of encroachment data. Data collected in the 1960s and 70s in North America (Hutchinson and Kennedy, 1966; Cooper, 1980) are still the main source of information on these manoeuvres (Pardillo-Mayora et al., 2010). In the development of the recently updated version ofthe Roadside Safety Analysis Programme (RSAP) – a computer programme for performing cost benefit analyses on roadside design developed under NCHRP Project 22–27 – the Cooper data was re-analyzed to attempt to resolve some of its longstanding problems (Ray et al., 2012).

A second approach is the development of generalized linear regression models fitted to cross-sectional data, to estimate ROR crash frequencies using exposure and relevant highway and roadside variables as covariates. The frequency of crashes in a given highway segment is treated as a random variable which takes discrete integer non-negative values distributed following a Poisson distribution. A generalization of the Poisson form that allows the variance of the model to be over-dispersed results in the Negative Binomial (NB) model. Lee and Mannering (2002) and Geedipally

that ROR crash frequency depends on the actual state of the system (flow, speed and density) in time and space (Hauer, 1997).

^∗ Corresponding author. Tel.: +351 218 443 970; fax: +351 218 443 029. E-mail addresses: croque@lnec.pt, carlos.alm.roque@gmail.com (C. Roque).

and Lord (2010) used Poisson and NB regression models to develop RORcrashpredictionmodels.Inbothcases, crashpredictionmodels use the following functional form:

$$\lambda_i = \beta_0 Q_i^{\beta_1} e^{\sum_i \gamma_i x_i} \tag{1}$$

The mean number of the crashes per unit of time for segment i (ⁱ) is a function of traffic flow, Qi, and a set of risk factors, xi. In Eq. (1) the effect of traffic flow on crashes is modelled in terms of an elasticity, which is a power, ˇ1, to which traffic volume is raised. The effects of various risk factors that influence the probability of crashes, given exposure, is generally modelled as an exponential function, that is as e (the base of natural logarithms) raised to a sum of the product of coefficients, i, and values of the variables, xi, denoting risk factors. In its most basic form, traffic flow is the sole explanatory variable in crash prediction models. For road segments, several crash prediction models also include segmentlength as covariate which in some cases is assumed as an offset, with ⁱ = 1. Depending on the road network characteristics, there are cases in which this simplification is not appropriate, both from the operational and purely statistical points of view (Mountain et al., 1996). According to Reurings et al. (2005), the segment length, the access density, the carriageway width and the shoulder width are a desirable minimum list of risk factors for road segments.

Eq. (1) is the conventional functional form of crash prediction models (Reurings et al., 2005). However, as noted by Hauer (2004), the crash phenomenon does not necessarily have to follow a simple monotonous, mathematical function. If the model functional form is not appropriate, the regression coefficients obtained are inaccurate, and in several cases the estimated values are ambiguous. One way to mitigate this problem is to set the traffic flow interval (lower and upper bounds) for which the modelled equations are valid; an alternative way is to fit equations with a different functional form, with overall shapes that are in better agreement with road operation characteristics. However, ADT is the average of the annual distribution of traffic volumes on selected road segments. Therefore, it is just a moment statistic of a highly non-uniform distribution, where seasonal, monthly, daily, hourly periodic trends may be detected. From traffic census, it is known that the ratios between daylight and nighttimes traffic or between winter and summertime traffic are notthe same for all segments on a road network. For that reason, it is still open to debate whether ADT is an appropriate macroscopic variable to solely represent exposure and the underlying crash mechanisms directly related to traffic volume, given the increasing availability of automatically collected data that may be used to calculate complimentary traffic distribution statistics. Nevertheless, one must acknowledge that the mentioned traffic time trends show, at least partially, a scale factor depending on the ADT value.

Following a comprehensive and systematic bibliographic search, few studies addressing the issue of crash prediction model functional form selection were identified. So far, little interest has been dedicated to study alternative functional forms to express the relationship between specific crash types and traffic flow. Reurings and Janssen (2007) compared models using the functional form expressed in Eq.(1) with models where Annual Average Daily Traffic (AADT) could be considered as a property of the carriageway under consideration and hence as a sort of continuous dummy-variable. They concluded that not only these last models did not have the desired structure but also that adding the variable AADT/1000 was indeed an improvement of the models for urban carriageways. Kononov et al. (2011) related traffic flow parameters, such as speed and density, to the choice of the functional form of crash prediction models. It compared models for urban freeways developed with sigmoid and exponential functional forms. Neural networks (NN) were used to explore the underlying relationship between accidents and other variables for urban freeway segments. The results were then compared with models calibrated by using these same data with generalized linear modelling and an NB error structure. The functional form generated through the training ofNNs suggests that a sigmoidmay be a reasonable approximation of the operational characteristics of crash occurrence on urban freeways.

The objective of this paper is to document the application of NB generalized linear models with different functional forms in the analysis of ROR crash data, exploring the underlying relationships between ROR crashes and traffic flow for interurban road segments.

The study objective was accomplished using observed and simulated datasets. The models were applied to single and dual carriageway road datasets. Subsequently, the models were applied to simulated datasets to show their general performance (Geedipally et al., 2012).

2. Methodology

This section describes the probabilistic structure of the negative binomial models, the functional form used for linking ROR crashes to covariates, the procedure employed for estimating the confidence intervals, and characteristics of the Monte Carlo simulation study.

2.1. Negative binomial models

In applying Poisson regression to crash frequency analysis, let yij be the number of ROR crashes on highway element i during period j. The Poisson model is (Washington et al., 2011):

$$p(y{ij}) = \frac{\exp(-\lambda{ij})\lambda{ij}^{y{ij}}}{y_{ij}!},\tag{2}$$

where P(yij) is the probability of y crashes occurring on highway element i during time period j and ij is the expected value of yij:

$$E(y{ij}) = \lambda{ij} = \exp(\beta X_{ij}), \tag{3}$$

for a roadway section i in time period j, ˇ is the vector of parameters to be estimated and Xij is a vector of explanatory variables describing roadway section geometric and environmental characteristics, as well as other relevant features such as traffic, that may affect crash frequency.

A feature of the Poisson distribution refers to the equality between the counts expected value and its variance. However, it is not always possible to assume that ij is constant. On the one hand, the decreasing trend in time of accident risk, as observed in many countries, weakens the validity of the hypothesis of constancy in time of probability of occurrence. On the other hand,there are unknown factors that may contribute to crash occurrence as well as factors which, although known, are quantified with measurement errors, in both cases justifying that the individual risks on each entity in a homogeneous group of entities are not identical. Thus,the ratio ofthe variance to the expected value differs from one, i.e., overdispersion or subdispersion are observed (Roque and Cardoso, 2013).

The negative binomial model is an extension of the Poisson regression model that accommodates data overdispersion. The negative binomial model is derived by rewriting Poisson parameter for each observation i at a given time interval j as (Washington et al., 2011):

$$\lambda{ij} = \exp(\beta X{ij} + \varepsilon_{ij}) \tag{4}$$

Fig. 1. The curve of the function ⁱ ⁼ ^ˇ0ADT^ˇ¹ ⁱ ^L^ˇ² ⁱ when ˇ¹ is below 1.

where exp (εij) is a random error term that follows a Gamma probability distribution with mean 1 and variance ˛. This term allows the variance to differ from the mean as below:

$$VAR\left[y{ij}\right] = E\left[y{ij}\right]\left[1 + \alpha E\left[y{ij}\right]\right] = E\left[y{ij}\right] + \alpha E\left[y_{ij}\right]^{2}$$ (5)

2.2. Functional forms

Two functional forms were considered. In both cases, the Average Daily Traffic (ADT) and the segmentlength (L) are the covariates considered. Based on crash data considered, the segment length was assumed to have a non-linear relation with the crash occurrence. Thus, the segment length is considered to be a covariate rather than considered as an offset.

The first statistical model considered in this study, hereinafter called Functional Form 1 (FF1), is similar to the one referred in Eq. (1), except that there are only two independent variables, namely Average Daily Traffic (ADT) and segment length (L). Although such model suffers from omitted variables bias (because many non-flow related factors are known to affect crash frequency, which are not considered), the empirical assessment carried out in this study still provides useful information for the development of improved ROR crash prediction models.

The mean of the ROR crashes per reference time period for segment i (ⁱ) can be calculated by:

$$\lambda_i = \beta_0 A D T_i^{\beta_1} L_i^{\beta_2} \tag{6}$$

where, ADTi =Average Daily Traffic of segment i, Li = length (in km) of segment i, ˇ⁰ = intercept (to be estimated), ˇ¹ = coefficient associated with ADT (to be estimated), ˇ² = coefficient associated with L (to be estimated).

It is usually assumed that the total number of crashes increases at a decreasing rate as the traffic flow increases (see Fig. 1). This relationship is characterized in predictive models by a coefficient for the traffic flow parameter (ˇ1) below 1.

In the second statistical model considered, the relationship between the number of ROR crashes, ⁱ, and ADT, is described by a modification of the Ricker (1954) model.

The Ricker (1954) model was initially used in ecological research, where it is important to estimate the survival probabilities of plants and animals and where population densities usually do not show an indefinitely increasing trend but show local maxima in the explanatory variable space.

Applied to crash prediction, the Ricker model is expressed by:

$$\lambda_i = \beta_0 A D T_i e^{-\beta_1 A D T_i} \tag{7}$$

Graphically this is a curve that passes through the origin that presents a maximum and decreases asymptotically to zero. The

Fig. 2. The curve of the function ⁱ ⁼ ^ˇ0ADT^ˇ¹ ⁱ ^e^ˇ² ADT ¹⁰,⁰⁰⁰ L^ˇ³ ⁱ when ˇ¹ is larger than 1 and ˇ² is negative.

Ricker model is not monotonously increasing. Consequently there are two values of traffic flow corresponding to each value of the number of crashes between zero and the maximum.

The modified Ricker model in the second functional form – hereinafter referred to as Functional Form 2 (FF2) – is calculated by adding to the referred model a parameter, ˇ1, and considering the segment length (L) as a covariate:

$$\lambda_i = \beta_0 A D T_i^{\beta_1} e^{\beta_2} \frac{A D T}{10,000} L_i^{\beta_3}$$ (8)

where, ⁱ = Mean number of ROR crashes in segment i, ADTi =Average Daily Traffic of segment i, Li = length (in km) of segment i, ˇ⁰ = intercept (to be estimated), ˇ¹ = coefficient associated with ADT (to be estimated), ˇ² = coefficient also associated with ADT (to be estimated), ˇ³ = coefficient associated with L (to be estimated).

In this case ADT/10,000 is used to avoid possible overflow during the estimation process. In accordance with empirical evidence and logical considerations described above, ˇ¹ should be larger than 1 and ˇ² should be negative (see Fig. 2).

With this functional form, plotting the number of ROR crashes against ADT, for a given road segment, gives an inverted "U-Shaped" curve. The ROR crash frequency increases with ADT up to a maximum; at ADT levels above the one corresponding to this maximum number of ROR crashes, the number of ROR crashes will decrease with increasing ADT.

2.3. Goodness-of-fit statistics

As pointed out by Mitra and Washington (2007) and Lord and Park (2008) it is recommended to use more than one goodnessof-fit measure for evaluating models. In this analysis, different methods were used for evaluating the goodness-of-fit and predictive performance of the fitted models, as described bellow:

(1) Akaike Information Criterion (AIC), which measures the goodness of fit of an estimated statistical model and is defined as (Akaike, 1974):

$$AIC = -2\log L + 2p \tag{9}$$

where L is themaximized value ofthe likelihood functionfor the estimated model, and p is the number of parameters in the statistical model. The AIC methodology attempts to find the model that best explains the data with a minimum of free parameters and thus it penalizes models with a large number of parameters. The model with the lowest AIC is considered to be the best model among all available models.

(2) Sum of model deviances (G2), being the model with the lowest G2, the one with the best fit. If G² is equal to zero, the model has a perfect fit (Washington et al., 2011). This statistic is given as:

$$G^{2} = 2\sum{i=1}^{n} y{i} \ln \left( \frac{y{i}}{\hat{\lambda}{i}} \right)$$ (10)

where yi is the observed number and -^ˆⁱ is the expected number of crashes.

(3) Mean Absolute Deviation (MAD), which gives a measure of the average magnitude of variability of prediction (Washington et al., 2011). Smaller values are preferred to larger values. It is computed using the following equation:

$$MAD = \frac{\sum_{i=1}^{n} \left| \varepsilon_i \right|}{n} \tag{11}$$

where εⁱ corresponds to the difference between the expected and the observed number of crashes.

(4) The proportion of systematic variation explained by the model may be quantified using the Elvik's index. There is systematic variation in the number of crashes whenever the variance exceeds the mean, meaning that overdispersion exists. The amount of overdispersion can be described by:

$$Var(x) = \lambda \times (1 + \mu \lambda) \tag{12}$$

where x is the observed number of accidents; is the average number of crashes; and is the overdispersion parameter. Eq. (11) can also be written as:

$$\mu = \frac{Var(x)/y - 1}{\lambda} \tag{13}$$

Elvik's index is calculated by comparing the over-dispersion parameter in the original dataset (0) with the over-dispersion parameter of a fitted model (m) (Elvik, 2007). The proportion of systematic variation explained by the model is calculated by the ratio between the two over-dispersion parameters previously determined, being Elvik's indexes closer to one associated with better-fitted models.

The methods described above evaluate the overall fit and performance of the models. However, it is important that the model predicts well fitted values for all practical values of a variable, not only a global good fit. The CURE (CUmulative REsiduals) method, as described in Hauer and Bamfo (1997), is particularly useful in this circumstance. CURE is a graphical method that was also used for comparing the two models performance.

This method consists of plotting the CURE for each independent variable of interest (usually ADT). The difference between the number of observed and predicted crashes for an entity (in this case, road section) and time period is called a 'residual'. By examining residuals one should be able to learn how well the functional form has been chosen and how well the model fits the data. To generate a CURE plot, road sections are sorted in ascending order by ADT. Because of the random nature of crash counts, the CURE line represents a so-called random walk. For a model that fits well in all ranges ofthe variable,the CURE plot should move up and down and oscillates around zero. If the CURE value progressively increases within a range of values of the independent variable, it means that, within that range,the model predicts fewer crashes than have been observed. In opposition, a decreasing CURE line indicates that, in that range, fewer crashes have been observed than are predicted by the model. If the model fit is poor, it can possibly be improved by a suitable modification of the functional form (Hauer, 2004)

2.4. Estimation of confidence intervals

Once goodness of fit is established for each model, it is of interest to provide confidence intervals and prediction intervals for models parameters and covariates, respectively. Wood (2005) has proposed a method for estimating the confidence intervals for the mean response (-), for the gamma mean (m), and the predicted response (y) at a new segment having similar characteristics as the segments used in the datasetfrom which the model was developed.

The following table gives the equations for calculating the 95% confidence and prediction intervals for NB models. In this table, ⁱ is the logarithm of the estimated mean response ⁱ for segment i, while is the inverse dispersion parameter estimated during the fitting process.

The primary purpose of this paper is to evaluate whether or not there will be any significant difference in modelling ROR crashes with different functional forms. To examine this difference, confidence intervals for the Poisson mean, the gamma mean and the predicted response were calculated for FF1 and FF2, in addition to the methods used for evaluating the goodness-of-fit and predictive performance of the models.

2.5. Modelling

The modelling procedure was accomplished using a 3-stage process:

(1) In the initial step, the two models were estimated using R (version 2.15.0) (R Development Core Team, 2011). "MASS" (Venables and Ripley, 2002) R package was used.
(2) The confidence interval for Poisson mean ⁱ, gamma mean m, and the prediction interval y at a new road section with the same traffic flow characteristics was then calculated using Table 1 equations for FF1 and FF2 models. Var(i) is calculated using XI−1XT where I ⁻¹ represents the variance–covariance matrices and X is the matrix containing the observed values. The variance–covariance matrices were calculated by R (R Development Core Team, 2011).
(3) The differences between confidence interval widths calculated by the two models were then estimated.

2.6. Simulated crash data

A Monte Carlo simulation study was used to verify the results produced using the empirical data. Using simulation, we can assess the general performance of the proposed models (Geedipally et al., 2012).

The simulation design was accomplished using the simulation protocol used by Geedipally and Lord (2010):

(1) The traffic flow and segment length from the dual carriageway dataset are used as the independent variables.
(2) Generate the "true" mean for ROR injury crashes at each segment using the two functional forms considered (FF1 and FF2). The parameters are defined prior to the simulation. These parameters are directly taken from the dual carriageway NB models in Tables 4 and 5 forRORinjury crashes on two-lane unidirectional sections of the Portuguese motorway net. Similarly, using the defined parameters, generate the "true" confidence intervals with the dual carriageway NB models.
(3) Simulate ROR injury crash counts for each site and for each functional form. Following the simulation, the parameters and confidence intervals are re-estimated by the two models (as described above in the modelling procedure).
(4) Steps (2) and (3) are repeated for 30 times (as in Lord, 2006) to obtain statistically reliable estimates. The average values for the

Table 1 95% confidence and prediction intervals (Wood, 2005).

Parameter	Intervals
	ˆ ˆie1.96 Var( ˆi) i , e1.96√Var( ˆi)
m	ˆi − 1.96 ˆi + 1.96 ˆ2 ˆ2 ˆ2 ˆ2 i Var( ˆi) + Var( ˆi)+ i Var( ˆi) + Var( ˆi)+ ˆ 2 ˆ 2 max 0, i i i i ,
y	19 ˆ2 ˆ2 Var( ˆi)+ ˆi + √ i Var( ˆi) + ˆ 2 ˆi 0, i i +

parameters and confidence intervals are finally used. For each functional form, the average confidence intervals obtained are then compared with the true intervals.

3. Data description

This section describes the characteristics of the single and dual carriageway datasets. In both cases, data characteristics and summary statistics for the ROR injury crash data are presented. Only single-vehicle ROR crashes involving roadside features were used in this study. Roadside features include impact attenuators, bridge parapet ends, bridge rails, guardrail fences, guardrail ends, median barriers, highway traffic sign posts, overhead sign supports, utility poles, culverts, curbs, ditches, embankments, trees, and other fixed objects. Multi-vehicle crashes were excluded: head-on and sideswipe opposite crashes are not included in ROR crashes on single carriageway roads.

The first dataset contained ROR injury crash and traffic data collected in a four-year period (2002–2005) at 608 rural single carriageway segments in Portugal. This dataset comprises mostly rural single carriageway segments situated in various regions across Portugal. The general speed limit on these segments is 90 km/h; nevertheless, on some short sections lower speed limits are posted (down to 60 km/in some cases). Access to trunk roads is only possible through grade separated interchanges; on other road categories there are signalized or stop-controlled at grade intersections. Cross section characteristics are not uniform, depending on the road category: carriageway widths vary between 6 and 9 metres; and shoulders of up to 2.5 m may be paved or non-paved. Table 2 provides relevant descriptive statistics for key explanatory variables and the crash data for single carriageway data.

The second dataset contained run-off-the-road injury crashes that occurred on two-lane unidirectional sections ofthe Portuguese freeway network in a four-year period (2007–2010). This dataset comprises dual carriageway motorway segments situated in various regions across Portugal. All segments have full access control, two lanes per carriageway and paved shoulders (with widths of less than 2.5 m and 4.0 m for left and right shoulders, respectively). Access to and from the motorway is only possible through interchange ramps. Table 3 provides relevant descriptive statistics for key explanatory variables and the crash data for dual carriageways.

4. Modelling results

This section describes the results of the comparison analysis between the two functional forms crash models using empirical

Table 2 Descriptive statistics of independent variables and crash data for single carriageway road segments.

Variable	Minimum	Maximum	Mean (Std. dev.)	Total
Segment length (km)	0.17	62.67	13.47 (8.70)	8189.21
ADT (vehicles/day)	283	47900	6927.82 (6700.03)	–
ROR injury crashes	0	129	11.07(11.54)	6733

Table 3 Descriptive statistics of explanatory variables and crash data for unidirectional dual carriageway road segments.

Variable	Minimum	Maximum	Mean (Std. dev.)	Total
Segment length (km)	0.32	33.28	6.56 (5.88)	3996.00
ADT (vehicles/day)	1056	78824	10689.61 (9901.45)	–
ROR injury crashes	0	45	5.00 (5.79)	3047

Table 4 Estimates of FF1 model coefficients.

Type of carriageway	ln(ˇ0)	ˇ1	ˇ2	˛
Single	−4.5555	0.5413	0.9261	2.6060
Dual	−6.4932	0.7246	0.8745	2.2791

Table 5 Estimates of FF2 model coefficients.

Type of carriageway	ln(ˇ0)	ˇ1	ˇ2	ˇ3	˛
Single	−6.4685	0.8112	−0.4537	0.9046	2.7362
Dual	−9.9360	1.1612	−0.4096	0.8472	2.4351

data as well as the output of the Monte Carlo simulation study. The confidence intervals were calculated for injury crashes for the single and dual carriageway empirical data and for the simulation study.

4.1. Empirical data

Table 4 contains the parameter estimates of FF1 models with their associated standard errors for single and dual carriageway road sections. Likewise, Table 5 provides the same parameter estimates of FF2 models.

With FF1, the coefficient (ˇ1) shows that the ROR injury crashes increase at a decreasing rate as traffic flow increases, whereas with FF2, at first ROR injury crashes increase with ADT, until reaching a maximum and thereafter decreases with a diminishing rate. The dispersion parameter (˛) reveals that the ROR injury crashes when modelled with FF2 are slightly more dispersed than the same crashes modelled with FF1.

Table 6 gives the goodness-of-fit statistics for each type of model. This table shows that models with FF2 are found to provide a better statistical fit than models with FF1.

Plots of the models predictions and CURE plots are provided in Figs. 3 and 4 for single and dual carriageways. In this example,

Table 6 Goodness-of-fit statistics.

Type of carriageway	Functional form	AIC	G2	MAD	Elvik's index
Single	FF1	3799.4	3028.415	5.498	0.681
	FF2	3776.9	2923.452	5.406	0.696
Dual	FF1	2943.0	1732.323	2.661	0.572
	FF2	2923.4	1674.070	2.625	0.600

Fig. 3. Comparison of crash prediction models for ROR injury crashes on single carriageways with (left) FF1 and (right) FF2.

a 1 km section length (between interchanges) was considered in model prediction plots; overall, different values for segment length generate shapes similar to the ones presented in this document. In each figure two graphs (on the left sides) correspond to the FF1; and the two graphs on the right sides correspond to the FF2 models.

The CURE plots for the crash prediction models, generated with FF2, consistently showedamodelfit superior to that ofthe FF1models. The model with FF2 was significantly less biased throughout the entire range of ADT; it showed a moderately amplified random walk without excursions outside of the 2 limits. In contrast to the model withFF2, FF1models showedmuchmore bias and frequently departed outside of the 2 boundaries.

For both the single and dual carriageway road segments data, the models with FF2 provided better estimates of expected crash frequencies at different levels of ADT in both data sets. It is worth to note that the influence of access control in crash frequency prediction models for single carriageways was also analyzed. In both cases the modified Ricker model shows better quality fit, similarly as with the aggregate models presented above.

Fig. 5 shows the confidence intervals for the Poisson mean '-This figure illustrates that the width of the intervals for Poisson mean '-' for ROR injury crashes on 1 kilometre segments of single and dual carriageways (i.e., the distance between the mean response and the upper confidence interval boundary) predicted by the FF1 models is slightly wider than the width predicted by FF2 models, for flows more than 20,000 vehicles/day and 40,000 vehicles/day, for single and dual carriageways, respectively.

Fig. 6 shows the confidence intervals for the gamma mean 'm' and for the predictive response 'y'.

The figure illustrates that the confidence interval for the gamma mean 'm' for ROR injury crashes in single carriageway road segments predicted by FF2 model is clearly narrower than the FF1 model for flows greater than 15,000 vehicles/day. For dual carriageway road segments, the interval predicted with the FF1 model is much wider than the one predicted with the FF2 model for flows above 30,000 vehicles/day.

Fig. 6 also shows the confidence intervals for the predictive response 'y' for FF1 and FF2 models. The figure illustrates that predictions with FF2 models have smaller confidence intervals than those models with FF1 for more than half of the time for single carriageways and, atleast,two thirds ofthe time for dual carriageways. According to Lord (2008), the predicted response and its associated variance are particularly important values to be computed, since in most cases, analysts will apply the model to an observation (or site) that was not used for developing the statistical models.

In sum, the comparison analysis showed that using two FF1 models to predict ROR injury crashes, for single and dual carriageway road segments, most often provides larger confidence intervals than using models with FF2, especially for the gamma mean values and predictive response.

Fig. 4. Comparison of crash prediction models for ROR injury crashes on dual carriageways with (left) FF1 and (right) FF2.

4.2. Simulation

Fig. 7 shows the comparison between the true crash means and the crash means predicted by FF1 and FF2 models. Values predicted with FF2 models are slightly closer to the true crash mean than those with FF1 models.

Fig. 8 shows that 95% confidence intervals predicted with FF2 model predicts are slightly closer to true intervals for Poisson mean (-), for ADT higher than 40,000 vehicles/day.

Fig. 9 illustrates the comparative analysis of the 95% confidence intervals for upper boundaries for the Gamma mean (m) and predictive response (y). Intervals predicted with the Gamma mean FF2 model are closer to true intervals. Nevertheless, for predictive responses, FF1 model predicted intervals are closer to true intervals than those from the FF2 model, for flows in excess of 40,000 vehicles/day.

The results of the analysis with the empirical and simulated data show that FF2 model has higher variance of mean and predicted response. However, the CURE plots for the crash prediction models consistently showed the models with FF2 provided better estimates of expected crash frequencies at different levels of ADT. Accordingly, the FF2 models confidence intervals will be narrower (with a very significant difference in case of the gamma mean values and predictive response). Although FF1 confidence intervals for predictive response (y) are closer to the "true intervals" for higher ADT values, the crash means predicted by a model with FF2 are closer to the true values and confidence intervals for Gamma mean (m) with this functional form are also closer to the "true intervals". Thus, this study supports the idea that ROR injury crashes should be modelled with a functional form similar to FF2.

5. Discussion

The results of the analysis presented above raise two important issues that deserve further discussion.

First, crash prediction models traditionally relate crash occurrence to exposure measures such as AADT, or ADT. Exposure may be defined as the number of traffic events in which there is a reasonable chain of event that could lead to a collision between road users or a collision with dangerous obstacles. Currently used, volume-based exposure measures, suffer many limitations. In case of aggregate measurements of traffic volume, such as ADT or AADT,temporal variations in demand (seasonal, monthly, daily and hourly) are not considered. Actually, a road with intense traffic flow during peak periodsmight have a differentRORcrash potentialthan other with the same ADT but with traffic flow more evenly spread out during daytime. Also, the distribution of weekday to weekend traffic volume may vary from one location to another or from daytime to night-time. To the extent that the actual hourly volume is an importantfactor in explaining the number of crashes,the hourly volume can accurately account for this effect in a way that ADT or other aggregate exposure measures cannot (Qin et al., 2006)

Furthermore, the fact that not all vehicles are interacting unsafely is not explicitly taken into account (El-Basyouny and

Fig. 5. 95-Percentile confidence intervals for the Poisson mean (-) for a 1 km segment of single and dual carriageways.

Sayed, 2013). In fact, if there are varying degrees of interaction between vehicles within a traffic stream, then, for example, the occurrence of crashes involving a single vehicle may have a different likelihood from those involving multiple vehicles.

These model limitations, due primarily to the cost of collecting and combining more detailed traffic and crash data using traditional systems, have been an obstacle to the development of detailed statistical models of the relationship between traffic flow and ROR crash frequency. Recent developments in traffic data collection and management systems allow for the development and test of other traffic distribution statistics that may complement ADT. Resorting to very small time periods, however, may generate new data consistency problems, as high accuracy in accident occurrence time reporting and full synchronization of traffic data collection stations and accident reporting teams are required.

Mensah and Hauer (1998) investigated what problems arise when data is in a form such that the causal link between the dependent variable (the effect) and the independent variables (causes) is indirect in the sense that the effect is observed over a "long" period of time during which the cause has assumed widely different values, but for which period only the average value of the cause is known. They concluded that by using additional information about the distribution of the hourly flows on such road segments over the year, the bias due to traffic flow averaging may be reduced. They also analyzed the so called function averaging, a problem that arises when a single crash prediction model is estimated for data to which two or more distinct models apply. They recommended to avoid this problem by, if possible, attempting to fit separate models to say, daytime and nighttime conditions, instead of a single model. Also Garbarino et al. (2001) highlighted the importance of daytime and nighttime conditions by demonstrating that the distribution of road crashes may be correlated with the sleep propensity curve.

As referred by Kononov et al. (2012) it would also be extremely valuable to know how safety varies with the volume-to-capacity (V/C) ratio and what V/C ratio provide the lower accident rates, in particular to more thoroughly assess the safety implications of highway improvements designed to increase capacity. In a study on freeways by Lord et al. (2005b) it was found that as density and V/C ratio increase, the number of single-vehicle crashes decreases.

Garber and Subramanyan (2001) related crashes to lane occupancy and concluded that peak crash rates do not occur during peak flows. In the particular case of ROR crashes this conclusion is intuitive. Kononov et al. (2011) observed that on uncongested freeways the number of crashes increases moderately with an increase in traffic; though, once some critical traffic density is reached, the number of crashes begins to increase at a much higher rate with an increase in traffic.

The second issue is related with the development of crash prediction models with a functional form that consider the fundamental relation between the ROR crash process and traffic flow. Miaou and Lord (2003) pointed out the need to develop new approaches for modelling of vehicle crashes on roadway networks that are based on the logic (e.g., reason, consistency, and coherency), flexibility, extensibility, and interpretability of the functional form.

Fig. 6. 95-Percentile confidence intervals for the Gamma mean (m) and for the predictive response (y) for a 1 km segment of single and dual carriageways.

It is very important to look at the relationship between the variables and whether the distribution or model is logically or theoretically sound, the so-called "goodness-of logic" (Miaou and Lord, 2003). FF2 follows the "goodness-of logic", starting by the logic of "no traffic flows, no crashes" and allowing the relationship between ROR crashes and traffic flows to be nonmonotonic. In fact, the number of ROR crashes increases, peaks and decreases as ADT increases, which implies that, on one hand, similar ROR crash frequencies may occur under very different traffic flow conditions and, on the other, that in some cases the number of ROR crashes may increase when ADT is reduced.

As shown,there is meritin improving the functional form of ADT based crash frequency prediction models, even though ADT may be regarded as too simple an explanatory variable to represent traffic distribution time series.

Fig. 7. Comparison of crash means for a 1 km segment of dual carriageways with (left) FF1 and (right) FF2.

Fig. 8. Comparison of 95-percentile confidence intervals for dual carriageway one kilometre segments. (a) Lower boundary for Poisson mean (-). (b) Upper boundary for Poisson mean (-

6. Summary and conclusions

This paper presents some evidence that the functional form of a ROR injury crash prediction models influences the quality of crash predictions as a function of traffic flow. For that purpose, this paper documented a research study that examined the potential differences between two different functional forms for the referred models.

To accomplish the comparison analysis, ROR injury crash data were collected on single and dual carriageway road sections in Portugal. Then, models were estimated for the two functional forms along with the confidence intervals for the Poisson mean (-), gamma mean (m), and predictive response (y). A simulation study was also conducted to complement the comparative analysis.

The following results were obtained from the analysis:

(1) CURE plots of models with FF2 consistently showed better quality model fit when compared with models with FF1.
(2) There is a clear difference in the prediction of confidence intervals for the gamma mean, and predictive response between models with FF2 and models with FF1. Globally, a model with FF1 provides wider confidence intervals than a model with FF2.
(3) Besides narrower confidence intervals were observed with models with FF2, the crash means and the confidence intervals for Gamma mean (m) predicted with FF2 are closer to the "true crash means" and "true intervals", respectively.
(4) This paper supports the idea that modelling run-off-the-road crashes with a modified Ricker model functional form that shapes the curve between crash frequency and traffic flow in a concave way, shows better quality model fit than a functional form with the assumption that model component for traffic flow continuously increases with traffic flow. This functionalformoffers a reasonably good estimate ofthe relationship between safety and exposure in a ROR crash.

As future work, there are some interesting topics to focus on and explore, namely the collection of hourly traffic volumes and ROR crash data in combination, allowing for the development of crash prediction models to assess how ROR crashes varies with the volume-to-capacity ratio. Also, there is merit in carrying out similar analyses using data collected on other types of crashes, such as multi-vehicle crashes, and other functional forms. Finally, to identify and test readily available AADT parameters (supplementary to the meanADT)that may be used to representthe operational mechanisms that explain improved functional forms for the relationship between traffic and crash frequencies.

Fig. 9. Comparison of 95-percentile confidence intervals for a 1 km section of dual carriageways. (a) Upper boundary for the Gamma mean (m). (b) Upper boundary for predictive response (y).

Acknowledgements

The authors gratefully acknowledge the scholarship SFRH/BD/82228/2011 by the Portuguese Science and Technology Foundation Agency (FCT-Fundac¸ ão para a Ciência e a Tecnologia). The authors express their gratitude to two anonymous reviewers for their valuable remarks that helped to improve the final version of this paper.

Appendix A. Supplementary data

Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.aap.2013.10.034. These data include GoogleMapsTM pinpoints from segment examples of the datasets described in this article.