# Bayesian Calibration of the Intelligent Driver Model

Chengyuan Zhang, and Lijun Sun (Corresponding author: Lijun Sun)C. Zhang and
L. Sun are with the Department of Civil Engineering, McGill University,
Montreal, QC H3A 0C3, Canada. E-mail<EMAIL_ADDRESS>(C. Zhang),
<EMAIL_ADDRESS>(L. Sun).The authors would like to thank the McGill
Engineering Doctoral Awards (MEDA), the Institute for Data Valorisation
(IVADO), the Interuniversity Research Centre on Enterprise Networks, Logistics
and Transportation (CIRRELT), Fonds de recherche du Québec – Nature et
technologies (FRQNT), and the Natural Sciences and Engineering Research
Council (NSERC) of Canada for providing scholarships and funding to support
this study.

###### Abstract

Accurate calibration of car-following models is essential for investigating
microscopic human driving behaviors. This work proposes a memory-augmented
Bayesian calibration approach, which leverages the Bayesian inference and
stochastic processes (i.e., Gaussian processes) to calibrate an unbiased car-
following model while extracting the serial correlations of residual. This
calibration approach is applied to the intelligent driver model (IDM) and
develops a novel model named MA-IDM. To evaluate the effectiveness of the
developed approach, three models with different hierarchies (i.e., pooled,
hierarchical, and unpooled) are tested. Experiments demonstrate that the MA-
IDM can estimate the noise level of unrelated errors by decoupling the serial
correlation of residuals. Furthermore, a stochastic simulation method is also
developed based on our Bayesian calibration approach, which can obtain
unbiased posterior motion states and generate anthropomorphic driving
behaviors. Simulation results indicate that the MA-IDM outperforms Bayesian
IDM in simulation accuracy and uncertainty quantification. With this Bayesian
approach, we can generate enormous but nonidentical driving behaviors by
sampling from the posteriors, which can help develop a realistic traffic
simulator.

###### Index Terms:

car-following, serial correlation, Bayesian inference, hierarchical model,
Gaussian processes

##  I Introduction

Microscopic car-following models are powerful tools to study and simulate
human driving behaviors in traffic flows at the trajectory level. They reveal
the mechanisms of complex interactions between a subject vehicle and its
leading vehicle. These interactions are the essential factors that affect the
dynamics of traffic flow and create diverse macroscopic traffic phenomena.
However, no model is perfect [1], and car-following models are no exception.
Experiences indicate that we can never calibrate a zero-error model that
perfectly fits the data [2]. Nevertheless, some car-following models are still
valid, with which one can adopt various calibrations (or parameter
identification) methods to approximate reality to different extents. Given
that the performance and fidelity of microscopic car-following models heavily
depend on accurate calibration, improving the reliability of calibration
methods is a critical research question.

In recent studies, probabilistic calibration methods have become an emerging
and promising approach with a solid statistical foundation. There are
essentially two probabilistic approaches for model calibration. The first
approach is maximum likelihood estimation (MLE): we often solve MLE as an
optimization problem and finally obtain a point estimate for model parameters;
see e.g., [3]. The second approach is Bayesian inference, which allows us to
exploit the full posterior distribution for model parameters through Markov
chain Monte Carlo (MCMC); see e.g., [4, 5]. In general, the Bayesian inference
approach is advantageous in two aspects: (1) we can perform uncertainty
quantification based on the full posterior distribution. This is particularly
important for a simulation model as we are often interested in the
distribution and uncertainty of the simulation results. (2) Bayesian inference
offers a hierarchical modeling scheme that allows us to learn parameter
distributions at both population and individual levels. A unique advantage of
the Bayesian hierarchical model is that it allows us to learn parameters for
each individual driver/sample/trajectory, while the population distribution
will prevent overfitting by imposing certain dependencies on the parameters.
This is particularly important for calibrating car-following models since
samples are often short trajectories collected from a large number of drivers
with diverse vehicle configurations/dynamics. Considering the substantial
differences in personal driving styles and vehicle configurations/dynamics, it
is critical to learn the population distribution that can generalize the
parameter distribution from short trajectory data. In addition, the learned
population distribution can help create diverse driving behaviors to create
more realistic simulations with driver heterogeneity. Note that in the
following of this paper we refer to “driver behavior” as the characteristics
resulting from both the driver and the vehicle.

A fundamental question in performing probabilistic calibration is how to
define the probabilistic model and data generation process. Most existing car-
following models seek parsimonious structures by simply taking observations
from the most recent (i.e., only one) step as input to generate
acceleration/speed as output for the current step. However, given the physical
inertia and delay in reaction, we should expect the output acceleration/speed
to be correlated with and determined by multiple steps of past observations
[6]. For instance, Wang et al. [6] show that the best performed deep learning
model essentially takes states/observations from the most recent $\sim$10 sec
as input. As a result, a parsimonious car-following model that only takes one-
step past observation/input will inevitably show temporally correlated
residuals. Nevertheless, most existing probabilistic calibration methods are
developed based on a simple assumption—the residuals are independent and
identically distributed (i.i.d.), and as a result, ignoring the
autocorrelation in the residuals will lead to biased calibration [7, 8]. For
example, Fig. 1 shows the residual process in acceleration ($a$), speed ($v$),
and gap ($s$) when calibrating an IDM model with i.i.d. noise, and we can
clearly see that the residuals have strong autocorrelation. To address the
serial correlation problem, two classes of methods are developed in the
literature. One approach is to directly process the time series data to
eliminate serial correlations. For instance, Hoogendoorn and Hoogendoorn [3]
adopted a difference transformation (see [9]) to eliminate the serial
correlation. However, they used empirical correlation coefficients to perform
the transformation instead of jointly learning the IDM parameters and the
correlations. Another approach is to explicitly model the serial correlations
(e.g., by stochastic processes). For example, Treiber and Kesting [10]
introduced Wiener processes to model the temporal correlation of residuals.

Figure 1: The residual trends of the calibrated IDM have a strong serial
correlation. The ground truth data refer to the unpooled model in Fig. 8.

In this paper, we propose a novel probabilistic calibration approach that
leverages the advantages of Bayesian methods and the power of stochastic
processes to calibrate the car-following model unbiasedly. We apply this
approach to calibrate the IDM and develop a memory-augmented Bayesian IDM (MA-
IDM), which models the serially-correlated residuals as zero-mean Gaussian
processes (GP). Taking advantage of the Bayesian methods, we can jointly learn
the model parameters and the GP hyperparameters by MCMC. In this work, the MA-
IDM is presented in three forms, i.e., pooled model, hierarchical model, and
unpooled model. We conduct numerical experiments on the Highway Drone (HighD)
Dataset [11], and the results demonstrate that our hierarchical model not only
learns disparate driving styles at the population level but also depicts
various driving behaviors at the individual level. Additionally, a stochastic
simulation method is developed to obtain the posterior motion states in terms
of acceleration, speed, and location. In stochastic simulations, one can
generate enormous but nonidentical sets of model parameters from the learned
posterior distributions, which help characterize anthropomorphic driving
behaviors.

The overall contributions of this work are threefold:

  1. 1.

We develop a novel Bayesian calibration approach to learn unbiased parameters
and their full posterior distribution. We introduce GP to characterize the
autocorrelation in residuals. This approach is applied to calibrate IDM, and
results in a novel car-following model, i.e., the MA-IDM.

  2. 2.

We implement the MA-IDM with three hierarchies. Especially with a hierarchical
MA-IDM, one can obtain diverse driving styles at the population level and
disparate driving behaviors at the individual level by sampling from the
posterior distributions of the well-calibrated hierarchical model. Therefore,
we can generate enormous drivers with heterogeneous driving behaviors/styles
governed by the same population distribution. Therefore, our model can help
create simulations with driver/car heterogeneity.

  3. 3.

We introduce an unbiased stochastic simulator, which is inspired by the
corresponding generative process of our Bayesian calibration approach. As a
result, the simulator can produce more realistic results than those with
homogeneous parameters or random parameters.

The structure of the remaining contents is organized as follows. First,
Section II introduces related literature about calibrating car-following
models and modeling serial correlations. Section III emphasizes some
preliminaries and formulates the car-following problem based on IDM. Then in
Section IV, we propose a novel Bayesian calibration method and develop several
novel car-following models, i.e., the Bayesian-IDM (B-IDM) and the MA-IDM,
with different hierarchies. Next, we conduct extensive experiments and
simulations, then thoroughly analyze the results in Section V and Section VI.
Finally, Section VIII wraps up this work and discusses several potential
directions worth further exploration.

##  II Related Work

###  II-A Car-Following Model Calibration

The proper choice of key parameters in a car-following model can help depict
and replicate complex driving behaviors. As in most cases we never know the
values of these parameters, we often collect observable measurements from
field car-following data in diverse scenarios and fit the model to the
observed data by adjusting the parameters to optimize certain objective
functions. The overall process is known as calibration [12]. We refer readers
to [2] about calibration of car-following models.

Various calibration methods have been studied in the literature. The genetic
algorithm (GA) is one of the most typical and traditional tools as a heuristic
evolutionary algorithm [2, 13]. For instance, Punzo et al. [2] conducted
extensive experiments to study $29$ goodness-of-fit functions (GoF) with the
combinations from $3$ measures of performance (MoP) and provided systematic
guidelines on GA-based calibration. However, GA is sensitive to different
choices of GoF and MoP; it is not only computationally expensive but also
requires extensive efforts to search for the best combination of GoF and MoP.
Besides, the best combination of one dataset may not be suitable for another.
Maximum likelihood (MLE) is a powerful approach. Hoogendoorn et al. [14, 3]
first proposed to calibrate car-following models based on MLE. Treiber and
Kesting [7] thoroughly investigated the MLE approach on the calibration and
validation. MLE seeks to estimate a single “best” value of an unknown
quantity, which relies on optimization algorithms. The Bayesian method, as an
alternative approach, aims to estimate the distribution of an unknown quantity
(i.e., posterior), and often relies on sampling-based algorithms. Rahman et
al. [4] calibrated car-following models based on a Bayesian approach. The
experiments showed that the Bayesian approach provided much better results
than the deterministic optimization algorithm. With the Bayesian approach,
Abodo et al. [5] developed a hierarchical model to calibrate Intelligent
Driver Model (IDM), which has also achieved promising results.

