---
category: literaturenote
citekey: zengapplyingbayesianmultivariatespatio2021
title: Applying a Bayesian multivariate spatio-temporal interaction model based approach to rank sites with promise using severity-weighted decision parameters
authors: "Zeng, Qiang; Xu, Pengpeng; Wang, Xuesong; Wen, Huiying; Hao, Wei"
year: 2021
date: 2021-00-00 2021
doi: 10.1016/j.aap.2021.106190
publication: "Accident Analysis & Prevention"
url: "https://doi.org/10.1016/j.aap.2021.106190"
zotero_key: RP3SSSDN
zotero_storage: V9T4FQBM
collections: transpordi_planeerimine
folder: 001_artiklid
firstAuthor: "Zeng, Qiang"
status: converted
---

Contents lists available at [ScienceDirect](www.sciencedirect.com/science/journal/00014575)
# Accident Analysis and Prevention
journal homepage: [www.elsevier.com/locate/aap](https://www.elsevier.com/locate/aap)

# Applying a Bayesian multivariate spatio-temporal interaction model based approach to rank sites with promise using severity-weighted decision parameters

Qiang Zeng a,b , Pengpeng Xu c , Xuesong Wang d , Huiying Wen a,b , Wei Hao e, \*
- a *School of Civil Engineering and Transportation, South China University of Technology, Guangzhou, 510641, PR China*
- b *Jiangsu Province Collaborative Innovation Center of Modern Urban Traffic Technologies, Southeast University Road #2, Nanjing, 211189, PR China*
- c *Department of Civil Engineering, The University of Hong Kong, Hong Kong, China*
- d *School of Transportation Engineering, Tongji University, Shanghai, 201804, PR China*
- e *School of Traffic and Transportation, Changsha University of Science and Technology, Changsha, 410114, PR China*
### ARTICLE INFO
#### *Keywords:* Multivariate spatio-temporal interaction model Site with promise Bayesian ranking Ranking criterion Crash severity
#### ABSTRACT
Ranking sites with promise is an essential step for cost-effective engineering improvement on roadway traffic safety. This study proposes a Bayesian multivariate spatio-temporal interaction model based approach for ranking sites. The severity-weighted crash frequency and crash rate are used as the decision parameters. The posterior expected rank and posterior mean of the decision parameters are adopted as the statistical criteria. The proposed approach is applied to rank road segments on Kaiyang Freeway in China, which is conducted via programming in the freeware WinBUGS. The results of Bayesian estimation and assessment indicate that incorporating spatio-temporal correlations and interactions into the crash frequency model significantly improves the overall goodness-of-fit performance and affects the identified crash-contributing factors and the estimated safety effects for each severity level. With respect to the ranking results, significant differences are found between those generated from the proposed approach and those generated from the naïve ranking approach and a Bayesian approach based on the multivariate Poisson-lognormal model. Besides, the ranks under the posterior mean criterion are found generally consistent with those under the posterior expected rank criterion.
# **1. Introduction**
Identifying and ranking sites with promise (SWiPs, also known as "crash hotspots" or "crash blackspots") for cost-effective engineering safety improvement on a roadway network, is the first step of roadway safety management process proposed by *Highway Safety Manual* ([AASHTO, 2010](#page-9-0)). The network screening for SWiPs requires that precise ranks of crash risks be estimated at sites of interest which could range from roadway segments ([Chen et al., 2020; Medury and Grembek, 2016](#page-9-0); [Wang et al., 2011, 2021](#page-9-0)), intersections [\(Huang et al., 2009](#page-9-0)) and ramps, to corridors and traffic analysis zones ([Dong et al., 2016;](#page-9-0) [Hao et al.,](#page-9-0) [2020; Huang et al., 2016](#page-9-0)). Due to the differences in the purpose of safety improvement programs, crash risks might be assessed in crash frequencies [\(Dong et al., 2016; Huang et al., 2009\)](#page-9-0), crash rates [\(Miaou and](#page-9-0) [Song, 2005\)](#page-9-0), crash cost ([Aguero-Valverde and Jovanis, 2009](#page-9-0); [Miaou and](#page-9-0)
[Song, 2005\)](#page-9-0), or their excess values [\(Aguero-Valverde and Jovanis, 2009](#page-9-0); [Huang et al., 2009](#page-9-0)).
In past decades, a large number of research efforts have been devoted to develop appropriate methods to rank SWiPs, which can be categorized into three groups: naïve or raw crash risk based approaches ([Deacon et al., 1975\)](#page-9-0), sliding window based approaches [\(Medury and](#page-9-0) [Grembek, 2016](#page-9-0)), and statistical model based approaches ([Miaou and](#page-9-0) [Song, 2005; Wang et al., 2011\)](#page-9-0). The raw crash risks of individual sites calculated without using any data from other sites, may suffer from two well-known statistical problems: the regression-to-the-mean problem ([Huang et al., 2009\)](#page-9-0) and the small area estimation problem [\(Miaou and](#page-9-0) [Song, 2005\)](#page-9-0). The former refers to the tendency that a lower crash risk will probably be observed in a period after a relatively high crash risk is observed. The latter stems from parameter estimation for a small geographical area or a small socio-demographic group in a large area.
*E-mail addresses:* [zengqiang@scut.edu.cn](mailto:zengqiang@scut.edu.cn) (Q. Zeng), [xh150723@163.com](mailto:xh150723@163.com) (P. Xu), [wangxs@tongji.edu.cn](mailto:wangxs@tongji.edu.cn) (X. Wang), [hywen@scut.edu.cn](mailto:hywen@scut.edu.cn) (H. Wen), [haowei@](mailto:haowei@csust.edu.cn) [csust.edu.cn](mailto:haowei@csust.edu.cn) (W. Hao).
0001-4575/© 2021 Published by Elsevier Ltd.
\* Corresponding author.
Specifically, through simulations, [Miaou and Song \(2005\)](#page-9-0) found that extraordinary high and extraordinary low raw crash rates tend to be observed on short segments with equal true crash risks. Technically, the naïve ranking methods are not able to examine the crash dispersion ([Elvik, 2007](#page-9-0)). In addition, without the ability of estimating normal expected crash frequency/rate/cost, they cannot use their excesses as the measure of crash risks [\(Huang et al., 2009](#page-9-0)).
A sliding window based approach can be a scan of road segments with specified sizes of window and moving steps ([Medury and Grembek,](#page-9-0) [2016\)](#page-9-0), or a grid-based scan with a pre-determined grid size to aggregate and generate averages for crash numbers, vehicle miles traveled, and other covariates. The simplest sliding window based approach would independently calculate crash frequency/rate for each window or grid, which also has the limitations of naïve ranking methods. Even though the estimates could be adjusted with reference to background population, a group with similar characteristics, or the estimation results of Poisson-based models, the ranking results are sensitive to the selection of the window/grid size and moving step size which are very difficult to optimize in engineering practice [\(Miaou and Song, 2005](#page-9-0)).
Compared with the approaches of the two types, statistical model based approaches have owned a significant advantage to increase the precision of estimated crash risks by pooling strength from similar sites based on available auxiliary variables, especially in the analysis with small sample sizes ([Dong et al., 2016](#page-9-0)). The Empirical Bayes (EB) method is a classical model based approach for ranking SWiPs and has been found more reliable than naïve ranking methods ([Elvik, 2007](#page-9-0)). However, it is often criticized for ignoring the uncertainty associated with the estimates of parameters in the safety performance function, which probably results in overestimation of the precision [\(Goldstein, 2003](#page-9-0); [Rao, 2003\)](#page-9-0).
As a good alternative to the EB method, the full Bayesian approach is one of the state-of-the-art techniques for site ranking ([Huang et al.,](#page-9-0) [2009\)](#page-9-0). An essential characteristic of the full Bayesian treatment is that the posterior distributions in Bayesian inference are able to quantify the uncertainties specific to model parameters and crash risks ([Miaou and](#page-9-0) [Song, 2005\)](#page-9-0). This ability may make that the ranking results based on the full Bayesian method are significantly different from those based on the EB method, even if they have the same model structure. Besides, it is more flexible to incorporate some of the decision-making functions (e.g., a pre-determined loss function) in the model of the full Bayesian method. Another distinct advantage of the full Bayesian method lies in its capability to explicitly represent the hierarchical structure of crash data ([Huang et al., 2009\)](#page-9-0). Using probabilistic and functional forms, unobserved site-specific and spatio-temporal heterogeneities can be taken into more realistic consideration in the specification of hierarchical models [\(Dong et al., 2016\)](#page-9-0).
The merit of accommodating spatial and temporal correlations under the hierarchical modeling framework has been demonstrated by quantities of studies on the analysis of crash frequency [\(Aguero-Valverde and](#page-9-0) [Jovanis, 2006](#page-9-0); [Li et al., 2019](#page-9-0); [Wen et al., 2019b; Zeng et al., 2020](#page-9-0)), crash rate ([Zeng et al., 2017a](#page-9-0), [2019b](#page-9-0)), and crash severity [\(Meng et al., 2017](#page-9-0); [Xu et al., 2016](#page-9-0); [Zeng et al., 2019a](#page-9-0); [Zhai et al., 2019](#page-9-0)), as well as SWiPs identification [\(Dong et al., 2016](#page-9-0); [Huang et al., 2009;](#page-9-0) [Miaou and Song,](#page-9-0) [2005\)](#page-9-0). For instance, [Miaou and Song \(2005\)](#page-9-0) compared the performance of Bayesian hierarchical models with and without spatial correlation in analyzing crash frequency and ranking SWiPs, and concluded that the inclusion of a spatial term in the model could remove most over-dispersion, significantly improve model fit, and affect the ranking results. [Huang et al. \(2009\)](#page-9-0) conducted an empirical analysis of identifying crash hotspots in Singapore using a ten-year crash dataset. They found that the full Bayesian approach using hierarchical model with temporal correlation significantly outperforms the naïve ranking and EB approaches in the correct identification of crash hotspots. [Dong et al.](#page-9-0) [\(2016\)](#page-9-0) proposed a Bayesian spatio-temporal interaction model based approach for identifying macroscopic hotspots (i.e., the hazardous traffic analysis zones), and compared it with a Poisson-lognormal model and a Bayesian spatial and temporal (without interaction) model. The results showed that the proposed model yields the best performance in both goodness-of-fit and hotspots identification, which highlights the importance of accounting for spatio-temporal correlations and interactions. More recently, [Ziakopoulos and Yannis \(2020\)](#page-9-0) reviewed spatial analyses in the field of roadway safety research, and indicated that spatial approaches can bring about more accurate identification of crash hotspots. Besides, [Gill et al. \(2019\)](#page-9-0) ranked intersections using the expected severity-weighted crash frequency estimated by Bayesian multivariate spatial models. However, only the last study considers crash severity and associated cost in ranking, which are the primary concern in many safety improvement programs. Moreover, [Hauer et al.](#page-9-0) [\(2004\)](#page-9-0) suggested that the use of severity-weighted crash frequency as the decision parameter for SWiPs ranking can lead to the most cost-effective projects.
When crash frequencies/rates at different severity levels are analyzed, multivariate models have often been proposed ([Aguer](#page-9-0)[o-Valverde and Jovanis, 2009; Zeng et al., 2017b](#page-9-0)), because they are able to account for the underlying correlation across crash severities due to the shared safety effects of unobserved factors. Among them, the Bayesian multivariate spatio-temporal interaction (MSTI) models are the state-of-the-art methods ([Cheng et al., 2018](#page-9-0); [Ma et al., 2017; Zeng](#page-9-0) [et al., 2019b](#page-9-0)), that can not only capture the respective spatial and temporal effects and their interactive effects of crash frequency/rate at each severity level, but also the correlations among them in these structured random effects. The empirical analyses indicate that they are superior to other alternatives in terms of model fit and predictive performance. As a consequence, the Bayesian MSTI models may also present superior performance in ranking SWiPs, especially when severity-weighted crash risks are employed as the decision parameters.
In the light of previous studies, this research aims to apply a Bayesian MSTI model to rank SWiPs, using the severity-weighted crash frequency and crash rate as the decision parameters. A crash dataset collected from Kaiyang Freeway, China in 2014 is used for the model development and site ranking. The ranking results based on the proposed model are compared with the counterparts based on the naïve ranking method and a multivariate model without spatio-temporal effects.
# **2. Methodology**
In this section, we first present the formulations of a multivariate Poisson-lognormal (MPLN) model and a MSTI model for analyzing crash frequency by severity in the Bayesian context (Section 2.1). We then introduce three criteria for assessing the performance of the two models in terms of goodness-of-fit (Section [2.2\)](#page-3-0). After that, two severityweighted decision parameters (Section [2.3\)](#page-3-0) and two ranking criteria (Section [2.4](#page-3-0)) are illustrated on the basis of the Bayesian models respectively.
#### *2.1. Model specification*
#### *2.1.1. Multivariate Poisson-lognormal model*
The MPLN regression is one of the widely used statistical methods for simultaneously modeling crash frequencies at various severity levels ([Aguero-Valverde and Jovanis, 2009;](#page-9-0) [Ma et al., 2008](#page-9-0)), since it can account for the common over-dispersion in crash data and correlation across crash severities. In the model, the crash count at roadway site *i* (*i* = 1*,* 2*,*…*, I*) during period *t* (*t* = 1*,* 2*,* …*, T*) at crash severity level *k* (*k* = 1*,* 2*,*…*,K*), Y*i,t,k*, when conditional on its mean *λi,t,k*, is assumed to follow a Poisson distribution:
$$P(Y_{i,l,k} = l) = \frac{(\lambda_{i,l,k})^l}{l!} exp(-\lambda_{i,l,k}), \quad l = 0, 1, 2, ....$$
(1)
As in the previous analysis of crash frequencies on roadway segments ([Wen et al., 2019a\)](#page-9-0), the mean of the Poisson distribution is structured as1 :
$$\ln \lambda_{i,l,k} = \beta_{k,0} + \sum_{m=1}^{M} \beta_{k,m} x_{i,l,m} + \beta_{k,M+1} \ln \epsilon_{i,l} + \theta_{i,k}$$
(2)
in which *xi,t,*1*, xi,t,*2*,* …*, xi,t,M*are a set of risk factors observed at roadway site *i* during period *t* and *βk,*1*, βk,*2*,* …*, βk,M*are their respective parameters specific to severity level *k*. *βk,*0 is a severity-specific constant term. *ei,t*is the total amount of vehicle kilometers traveled (VKT) at roadway site *i* during period *t* and *αk*is its coefficient corresponding to severity level *k*.
The random term *θi, k*which represents the unobserved heterogeneous effect specific to site *i* and severity level *k*, is able to handle the over-dispersion issue. To capture the correlation across crash severities, it is assumed to follow a multi-normal distribution with zero means:
$$\Theta_i \sim N_K(0, \Omega)$$
(3)
$$\mathbf{\Theta}_{i} = (\theta_{i,1}, \theta_{i,2}, \dots, \theta_{i,K})' \tag{4}$$
$$\mathbf{\Omega} = \begin{pmatrix} \delta_{1,1} & \delta_{1,2} & \cdots & \delta_{1,K} \\ \delta_{2,1} & \delta_{2,2} & \cdots & \delta_{2,K} \\ \cdots & \cdots & \cdots \\ \delta_{K,1} & \delta_{K,2} & \cdots & \delta_{K,K} \end{pmatrix}$$
(5)
In the specification of the variance-covariance matrix **Ω**, the diagonal element *δk,k* is the variance of the random term *θi,k*, while the offdiagonal element *δk*1*,k*2 (*k*1 ∕= *k*2) is the covariance between *θi, k*1 and *θi, k*2 .
# *2.1.2. Multivariate spatio-temporal interaction model*
Crash data are usually associated with certain location and time dimensions. Spatial correlations may exist among adjacent road sites due to spatial confounding factors ([Aguero-Valverde and Jovanis, 2006\)](#page-9-0). In a similar vein, unobserved time-dependent factors and those with unspecified time-dependent effects may lead to temporal correlations among successive periods [\(Huang et al., 2009\)](#page-9-0). The spatial or temporal effects may vary across the observed temporal or spatial units, which result in spatio-temporal interaction in crash data [\(Dong et al., 2016](#page-9-0)). Moreover, the spatial and temporal effects and their interactive effects may be correlated among crash frequencies at different severity levels. The MSTI model is an appropriate method to capture all these structured random effects ([Cheng et al., 2018](#page-9-0)). Specifically, a random term *φi,k*with multivariate conditional autoregressive (CAR) priors and a random term *εi,t,k*with a lag-1 dependence are incorporated in the link function to accommodate the spatial and temporal correlations respectively; while the spatio-temporal interaction is formulated as the product of a spatial term and the time trend, i.e., *tϕi,k*. The link function is expressed as:
$$\ln \lambda_{i,t,k} = \beta_{k,0} + \sum_{m=1}^{M} \beta_{k,m} x_{i,t,m} + \beta_{k,M+1} \ln e_{i,t} + \theta_{i,k} + \varphi_{i,k} + \epsilon_{i,t,k} + t \phi_{i,k}$$
(6)
The spatial term *φi,k*and the spatial component *ϕi,k*of the spatiotemporal interaction are specified as [\(Zeng et al., 2019b](#page-9-0)):
$$\Psi_i \sim N_K \left( \overline{\Psi}_i, \Omega_i / \sum_{j \neq i} \omega_{ij} \right)$$
(7)
$$\Phi_i \sim N_K \left( \overline{\Phi}_i, \Omega_a / \sum_{i \neq i} \omega_{i,j} \right)$$
(8)
with
$$\mathbf{\psi}_i = (\varphi_{i,1}, \varphi_{i,2}, \dots, \varphi_{i,K})$$
(9)
$$\overline{\Psi}_{i} = \left(\sum_{i \neq i} \varphi_{i,1} \omega_{i,j}, \sum_{j \neq i} \varphi_{i,2} \omega_{i,j}, \dots, \sum_{j \neq i} \varphi_{i,K} \omega_{i,j}\right) \cdot / \sum_{i \neq i} \omega_{i,j}$$
(10)
$$\Omega_{i} = \begin{pmatrix}
\delta_{1,1}^{i} & \delta_{1,2}^{i} & \cdots & \delta_{1,K}^{i} \\
\delta_{2,1}^{i} & \delta_{2,2}^{i} & \cdots & \delta_{2,K}^{i} \\
\vdots & \vdots & \vdots & \vdots \\
\delta_{K,1}^{i} & \delta_{K,2}^{i} & \cdots & \delta_{K,K}^{i}
\end{pmatrix}$$
(11)
$$\mathbf{\Phi}_i = \left(\phi_{i,1}, \phi_{i,2}, \dots, \phi_{i,K}\right)' \tag{12}$$
$$\overline{\Phi}_{i} = \left(\sum_{j \neq i} \phi_{i,1} \omega_{i,j}, \sum_{j \neq i} \phi_{i,2} \omega_{i,j}, \dots, \sum_{j \neq i} \phi_{i,K} \omega_{i,j}\right) \cdot / \sum_{j \neq i} \omega_{i,j}$$
(13)
and
$$\Omega_{u} = \begin{pmatrix}
\delta_{1,1}^{u} & \delta_{1,2}^{u} & \cdots & \delta_{1,K}^{u} \\
\delta_{2,1}^{g} & \delta_{2,2}^{g} & \cdots & \delta_{2,K}^{g} \\
\cdots & \cdots & \cdots & \cdots \\
\delta_{K,1}^{u} & \delta_{K,2}^{u} & \cdots & \delta_{K,K}^{u}
\end{pmatrix}$$
(14)
In the above equations, *ωi,j*is the adjacent weight between roadway sites *i* and. Herein, the adjacency-based first order proximity structure is used. That is, if sites *i* and *j* are spatially adjacent, *ωi,j* = 1; otherwise, *ωi,j* = 0. **Ω***s*and **Ω***a*are the variance-covariance matrices for *φi,k*and *ϕi,k* respectively. In the matrices, the diagonal elements represent their respective variance and the off-diagonal elements represent their respective covariance of different severity levels.
An independent form for the temporal effects of different severity levels suggested by [Zeng et al. \(2018\)](#page-9-0) is adopted in this study. Specifically,
$$\varepsilon_{i,1,k} \sim N\left(0, \frac{\sigma_{\varepsilon_k}^2}{1 - \gamma_k^2}\right)$$
(15)
$$\varepsilon_{i,t,k} \sim N\left(\gamma_k \varepsilon_{i,t-1,k}, \ \sigma_{\varepsilon_k}^2\right)$$
(16)
where *γk*and *σεk*are the autocorrelation coefficient and standard deviation of the temporal term for injury severity level *k*, respectively.
To highlight the importance of spatial correlation, temporal correlation, and spatio-temporal interaction, the proportions of extra-Poisson variations respectively explained by them for each severity level *k* are calculated as:
$$\alpha_{s,k} = \frac{sd(\mathbf{\varphi}_k)}{sd(\mathbf{\varphi}_k) + sd(\mathbf{E}_k) + sd(\mathbf{E}_k) + sd(\mathbf{E}_k)}$$
(17)
$$\alpha_{p,k} = \frac{sd(\mathbf{e}_k)}{sd(\mathbf{\phi}_k) + sd(\mathbf{e}_k) + sd(\mathbf{A}_k) + sd(\mathbf{\theta}_k)}$$
(18)
$$\alpha_{a,k} = \frac{sd(\mathbf{A}_k)}{sd(\mathbf{\varphi}_k) + sd(\mathbf{E}_k) + sd(\mathbf{A}_k) + sd(\mathbf{\theta}_k)}$$
(19)
with
$$\mathbf{\phi}_k = (\varphi_{1,k}, \varphi_{2,k}, \dots, \varphi_{I,k})$$
(20)
$$\mathbf{\varepsilon}_{k} = (\varepsilon_{1,1,k}, \varepsilon_{1,2,k}, ..., \varepsilon_{I,T,k}) \tag{21}$$
1 Although the specifications of model structures and decision parameters (in Section [2.3](#page-3-0)) are specific to roadway segments (i.e., the type of sites in the empirical analysis), they can be easily extended to other types, such as intersections.
$$\mathbf{A}_{k} = (\phi_{1,k}, 2\phi_{1,k}, ..., T\phi_{I,k})$$
(22)
$$\mathbf{\theta} = (\theta_{1,k}, \theta_{2,k}, \dots, \theta_{I,k}) \tag{23}$$
# *2.2. Model assessment*
The goodness-of-fit of the multivariate models are compared based on the measures of three common assessment criteria [\(Dong et al., 2016](#page-9-0); [Zeng et al., 2019b\)](#page-9-0), that is, deviance information criterion (DIC), mean absolute deviance (MAD), and mean squared prediction error (MSPE). DIC is deemed as a Bayesian equivalent of Akaike information criterion and Bayesian information criterion. It provides a hybrid measure of model complexity and fitting and is defined as ([Spiegelhalter et al.,](#page-9-0) [2002\)](#page-9-0),
$$DIC = \overline{D} + pD \tag{24}$$
where *D* is the posterior mean deviance that can be used as a fitness measure of the model, and *pD* is the effective number of parameters used to measure model complexity. Generally, a model with a lower DIC value is preferred. Particularly, difference greater than 10 indicates that the model is superior ([Spiegelhalter et al., 2005\)](#page-9-0).
With respect to all observations, the MAD and MSPE are defined as follows:
$$MAD = \frac{1}{ITK} \sum_{i=1}^{I} \sum_{k=1}^{T} \sum_{k=1}^{K} |Y_{i,t,k} - \lambda_{i,t,k}|$$
(25)
$$MSPE = \frac{1}{ITK} \sum_{i=1}^{I} \sum_{t=1}^{T} \sum_{k=1}^{K} (Y_{i,t,k} - \lambda_{i,t,k})^{2}$$
(26)
In addition, for the evaluation of model fitting in the crash frequencies at each severity level, the severity-specific MAD and MSPE are also computed:
$$MAD_k = \frac{1}{IT} \sum_{i=1}^{I} \sum_{t=1}^{T} |Y_{i,t,k} - \lambda_{i,t,k}|$$
(27)
$$MSPE_k = \frac{1}{IT} \sum_{i=1}^{I} \sum_{i=1}^{T} (Y_{i,i,k} - \lambda_{i,i,k})^2$$
(28)
### *2.3. Decision parameters*
The decision parameter establishes a link between the calibrated model based on the historical crash data and the SWiPs to be identified. The context under which the rank is to be used (especially the range of treatments to be performed) should be taken into consideration, when choosing an appropriate decision parameter. In the previous studies ([Gill et al., 2019;](#page-9-0) [Huang et al., 2009](#page-9-0); [Miaou and Song, 2005](#page-9-0)), crash frequency and crash rate are the most widely used decision parameters, which provide absolute and relative criteria on assessing roadway safety performance, respectively. Moreover, as mentioned earlier, crash severity is also of great interest in many safety improvement programmes. Therefore, the severity-weighted crash frequency and crash rate are used herein and are formulated in the following equations respectively:
$$P_i^f = \sum_{k=1}^{T} \sum_{k=1}^{K} c_k \lambda_{i,t,k}$$
(29)
$$\mathbf{P}_{i}^{r} = \left(\sum_{t=1}^{T} \sum_{k=1}^{K} c_{k} \lambda_{i,t,k}\right) / \sum_{t=1}^{T} e_{i,t}$$
(30)
Noticeably, the formulations are specific to the Bayesian modelsbased approaches. For the naïve ranking approach, the *λi,t,k*at Eqs. (29) and (30) should be replaced by Y*i,t,k*.
# *2.4. Ranking criteria*
Given a decision parameter, we can rank sites according to an absolute standard (i.e., selecting sites based on whether they exceed a predetermined threshold in risk levels) or a relative standard (i.e., selecting a pre-determined number of sites based on their relative risk levels). As in most previous studies [\(Dong et al., 2016; Huang et al., 2009; Miaou](#page-9-0) [and Song, 2005\)](#page-9-0), the relative standard is used in this research.
In the context of Bayesian inference, ranking by posterior mean of the decision parameter which denotes its expected value taken over its posterior distribution given all observed data, is one of the most popular criteria [\(Aguero-Valverde and Jovanis, 2009;](#page-9-0) [Dong et al., 2016; Huang](#page-9-0) [et al., 2009](#page-9-0)). If the rank order is the main interest, however, ranking by the posterior mean may perform poorly, especially when there are many overlappings in the Bayesian credible intervals, say at the 95 % level, of the decision parameter for the observed sites. Alternatively, ranking by posterior expected rank of the decision parameter may be better ([Shen](#page-9-0) [and Louis, 1998\)](#page-9-0). This is due to the fact that posterior expected rank takes both the mean and variance of the posterior distribution of the decision parameter into account ([Goldstein and Spiegelhalter, 1996](#page-9-0)). Specifically, it is defined as the expected value of the true rank order of the decision parameter taken over its posterior distribution given all observed data. [Miaou and Song \(2005\)](#page-9-0) have proposed this statistical criterion for ranking crash hotspots. In this study, the ranking results based on the two criteria are compared.
#### **3. Data**
A comprehensive database from Kaiyang Freeway in Guangdong Province of China in the entirety of 2014 is used for the model development and SWiP identification. The database is assembled with information from four different data sources: crash incidence, roadway design, traffic flow, and weather condition. The 125-kilometer freeway is split into 154 homogeneous segments. Please refer to [Zeng et al.](#page-9-0) [\(2017c\)](#page-9-0) for the detailed segmentation method.
A total of 692 traffic crashes on the freeway in 2014 are documented in the Highway Maintenance and Administration Management Platform maintained by the Guangdong Transportation Group. According to the criteria defined by the Ministry of Public Security in China, the crash severities are divided into four levels: slight crash, medium crash, severe crash, and very severe crash. Please see [Zeng et al. \(2019a\)](#page-9-0) for the detailed definitions of each severity level which take injury severity, the amount of property damage and casualties into consideration. These levels can be viewed as a more comprehensive metric for crash severity than the common one (i.e., the most severe injury sustained by road users involved in the crash). Due to the rareness of severe and very severe crashes, they are combined with medium crashes into one level: medium and severe (MS) crash. Given their definitions and proportions in the dataset, the cost ratios for slight crash and MS crash are set roughly 1:18. To reveal the impacts of the seasonal variations of traffic and weather attributes on crash occurrence, the crashes are aggregated by roadway segment and crash severity at the seasonal level ([Shankar](#page-9-0) [et al., 2004](#page-9-0)). Four seasons are categorized in Guangdong Province: (1) spring—from February to April; (2) summer—from May to July; (3) autumn—from August to October; and (4) winter—November, December, and January.
We obtain the information on roadway design from the horizontal and longitudinal profile, provided by Guangdong Province Communication Planning and Design Institute Co., Ltd. Contributing to the design standards specified in the [Design Specification for Highway Alignment](#page-9-0) [\(2006\),](#page-9-0) many roadway-specific attributes are constant along the freeway. Consequently, vertical grade, horizontal curvature, and the presence of bridges and ramps, are used as roadway factors.
The traffic data are collected from the Guangdong Freeway Network Toll System which groups motor vehicles into five classes, based on their head height, axis number, wheel number, and wheelbase. The detailed
vehicle categorization criteria can be found in [Zeng et al. \(2017c\).](#page-9-0) The toll system records daily traffic volume of each vehicle class. We calculate the normalized seasonal traffic volume (NSTV) as the weighted sum of the daily traffic volumes of the five classes within a certain season. The weights of the five vehicle classes are 1, 1.5, 2, 3, and 3.5, respectively, upon the recommendation of the Guangdong Transportation Department. The VKT is calculated as the product of segment length and NSTV. We also compute the proportion of each vehicle class to represent the traffic composition. Due to the perfect collinearity of the vehicle-class proportions, that of vehicles in class 1, *Veh\_1*, is set as the reference case. The proportions of vehicles in classes 2 and 5 are also excluded, due to their mutual correlation and significant correlation to the proportion of vehicles in class 3.
The daily data on precipitation and wind speed from three countylevel weather stations along the freeway, are extracted from the Meteorological Information Management System which is managed by the Guangdong Climate Center. We assign the roadway segments to the nearest weather stations according to the Euclidean distances between the segments and weather stations. The seasonal average daily precipitation and wind speed are calculated and used as covariates in this study. Table 1 summarizes the definitions and descriptive statistics of the variables used in the model development and SWiPs ranking.
### **4. Model estimation and result analysis**
# *4.1. Model estimation and diagnostics*
The full Bayesian analysis is conducted in WinBUGS ([Lunn et al.,](#page-9-0) [2000\)](#page-9-0), which adopts the Metropolis-Hastings algorithm to sample from the unnormalized posterior distribution, producing Markov-chain-Monte-Carlo (MCMC) runs for each parameter of interest. A specification of the prior distributions of the parameters is required to obtain the Bayesian estimates. In the absence of sufficient knowledge, we use uninformative prior distributions for these parameters. To be specific, a diffused normal distribution, *Normal*(0*,* 104), is used as the priors of the regression coefficients *βk,m*(*k* = 1*,* 2; *m* = 0*,* 1*,* …*,* 9), and the severity-specific autocorrelation coefficients *γk*(*k* = 1*,* 2) for the temporal correlation. A Wishart distribution, *W* (( 1 0 0 1 ) *,* 2 ) , is used as the priors of the inverse of the variance-covariance matrices (i.e., **Ω**− 1, **Ω**− 1 *s* , and **Ω**− 1 *a* ). A gamma distribution *gamma*(0*.*001*,* 0*.*001) is used as the priors of the precisions of temporal effects, 1*/σ*2 *εk* (*k* = 1*,* 2) ([Zeng](#page-9-0)
**Table 1** Descriptive Statistics of Variables in the Model.
| Variables | Description | Mean | S.D. | Min. | Max. |
|---------------|--------------------------------------------------------------------|------|------|------|-------|
| Slight crash | Seasonal slight crash count | 0.58 | 0.90 | 0 | 4 |
| MS crash | Seasonal medium and severe
crash count | 0.55 | 0.85 | 0 | 6 |
| VKT | Seasonal vehicle kilometers
traveled (103
km•pcua
) | 4046 | 1541 | 893 | 12209 |
| Veh_3 | The proportion of vehicles in
class 3 (%) | 21.8 | 0.97 | 19.6 | 23.4 |
| Veh_4 | The proportion of vehicles in
class 4 (%) | 6.58 | 0.35 | 5.90 | 7.37 |
| Curvature | The horizontal curvature of the
freeway segment (0.1 km− 1
) | 1.77 | 1.27 | 0 | 4.35 |
| Gradient | The vertical grade of the
freeway segment (%) | 0.74 | 0.57 | 0 | 2.91 |
| Bridge | A part of a bridge: Yes = 1, No
= 0 | 0.50 | 0.50 | 0 | 1 |
| Ramp | Presence of ramp: Yes = 1, No
= 0 | 0.21 | 0.41 | 0 | 1 |
| Wind speed | Average seasonal wind speed
(m/s) | 2.49 | 0.85 | 1.67 | 3.92 |
| Precipitation | Seasonal average daily
precipitation (m/s) | 4.47 | 2.75 | 0.64 | 8.94 |
a pcu: passenger car units.
[et al., 2019b\)](#page-9-0). The posterior expected rank of the decision parameters are obtained via the *rank* function available in WinBUGS. With reference to the convergence speed and the size of dataset, a chain of 50,000 MCMC simulation iterations is set up for each model, and the first 30, 000 iterations are excluded as a burn-in. The MCMC trace plots for the model parameters and the Gelman-Rubin convergence statistics are used to monitor the simulations convergence.
The results of model estimation and comparison are shown in Table 2., where only the factors with statistically significant (at least at the 90 % credibility level) effects on slight or MS crash frequency are included. According to the results, we can find that the MSTI model yields a smaller DIC value (= 2363) than that of the MPLN model (= 2381). It indicates that the former fits the crash data better than the latter, which is further demonstrated by the lower values of MAD and MSPE for each crash severity levels and all observations of the MSTI model. These results are reasonable, because there are considerable spatio-temporal correlations and interactions for both slight and MS crashes in the dataset, as suggested by the Bayesian estimates of the standard deviations of spatial terms (*sd*(**φ**1) = 0.33, *sd*(**φ**2) = 0.16), temporal terms (*sd*(**ε**1) = 0.37, *sd*(**ε**2) = 0.32), and spatio-temporal interactive terms (*sd*(**A**1) = 0.24, *sd*(**A***k*) = 0.24). The results conform to a large number of previous studies which have shown that explicitly accounting for spatial correlation, temporal correlation, and spatiotemporal interaction can mitigate the model misspecification and improve model fit ([Aguero-Valverde and Jovanis, 2006;](#page-9-0) [Dong et al.,](#page-9-0) [2016; Wen et al., 2019b; Zeng et al., 2019b](#page-9-0)). Besides, the Bayesian estimates of α*s, k*, α*p,k*, and α*a,k*(*k* = 1*,* 2) indicate that spatial correlations, temporal correlations and spatio-temporal interactions account for 28 % (20 %), 31 % (40 %), and 21 % (19 %) of extra-variability in the count distributions of slight crashes (MS crashes), respectively.
Comparison of the coefficient estimates in the two models shows that the effect of *Wind speed* on slight crashes is significant only in the MPLN model, although the safety effects of other significant factors are generally consistent in the two models. These results may attribute that omitting spatio-temporal correlations or interactions might
**Table 2** Parameter estimation and model comparison results.
| | Multivariate Poisson-lognormal
model | | Multivariate spatio-temporal
interaction model | | |
|---------------|-----------------------------------------|----------------------|---------------------------------------------------|----------------------|--|
| | Slight crash | MS crash | Slight crash | MS crash | |
| Constant | − 7.81 (2.22)
** | − 11.39 (2.73)
** | − 6.55 (2.19)
** | − 11.54 (1.46)
** | |
| VKT | 1.03 (0.17)** | 1.23 (0.23)** | 1.06 (0.19)** | 1.30 (0.13)** | |
| Veh_3 | − 0.17 (0.05)
** | — | − 0.13 (0.05)
** | — | |
| Curvature | — | 0.12 (0.06)** | — | 0.13 (0.05)** | |
| Gradient | 0.26 (0.11)** | 0.29 (0.11)** | 0.24 (0.12)* | 0.31 (0.11)** | |
| Wind speed | 0.42 (0.07)** | — | — | — | |
| Precipitation | — | 0.08 (0.02)** | — | 0.08 (0.02)** | |
| sd(θk) | 0.34 (0.12)** | 0.27 (0.10)** | 0.23 (0.12)** | 0.17 (0.10)** | |
| sd(φk) | — | — | 0.33 (0.12)** | 0.16 (0.08)** | |
| sd(εk) | — | — | 0.37 (0.17)** | 0.32 (0.14)** | |
| sd(Ak) | — | — | 0.24 (0.08)** | 0.15 (0.06)** | |
| αs, k | — | — | 0.28 (0.10)** | 0.20 (0.10)** | |
| αp,k | — | — | 0.31 (0.13)** | 0.40 (0.15)** | |
| αa,k | — | — | 0.21 (0.06)** | 0.19 (0.07)** | |
| MADk | 0.60 | 0.59 | 0.56 | 0.57 | |
| MSPEk | 0.65 | 0.59 | 0.59 | 0.55 | |
| MAD | 0.60 | | 0.57 | | |
| MSPE | 0.62 | | 0.57 | | |
| D | 2314 | | 2227 | | |
| pD | 67 | | 136 | | |
| DIC | 2381 | | 2363 | | |
| | | | | | |
\* Significant at the 90 % credibility level. \*\* Significant at the 95 % credibility level.
underestimate the variance of coefficients for certain factors, which have been founded in many previous studies [\(Aguero-Valverde and](#page-9-0) [Jovanis, 2006](#page-9-0); [Ahmed et al., 2011](#page-9-0)). Besides, some factors (e.g., *Veh\_3*, *Curvature*, and *Precipitation*) have significant effects on either slight crashes or MS crashes. Overexposure of a specific crash severity may uncover some undetected deficiencies and provide a deep insight into the safety performance of roadway entities ([Ye et al., 2013; Zeng et al.,](#page-9-0) [2017b\)](#page-9-0).
Specifically, the estimation results of the MSTI model indicate that *VKT* is significantly positively associated with the crash frequencies at the two severity levels. The estimated coefficients are near 1.0, which imply that there may be linear relationships between vehicle kilometers traveled and slight and MS crash frequencies. More slight and MS crashes tend to occur roadway segments with higher vertical grades. In addition, slight crashes are more likely to occur in the traffic conditions with less vehicles in class 3; MS crashes are more likely to occur on freeway segments with greater horizontal curvature and seasons with higher precipitation. These findings are generally consistent and confirm to engineering intuitions and the counterparts in previous studies [\(Wen et al., 2019a;](#page-9-0) [Zeng et al., 2017c](#page-9-0)).
## *4.2. Comparison of ranking results*
The ranking results of the approaches based on Bayesian multivariate models and naïve approach are compared with severity-weighted crash frequency and rate used as the decision parameters. As described above, if the ranks are the primary interest, ranking by the posterior expected rank of decision parameters is recommended for Bayesian models-based approaches. Meanwhile, the results of model assessment indicate that the MSTI model outperforms the MPLN model. Therefore, ranking results by the posterior expected rank based on the MSTI model are used as baselines in the comparisons. In particular, Figs. 1(a) and [2](#page-6-0) (a) show the


**Fig. 1.** Comparison of severity-weighted crash frequency ranking results based on the alternative approaches.


**Fig. 2.** Comparison of severity-weighted crash rate ranking results based on the alternative approaches.
comparisons of ranking results for all roadway segments, using the severity-weighted crash frequency and rate, respectively. To get a closer view of how the top-ranked sites compared, [Figs. 1\(](#page-5-0)b) and 2 (b) zoom in the comparisons on the top 15 (about 10 %) road segments (based on the PER criterion and the MSTI model).
From these figures, one can observe that the ranks produced by the naïve ranking approach are significantly different from the counterparts generated by the Bayesian MSTI model, for both decision parameters. For instance, under the severity-weighted crash frequency, road segment No. 49 ranks 9 by the Bayesian MSTI model, but ranks 40 by the naïve ranking approach; under the severity-weighted crash rate, road segment No. 23 ranks 10 by the Bayesian MSTI model, but ranks 148 by the naïve ranking approach. The substantial differences are probably due to the big disparity in cost ratios that give considerably more weights to severe crashes than slight and medium crashes in the decision parameters, as well as the highly stochastic and sporadic characteristics of severe crashes [\(Miaou and Song, 2005\)](#page-9-0). In other words, the "small area" issue of the traffic crashes estimation becomes more conspicuous because smaller domains are used (in this case, by crash severity) and considerably greater weights are given to the rarest events in the decision parameters.
The figures show that the ranks produced by the two Bayesian multivariate models are more consistent than the ranks produced by the naïve ranking approach, especially for the top-ranked sites. Nonetheless, there are also significant differences in the ranks of certain roadway segments. For instance, under the severity-weighted crash frequency, road segment No. 90 ranks 94 by the Bayesian MSTI model, but ranks 66 by the MPLN model; under the severity-weighted crash rate, road segment No. 99 ranks 104 by the Bayesian MSTI model, but ranks 46 by the MPLN model. The findings are generally consistent to those in the
previous studies ([Dong et al., 2016; Miaou and Song, 2005](#page-9-0)), which have found that accommodating spatio-temporal correlation and interaction in Bayesian models could significantly affect the ranking results, due to their abilities of increasing precision of parameter estimation by borrowing strengths across adjacent sites and successive periods.
Based on the MSTI model, the ranking results under the two statistical criteria——posterior mean and posterior expected rank, are also compared. The comparisons for severity-weighted crash frequency and rate are shown in Figs. 3 and 4 respectively, where Figs. 3(a) and [4](#page-8-0) (a) correspond to all road segments and Figs. 3(b) and [4](#page-8-0) (b) focus on the topranked road segments. According to these figures, we can see that the ranks produced under the posterior mean criterion are very consistent with the counterparts produced under the posterior expected rank criterion, especially for the top 15 segments. Specifically, ten and six rankpairs of the top-ranked sites are identical under the two criteria for severity-weighted crash frequency and rates respectively. In addition, for each decision parameter, only one top-ranked site under the posterior mean criterion are not among the top 15 sites under the posterior expected rank criterion.
## **5. Conclusions and future research**
This paper advocates a Bayesian approach based on the MSTI model for ranking SWiPs, with severity-weighted crash frequency and rate used as decision parameters. The MSTI model is developed under the Bayesian hierarchical modeling framework, which can account for the spatial correlations across adjacent sites, the temporal correlations across successive periods, the spatio-temporal interactive effects, and the correlations across crash severities in the random effects. A one-year dataset (including crash, traffic, road design, and weather condition information on Kaiyang Freeway in 2014) was collected to demonstrate the approach, which is implemented and assessed by Bayesian methods via the WinBUGS program.
The results of Bayesian estimation indicate that there are significant spatio-temporal correlations and interactions in the crash data. Accounting for these structured random effects significantly affects the estimated safety effects of the observed risk factors and improves the goodness-of-fit performance for each severity level, as compared with the counterparts of the MPLN model. Moreover, the ranking results


**Fig. 3.** Comparison of severity-weighted crash frequency ranking results between PER and PM.


**Fig. 4.** Comparison of severity-weighted crash rate ranking results between PER and PM.
generated from the MSTI model are found different from those from the MPLN model and the naïve ranking approach, which highlights the importance of considering spatio-temporal correlations and interactions in ranking sites. The ranking results under the posterior mean and posterior expected rank of the decision parameters are also compared and found generally consistent, especially for the top-ranked sites.
Although this study has dealt with a number of aspects in applying the proposed model based approach for ranking SWiPs, several limitations of this study and directions for future research are worth mentioning. First, only one-year crash dataset is obtained for the approach demonstration. While the differences in ranking results have been found, the dataset cannot be used to evaluate the performance of the alternative models, as in the previous studies ([Cheng and Wash](#page-9-0)[ington, 2008](#page-9-0); [Dong et al., 2016;](#page-9-0) [Huang et al., 2009\)](#page-9-0). Crash data collected for more years is needed to further justify the strength of the Bayesian MSTI model in sites ranking. Second, the correlations across crash severities in the random effects are found statistically insignificant. A plausible cause is that their correlations mainly stem from the safety effects of the observed factors. Nevertheless, the MSTI model is found to outperform its univariate version (the results are not shown in the text). Further investigation using other field datasets (especially those at macro level) or simulated data is full merited. Besides, it is worth incorporating spatio-temporal correlation and interaction into deep learning techniques and applying to ranking SWiPs in future research ([Bao et al., 2019](#page-9-0); [Cai et al., 2019](#page-9-0)).
#### **Author statement**
Qiang Zeng: Conceptualization, Methodology, Writing -original draft preparation; Pengpeng Xu: Methodology, Validation; Xuesong Wang: Funding acquisition, Writing -review & editing; Huiying Wen: Resources, Supervision; Wei Hao: Funding acquisition,Supervision.
#### **Declaration of Competing Interest**
The authors report no declarations of interest.
#### **Acknowledgments**
This work was jointly supported by the International Science &
Technology Cooperation Program of China (No. 2017YFE0134500); the National Natural Science Foundation of China (Nos. 71801095 and 51808057); Hunan Provincial Natural Science Foundation of China (No. 2019JJ30026); and the Fundamental Research Funds for the Central Universities (No. 2020ZYGXZR007).
#### **References**
- [AASHTO, 2010. Highway Safety Manual, 1st edition](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0005).
- [Aguero-Valverde, J., Jovanis, P.P., 2006. Spatial analysis of fatal and injury crashes in](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0010) [Pennsylvania. Accid. Anal. Prev. 38 \(3\), 618](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0010)–625.
- [Aguero-Valverde, J., Jovanis, P.P., 2009. Bayesian multivariate Poisson lognormal](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0015) [models for crash severity modeling and site ranking. Transp. Res. Rec. 2136, 82](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0015)–91.
- [Ahmed, M., Huang, H., Abdel-Aty, M., Guevara, B., 2011. Exploring a Bayesian](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0020) [hierarchical approach for developing safety performance functions for a](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0020) [mountainous freeway. Accid. Anal. Prev. 43 \(4\), 1581](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0020)–1589.
- [Bao, J., Liu, P., Ukkusuri, S.V., 2019. A spatiotemporal deep learning approach for](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0025) [citywide short-term crash risk prediction with multi-source data. Accid. Anal. Prev.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0025) [122, 239](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0025)–254.
- [Cai, Q., Abdel-Aty, M., Sun, Y., Lee, J., Yuan, J., 2019. Applying a deep learning](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0030) [approach for transportation safety planning by using high-resolution transportation](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0030) [and land use data. Transp. Res. Part A: Policy Pract. 127, 71](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0030)–85.
- [Chen, Y., Wu, L., Huang, Z., 2020. Assessing quality of crash modification factors](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0035) [estimated by empirical Bayes before-after methods. J. Cent. South Univ. 27 \(8\),](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0035) [2259](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0035)–2268.
- [Cheng, W., Washington, S., 2008. New criteria for evaluating methods of identifying hot](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0040) [spots. Transp. Res. Rec. 2083, 76](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0040)–85.
- [Cheng, W., Gill, G.S., Ensch, J.L., Kwong, J., Jia, X., 2018. Multimodal crash frequency](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0045) [modeling: multivariate space-time models with alternate spatiotemporal](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0045) [interactions. Accid. Anal. Prev. 113, 159](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0045)–170.
- [Deacon, J.A., Zegeer, C.V., Deen, R.C., 1975. Identification of hazardous rural highway](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0050) [locations. Transp. Res. Rec. 543, 16](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0050)–33.
- [Design Specification for Highway Alignment, 2006. Ministry of Transport of the People](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0055)'s [Republic of China \(JTG D20-2006\)](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0055).
- [Dong, N., Huang, H., Lee, J., Gao, M., Abdel-Aty, M., 2016. Macroscopic hotspots](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0060) [identification: a Bayesian spatio-temporal interaction approach. Accid. Anal. Prev.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0060) [92, 256](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0060)–264.
- [Elvik, R., 2007. State-of-the-art Approaches to Road Accident Black Spot Management](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0065) [and Safety Analysis of Road Networks. Institute of Transport Economics, Oslo](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0065).
- [Gill, G.S., Cheng, W., Zhou, J., Park, V.S., 2019. Comparative analysis of cost-weighted](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0070) [site ranking using alternate distance-based neighboring structures for spatial crash](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0070) [frequency modeling. J. Transp. Saf. Secur. 11 \(1\), 79](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0070)–101.
- [Goldstein, H., 2003. Multilevel Statistical Models, third ed. Hodder Arnold, London.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0075) [Goldstein, H., Spiegelhalter, D.J., 1996. League tables and their limitations: statistical](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0080) [issues in comparisons of institutional performance. J. R. Stat. Soc. Ser. A Stat. Soc.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0080) [159 \(3\), 385](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0080)–409.
- Hao, W., Liu, L., Yang, X., Li, Y., Byon, Y.J., 2020. Reducing CACC platoon disturbances caused by state jitters by combining two stages driving state recognition with multiple Platoons' strategies and risk prediction. IEEE Trans. Intell. Transp. Syst. [https://doi.org/10.1109/TITS.2020.3033436.](https://doi.org/10.1109/TITS.2020.3033436) In press.
- [Hauer, E., Allery, B.K., Kononov, J., Griffith, M.S., 2004. How best to rank sites with](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0090) [promise. Transp. Res. Rec. 1897, 48](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0090)–54.
- [Huang, H., Chin, H.C., Haque, M.M., 2009. Empirical evaluation of alternative](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0095) [approaches in identifying crash hot spots: naive ranking, empirical Bayes, full Bayes](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0095) [methods. Transp. Res. Rec. 2103, 32](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0095)–41.
- [Huang, H., Song, B., Xu, P., Zeng, Q., Lee, J., Abdel-Aty, M., 2016. Macro and micro](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0100) [models for zonal crash prediction with application in hot zones identification.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0100) [J. Transp. Geogr. 54, 248](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0100)–256.
- [Li, Z., Chen, X., Ci, Y., Chen, C., Zhang, G., 2019. A hierarchical Bayesian spatiotemporal](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0105) [random parameters approach for alcohol/drug impaired-driving crash frequency](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0105) [analysis. Anal. Methods Accid. Res. 21, 44](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0105)–61.
- [Lunn, D.J., Thomas, A., Best, N., Spiegelhalter, D., 2000. WinBUGS-a Bayesian modelling](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0110) [framework: concepts, structure, and extensibility. Stat. Comput. 10 \(4\), 325](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0110)–337.
- [Ma, J., Kockelman, K.M., Damien, P., 2008. A multivariate Poisson-lognormal regression](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0115) [model for prediction of crash counts by severity, using Bayesian methods. Accid.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0115) [Anal. Prev. 40 \(3\), 964](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0115)–975.
- [Ma, X., Chen, S., Chen, F., 2017. Multivariate space-time modeling of crash frequencies](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0120) [by injury severity levels. Anal. Methods Accid. Res. 15, 29](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0120)–40.
- [Medury, A., Grembek, O., 2016. Dynamic programming-based hot spot identification](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0125) [approach for pedestrian crashes. Accid. Anal. Prev. 93, 198](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0125)–206.
- [Meng, F., Xu, P., Wong, S.C., Huang, H., Li, Y.C., 2017. Occupant-level injury severity](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0130) [analyses for taxis in Hong Kong: a Bayesian space-time logistic model. Accid. Anal.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0130) [Prev. 108, 297](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0130)–307.
- [Miaou, S.P., Song, J.J., 2005. Bayesian ranking of sites for engineering safety](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0135) [improvements: decision parameter, treatability concept, statistical criterion, and](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0135) [spatial dependence. Accid. Anal. Prev. 37 \(4\), 699](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0135)–720.
- [Rao, J.N.K., 2003. Small Area Estimation. Wiley, New York](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0140).
- [Shankar, V.N., Chayanan, S., Sittikariya, S., Shyu, M.B., Juvva, N.K., Milton, J.C., 2004.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0145) [Marginal impacts of design, traffic, weather, and related interactions on roadside](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0145) [crashes. Transp. Res. Rec. 1897, 156](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0145)–163.
- [Shen, W., Louis, T.A., 1998. Triple-goal estimates in two-stage hierarchical models. J. R.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0150) [Stat. Soc.: Ser. B Stat. Methodol. 60 \(2\), 455](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0150)–471.
- [Spiegelhalter, D.J., Best, N.G., Carlin, B.P., Van Der Linde, A., 2002. Bayesian measures](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0155) [of model complexity and fit. J. R. Stat. Soc.: Series B Stat. Methodol. 64 \(4\),](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0155) 583–[639](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0155).
- [Spiegelhalter, D.J., Thomas, A., Best, N., Lunn, D., 2005. WinBUGS User Manual. MRC](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0160) [Biostatistics Unit, Cambridge, United Kingdom](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0160).
- [Wang, C., Quddus, M.A., Ison, S.G., 2011. Predicting accident frequency at their severity](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0165) [levels and its application in site ranking using a two-stage mixed multivariate model.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0165) [Accid. Anal. Prev. 43 \(6\), 1979](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0165)–1990.
- [Wang, Z., Yu, J., Hao, W., Xiang, J., 2021. Joint optimization of running route and and](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0170) [scheduling for the mixed demand responsive feeder transit witn time-dependent](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0170) [travel times. IEEE Trans. Intell. Transp. Syst. 22 \(4\), 2498](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0170)–2509.
- [Wen, H., Zhang, X., Zeng, Q., Lee, J., Yuan, Q., 2019a. Investigating spatial](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0175) [autocorrelation and spillover effects in freeway crash-frequency data. Int. J. Environ.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0175) [Res. Public Health 16 \(2\), 219](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0175).
- [Wen, H., Zhang, X., Zeng, Q., Sze, N.N., 2019b. Bayesian spatial-temporal model for the](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0180) [main and interaction effects of roadway and weather characteristics on freeway](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0180) [crash incidence. Accid. Anal. Prev. 132, 105249.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0180)
- [Xu, X., Xie, S., Wong, S.C., Xu, P., Huang, H., Pei, X., 2016. Severity of pedestrian injuries](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0185) [due to traffic crashes at signalized intersections in Hong Kong: a Bayesian spatial](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0185) [logit model. J. Adv. Transp. 50 \(8\), 2015](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0185)–2028.
- [Ye, X., Pendyala, R.M., Shankar, V., Konduri, K.C., 2013. A simultaneous equations](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0190) [model of crash frequency by severity level for freeway sections. Accid. Anal. Prev.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0190) [57, 140](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0190)–149.
- [Zeng, Q., Wen, H., Huang, H., Abdel-Aty, M., 2017a. A Bayesian spatial random](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0195) [parameters Tobit model for analyzing crash rates on roadway segments. Accid. Anal.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0195) [Prev. 100, 37](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0195)–43.
- [Zeng, Q., Wen, H., Huang, H., Pei, X., Wong, S.C., 2017b. A multivariate random](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0200)[parameters Tobit model for analyzing highway crash rates by injury severity. Accid.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0200) [Anal. Prev. 99, 184](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0200)–191.
- Zeng, Q., Sun, J., Wen, H., 2017c. Bayesian hierarchical modeling monthly crash counts on freeway segments with temporal correlation. J. Adv. Transp., 5391054 [https://](https://doi.org/10.1155/2017/5391054) [doi.org/10.1155/2017/5391054.](https://doi.org/10.1155/2017/5391054)
- [Zeng, Q., Wen, H., Huang, H., Pei, X., Wong, S.C., 2018. Incorporating temporal](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0210) [correlation into a multivariate random parameters Tobit model for modeling crash](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0210) [rate by injury severity. Transp. A: Transp. Sci. 14 \(3\), 177](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0210)–191.
- [Zeng, Q., Gu, W., Zhang, X., Wen, H., Lee, J., Hao, W., 2019a. Analyzing freeway crash](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0215) [severity using a Bayesian spatial generalized ordered logit model with conditional](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0215) [autoregressive priors. Accid. Anal. Prev. 127, 87](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0215)–95.
- [Zeng, Q., Guo, Q., Wong, S.C., Wen, H., Huang, H., Pei, X., 2019b. Jointly modeling area](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0220)[level crash rates by severity: a Bayesian multivariate random-parameters spatio](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0220)[temporal Tobit regression. Transp. A: Transp. Sci. 15 \(2\), 1867](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0220)–1884.
- Zeng, Q., Wang, X., Wen, H., Yuan, Q., 2020. An empirical investigation of the factors contributing to local-vehicle and non-local-vehicle crashes on freeway. J. Transp. Saf. Secur., 1779422 . In press.
- [Zhai, X., Huang, H., Sze, N.N., Song, Z., Hon, K.K., 2019. Diagnostic analysis of the](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0230) [effects of weather condition on pedestrian crash severity. Accid. Anal. Prev. 122,](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0230) 318–[324](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0230).
- [Ziakopoulos, A., Yannis, G., 2020. A review of spatial approaches in road safety. Accid.](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0235) [Anal. Prev. 135, 105323](http://refhub.elsevier.com/S0001-4575(21)00221-9/sbref0235).