###  II-B Serial Correlation Modeling

Regression plays a vital role in model calibration. In the literature, many
regression models are developed heavily based on an assumption—the residuals
are i.i.d. random variables. However, time series are successive observations
with serial correlations, which violates this i.i.d. assumption. For instance,
the i.i.d. assumption in the least squares regression methods was challenged
by Durbin and Watson [15], in which they tested the errors for independence.
Here, we argue that the calibration methods based on i.i.d. assumption would
result in biased parameter identification. Generally, as mentioned in [12],
the observations are usually regarded as the _true process_ with an i.i.d.
_observation error_ ; while in calibration methods, observations are modeled
as the combination of a _calibrated model_ and an independent _model
inadequacy function_ with an i.i.d. _observation error_. In calibrations
involving time series data, the _model inadequacy function_ is mainly used to
capture the serial correlations.

There are generally two ways to consider the _model inadequacy_ part by
investigating the serial correlations: (1) by directly processing the
nonstationary data and eliminating serial correlations, such that one can
safely ignore the _model inadequacy function_ and obtain stationary time
series; or (2) by explicitly modeling the serial correlations based on
specific _model inadequacy functions_. For instance, Hoogendoorn [3] performed
a differencing transformation to eliminate serial correlations, which then
showed no significant difference between autocorrelation coefficients and
zeros in the previously mentioned Durbin–Watson test [15]. However, the
information conveyed by the serial correlations is directly discarded, which
prevents us from modeling the generative processes of observations. Another
way is to explicitly develop the formation of serial correlations based on
further assumptions. For example, state-space models (e.g., Kalman filter
[16]) regard the residual trends as a filtering problem based on Markovian
assumptions. Dynamic regression models leverage linear regression and
autoregressive integrated moving average (ARIMA) model into a single powerful
regression model for forecasting time series data [17]. Stochastic processes
(e.g., Gaussian processes [18]) are also helpful tools for analyzing time
series, which assume time series as manifestations of stochastic processes.
These assumptions benefit the methods by deriving concise solutions but
constraining some capability. For the main scope of this paper, i.e.,
calibration and simulation of car-following models, we use Gaussian processes
to model the serially correlated residuals (i.e., for the _model inadequacy_
part). Gaussian process provides a solid statistical solution to learn
unbiased model parameters, and more importantly, it allows us to understand
the temporal effect in driving behavior through the lengthscale
hyperparameters $l$, which partially explains the memory effect of human
driving behaviors.

##  III Preliminaries and Problem Formulation

###  III-A IDM and Probabilistic IDM

The traditional IDM [19] is a continuous nonlinear function
$f:\mathbb{R}^{3}\mapsto\mathbb{R}$ which maps the gap, the speed (not
velocity111Note that the term ‘velocity’ represents the directional speed,
which is usually denoted as a vector. In the car-following case, motion
direction is longitudinal in the vehicle coordinate system at default.
Therefore, we use the term ‘speed’ in this paper.), and the speed difference
(approach rate) to acceleration $a_{\mathrm{IDM}}$ at a certain timestamp.
Here, we denote $s$ as the gap between the following vehicle and the leading
vehicle, $v$ as the speed of the following vehicle, and $\Delta v=-ds/dt$
(note the minus sign) as the speed difference. The physical meaning of these
notations are illustrated in Fig. 2, where $v_{l}$ denotes the speed of the
leading vehicle. With these notations, IDM is defined as

$\displaystyle\begin{split}&a_{\mathrm{IDM}}=f(s,v,\Delta
v)\triangleq\alpha\,\left(1-\left({\frac{v}{v_{0}}}\right)^{\delta}-\left({\frac{s^{\ast}(v,\Delta
v)}{s}}\right)^{2}\right),\end{split}$ (1)
$\displaystyle\begin{split}s^{\ast}(v,\Delta
v)=s_{0}+s_{1}\,\sqrt{\frac{v}{v_{0}}}+v\,T+{\frac{v\,\Delta
v}{2\,{\sqrt{\alpha\,\beta}}}},\end{split}$ (2)

where, $v_{0},s_{0},T,\alpha,\beta,$ and $\delta$ are model parameters with
the following meaning: the desired speed $v_{0}$ is the free-flow speed; the
jam spacing $s_{0}$ denotes a minimum gap distance from the leading vehicle;
the safe time headway $T$ represents the minimum time interval between the
following vehicle and the leading vehicle; the acceleration $\alpha$ and the
comfortable braking deceleration $\beta$ are the maximum vehicle acceleration
and the desired deceleration to keep safe, respectively; the exponential
coefficient $\delta$ is a constant, usually set as $4$ at default. In IDM, the
deceleration is controlled by the desired minimum gap $s^{\ast}$ in Eq. (2),
in which we set $s_{1}=0$ following [19] to obtain a model with interpretable
and easily measurable parameters.

Figure 2: Physical settings of a car-following scenario.

To make the notations concise and compact, we define vector
$\boldsymbol{\theta}=[v_{0},s_{0},T,\alpha,\beta]\in\mathbb{R}^{5}$
representing the IDM parameters that we are interested in. For a certain
vehicle $d$, we have the IDM decisions
$a_{\mathrm{IDM}\,d}^{(t)}=f(s_{d}^{(t)},v_{d}^{(t)},\Delta v_{d}^{(t)})$,
where the subscript $d$ represents the index for each driver and the
superscript $(t)$ indicates the timestamp. Compactly, we denote the inputs at
time $t$ as a vector $\boldsymbol{i}_{d}^{(t)}=[s_{d}^{(t)},v_{d}^{(t)},\Delta
v_{d}^{(t)}],\forall t\in\\{t_{0},\dots,t_{0}+(T-1)\Delta t\\}$. Here, we
adopt the scheme in [7, 20] with a step of $\Delta t$ to update vehicle speed
and position as following:

$\displaystyle v^{(t+\Delta t)}$ $\displaystyle=v^{(t)}+a^{(t)}\Delta t,$ (3a)
$\displaystyle x^{(t+\Delta t)}$ $\displaystyle=x^{(t)}+\frac{v^{(t+\Delta
t)}+v^{(t)}}{2}\Delta t.$ (3b)

To this end, one can define a concise function $\mathcal{F}_{\text{IDM}}$ that
updates the following vehicle’s speed at $t+\Delta t$ based on Eq. (3a),
written as

$\displaystyle\begin{split}v_{d}^{(t+\Delta
t)}|\boldsymbol{i}_{d}^{(t)},\boldsymbol{\theta}_{d}=\mathcal{F}_{\text{IDM}}(\boldsymbol{i}_{d}^{(t)};\boldsymbol{\theta}_{d})\triangleq
v_{d}^{(t)}+a_{\mathrm{IDM}\,d}^{(t)}\Delta t.\end{split}$ (4)

However, imperfect driving and irregularities in driving behaviors result in
erratic components of the driver’s action, as stated by Treiber and Kesting
[10] in Chapter 12.3.2. It can be modeled by adding to $a_{\text{IDM}}$ some
acceleration noise of standard deviation $\sigma_{\epsilon}$. In such a
setting, we consider $a_{\text{IDM}}$ as a rational behavior model, while the
random term as the imperfect driving behaviors that cannot be modeled by IDM.
By taking the random term with _i.i.d_ assumption into consideration, one can
develop a stochastic version of IDM, i.e., the probabilistic IDM [21, 22],
written as

$a_{d}^{(t)}|\boldsymbol{i}_{d}^{(t)},\boldsymbol{\theta}_{d}\,\stackon{\sim}{\tiny{\textit{i.i.d.}}}\,\mathcal{N}(a_{\mathrm{IDM}\,d}^{(t)},\sigma_{\epsilon}^{2}),$
(5)

where $a_{\mathrm{IDM}}$ and $\sigma^{2}_{\epsilon}$ represent the mean and
variance. Given the dynamic updating in Eq. (4), we have ${v}_{d}^{(t+\Delta
t)}=v_{d}^{(t)}+a_{d}^{(t)}\Delta t+\epsilon_{v},\
\epsilon_{v}\,\stackon{$\sim$}{\tiny{\emph{i.i.d.}}}\,\mathcal{N}(0,(\sigma_{\epsilon}\Delta
t)^{2})$. Thus,

${v}_{d}^{(t+\Delta
t)}|\boldsymbol{i}_{d}^{(t)},\boldsymbol{\theta}_{d}\,\stackon{\sim}{\tiny{\textit{i.i.d.}}}\,\mathcal{N}(\mathcal{F}_{\text{IDM}}(\boldsymbol{i}_{d}^{(t)};\boldsymbol{\theta}_{d}),(\sigma_{\epsilon}\Delta
t)^{2}).$ (6)

observation individual level (driving behaviors) population level (driving
styles) $\boldsymbol{\theta}$ (a) Pooled model. $y_{1}$$y_{2}$ $\cdots$
$y_{D}$$\boldsymbol{\theta}$ (b) Hierarchical model.
$\boldsymbol{\theta}_{1}$$y_{1}$$\boldsymbol{\theta}_{2}$$y_{2}$ $\cdots$
$\cdots$ $\boldsymbol{\theta}_{D}$$y_{D}$$\boldsymbol{\theta}_{1}$ (c)
Unpooled model. $y_{1}$$\boldsymbol{\theta}_{2}$$y_{2}$ $\cdots$ $\cdots$
$\boldsymbol{\theta}_{D}$$y_{D}$

Figure 3: Pooled model, hierarchical model, and unpooled model.

The calibration of car-following model mainly concerns estimating the model
parameter $\boldsymbol{\theta}$ to learn the best mapping from
$\boldsymbol{i}_{d}^{(t)}$ to a measure of performance [2], e.g.,
${v}_{d}^{(t+\Delta t)}$ here. Various approaches such as utility-based
optimization methods and maximum-likelihood techniques are usually adopted
[7].

###  III-B Pooled, Hierarchical, and Unpooled Model

Here we introduce three different hierarchies of the model, i.e., pooled,
hierarchical, and unpooled models, in which hierarchical models are often
regarded in terms of “partial pooling,” compromising between pooled and
unpooled models [23]. As shown in Fig. 3, we illustrate the pooling technique
and explicitly compare the differences among pooled, hierarchical, and
unpooled models. With a pooled model (referred to as complete pooling), all of
the drivers $\forall d\in\\{1,\dots,D\\}$ is assumed to have the same driving
behaviors (characterized by the parameter set $\boldsymbol{\theta}$). While
with an unpooled model (referred to as no pooling), each driver $d$ is assumed
to have a completely unrelated parameter set $\boldsymbol{\theta}_{d}$. But
with a hierarchical model (referred to as partial pooling), each driver is
assumed to have a different parameter set $\boldsymbol{\theta}_{d}$, while the
data for all the observed drivers $d\in\\{1,\dots,D\\}$ inform the parameter
estimates for each driver $d$. In general, hierarchical models consider a
population model of drivers as the prior for each specific driver. And from
this population model perspective, the unpooled model corresponds to infinite
population variance, whereas pooled model corresponds to zero population
variance.

${v}_{d}^{(t+1)}$$\boldsymbol{i}_{d}^{(t)}$$\boldsymbol{\theta}$$\sigma_{\epsilon}$
Time $t\in\\{1,\dots,T\\}$ Drivers $d\in\\{1,\dots,D\\}$
$\mathcal{F}_{\text{IDM}}$

(a) Pooled B-IDM.

${v}_{d}^{(t+1)}$$\boldsymbol{i}_{d}^{(t)}$$\boldsymbol{\theta}_{d}$$\boldsymbol{\Sigma}$$\boldsymbol{\theta}$$\sigma_{\epsilon}$$\boldsymbol{\sigma}_{0}$
Drivers $d\in\\{1,\dots,D\\}$ Time $t\in\\{1,\dots,T\\}$
$\mathcal{F}_{\text{IDM}}$

(b) Hierarchical B-IDM.

${v}_{d}^{(t+1)}$$\boldsymbol{i}_{d}^{(t)}$$\boldsymbol{\theta}_{d}$$\sigma_{\epsilon}$
Time $t\in\\{1,\dots,T\\}$ $\mathcal{F}_{\text{IDM}}$

(c) Unpooled B-IDM.

$\boldsymbol{v}_{d}$$\boldsymbol{i}_{d}$$\boldsymbol{\theta}$$\ell$$\sigma_{k}$$\sigma_{\epsilon}$
Drivers $d\in\\{1,\dots,D\\}$ $\mathcal{F}_{\text{IDM}}$ $GP$

(d) Pooled MA-IDM.

$\boldsymbol{v}_{d}$$\boldsymbol{i}_{d}$$\boldsymbol{\theta}_{d}$$\boldsymbol{\theta}$$\boldsymbol{\Sigma}$$\boldsymbol{\sigma}_{0}$$\ell$$\sigma_{k}$$\sigma_{\epsilon}$
Drivers $d\in\\{1,\dots,D\\}$ $GP$ $\mathcal{F}_{\text{IDM}}$

(e) Hierarchical MA-IDM.

$\boldsymbol{v}_{d}$$\boldsymbol{i}_{d}$$\boldsymbol{\theta}_{d}$$\ell$$\sigma_{k}$$\sigma_{\epsilon}$
$\mathcal{F}_{\text{IDM}}$ $GP$

(f) Unpooled MA-IDM.

Figure 4: Probabilistic graphical models.

##  IV Bayesian Calibration of Car Following Model

In this section, we develop a novel Bayesian calibration method and apply it
to calibrate IDM. Specifically, we present the formulations of four models,
i.e., the pooled and the hierarchical models of the B-IDM and the MA-IDM,
respectively. The corresponding probabilistic graphical models are shown in
Fig. 4.

###  IV-A Pooled and Hierarchical Bayesian IDM (B-IDM)

Firstly, we introduce the form of the pooled B-IDM:
$\forall(t,d)\in\\{(t,d)\\}_{t=t_{0},d=1}^{t_{0}+(T-1)\Delta t,D}$, we have

$\displaystyle\ln(\boldsymbol{\theta})$
$\displaystyle\sim\mathcal{N}(\boldsymbol{\mu}_{0},\boldsymbol{\Sigma}_{0})\in\mathbb{R}^{5},$
(7a) $\displaystyle\ln(\sigma_{\epsilon})$
$\displaystyle\sim\mathcal{N}(\mu_{\epsilon},\sigma_{1})\in\mathbb{R},$ (7b)
$\displaystyle{v}_{d}^{(t+\Delta
t)}|\boldsymbol{i}_{d}^{(t)},\boldsymbol{\theta}\,$
$\displaystyle\stackon{$\sim$}{\tiny{\textit{i.i.d.}}}\,\mathcal{N}(\mathcal{F}_{\text{IDM}}(\boldsymbol{i}_{d}^{(t)};\boldsymbol{\theta}),(\sigma_{\epsilon}\Delta
t)^{2})\in\mathbb{R}.$ (7c)

In the previous literature [5], some samples that drawn from the Markov chains
are unrealistically with negative values. To avoid this, here we assume
$\boldsymbol{\theta}$ and $\sigma_{\epsilon}$ are log-normal distributed. With
such a setting, one can ensure the corresponding parameters fall in the set of
positive real number $\mathbb{R}=\left\\{x\in\mathbb{R}\mid x>0\right\\}$.

However, the pooled model can only capture the general driving styles of a
“population driver.” To characterize and model the heterogeneity in driving
behavior, one can learn different parameters for each sample of the data,
which translates to introducing heterogeneity at the level of the driver
population [10]. Here we develop a hierarchical B-IDM that explicitly models
the general driving behaviors at the population level, then depicts
heterogeneous behaviors at the individual level, described as follows:
$\forall(t,d)\in\\{(t,d)\\}_{t=t_{0},d=1}^{t_{0}+(T-1)\Delta t,D}$, we have

$\displaystyle\boldsymbol{\sigma}_{0}$
$\displaystyle\sim\mathrm{Exp}(\lambda)\in\mathbb{R}^{5},$ (8a)
$\displaystyle\boldsymbol{\Sigma}$
$\displaystyle\sim\mathrm{LKJCholeskyCov}(\eta,\boldsymbol{\sigma}_{0})\in\mathbb{R}^{5\times
5},$ (8b) $\displaystyle\ln(\boldsymbol{\theta})$
$\displaystyle\sim\mathcal{N}(\boldsymbol{\mu}_{0},\boldsymbol{\Sigma})\in\mathbb{R}^{5},$
(8c) $\displaystyle\ln(\boldsymbol{\theta}_{d})$
$\displaystyle\sim\mathcal{N}(\ln(\boldsymbol{\theta}),\boldsymbol{\Sigma})\in\mathbb{R}^{5},$
(8d) $\displaystyle\ln(\sigma_{\epsilon})$
$\displaystyle\sim\mathcal{N}(\mu_{\epsilon},\sigma_{1})\in\mathbb{R},$ (8e)
$\displaystyle{v}_{d}^{(t+\Delta
t)}|\boldsymbol{i}_{d}^{(t)},\boldsymbol{\theta}_{d}\,$
$\displaystyle\stackon{$\sim$}{\tiny{\textit{i.i.d.}}}\,\mathcal{N}(\mathcal{F}_{\text{IDM}}(\boldsymbol{i}_{d}^{(t)};\boldsymbol{\theta}_{d}),(\sigma_{\epsilon}\Delta
t)^{2})\in\mathbb{R}.$ (8f)

where LKJ stands for the LKJ distribution developed by Lewandowski, Kurowicka,
and Joe [24]. And LKJCholeskyCov is a distribution over Cholesky decomposed
covariance matrices, such that the underlying correlation matrices follow an
LKJ distribution and the standard deviations follow an exponential
distribution, as specified by PyMC [25]. In practice, LKJCholeskyCov
distribution is a more robust prior for covariance matrix than the inverse-
Wishart distribution [26].

###  IV-B Pooled and Hierarchical Memory-Augmented IDM (MA-IDM)

For daily driving tasks, temporally-consistent inertia of actions/motions is
one of the vital driving profiles for human drivers. This phenomenon is termed
as _driving persistence_ in [10]. For instance, if drivers slam on their
brakes too hard at a certain timestamp, they will likely persist in a vast
deceleration in the next seconds. However, in most stochastic models, as well
as B-IDM introduced above, the persistence of the acceleration noise is
ignored, and its time dependence is wrongly modeled by white noise, such as in
Eq. (5). To make this point clear, we reformulate Eq. (5) as

$a_{d}^{(t)}=a_{\mathrm{IDM}\,d}^{(t)}+\epsilon_{t},\
\epsilon_{t}\,\stackon{\sim}{\tiny{\emph{i.i.d.}}}\,\mathcal{N}(0,\sigma_{\epsilon}^{2}).$
(9)

However, the developed models and the calibrated parameters are only valid if
the residuals conform to the i.i.d. assumption, and the temporal correlation
in residuals leads to the biased estimation of model parameters. In the
following, we will derive an unbiased Bayesian model in the pooled and
hierarchical form by modeling the time-dependent stochastic residual term with
Gaussian processes.

Given the nature of time-series autocorrelation, we model the $\epsilon_{t}$
in Eq. (9) by Gaussian processes, such that

$\displaystyle
a_{d}^{(t)}=a_{\mathrm{IDM}\,d}^{(t)}+a_{\mathrm{GP}\,d}^{(t)}+\epsilon$ (10)

with $\epsilon$ as time-independent i.i.d. noise, and
$a_{\mathrm{GP}\,d}^{(t)}\sim\mathcal{GP}(0,k(t,t^{\prime}))$, where $k$ is
the kernel function. Specifically, if we evaluate the continuous function at
several discrete timestamps $t=t_{0},\dots,t_{0}+(T-1)\Delta t$ with step
$\Delta t$ to obtain
$\boldsymbol{a}_{\mathrm{GP}\,d}=\\{a_{\mathrm{GP}\,d}^{(t)}\\}$, from the
definition of Gaussian processes, the prior distribution
$p(\boldsymbol{a}_{\mathrm{GP}\,d})$ is usually set as a zero-mean Gaussian
with a Gram matrix $\boldsymbol{K}$ being the covariance matrix, written as

$\displaystyle\boldsymbol{a}_{\mathrm{GP}\,d}\sim\mathcal{N}(\boldsymbol{0},\boldsymbol{K})\in\mathbb{R}^{T},$
(11)

where $\boldsymbol{K}\in\mathbb{R}^{T\times T}$ is the Gram matrix

$\displaystyle\boldsymbol{K}=\left[\begin{array}[]{ccc}k(\Delta t,\Delta
t)&\cdots&k(\Delta t,T\Delta t)\\\ \vdots&\ddots&\vdots\\\ k(T\Delta t,\Delta
t)&\cdots&k(T\Delta t,T\Delta t)\end{array}\right]\in\mathbb{R}^{T\times T}.$
(15)

Here we specify that $k$ is the squared exponential kernel
$k(t,t^{\prime})=\sigma_{k}^{2}\,\mathrm{exp}\left(-\frac{(t-t^{\prime})^{2}}{2\ell^{2}}\right)$
with $\ell$ being the lengthscale and $\sigma_{k}^{2}$ being the noise
variance. Notice that we reformulated $\boldsymbol{K}$ for concise notations
with the property $k(t_{0},t)\equiv k(t_{0}+\Delta t,t+\Delta t)\equiv
k(\Delta t,t-t_{0}+\Delta t)$.

In a compact form, by combining $\boldsymbol{a}_{\mathrm{IDM}\,d}$,
$\boldsymbol{a}_{\mathrm{GP}\,d}$ and $\epsilon$, we have

$\displaystyle\boldsymbol{a}_{d}|\boldsymbol{i}_{d},\boldsymbol{\theta}_{d}\sim\mathcal{N}(\boldsymbol{a}_{\mathrm{IDM}\,d},\boldsymbol{K}+\sigma_{\epsilon}^{2}\boldsymbol{I})\in\mathbb{R}^{T},$
(16)

where $\boldsymbol{I}\in\mathbb{R}^{T\times T}$.

Upon the above discussion, we absorb the $\boldsymbol{a}_{d}$ from Eq. (16)
into the B-IDM, one can obtain the pooled and the hierarchical MA-IDM which
will be discussed in the following.

For the pooled MA-IDM, $\forall d\in\\{1,\dots,D\\}$, we have

$\displaystyle\ln(\boldsymbol{\theta})$
$\displaystyle\sim\mathcal{N}(\boldsymbol{\mu}_{0},\boldsymbol{\Sigma}_{0})\in\mathbb{R}^{5},$
(17a) $\displaystyle\ln(\sigma_{\epsilon})$
$\displaystyle\sim\mathcal{N}(\mu_{\epsilon},\sigma_{1})\in\mathbb{R},$ (17b)
$\displaystyle\ln(\sigma_{k})$
$\displaystyle\sim\mathcal{N}(\mu_{k},\sigma_{2})\in\mathbb{R},$ (17c)
$\displaystyle\ln(\ell)$
$\displaystyle\sim\mathcal{N}(\mu_{\ell},\sigma_{\ell})\in\mathbb{R},$ (17d)
$\displaystyle\boldsymbol{v}_{d}|\boldsymbol{i}_{d},\boldsymbol{\theta}$
$\displaystyle\sim\mathcal{N}\left(\mathcal{F}_{\text{IDM}}(\boldsymbol{i}_{d};\boldsymbol{\theta}),(\boldsymbol{K}+\sigma_{\epsilon}^{2}\boldsymbol{I})\Delta
t^{2}\right)\in\mathbb{R}^{T}.$ (17e)

Note that different from Eq. (7c), we have a multivariate normal distribution
in Eq. (17e).

Similarly, we can develop the hierarchical MA-IDM as

$\displaystyle\boldsymbol{\sigma}_{0}$
$\displaystyle\sim\mathrm{Exp}(\lambda)\in\mathbb{R}^{5},$ (18a)
$\displaystyle\boldsymbol{\Sigma}$
$\displaystyle\sim\mathrm{LKJCholeskyCov}(\eta,\boldsymbol{\sigma}_{0})\in\mathbb{R}^{5\times
5},$ (18b) $\displaystyle\ln(\boldsymbol{\theta})$
$\displaystyle\sim\mathcal{N}(\boldsymbol{\mu}_{0},\boldsymbol{\Sigma})\in\mathbb{R}^{5},$
(18c) $\displaystyle\ln(\boldsymbol{\theta}_{d})$
$\displaystyle\sim\mathcal{N}(\ln(\boldsymbol{\theta}),\boldsymbol{\Sigma})\in\mathbb{R}^{5},$
(18d) $\displaystyle\ln(\sigma_{\epsilon})$
$\displaystyle\sim\mathcal{N}(\mu_{\epsilon},\sigma_{1})\in\mathbb{R},$ (18e)
$\displaystyle\ln(\sigma_{k})$
$\displaystyle\sim\mathcal{N}(\mu_{k},\sigma_{2})\in\mathbb{R},$ (18f)
$\displaystyle\ln(\ell)$
$\displaystyle\sim\mathcal{N}(\mu_{\ell},\sigma_{\ell})\in\mathbb{R},$ (18g)
$\displaystyle\boldsymbol{v}_{d}|\boldsymbol{i}_{d},\boldsymbol{\theta}_{d}\,$
$\displaystyle\stackon{$\sim$}{\tiny{\textit{i.i.d.}}}\,\mathcal{N}\left(\mathcal{F}_{\text{IDM}}(\boldsymbol{i}_{d};\boldsymbol{\theta}_{d}),(\boldsymbol{K}+\sigma_{\epsilon}^{2}\boldsymbol{I})\Delta
t^{2}\right)\in\mathbb{R}^{T}.$ (18h)

##  V Experiments

In this section, we access the performance of our proposed car-following
models by calibrating on naturalistic human-driving car-following data.

###  V-A Experimental Settings

####  V-A1 Dataset

The behavior of car-following models can be affected and contaminated by the
noise in empirical data [27]. Therefore, selecting an appropriate dataset can
keep the model performances away from data quality impact at the most while
relieving us from the redundant work of filtering the data noise. In this
paper, we use the HighD dataset [11]—a high-resolution trajectory data
collected using drones. Compared with the widely-used NGSIM dataset [28],
HighD benefits from advanced computer vision techniques and more reliable
data-capture methods.

The HighD dataset has $60$ video recordings, logged with the sampling
frequency of $25$ Hz on several German highway sections with a length of $420$
m. In each recording, the trajectories, velocities, and accelerations of two
types of vehicles (car and truck) are measured and estimated. We follow the
same data processing procedures as in [29] to transform the data into a new
coordinate system.

(a) Trajectory visualization. (b) Empirical histogram.

Figure 5: Visualizations of the $16$ selected leader-follower pairs.

####  V-A2 Car-Following Data Extraction

According to [30], long trajectories are preferred for robust estimation of
IDM parameters. Thus, we first filter out the data with a car-following
duration lasting less than a certain threshold $t_{0}$ to extract leader-
follower pairs. To obtain informative stop-and-go data with frequent
interactions, we set $t_{0}=50\,\mathrm{s}$. Then we randomly select $8$
leader-follower pairs from these interactive data for each type of vehicle,
respectively. Fig. 5 shows examples of the trajectories used for model
training. It should be noted that here we only consider the heterogeneity of
the follower; if interested, readers can also investigate different
combinations of more types of vehicles for both the leader and the follower,
e.g., as discussed in [31].

####  V-A3 MCMC Settings

Here we assume that the internal states of each driver are consistent, such
that human drivers do not change their driving styles in the region of
interest (ROI). In this work, we utilize PyMC [25] to estimate the posterior
distributions. Specifically, we adopt one kind of Hamiltonian Monte Carlo
(HMC) method, i.e., the No-U-Turn Sampler (NUTS) [32], and we set the burn-in
steps as $20000$ to make sure the detailed balance condition holds for the
Markov chains. The IDM priors are set by mainly following the recommended
values $\boldsymbol{\theta}_{\text{rec}}=[33.3,2.0,1.6,0.73,1.67]$ in [30],
e.g., the population prior $\boldsymbol{\mu}_{0}$ are set as
$\ln(\boldsymbol{\theta}_{\text{rec}})$. As for other hyperparameters, we set
$\lambda=100$ for the exponential distribution; $\eta=2$ for the
LKJCholeskyCov distribution; and the hyperparameters for the log-normal
distributions of $\sigma_{\epsilon}$, $\sigma_{k}$, and $\ell$ are flexible,
one can adjust them to ensure that they are positive values falling in
reasonable ranges.

####  V-A4 Reproducibility

All codes for reproducing the experiment results reported in this paper will
be available at https://github.com/Chengyuan-Zhang/IDM_Bayesian_Calibration.

(a) Hierarchical B-IDM posteriors of car $\\#355$. (b) Hierarchical MA-IDM
posteriors of car $\\#355$.

Figure 6: The posterior distributions and the correlation matrices of the
parameters of a car driver $\\#355$ in the hierarchical B-IDM and MA-IDM.

(a) Hierarchical B-IDM posteriors of truck $\\#211$. (b) Hierarchical MA-IDM
posteriors of truck $\\#211$.

Figure 7: The posterior distributions and the correlation matrices of the parameters of two representative driver $\\#211$ in the hierarchical B-IDM and MA-IDM. TABLE I: Posterior Expectations of Model Parameters. | $\boldsymbol{\theta}=[v_{0},s_{0},T,\alpha,\beta]$ | $\sigma_{\epsilon}$ | $\sigma_{k}$ | $\ell$
---|---|---|---|---
Pooled B-IDM (car) | $[14.552,2.906,0.469,0.622,1.267]$ | $0.517$ | $/$ | $/$
Pooled B-IDM (truck) | $[14.844,5.305,0.229,0.568,1.205]$ | $0.487$ | $/$ | $/$
Hierarchical B-IDM (population) | $[27.142,3.647,0.810,0.424,2.065]$ | $0.334$ | $/$ | $/$
Hierarchical B-IDM (car $\\#355$) | $[29.340,3.173,0.709,0.917,1.319]$ | $0.334$ | $/$ | $/$
Hierarchical B-IDM (truck $\\#211$) | $[31.013,2.800,1.680,0.287,1.246]$ | $0.334$ | $/$ | $/$
Unpooled B-IDM (car $\\#355$) | $[22.633,3.190,0.697,0.935,1.211]$ | $0.410$ | $/$ | $/$
Unpooled B-IDM (truck $\\#211$) | $[32.580,1.467,2.050,0.325,1.555]$ | $0.260$ | $/$ | $/$
Pooled MA-IDM (car) | $[13.629,2.579,0.422,0.531,1.754]$ | $0.105$ | $0.250$ |  $32.778$ ($1.311\,\mathrm{s}$)
Pooled MA-IDM (truck) | $[14.134,5.101,0.275,0.539,1.758]$ | $0.102$ | $0.242$ |  $37.226$ ($1.489\,\mathrm{s}$)
Hierarchical MA-IDM (population) | $[20.617,3.672,0.710,0.693,2.657]$ | $0.103$ | $0.212$ |  $33.141$ ($1.326\,\mathrm{s}$)
Hierarchical MA-IDM (car $\\#355$) | $[24.514,3.130,0.732,0.982,1.709]$ | $0.103$ | $0.212$ |  $33.141$ ($1.326\,\mathrm{s}$)
Hierarchical MA-IDM (truck $\\#211$) | $[29.938,2.013,1.907,0.397,1.788]$ | $0.103$ | $0.212$ |  $33.141$ ($1.326\,\mathrm{s}$)
Unpooled MA-IDM (car $\\#355$) | $[18.334,3.182,0.699,0.969,1.400]$ | $0.103$ | $0.164$ |  $25.542$ ($1.022\,\mathrm{s}$)
Unpooled MA-IDM (truck $\\#211$) | $[27.059,1.294,2.086,0.408,1.942]$ | $0.102$ | $0.137$ |  $32.938$ ($1.318\,\mathrm{s}$)
* Recommendation values [19]: $\boldsymbol{\theta}_{\text{rec}}=[33.3,2.0,1.6,0.73,1.67]$.

###  V-B Calibration Results and Analysis

We compare the calibration results of the pooled, hierarchical, and unpooled
models for the B-IDM and the MA-IDM in the following. The most straightforward
results are the calibrated IDM parameters in the form of posterior
distributions, for which we utilized _corner.py_ [33] for visualizations. Here
we use the hierarchical models as examples to illustrate our results. We
visualized the posterior distributions of $\boldsymbol{\theta}_{d}$ (i.e., at
the individual level) for two drivers in Fig. 6 and Fig. 7, respectively. The
posterior distributions of the B-IDM seem pretty similar to that of the MA-
IDM. To make it clearer to compare, we summarize the posterior expectations of
each parameter in Table I.

####  V-B1 Correlation matrices

From the upper right part of each subfigure in Fig. 6 and Fig. 7, the
correlation matrices indicate several strong correlations among IDM
parameters. The strong positive correlations exist in the pairs of
$(T,\alpha)$, $(T,\beta)$, and $(\alpha,\beta)$; while the strong negative
correlations exist in the pairs of $(v_{0},s_{0})$, $(v_{0},\alpha)$,
$(s_{0},T)$, $(s_{0},\alpha)$, and $(s_{0},\beta)$. These correlations
basically consistent with the results in the literature [7], in which the
calibrations were based on the likelihood but not the posterior. Also, the
correlation matrices of $\boldsymbol{\theta}_{d}$ at the individual level are
also consistent by using the hierarchical models and the unpooled models.

####  V-B2 Analysis of the calibrated IDM parameters

We verify that all of the calibrated parameters fall in the recommendation
ranges according to [34]. And the calibration results in Table I demonstrate
the reliability of our methods to calibrate plausible IDM parameters. In the
following, we make some comparisons among several pairs of calibration
results.

Firstly, one obvious conclusion is that cars tend to have larger maximum
vehicle acceleration $\alpha$ than trucks. Recall that the $\alpha$ in Eq. (1)
is an amplification coefficient that acts on the entire function. This finding
aligns with the truth that trucks are heavier, and truck drivers often use
eco-friendly driving strategies to avoid aggressive accelerations.

Secondly, the results of the pooled model indicate the general driving
behavior of a ‘population driver,’ i.e., the model learns a set of IDM
parameters with a moderate driving policy, which expects to explain all
behaviors in the data. However, as we noticed that human driving behaviors
vary from one driver to another because of disparate driving situations, the
capability of the pooled model is limited in explaining heterogeneous data.
Consequently, the pooled model considers the diversity as errors to compensate
for the uncertain behaviors and assign a significant variance
$\sigma_{\epsilon}^{2}$ to the residual term. The unpooled model, on the other
hand, is an example with smaller residual variance $\sigma_{\epsilon}^{2}$,
since it creates and fits a separable model for each specific driver at the
risk of overfitting.

Thirdly, we obtain mostly similar identified parameter values between
hierarchical and unpooled models at the individual driver’s level. The results
demonstrate the consistency of our proposed models with different hierarchies.
Also, it strongly supports the claim that the hierarchical models can
simultaneously identify style-specific IDM parameters for diverse drivers.

####  V-B3 Analysis of GP parameters

Recall that the RBF kernel function $k$ has two essential hyperparameters, the
lengthscale $\ell$ and the noise variance $\sigma_{k}^{2}$. One can refer to
Table I for the identified parameters of the MA-IDM.

First, we focus on analyzing the lengthscale $\ell$. Theoretically, for
$f\sim\mathcal{GP}(0,k)$, the correlation between $f(x)$ and $f(y)$ is
$k(x,y)$, e.g., when $\|x-y\|=\ell$, the correlation is $\exp(-1/2)\approx
0.61$. In the setting of car-following scenarios, $\ell$ can be interpreted as
an action delay in the stimulus-response mechanism. From the results, we
conclude that in our case, the car-following motion delay has: a decent amount
of correlation within $\ell\approx 1.33\,\mathrm{s}$; and is only slightly
dependent within $2\ell\approx 2.66\,\mathrm{s}$; but essentially independent
out of $3\ell\approx 4.00\,\mathrm{s}$. These results are also consistent with
the literature; a driving reaction time of $1.5$ seconds was concluded in
[35]. Additionally, all pooled, hierarchical, and unpooled models show that
truck drivers have longer reaction delays $\ell$ than car drivers. We also
noticed that the unpooled models have smaller $\ell$ than the pooled models,
which is reasonable—the pooled models are less accurate and have to compensate
for the uncertainty by increasing $\ell$ to obtain more likelihood of
explaining the data.

Besides, experiments demonstrate that the variance $\sigma_{\epsilon}^{2}$ of
the time-independent i.i.d. noise $\epsilon$ could be decreased at a
significance level by introducing GP. One should note that the
$\sigma_{\epsilon}^{2}$ of the B-IDM varies as we use disparate models or
calibrate on different data. However, we would expect the random i.i.d. noise
to be at almost the same level for the same dataset. On the other side, for
the MA-IDM, the $\sigma_{\epsilon}^{2}$ stabilized approximately to
$0.1\,\mathrm{m/s^{2}}$, which demonstrates that the errors captured by
$\epsilon$ are i.i.d. In addition, the overall residual noise level (including
the parts explained by GP) is significantly decreased, which reflects that the
calibrated IDM can better fit the data. Thus, we can identify the
effectiveness of the temporally correlated residual modeling in the MA-IDM.

###  V-C Data Completeness and Parameter Orthogonality

Several questions are worth considering before calibrating the model: What
scenarios may the data contain to train a well-calibrated model? There are
several parameters in the IDM formulation; are these parameters equally
important? Do we need to calibrate all the parameters [30]? The previous
literature has reported that parts of the parameters could not be determined
directly from the data due to insufficient information [3]. Researchers also
developed a reduced IDM that fixes some unimportant parameters [30]. Also, the
identifiability of car-following systems is discussed in [36]. We analyze
these problems from data completeness and parameter orthogonality perspectives
in this part.

Parameter orthogonality of IDM was first discussed in [7], in which the
parameter orthogonality is illustrated as “different identifiable driving
situations are represented by disparate parameters, ideally one for each
situation.” We notice that in Eq. (1) the two parameters $v_{0}$ and $\delta$
are correlated with each other and orthogonal to the remaining parameters in
$\boldsymbol{\theta}$. Specifically, they are shown in the same term
$(\frac{v}{v_{0}})^{\delta}$ but nowhere else. Given a certain $v$ and a
specific quantity of this term, one can find infinite groups of the
combination of $v_{0}$ and $\delta$. Experiments show that if we jointly
calibrate $v_{0}$ and $\delta$, they would both converge to their priors.
Therefore, for the sake of avoiding fractional exponential order and
numerically unstable issues, we followed the conventional settings and fix
$\delta=4$ to obtain a parameter orthogonal model.

Then we verify the completeness of the car-following data used in this work.
The IDM parameters need to be learned from corresponding data collected in
specific scenarios. For instance, $v_{0}$ is learned from the free-flow data;
$s_{0}$ and $T$ are learned from the car-following data; $\alpha$ is learned
from the freely accelerating data, and $\beta$ is learned from the approaching
(with braking) data. However, since we only selected the interactive stop-and-
go data, we don’t expect the calibrated model to identify $v_{0}$ well. The
experiments are consistent with our expectation—from Table I, we can see the
identified $v_{0}$ varies in a large range for disparate models; and from Fig.
6 and Fig. 7, the posterior distributions of $v_{0}$ are noisy and tend to
larger $v_{0}$. However, experiments indicate that this heavily depends on
whether the prior is strong enough—if we set a noninformative prior, the
posterior would be unrealistically high; and if we set a very strong prior, it
would be almost the same as our prior. Another interesting fact is that if we
set a weak prior with a very large $v_{0}$, the model would still converge but
with an unrealistic high posterior value of $v_{0}$; the other parameters,
however, are slightly changed.

##  VI Simulation

Based on the Bayesian calibration approach, in this section we propose a
stochastic simulation method corresponding to the generative processes of the
proposed Bayesian models, which is different from the traditional
deterministic car-following simulations. Then, the quantitative analysis is
demonstrated for comparison.

TABLE II: Evaluations of the deterministic and stochastic simulation results with different proposed models. Data $\backslash$ Model | Pooled model | Hierarchical model | Unpooled model
---|---|---|---
B-IDM | MA-IDM | B-IDM | MA-IDM | B-IDM | MA-IDM
Car $\\#355$ | $e(a)\times 10$ | $3.93\pm 0.02$ | $4.13\pm 0.02$ | $3.37\pm 0.05$ | $\mathbf{3.30}\pm 0.06$ | $3.37\pm 0.03$ | $3.33\pm 0.04$
$e(v)\times 10$ | $8.57\pm 0.17$ | $10.98\pm 0.26$ | $4.82\pm 0.12$ | $\mathbf{4.69}\pm 0.13$ | $4.77\pm 0.06$ | $4.74\pm 0.10$
$e(s)$ | $3.02\pm 0.09$ | $4.26\pm 0.20$ | $1.39\pm 0.07$ | $\mathbf{1.30}\pm 0.06$ | $1.37\pm 0.04$ | $1.32\pm 0.05$
${\mathrm{CRPS}}(a)$ | $0.23/0.54$ | $0.25/0.55$ | $\mathbf{0.19}/\mathbf{0.51}$ | $\mathbf{0.19}/0.52$ | $\mathbf{0.19}/\mathbf{0.51}$ | $0.20/\mathbf{0.51}$
${\mathrm{CRPS}}(v)$ | $0.59/0.58$ | $0.60/0.89$ | $0.32/0.36$ | $\mathbf{0.26}/0.27$ | $0.33/0.32$ | $0.27/\mathbf{0.26}$
${\mathrm{CRPS}}(s)$ | $2.21/4.74$ | $2.59/5.85$ | $0.89/3.32$ | $\mathbf{0.72}/2.79$ | $0.90/3.18$ | $0.76/\mathbf{2.73}$
Truck $\\#211$ | $e(a)\times 10$ | $2.69\pm 0.02$ | $2.54\pm 0.03$ | $2.14\pm 0.17$ | $\mathbf{1.93}\pm 0.09$ | $2.00\pm 0.06$ | $\mathbf{1.93}\pm 0.05$
$e(v)\times 10$ | $9.06\pm 0.10$ | $8.76\pm 0.13$ | $6.28\pm 1.53$ | $4.58\pm 0.57$ | $4.87\pm 0.39$ | $\mathbf{4.31}\pm 0.26$
$e(s)$ | $5.71\pm 0.06$ | $5.49\pm 0.09$ | $4.96\pm 1.95$ | $2.58\pm 0.84$ | $3.09\pm 0.50$ | $\mathbf{2.12}\pm 0.42$
${\mathrm{CRPS}}(a)$ | $0.15/0.07$ | $0.13/0.06$ | $0.12/0.07$ | $\mathbf{0.11}/0.06$ | $0.12/0.07$ | $\mathbf{0.11}/\mathbf{0.05}$
${\mathrm{CRPS}}(v)$ | $0.75/1.47$ | $0.56/0.92$ | $0.37/0.22$ | ${0.25/\mathbf{0.17}}$ | $0.33/0.08$ | $\mathbf{0.24/0.17}$
${\mathrm{CRPS}}(s)$ | $4.53/2.56$ | $3.64/1.29$ | $2.64/7.72$ | $1.31/3.91$ | $1.83/6.07$ | $\mathbf{1.14/3.41}$

(a) Predicted motion states with hierarchical models. (b) Predicted motion
states with unpooled models.

(c) Posterior motion states with hierarchical models. (d) Posterior motion
states with unpooled models.

Figure 8: The deterministic (top) and stochastic (bottom) simulation results
of a representative truck driver $\\#211$. The black lines indicate the ground
truth of human driving data; The yellow and fuchsia dotted lines are the
predicted motion states with the expectations of parameter posteriors; The red
and blue lines are the predicted motion states with the parameter samples
drawn from posteriors.

###  VI-A Deterministic and Stochastic Simulations

Deterministic simulation is a commonly used approach in trajectory
calibrations [7]. Typically, given the leading vehicle’s trajectory and the
initial motion state of the following vehicle, we can use the calibrated IDM
model to conduct simulations and obtain predicted motion states. Then, one can
compare the predicted motion states with the observations to further help with
either calibration or evaluation. In the previous literature, most motion
states updating methods are based on Eq. (3), which is deterministic. Such
deterministic simulations limit the simulator’s ability of reproducing
anthropomorphic driving behaviors. Furthermore, as we elaborated in Section
IV-B, only biased models would be obtained without considering serial
correlation and the deterministic simulations are no exception.

Given that the intrinsic nature of our Bayesian model is probabilistic
distribution, taking benefit from it, we propose a stochastic simulation
method. The key idea is to sample an action $a_{d}^{(t)}$ from Eq. (9) for the
B-IDM or Eq. (16) for the MA-IDM instead of using a deterministic
$a_{\mathrm{IDM}\,d}^{(t)}$. The stochasticity of such simulation is
determined by the noise level $\sigma_{\epsilon}$, calibrated from the data as
discussed in Section V-B3. According to the parameter identification results
in Table I, we control the noise level in simulations by fixing
$\sigma_{\epsilon}=0.3\,\mathrm{m/s^{2}}$ for the B-IDM and
$\sigma_{\epsilon}=0.1\,\mathrm{m/s^{2}}$ for the MA-IDM.

###  VI-B Parameters Generation and Error Analysis for Deterministic
Simulations

It was reported in the previous literature that car-following parameters
simply drawn from uncorrelated marginal distributions could yield unreliable
results in simulation and, consequently, inaccurate interpretation [37]. Here
we emphasize that one advantage of our Bayesian method compared with the
traditional optimization-based solution is obtaining the joint distributions
of parameters instead of point estimations. Therefore, rather than using the
expectation values from Table I, one can draw enormous samples from the joint
posterior distributions to generate anthropomorphic driving behaviors.
Specifically, from the viewpoint of a generative process, it is not difficult
to sample $\boldsymbol{\theta}_{d},d=1,\dots,D$ from the posteriors of
$\boldsymbol{\theta}$ (at the population level) to obtain disparate individual
driving styles; and sample style-specific IDM parameters
$\boldsymbol{\theta}_{d,i},i=1,\dots,N$ from the posteriors of
$\boldsymbol{\theta}_{d}$ (at the individual level) to obtain enormous but
nonidentical driving behaviors.

To evaluate the calibration and generation results, we use one representative
driver for each vehicle type, i.e., the car driver $\\#355$ and the truck
driver $\\#211$. We generate $N=5000$ sets of IDM parameters for each model.
Then, by providing the leading vehicle’s trajectories and starting the
followers with the same initial conditions, we simulate the trajectories of
the followers, whose behaviors are controlled by a calibrated IDM with a set
of sampled parameters [38]. The simulation step is therefore set as
$0.04\,\mathrm{sec}$ which corresponds to the $25$ Hz sampling frequency of
the leading vehicle’s data. To illustrate the simulated trajectories, we show
the results with parameter samples from $\boldsymbol{\theta}_{\\#211}$ in the
top row of Fig. 8. In the following, we compare the ground-truth trajectories
with the fully simulated trajectories via measuring the absolute root mean
square errors (RMSE) $e$ on their motion states (i.e., $a$, $v$, and $s$),
written as

$\displaystyle
e(x,\boldsymbol{\hat{f}},\boldsymbol{\theta}_{d,i})=\sqrt{\frac{1}{N}\sum_{i=1}^{N}\left(f(x_{i}|\boldsymbol{\theta}_{d,i})-\hat{f}_{i}\right)^{2}}.$
(19)

The simulation errors are reported in Table II in the form of (mean $\pm$
standard deviation), where the models with the lowest error are bolded. Recall
that the pooled models, acting as general behavioral models, were trained to
fit all the car-following data of several drivers. Apparently, given that
different drivers have distinct driving styles and have experienced dissimilar
driving situations, the data could reflect diverse operating characteristics.
Thus, a model that best fits a particular driver does not necessarily do so
for a different driver [34]. Therefore, it is reasonable that the pooled
models have limited performances. Besides, the simulation results indicate
that the errors of the hierarchical models are at the same level as the
unpooled models, which agrees with the findings in Section V-B2 that the
hierarchical models can identify disparate IDM parameters of diverse drivers.
Furthermore, although we can only obtain biased results in deterministic
simulations, the results of the hierarchical models and the unpooled models
still demonstrate that the MA-IDM has lower simulation errors than the B-IDM.
As shown in Fig. 8, the predicted motion states with the calibrated MA-IDM
tend to be closer to ground truth when compared with those driven by the
calibrated B-IDM. These evidences support that considering serial correlations
during parameter identification can reduce bias. However, it is worth noting
that even though the identified parameters are unbiased since we considered
the serial correlation during calibration, the predicted motion states are
still biased in deterministic simulations.

###  VI-C Unbiased Estimation and Bayesian Evaluations for Stochastic
Simulations

Benefiting from the Bayesian method, we can obtain posterior motion states in
stochastic simulations rather than the predicted motion states in
deterministic simulations. Here we propose a memory-augmented method by
considering the serial correlation in stochastic simulations to obtain
unbiased posterior motion states. Specifically, recall that we assumed
$a_{d}^{(t)}=a_{\mathrm{IDM}\,d}^{(t)}+a_{\mathrm{GP}\,d}^{(t)}+\epsilon$ in
Section IV-B, we simulate the motion of the following vehicle at a certain
timestamp by simultaneously taking several independent steps as follows. To
initialize, we draw several random functions $a_{\mathrm{GP}\,d},d=1,\dots,D$
from GP with the RBF kernel parameterized by
$\sigma_{k}=0.2\,\mathrm{m/s^{2}}$ and $\ell=33\,\mathrm{s}$ at the beginning
of the simulations. Then for each driver $d$ at timestamp $t$: (1) we obtain
the first term $a_{\mathrm{IDM}\,d}^{(t)}$ from IDM with the typical
procedures in deterministic simulations; (2) we evaluate the $d$-th drawn
function $a_{\mathrm{GP}\,d}^{(t)}$ at time $t$ as the second term; (3) we
sample a random noise $\epsilon$ from $\mathcal{N}(0,\sigma_{\epsilon}^{2})$.
The posterior motion states are shown in the bottom row of Fig. 8.

To evaluate the performances of stochastic simulations and quantify the
uncertainty of posterior motion states, we adopted continuous ranked
probability score (CRPS) [39], which can be written as

$\mathrm{CRPS}(y_{t})=\int_{-\infty}^{+\infty}\left(F(y)-\mathbbm{1}\\{y>y_{t}\\}\right)^{2}dy,$
(20)

where $y_{t}$ is the observation at time $t$, $F$ is the forecast cumulative
distribution function, and $\mathbbm{1}$ is the indicator function. To build
the evaluation metric, we first evaluate CRPS for $a$, $v$, and $s$ at each
timestamp $t$, then averaging over time by
$\overline{\mathrm{CRPS}}(y)=\frac{1}{T}\sum_{t=1}^{T}\mathrm{CRPS}(y_{t})$.
By noticing the gaps between the biased motion state predictions and the
ground truth from the top row of Fig. 8, we also evaluate CRPS at
$t_{0}=37\,\mathrm{s}$ and $t_{0}=45\,\mathrm{s}$ for car and truck data,
respectively. The Bayesian evaluation results are shown in the form of
($\overline{\mathrm{CRPS}}(y)/\mathrm{CRPS}(y_{t}|t=t_{0})$) in Table II.

From Fig. 8, an interesting but not surprising finding is the red lines in the
bottom row (i.e., the posterior motion states of the MA-IDM) perfectly fix the
gaps of predicted motion states shown in the top row. And the posterior motion
states in unbiased simulations tightly contain the ground truth in the
envelope. In addition, in stochastic simulations, although the B-IDM (the blue
lines) loosely and easily contains the ground truth of the acceleration curve
in its envelope, the simulated results for $v$ and $s$ are still biased. From
Table II, we can conclude that CRPS is lowered in the memory-augmented
simulations. Therefore, by introducing the memory-augmented method to model
the serial correlation, the unbiased stochastic simulations with the unbiased
parameters identified by the MA-IDM outperform the biased results in the
B-IDM.

(a) Simulation with homogeneous IDM parameters. (b) Simulation with
heterogeneous IDM parameters.

Figure 9: Time-space diagram of multi-vehicle car-following simulations. (a)
the parameters are taken as the recommendation values [19]. (b) the parameters
are sampled from the posteriors in the hierarchical MA-IDM.

Figure 10: Fundamental diagrams of ring-road simulation.

###  VI-D Multi-Vehicle Simulations in A Ring Scenario

The single-lane ring road network has been extensively studied to investigate
the microscopic car-following behaviors in a closed traffic system [40, 41].
Several key elements determine the dynamics on the ring road: The ring radius
is set as $128$ m, which results in approximately $804$ m of the
circumference; The initial speed of each vehicle is set as $11.6$ m/s; There
are a total of $37$ vehicles simulated for $3000$ sec; The simulation step is
set as $0.5$ sec.

Here we conduct stochastic simulations with two different settings to compare
the microscopic behaviors, as shown in Fig. 9. The left of Fig. 9 is simulated
with homogeneous IDM parameters taken as the recommended values listed in
Table I, and the noise level is controlled to the same level discussed in
Section VI-A. The right part shows the results of the memory-augmented
simulations, in which the GP settings are the same as in Section VI-C. As can
be seen, we can observe a recurring pattern from the simulations with fixed
IDM parameters, although a stochastic term is introduced. On the contrary, we
can obtain different random car-following behaviors in the heterogeneous
setting, which reflects the dynamic and diverse driving styles of different
drivers. Additionally, we demonstrated the speed-density and flow-density
fundamental diagrams of the simulation results, as shown in Fig. 10. These
diagrams indicate that the simulations with heterogeneous IDM parameters can
generate a wide variety of car-following behaviors, thus, are more realistic
and stochastic than those with homogeneous ones.

##  VII Future Scope and Discussions

###  VII-A Future Works

There are several aspects that can be further explored in future works. First,
our method is not limited to calibrating IDM; one can easily extend it to
other microscopic car-following models (even other regression problems, see
e.g., [42]) using a similar Bayesian scheme introduced in this work. Besides,
the MoP is not restricted to using speed, as mentioned by [2]. Although the
likelihood here is developed based on speed, the spacing or acceleration could
also be a good MoP. One can easily make minor adjustments based on our method
to calibrate car-following models using spacing or acceleration.

Additionally, human drivers have personalized styles, which are not life-long
determined but time-varying [43]; thus, we ask: “how to model the dynamic
time-varying IDM parameters for every single driver?” Following this concern,
we further emphasize that in reality the residual variance and lengthscale
should be scenario specific. For instance, one could freely take actions with
a significant variance in free flow, but most drivers take almost the same
actions (e.g., hard brake) in some risky scenarios. Thus, the GP
hyperparameters (e.g., noise variance and temporal lengthscale) should be
input-dependent, and the GP itself will become nonstationary. An interesting
research question is to quantify how the GP hyperparameters vary with the
input/scenario using non-stationary GP regression [44].

Furthermore, more realistic scenarios involving multiple factors can be
studied by extending our model into a more complex architecture. For instance,
human driving behaviors result from frequent transitions among various driving
styles [45, 43]. Thus, it is no longer reasonable to model the driving
behaviors using a single distribution; instead, a mixture model is suitable
[46]. In addition, as mentioned in [47], power dynamics of vehicle are not
uniform across the speed range; also, different vehicles (e.g., cars and
trucks in this work) have differences in power capabilities. These factors
lead to vehicle dynamics heterogeneity [48, 47], which can be modeled by
constructing deeper hierarchical models, e.g., add another layer in Fig.3 to
depict different vehicle dynamic characteristics.

###  VII-B Open Discussions

Here we post some open questions to be discussed with the community. Although
we have calibrated an unpooled model, we note that the calibration results
will not be so satisfactory since the segment is too short and contains
minimal information (e.g., with only congestion but without high-speed free
flow data, or vice versa.) So, we must reconsider when collecting data:

  * •

What kinds of data are preferable for calibrating an accurate car-following
model? Which maneuvers (e.g., free-flow, freely accelerating, approaching,
keeping distance, etc.) that are missing in the training data may
significantly impact the calibration performance? Is the model calibrated on
data with partial maneuvers meaningful? If so, how to evaluate the calibration
performance in the scenarios of the missing maneuvers?

These questions further lead us to the discussion about the data collection
methods. By leveraging different methods and sensors, a lot of microscopic
datasets are constructed. For instance, sensor data collected by stationary
loop detectors, extended floating-car data (xFCD) collected by GPS, trajectory
data collected by drone, etc. Limitations are obvious for each data collection
method: stationary loop detectors are mounted on the road, thus can only
obtain information at a specific location; xFCD care more about ego vehicle’s
motion states but only with limited information regarding other traffic
participants; drones collect data based computer vision approach, which can
only extract trajectory from a limited ROI222The pNEUMA dataset [49] is a
great attempt trying to extend the ROI by sharing information among a swarm of
drones.. In the near future of Internet of Things (IoT) for transportation, we
still wondering:

  * •

How to collect high-quality public-available datasets with rich information?
How to leverage different data collection methods with IoT for transportation?
Whose data should we use/trust when obtaining different values from different
sensors regarding the same variable?

Additionally, Gaussian processes still have some limitations. For instance, as
we mentioned in Section V-C—when data is limited, setting a strong prior would
constraint the posterior near to our prior; but setting a noninformative prior
would have unrealistic results. Also, the computational complexity has limited
its application to larger datasets.

  * •

How to select an appropriate prior in Bayesian methods? How to construct the
kernel function that suits the driving characteristics? How to compare and
incorporate different stochastic processes to modeling driving behaviors?

##  VIII Conclusion

Most literature developed their models upon the i.i.d. assumption of
residuals, which does not hold. This paper developed a novel car-following
model (i.e., the MA-IDM) by leveraging Bayesian methods and Gaussian
processes, which model the serial correlation of residuals. In such a way, it
makes the i.i.d. assumption valid and has the model calibration unbiased.
Extensive experiments compare the pooled, hierarchical, and unpooled models
for both of the B-IDM and the MA-IDM. Experiment results demonstrate the
effectiveness of the MA-IDM. With the well-calibrated MA-IDM, we also
developed an unbiased stochastic simulator, which generates enormous drivers
with nonidentical driving behaviors, and provides insight into creating a
realistic simulator.

## References

  * [1] I. D. Ofiteru and C. Picioreanu, “No model is perfect, but some are useful,” _Science_ , vol. 376, no. 6596, pp. 914–916, 2022.
  * [2] V. Punzo, Z. Zheng, and M. Montanino, “About calibration of car-following dynamics of automated and human-driven vehicles: Methodology, guidelines and codes,” _Transportation Research Part C: Emerging Technologies_ , vol. 128, p. 103165, 2021.
  * [3] S. Hoogendoorn and R. Hoogendoorn, “Calibration of microscopic traffic-flow models using multiple data sources,” _Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences_ , vol. 368, no. 1928, pp. 4497–4517, 2010.
  * [4] M. Rahman, M. Chowdhury, T. Khan, and P. Bhavsar, “Improving the efficacy of car-following models with a new stochastic parameter estimation and calibration method,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 16, no. 5, pp. 2687–2699, 2015.
  * [5] F. Abodo, A. Berthaume, S. Zitzow-Childs, and L. Bobadilla, “Strengthening the case for a bayesian approach to car-following model calibration and validation using probabilistic programming,” in _2019 IEEE Intelligent Transportation Systems Conference (ITSC)_. IEEE, 2019, pp. 4360–4367.
  * [6] X. Wang, R. Jiang, L. Li, Y. Lin, X. Zheng, and F.-Y. Wang, “Capturing car-following behaviors by deep learning,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 19, no. 3, pp. 910–920, 2017.
  * [7] M. Treiber and A. Kesting, “Microscopic calibration and validation of car-following models–a systematic approach,” _Procedia-Social and Behavioral Sciences_ , vol. 80, pp. 922–939, 2013.
  * [8] V. Punzo and M. Montanino, “A two-level probabilistic approach for validation of stochastic traffic simulations: impact of drivers’ heterogeneity models,” _Transportation research part C: emerging technologies_ , vol. 121, p. 102843, 2020.
  * [9] D. Cochrane and G. H. Orcutt, “Application of least squares regression to relationships containing auto-correlated error terms,” _Journal of the American statistical association_ , vol. 44, no. 245, pp. 32–61, 1949.
  * [10] M. Treiber and A. Kesting, “Traffic flow dynamics,” _Traffic Flow Dynamics: Data, Models and Simulation, Springer-Verlag Berlin Heidelberg_ , pp. 983–1000, 2013.
  * [11] R. Krajewski, J. Bock, L. Kloeker, and L. Eckstein, “The highd dataset: A drone dataset of naturalistic vehicle trajectories on german highways for validation of highly automated driving systems,” in _2018 21st International Conference on Intelligent Transportation Systems (ITSC)_ , 2018, pp. 2118–2125.
  * [12] M. C. Kennedy and A. O’Hagan, “Bayesian calibration of computer models,” _Journal of the Royal Statistical Society: Series B (Statistical Methodology)_ , vol. 63, no. 3, pp. 425–464, 2001.
  * [13] C. Chen, L. Li, J. Hu, and C. Geng, “Calibration of mitsim and idm car-following model based on ngsim trajectory datasets,” in _Proceedings of 2010 IEEE International Conference on Vehicular Electronics and Safety_. IEEE, 2010, pp. 48–53.
  * [14] S. P. Hoogendoorn and S. Ossen, “Parameter estimation and analysis of car-following models,” in _Transportation and Traffic Theory. Flow, Dynamics and Human Interaction. 16th International Symposium on Transportation and Traffic TheoryUniversity of Maryland, College Park_ , 2005.
  * [15] J. Durbin and G. S. Watson, “Testing for serial correlation in least squares regression: I,” _Biometrika_ , vol. 37, no. 3/4, pp. 409–428, 1950.
  * [16] G. Welch, G. Bishop, _et al._ , “An introduction to the kalman filter,” 1995\.
  * [17] R. J. Hyndman and G. Athanasopoulos, _Forecasting: principles and practice_. OTexts, 2018.
  * [18] C. E. Rasmussen, “Gaussian processes in machine learning,” in _Summer school on machine learning_. Springer, 2003, pp. 63–71.
  * [19] M. Treiber, A. Hennecke, and D. Helbing, “Congested traffic states in empirical observations and microscopic simulations,” _Physical review E_ , vol. 62, no. 2, p. 1805, 2000.
  * [20] M. Treiber, A. Kesting, and D. Helbing, “Delays, inaccuracies and anticipation in microscopic traffic models,” _Physica A: Statistical Mechanics and its Applications_ , vol. 360, no. 1, pp. 71–88, 2006.
  * [21] R. P. Bhattacharyya, R. Senanayake, K. Brown, and M. J. Kochenderfer, “Online parameter estimation for human driver behavior prediction,” in _2020 American Control Conference (ACC)_. IEEE, 2020, pp. 301–306.
  * [22] M. Treiber and A. Kesting, “The intelligent driver model with stochasticity-new insights into traffic flow oscillations,” _Transportation research procedia_ , vol. 23, pp. 174–187, 2017.
  * [23] A. Gelman, “Multilevel (hierarchical) modeling: what it can and cannot do,” _Technometrics_ , vol. 48, no. 3, pp. 432–435, 2006.
  * [24] D. Lewandowski, D. Kurowicka, and H. Joe, “Generating random correlation matrices based on vines and extended onion method,” _Journal of multivariate analysis_ , vol. 100, no. 9, pp. 1989–2001, 2009.
  * [25] J. Salvatier, T. V. Wiecki, and C. Fonnesbeck, “Probabilistic programming in python using pymc3,” _PeerJ Computer Science_ , vol. 2, p. e55, 2016.
  * [26] I. Alvarez, J. Niemi, and M. Simpson, “Bayesian inference for a covariance matrix,” _arXiv preprint arXiv:1408.4050_ , 2014.
  * [27] M. Montanino and V. Punzo, “Trajectory data reconstruction and simulation-based validation against macroscopic traffic patterns,” _Transportation Research Part B: Methodological_ , vol. 80, pp. 82–106, 2015\.
  * [28] V. Punzo, M. T. Borzacchiello, and B. Ciuffo, “On the assessment of vehicle trajectory data accuracy and application to the next generation simulation (ngsim) program data,” _Transportation Research Part C: Emerging Technologies_ , vol. 19, no. 6, pp. 1243–1262, 2011.
  * [29] C. Zhang, J. Zhu, W. Wang, and J. Xi, “Spatiotemporal learning of multivehicle interaction patterns in lane-change scenarios,” _IEEE Transactions on Intelligent Transportation Systems_ , 2021.
  * [30] V. Punzo, M. Montanino, and B. Ciuffo, “Do we really need to calibrate all the parameters? variance-based sensitivity analysis to simplify microscopic traffic flow models,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 16, no. 1, pp. 184–193, 2014.
  * [31] L. Liu, L. Zhu, and D. Yang, “Modeling and simulation of the car-truck heterogeneous traffic flow based on a nonlinear car-following model,” _Applied Mathematics and Computation_ , vol. 273, pp. 706–717, 2016.
  * [32] M. D. Hoffman, A. Gelman, _et al._ , “The no-u-turn sampler: adaptively setting path lengths in hamiltonian monte carlo.” _J. Mach. Learn. Res._ , vol. 15, no. 1, pp. 1593–1623, 2014.
  * [33] D. Foreman-Mackey, “corner.py: Scatterplot matrices in python,” _The Journal of Open Source Software_ , vol. 1, no. 2, p. 24, jun 2016. [Online]. Available: https://doi.org/10.21105/joss.00024
  * [34] A. Kesting and M. Treiber, “Calibrating car-following models by using trajectory data: Methodological study,” _Transportation Research Record_ , vol. 2088, no. 1, pp. 148–156, 2008.
  * [35] T. J. Triggs, W. G. Harris, _et al._ , “Reaction time of drivers to road stimuli,” 1982.
  * [36] Y. Wang, M. L. Delle Monache, and D. B. Work, “Identifiability of car-following dynamics,” _Physica D: Nonlinear Phenomena_ , vol. 430, p. 133090, 2022.
  * [37] J. Kim and H. S. Mahmassani, “Correlated parameters in driving behavior models: Car-following example and implications for traffic microsimulation,” _Transportation Research Record_ , vol. 2249, no. 1, pp. 62–77, 2011.
  * [38] B. Ciuffo, V. Punzo, and V. Torrieri, “Comparison of simulation-based and model-based calibrations of traffic-flow microsimulation models,” _Transportation Research Record_ , vol. 2088, no. 1, pp. 36–44, 2008.
  * [39] J. E. Matheson and R. L. Winkler, “Scoring rules for continuous probability distributions,” _Management science_ , vol. 22, no. 10, pp. 1087–1096, 1976\.
  * [40] M. Treiber and A. Kesting, “An open-source microscopic traffic simulator,” _IEEE Intelligent Transportation Systems Magazine_ , vol. 2, no. 3, pp. 6–13, 2010.
  * [41] Y. Sugiyama, M. Fukui, M. Kikuchi, K. Hasebe, A. Nakayama, K. Nishinari, S.-i. Tadaki, and S. Yukawa, “Traffic jams without bottlenecks—experimental evidence for the physical mechanism of the formation of a jam,” _New journal of physics_ , vol. 10, no. 3, p. 033001, 2008.
  * [42] Z. Cheng, X. Wang, X. Chen, M. Trepanier, and L. Sun, “Bayesian calibration of traffic flow fundamental diagrams using gaussian processes,” _arXiv preprint arXiv:2208.02799_ , 2022.
  * [43] J. Yarlagadda and D. S. Pawar, “Heterogeneity in the driver behavior: an exploratory study using real-time driving data,” _Journal of advanced transportation_ , vol. 2022, 2022.
  * [44] M. Heinonen, H. Mannerström, J. Rousu, S. Kaski, and H. Lähdesmäki, “Non-stationary gaussian process regression with hamiltonian monte carlo,” in _Artificial Intelligence and Statistics_. PMLR, 2016, pp. 732–740.
  * [45] W. Wang, J. Xi, and D. Zhao, “Driving style analysis using primitive driving patterns with bayesian nonparametric approaches,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 20, no. 8, pp. 2986–2998, 2018.
  * [46] ——, “Learning and inferring a driver’s braking action in car-following scenarios,” _IEEE Transactions on Vehicular Technology_ , vol. 67, no. 5, pp. 3887–3899, 2018.
  * [47] M. A. Makridis, A. Anesiadou, K. Mattas, G. Fontaras, and B. Ciuffo, “Characterising driver heterogeneity within stochastic traffic simulation,” _Transportmetrica B: Transport Dynamics_ , vol. 0, no. 0, pp. 1–19, 2022\.
  * [48] B. Ciuffo, M. Makridis, T. Toledo, and G. Fontaras, “Capability of current car-following models to reproduce vehicle free-flow acceleration dynamics,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 19, no. 11, pp. 3594–3603, 2018.
  * [49] E. Barmpounakis and N. Geroliminis, “On the new era of urban traffic monitoring with massive drone data: The pneuma large-scale field experiment,” _Transportation research part C: emerging technologies_ , vol. 111, pp. 50–71, 2020.

|  Chengyuan Zhang (Student Member, IEEE) is a Ph.D. student with the
Department of Civil Engineering, McGill University. He received his B.S.
degree in Automotive Engineering from Chongqing University, Chongqing, China,
in 2019. From 2019 to 2020, he worked as a visiting researcher with the
Department of Mechanical Engineering, UC Berkeley. His research interests are
Bayesian learning, macro/micro human driving behavior analysis, and multi-
agent interaction modeling in intelligent transportation systems.
---|---
|  Lijun Sun (Senior Member, IEEE) received the B.S. degree in civil
engineering from Tsinghua University, Beijing, China, in 2011, and the Ph.D.
degree in civil engineering (transportation) from the National the University
of Singapore in 2015. He is currently an Assistant Professor in the Department
of Civil Engineering, McGill University, Montreal, QC, Canada. His research
centers on intelligent transportation systems, machine learning,
spatiotemporal modeling, travel behavior, and agent-based simulation.
---|---