A phase-based smoothing method for accurate traffic speed estimation with floating car data

Transportation Research Part C 85 (2017) 644–663

Contents lists available at ScienceDirect

Transportation Research Part C journal homepage: www.elsevier.com/locate/trc

A phase-based smoothing method for accurate traffic speed estimation with floating car data☆

MARK

⁎

Felix Rempea,b, , Philipp Franecka, Ulrich Fastenratha, Klaus Bogenbergerb a b

Location-Based Services, Traffic Information and Routing, BMW Group, Petuelring 130, 80788 München, Germany Department of Traffic Engineering, Munich University of the Federal Armed Forces, Werner-Heisenberg-Weg 39, 85577 Neubiberg, Germany

AR TI CLE I NF O

AB S T R A CT

Keywords: Traffic state estimation Floating car data Three-phase traffic theory Generalized adaptive smoothing method Congestion patterns

In this paper, a novel freeway traffic speed estimation method based on probe data is presented. In contrast to other traffic speed estimators, it only requires velocity data from probes and does not depend on any additional data inputs such as density or flow information. In the first step the method determines the three traffic phases free flow, synchronized flow, and Wide Moving Jam (WMJ) described by Kerner et al. in space and time. Subsequently, reported data is processed with respect to the prevailing traffic phase in order to estimate traffic velocities. This two-step approach allows incorporating empirical features of phase fronts into the estimation procedure. For instance, downstream fronts of WMJs always propagate upstream with approximately constant velocity, and downstream fronts of synchronized flow phases usually stick to bottlenecks. The second step assures the validity of measured velocities is limited to the extent of its assigned phase. Effectively, velocity information in space-time can be estimated more distinctively and the result is therefore more accurate even if the input data density is low. The accuracy of the proposed Phase-Based Smoothing Method (PSM) is evaluated using real floating car data collected during two traffic congestions on the German freeway A99 and compared to the performance of the Generalized Adaptive Smoothing Method (GASM) as well as a naive algorithm. The quantitative and qualitative results show that the PSM reconstructs the congestion pattern more accurately than the other two. A subsequent analysis of the computational efficiency and sensitivity demonstrates its practical suitability.

1. Introduction For every type of traffic application as well as for traffic research, the best-possible estimate of current and historical traffic speed is fundamental. Despite the fact that there are usually only a limited number of sensors available that provide sparse measurements in space and time of the real traffic situation. The task of a traffic estimation method is to analyze the given measurements and provide an estimate of the desired traffic variable for every position on the road and every point in time with the main purpose of reconstructing and providing accurate estimates of prevailing real traffic conditions. For traffic estimation, the quick spread of Floating Car (FC) data in recent years provides great potentials; however, still faces new challenges. One great advantage is the possibility to obtain traffic data on every position on a road, whereas fixed sensors can only provide data for pre-determined locations. Besides the high spatial resolution, another advantage is the ever-decreasing cost of this technology, which in turn increases the expectation of developing availability of FC data in the upcoming years. The installation and

☆ ⁎

This article belongs to the Virtual Special Issue on “Future Traffic Management”. Corresponding author at: Location-Based Services, Traffic Information and Routing, BMW Group, Petuelring 130, 80788 München, Germany. E-mail address: [email protected] (F. Rempe).

http://dx.doi.org/10.1016/j.trc.2017.10.015 Received 24 October 2016; Received in revised form 11 October 2017; Accepted 15 October 2017 0968-090X/ © 2017 Published by Elsevier Ltd.

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

maintenance of fixed detectors, on the one hand is costly, especially when considering the vastness of networks that have to be covered. In spite of this, they still provide great benefits in the provision of macroscopic traffic speeds as well as density and flow values in short time intervals. In contrast, time intervals between two reported trajectories can vary drastically and FC data usually does not contain explicit density and flow information. Because of the sparse data and the under-determined traffic state coupled with high spatio-temporal traffic dynamics accurate traffic state estimation with FC data is a major challenge. For decades traffic state estimation is an ongoing topic of research. A comprehensive survey of current traffic estimation approach on highways is given in Seo et al. (2017). The main approaches that consider traffic dynamics can be classified into two categories. The first category comprises analytical flow models coupled with data assimilation techniques. First order models are usually based on the Lighthill-Whitham-Richards model (Lighthill and Whitham, 1955; Richards, 1956) and use Kalman filters (van Lint and Djukic, 2014) in order to match model expectation and observation. Some examples are given in Suzuki et al. (2003), where a Kalman filter is applied in order to estimate traffic conditions based on a mix of loop detector and FC data. Also, results presented in Yuan et al. (2012) and Duret and Yuan (2017) studied the benefits of a Lagrangian model compared to Eulerian approaches when applied to detector and FC data. In addition, higher order models have been proposed that account for more sophisticated traffic dynamics in Aw and Rascle (2000). Wang and Papageorgiou (2005) describe a second order model that estimates traffic conditions on a freeway in real-time. Computational issues are addressed by van Hinsbergen et al. (2012), where the authors demonstrate the development of a localized filter that performs real-time computations. All models mentioned above rely strongly on flow and density data. In contrast, Herrera and Bayen (2010) and Work et al. (2008, 2009) focus on probe data and develop models using a fundamental diagram to be able to estimate densities from probe velocities. Furthermore, Bekiaris-Liberis et al. (2016) propose a macroscopic model for traffic density estimation using a linear parameter-varying system that relies mainly on probe velocity measurements. Although in the latter mentioned approaches most of the information is obtained from probe data, the proposed models still require flow or density measurements at the boundaries. In practical applications, the need for additional flow or density information drastically limits the applicability of an approach on a large scale since it adds further complexity to data acquisition and processing. The second category of algorithms comprises of estimation methods that are based on empirical traffic theory. First, in Kerner et al. (2004) a model called ASDA/FOTO based on Kerner’s Three Phase traffic theory (Kerner, 2004, 2009) is introduced with the purpose of providing current and predictive traffic information. The approach analyzes velocity and flow data provided by detectors and reconstructs spatio-temporal regions of free flow, synchronized flow and Wide Moving Jams. While that approach is completely based on loop detector data, Palmer et al. (2011) study the reconstruction of phase regions with trajectory data exclusively. The phase transitions in space-time of individual vehicles are identified by means of velocity and time conditions and aggregated into phase objects using a clustering approach. The advantages of this approach include the exclusive use of FC data in order to estimate velocities and phases, nonetheless, velocities inside a traffic phase are estimated to a constant over space and time, which in turn limits the accuracy of an estimation. Additionally, the trajectories without phase transitions are discarded, which subsequently results in a loss of valuable information. For example, a trajectory in free flow state passing through an estimated synchronized flow region will not influence the phase since it does not contain phase transitions. Another well-known method is the Generalized Adaptive Smoothing Method (GASM). It is based on the observation that shockwaves in congested traffic propagate upstream and shockwaves in free traffic propagate downstream (Treiber and Helbing, 2003). Using two convolution processes, traffic data is smoothed in these directions and aggregated adaptively. The advantages of the GASM are that it can be applied to velocity data of different sources (Treiber et al., 2010; van Lint and Hoogendoorn, 2009; Jiang et al., 2017) that it allows for an efficient implementation (Schreiter et al., 2010) and that it proved to be significantly more accurate than isotropic smoothing (Treiber and Helbing, 2003; Rempe et al., 2016). On the other hand, it tends to propagate low velocities up- and downstream unconditionally although they might be part of stationary congestion upstream a bottleneck (Treiber et al., 2010). Thus, when it is applied to sparse probe data the estimated velocities in stationary congestion patterns lack accuracy. Inspired by the GASM and works by Kerner et al., a novel Phase-based Smoothing Method (PSM) is presented with the purpose of estimating traffic velocities using FC data with higher accuracy. The key idea is that the PSM distinguishes between the three states: free flow, synchronized flow and Wide Moving Jam (Kerner, 2004, 2009). The distinction between these phases allows one to account for common propagation characteristics of each phase. After identifying these phases, traffic measurements are smoothed accordingly depending on the phase they belong to, which makes it possible to incorporate prevailing traffic measurements into the estimation algorithm. This approach enables higher estimation accuracies, while keeping the efficiency and robustness of the GASM. In the following preliminary section a summary of the Three Phase traffic theory and a description of fundamental convolution operations used frequently in the PSM are given. Section 3 describes the PSM. The first step considers the estimation of phases in space-time. Step two explains how the final velocity estimate is obtained using the combination of phase information and raw data. Evaluation of the method is described in s 4 and 5. The trajectory data reported by a huge fleet of vehicles during a congestion on German freeway A99 is used to assess the performance of the PSM. Finally, the results are discussed and future work is proposed. 2. Preliminaries Two main concepts form the basis of the PSM. The first concept considered is Kerner’s Three Phase traffic theory. In order to motivate how his theory is applied to traffic estimation algorithms, it is summarized briefly and the most important empirical findings relevant for the development of the PSM are presented. The second is the use of smoothing operations. Since FC data are prone to noise due to GPS sampling and further processing and FC data of individual vehicles constitute only a sample of the macroscopic traffic speed, smoothing of data is inevitable in order to remove outliers and average out velocity differences among several trajectories. In discrete space, smoothing is a standard operation and does not require an explanation, though, for reasons of generality 645

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 1. The three phases in vehicle speed vs. vehicle density plot (left) and a contour plot depicting a typical congestion pattern consisting of a WMJ inducing a synchronized flow phase at a bottleneck (right) (Kerner, 2009).

and clarity the PSM is described in continuous space. Therefore, smoothing of trajectories in continuous space is outlined formally in the second part of the preliminary section. 2.1. Summary of the three phase traffic theory In the Three Phase traffic theory the three aforementioned phases are postulated: free flow, synchronized flow and Wide Moving Jams (see Fig. 1). The free flow state describes the state of traffic where flow and density of traffic is virtually linearly related. In this state, traffic demand is lower than the (free flow) capacity of the road. Additionally, there is a possibility of a high spread of velocities between vehicles on different lanes (Kerner, 2004, 2009). A transition from free flow state to synchronized flow state requires high vehicle densities and can happen spontaneously due to local perturbations such as lane changing maneuvers, or can be induced by moving traffic jams propagating upstream (Kerner, 2004, 2009). A transition from free flow state to congested state often infers a capacity drop (Schoenhof and Helbing, 2007; Laval, 2007; Hall and Agyemang-Duah, 1991; Daganzo, 1996). The consequence of a reduced maximal road capacity it that once a synchronized flow state originates and the demand does not change significantly, the synchronized flow state will persist. Typically, the downstream front of the phase sticks to a bottleneck such as an on- or off-ramp, the position of a lane closure, and/or the beginning of a construction site. The vehicles’ velocities in this phase are significantly reduced compared to the free flow phase. Additionally, the variance of velocities among vehicles on different lanes is small, which motivates the name synchronized flow. The emergence of WMJs is also a spontaneous process. Microscopically, it can be explained with the over-deceleration effect (Kerner, 2009). The effect describes a vehicle with relatively small time headway to a preceding, decelerating vehicle which tends to decelerate stronger than necessary. As following vehicles behave similarly a shockwave emerges that propagates upstream forming the upstream front of a WMJ phase. If vehicle density is high, the over-deceleration effect is more likely to occur. The mean velocity of vehicles after deceleration can decrease down to 0 km/h. The moments when vehicles start to accelerate again form the downstream front of a WMJ. Its velocity has been estimated by many researchers to an approximate value of −15 km/ h (Lighthill and Whitham, 1955; Richards, 1956; Kerner, 2004, 2009; Treiber and Helbing, 2003; Treiber et al., 2010; Schoenhof and Helbing, 2007). 2.2. Data representation and convolution Since trajectory data provides the foundation of traffic speed estimation in this article, data representation and information density are described in a general form. The trajectory of a vehicle c ∈ {1,…,Nc } can be described as a function x c (t ) ∈ [0,L] denoting the position of that vehicle along a road segment with length L ∈  +. Accordingly, vehicle velocity vc (t ) denotes the derivative of x c . Each vehicle that passes through space-time domain [0,L] × [0,T ] of a road with length L , observed for time period T , provides partial information about the domain. In order to model the space-times for which information is available, a simple car-following model with parameters vehicle length and time headway is assumed. x 0 denotes the length of the vehicle c and the minimal distance to the preceding vehicle in queueing traffic. T0 denotes the time headway to a preceding vehicle. Therefore, at time t , vehicle c occupies the space interval [x c (t ),x c (t ) + x 0 + T0 vc (t )] (see Fig. 2 left). Occupying a region in space-time means that in the interval only one vehicle can exist and, furthermore, it is assumed that in this interval, the vehicle’s velocity is a representation of the traffic velocity. Note that this holds only for single-lane roads. Let Ψ(t ,x ) be a function that indicates whether time t and position x is occupied by any observed vehicle:

1 if ∃ c: x c (t ) < x < x c (t ) + x 0 + T0 vc (t ), Ψ(t ,x ) = ⎧ ⎨ ⎩ 0 otherwise.

(1)

Ψ(t ,x ) describes a binary system which has the following properties: if Ψ(t ,x ) = 1 ,∀ t ∈ [0,T ],x ∈ [0,L], then the traffic density is high, 646

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 2. Illustration of the region Ψ occupied by a vehicle with velocity vc (left), convolution process ΓVFCD for point (t ,x ) given Ψ of several trajectories, respective velocity VFCD and convolution kernel Φ .

all vehicles maintain a time headway of T0 to their preceding vehicles, and all vehicle positions as well as velocities are known. Otherwise, there are space-times for which no velocity information is available. The reason can be that traffic density is low such that there are gaps between vehicles, or for a part of the vehicles no position and velocity data is given. Accordingly, VFCD (t ,x ) ∈  + denotes the velocity information reported by all vehicles, combined into a single two dimensional function. This velocity is only valid in space-time (t ,x ) which is occupied by any vehicle c∗, i.e. Ψ(t ,x ) = 1, and is set to the velocity of the closest vehicle upstream of (t ,x ) :

v c∗ (t ) if Ψ(t ,x ) = 1,where c∗ = argmin∀ c : x − x c (t ) ⩾ 0 (x −x c (t )) VFCD (t ,x ) = ⎧ ⎨ otherwise. ⎩− 1

(2)

2.3. Convolution of trajectory data with smoothing kernels In order to filter data, the two-dimensional convolution operation is commonly used. The function:

ΓVFCD (w,Φ,t,x) =

T

∫0 ∫0

L

̂ w ( t ,̂ x )· ̂ VFCD ( t ,̂ x )·Ψ( ̂ ̂ ̂ Φ(t − t ,̂ x −x )· t ,̂ x )̂ dx dt

(3)

represents a weighted continuous convolution, where the function ϕ (t ,x ) ∈  denotes the kernel and w (t ,x ) ∈  represents a spacetime-dependent weight of the input data VFCD (t ,x ) . Definitions of kernel functions and weightings are given in section 2.4. In order to use the convolution equation for smoothing operations, the results of ΓVFCD need to be normalized. The normalization term is similar to the aforementioned convolution but omits velocity input VFCD :

D (w,Φ,t,x) =

T

∫0 ∫0

L

̂ w ( t ,̂ x )·Ψ( ̂ ̂ .̂ Φ(t − t ,̂ x −x )· t ,̂ x )̂ dx dt

(4)

Furthermore, the kernel Φ(t ,x ) and weighting w (t ,x ) are chosen in such a way that D (w,Φ,t ,x ) represents an estimate of the local data density. The local density quantifies the amount of information available for the velocity estimation in (t ,x ) . The normalized convolution of weighted and velocity function Vc (t ,x ) with kernel Φ(t ,x ) is given by:

VVFCD (w,Φ,t ,x ) =

ΓVFCD (w,Φ,t ,x ) . D (w,Φ,t ,x )

(5)

Given, for example, a kernel function that returns values greater or equal to zero with its maximum at (0,0) , Eq. (5) describes a common smoothing process (see Fig. 2 right). Then, for space-time (t ,x ) a weighted average velocity of all nearby velocities valid in their respective occupied regions Ψ(t ,x ) is computed. The weights depend on the kernel function, distance to (t ,x ) , sampled data as well as the weighting function w (t ,x ) . Notice that Eq. (5) implies that the speeds of those vehicles which occupy more space-time (either because they are longer, or they driver at higher speeds) are weighted stronger than the speeds of those which occupy less space-time. This pertains the interpretation of the resulting macroscopic traffic speed. As explained thoroughly in Treiber and Kesting (2013), there exist several ways to aggregate individual vehicle speeds over time and space such as the (arithmetic or harmonic) time mean space, the space mean speed or other methods such as Edie’s (Edie, 1963). The average that the PSM resembles has similarities with all of these definitions, but differs in one important aspect: Individual vehicles are converted into homogeneous regions in space-time for which a velocity is available. During that process the information whether a region was occupied by one or several vehicles may get lost: the region may emerge from several short vehicles following each other with time headways T0 , or one long truck. Consequently, the smoothed resulting velocity does not correspond to the mean vehicle speed averaged over all reported vehicles but to the mean speed averaged over all occupied regions in space-time. In that sense, the presented way to average is just another (slightly different) way to define a macroscopic traffic speed. It is motivated by the representation of vehicle speed’s with validities Ψ , which can be seen as an approach to advance from a macroscopic to a mesoscopic data representation in time and space. 647

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

2.4. Traffic-motivated convolution kernels In general, modifying the convolution kernel Φ is seen in many applications, however, in this article it is used mainly for smoothing operations. For traffic related smoothing with characteristic wave velocity vdir , the GASM uses a simple and effective formulation kernel where

⎧ ⎛ ⎪ exp ⎜− τ ,σ Φvdir (t ,x ) = ⎝ ⎨ ⎪ exp − ⎩

(

t−

x vdir

τ t τ

−

x σ

−

x σ

⎞ if vdir ≠ 0 ⎟ ⎠

)

if vdir = 0.

(6)

It defines a main smoothing direction that is used to interpolate velocities in a shockwave characteristic direction (see (Treiber and Helbing, 2003) for more details). Its maximal value is located at (0,0) and function values decay exponentially, converging to zero. In Fig. 2 (right) a kernel with negative velocity vdir is depicted as a discrete contour plot (red shape). 3. A phase-based smoothing method 3.1. Overview The development of PSM uses a probabilistic approach. The advantage is that in addition to a final velocity estimate also a trustworthiness of the estimate can be determined. In the following discussion, the steps for calculating the final velocity are introduced in order to provide an overview and also define a few fundamental variables. Exact definitions are given in the following sections. Fig. 3 summarizes the workflow for processing raw trajectory data VFCD into a continuous velocity estimate VE (t ,x ) . Let Ω = {F ,S,J } be a set of the phases, free flow (F), synchronized flow (S) and Wide Moving Jam (J). Then, Pp (t ,x ) ∈ [0,1] with p ∈ Ω denotes the probability for the traffic at time t and at position x to be in phase p . As first step, raw data is convolved with different smoothing kernels. Several fuzzy phase criteria are defined and applied to resulting values. Respective criteria probabilities Ppi denote their degree of fulfillment. Next, for each phase, the criteria probabilities

Fig. 3. Flow-chart of the steps to process raw trajectory data VFCD into a continuous velocity estimate VE .

648

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

are aggregated into preliminary phase probabilities Pp ′. The probability PU ∈ [0,1] is estimated to establish the level of uncertainty in assigning (t ,x ) to any of the phases. The relations between phases are modelled in the following step, which results in the final phase probabilities Pp . Based on these and raw trajectory data, for each phase and each point in space-time, a velocity estimate VpH (t ,x ) is computed. Additionally, a fallback velocity VU (t ,x ) is assumed that serves as a best-guess velocity in case the uncertainty PU (t ,x ) is high. Finally, the resulting estimate VE (t ,x ) is determined by aggregating the probabilities Pp (t ,x ) , PU (t ,x ) and their respective velocity estimates VpH (t ,x ) and VU (t ,x ) . 3.2. Modelling phase probabilities As described in ection 2.1 each traffic phase has different characteristics. These empirical characteristics can be used to identify the most likely traffic phase p in space-time (t ,x ) . Let Pp1,Pp2,…,PpNk , Ppi ∈  be a number of Nk (p) ∈  ⧹ {0} criteria that (t ,x ) needs to fulfill in order to belong to phase p . Each criterion is modelled as a fuzzy decider Ppi ∈ [0,1]. From there the independent probability Pp′ (t ,x ) is determined as the product of all criteria probabilities for each phase individually: Nk (p)

Pp′ (t ,x ) =

∏

Ppi (t ,x ).

(7)

i

Ppi .

Consequently, the preliminary phase probability Pp′ is always lower or equal to the lowest criteria probability In other words all criteria need to be fulfilled to a certain degree to be able to assign a point in space-time to phase p . In applications where several phase criteria are proposed, it might be desired to accumulate criteria of the less rigid probabilities such that the failure of one or more criteria will be tolerated. This would require Eq. (7) to be adapted. The applied fuzzy logic for now only makes use of the logical ‘AND’ relation as the product of probabilities. An extension with a logical ‘OR’, modelled as a OR b = 1−(1−a)(1−b) = a + b−ab with a,b ∈ [0,1], allows defining more complex rules. Another potential candidate would be an ‘opinion pool’ described in Jacobs (1995) where several classifiers are combined with respect to their expected trustworthiness. For now, two classes of criteria are presented and later applied. The first is the velocity criterion Pp1 that uses velocity information of smoothed data for determining the phase probability. The second, density criterion Pp2 is applied ensuring that a phase hypothesis is supported by nearby data. Table 1 is provided below to give an overview of all probabilities and will be explained thoroughly in the following sections. Note that the previous description does not adhere to classical probability axioms since it causes a contradiction in this context. E.g. if one considers a hypothesis H with a probability of P (H ) that it is true, then the hypothesis being false will have probability 1−P (H ) . Now, assume for each phase p there is at least one fuzzy decider Ppi that is zero. This in turn causes all preliminary phase probabilities Pp′ to be also zero. In a classical interpretation of probabilities it would follow that traffic does not appear in one of the phases. This naturally contradicts the theory model since according to the theory, traffic must always be in one of the three phases. Therefore, with respect to the evidence theory (Shafer, 1976), the modelled phase probabilities Pp′ need to be interpreted as independent beliefs that a phase hypothesis is correct. Thus, each probability is an estimate of a probability. On the other hand, the probability 1−Pp′ needs to be interpreted as the probability that traffic is in any one of the three phases. The higher 1−Pp′ , the more uncertain is a classification as phase p . That principle of uncertainty is important to fuse data and provide a quality estimate of reconstructed velocities. As a consequence, probabilities in Table 1 do not sum up to one. The probabilistic approach allows the usage of further criteria instead of just the presented two. For example, if flow or density information is available, phase characteristic ranges of flow and density can be used to further differentiate between free flow, synchronized flow and WMJs. Another idea is to fuse heterogeneous traffic information. Traffic messages provided by public or commercial traffic feeds can be used as a prior to provide an expectation of a synchronized flow phase. Subsequently, real measured data can be used to refine jam fronts and traffic dynamics inside the phase. 3.2.1. Velocity criterion According to Kerner (2009), traffic velocity information is crucial in order to determine to which traffic phase space-time (t ,x ) probably belongs. The velocity information measured in phase-characteristic directions around (t ,x ) are used to determine the probability Pp1. Modern traffic theories state that traffic breakdown, which is the transition from free to congested flow, is a probabilistic event triggered by perturbations (Kerner, 2009; Schoenhof and Helbing, 2007) A traffic breakdown is usually connected to a capacity drop (Schoenhof and Helbing, 2007; Laval, 2007) and a significant drop in average velocities (Kerner, 2009; Schoenhof and Helbing, 2007). Therefore, velocity is a good candidate in order to distinguish between free flow and congested flow. Here, the fuzzy thresholds vFthres and vSthres are applied to differentiate between these two flow regimes. Table 1 Overview of probabilities computed during phase identification using a velocity and a density criterion. Velocity Criterion Pp1

Density Criterion Pp2

Preliminary Phase Probability

Final Phase Probability

Free Flow

PF1

PF2

PF′ = PF1 ·PF2

PF = PF′ ·(1−PJ′ )

Synchronized Flow

PS1

PS2

PS′ = PS1·PS2

PS = PS′ ·(1−PJ′ )

Wide Moving Jam

PJ1

PJ2

PJ′ = PJ1·PJ2

PJ = PJ′

Uncertain State

PU = (1−PF′ )(1−PS′ )(1−PJ′ )

649

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 4. Sigmoid functions parametrized for different phases.

The distinction between WMJ and synchronized flow, as congested regimes, is less obvious. In fact, the upper velocity of the WMJ is significantly lower when compared to the threshold to free flow. We set threshold vJthres to 30 km/h (Kerner et al., 2004). Also, velocities in this particular phase can decrease down to 0 km/h, therefore in such regions, no lower bound is required. Lastly, we will discuss the lower velocity bound of the synchronized flow phase. Kerner suggests that low velocities should be assigned to the WMJ phase. He argues that the type of congestion pattern changes from synchronized pattern, over general pattern to mega-jams (Kerner, 2009; Helbing et al., 2009) with stronger bottlenecks. According to his definition, mega-jams are WMJs with a great width, which have a very low average speed. In two-phase theories these patterns are usually called Homogeneous Congested Traffic (HCT) (Schoenhof and Helbing, 2007). In contrast to the ASDA/FOTO model (Kerner et al., 2004) which applies a lower velocity threshold for the synchronized flow, in this approach no threshold is set. Instead, a dominance of the WMJ phase is modeled (see ection 3.2.4) which assures that in the presence of low velocities the probability of a WMJ is increased. This turned out to be more robust than the definition of a lower threshold. A fuzzy decider function σ (v,v thres,λ ) (sigmoid function) is applied in order to translate a velocity v into a probability Pp1 ∈ [0,1]. The parameters are thresholds v thres and λ , which determines the strictness of the threshold, i.e. the higher λ , the higher the gradient of the transition region:

σ (v,v thres,λ ) = 1−

1 . 1 + exp(−λ (v−v thres ))

(8)

Fig. 4 depicts the results of the sigmoid functions parametrized for free flow (blue), synchronized flow (red), and WMJ (green) with respect to different velocities. In the following the velocity criteria are defined for the three phases. 3.2.1.1. Phases: free flow and synchronized flow. Probabilities PF1 (t ,x ) and PS1 (t ,x ) are computed as:

PF1 (t ,x ) = 1−σ (VF (t ,x ),vFthres,λF )

(9)

PS1 (t ,x ) = σ (VS (t ,x ),vSthres,λS )

(10)

vpthres

where and λp denote the parameters of the decider function with respect to the characteristic velocity of phase p . The velocities VF and VS are computed according to the normalized convolution process Eq. (5) as:

VF (t ,x ) = VVFCD (w 0,ΦF ,t ,x )

(11)

VS (t ,x ) = VVFCD (w 0,ΦS ,t ,x )

(12)

vpdir

with ΦF and ΦS denoting phase specific smoothing kernels, with characteristic velocities of and parameters τp and σp . The standard weighting w 0 is defined as w 0 (t ,x ) = 1. The standard weighting implies that all smoothed raw data have an equal significance for the determination of the phase. 3.2.1.2. Phase: wide moving jam. The velocity criterion for the J phase is more distinct. It is postulated that (t ,x ) can only be assigned to the J phase if, both, up and downstream of (t ,x ) low velocities are observed. That ensures the WMJs are not extrapolated far beyond measurements in order to be able to reduce wrongly estimated congested regions. The trick, in order to compute this condition efficiently, is to check for measurements up- and downstream of (t ,x ) smoothing data with differing kernels and compute the probabilities independently. Then, the product of both probabilities represents the need to fulfill both requirements. Therefore, PJ1 is 650

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

computed as:

PJ1 (t ,x ) = σ (V Ju (t ,x ),vSthres,λS )·σ (V Jd (t ,x ),vJthres,λJ )

(13)

where V Ju and V Jd denote the velocity estimates computed as:

V Ju (t ,x ) = VVFCD (w 0,ΦuJ ,t ,x ),

(14)

V Jd (t ,x ) = VVFCD (w 0,ΦdJ ,t ,x ).

(15)

Different kernel functions ΦuJ and ΦdJ defined as: ,σ Φ τdir (t ,x ) if x ⩽ 0 ΦuJ (t ,x ) = ⎧ vJ ⎨0 otherwise ⎩

(16)

,σ Φ τdir (t ,x ) if x ⩾ 0 ΦdJ (t ,x ) = ⎧ vJ ⎨0 otherwise ⎩

(17)

ΦdJ

ΦuJ

are applied. is defined in such a way that it considers upstream data of (t ,x ) , where only considers downstream data. In that way, the data is smoothed in different directions. The combination of the sigmoid functions of the velocities V Ju and V Jd ensures that space-time (t ,x ) is only estimated as a WMJ if upstream and downstream data supports the phase hypothesis. Note that Eq. (13) uses different velocity thresholds vJthres and vSthres for the sigmoid functions. The difference stems from the expectation that a WMJ, which is detected by downstream velocities below vJthres , will propagate upstream as long as the upstream traffic is in a state of critical flowdensity (Kerner, 2009). Since no density or flow data is available, this state is assumed to be the congested region with the velocity threshold vSthres . By requiring both criteria to be fulfilled it is assured that the J phase is only reconstructed for low velocity measurements but never extrapolated. It is quite possible that the moving jam emerged earlier than the time from when the first equipped vehicle perceived it and propagated further upstream than the last equipped vehicle passing through the moving jam. Unfortunately, sparse data does not allow one to recognize exactly when a WMJ emerged and when it dissolved. Extrapolating a shockwave upstream or downstream means to risk overestimating. Thus, this approach can be described as the cautious way aimed at minimizing wrongly estimated low velocities. 3.2.2. Density criterion The second criterion, called the density criterion Pp2 (t ,x ) uses data density D Eq. (4), in order to quantify how well a phase hypothesis is supported by nearby data. The necessity for this criterion stems from varying data density that comes along with FCD. Data density D (w 0,Φp,t ,x ) is computed for each phase p using the respective kernel Φp . Since weighting w 0 and the applied kernels are greater than zero, D (w 0,Φp ,t,x) is always greater or equal to zero. Note that it is not normalized such that its values often exceed the value one. In order to translate density into the probability Pp2 , its values are limited to an upper bound of one:

Pp2 (t ,x ) = min(1,D (w 0,Φp,t ,x )).

(18)

This forces the density criteria Pp2 to equal one if data is nearby, and otherwise converge to zero when the distance between (t ,x ) and the measurements grows larger. The validity of a measurement in space and time can be parametrized by adapting the kernel function or by modifying the weighting w 0 . Note that the applied minimum operator is one way to translate densities into probabilities and which is chosen due to its simplicity. Other operators with smoother transitions such as sigmoid functions may also be applied. 3.2.3. Probability of uncertain state After computing the aforementioned probabilities Pp ′ as the product of all phase criteria Pp1 and Pp2 , uncertainty PU (t ,x ) is determined. PU models the degree of uncertainty in assigning (t ,x ) to any of the phases:

PU (t ,x ) = Πp ∈ Ω (1−Pp′ (t ,x )).

(19)

For space-time (t ,x ) with a high uncertainty a best-guess velocity VU needs to be assumed. In this paper, the speed limit serves as the best guess, though other kinds of velocities, e.g. historical averages, can be considered. Accordingly, a probability Q (t ,x ) ∈ [0,1] can be estimated that estimates the certainty of the result. It can be interpreted as a quality measure and is computed as: (20)

Q (t ,x ) = 1−PU (t ,x ).

3.2.4. Dominance of the J phase Up until now, the independent phase probabilities PP′ have been determined. Since Pp′ are independent, particular region characteristics could occur, especially in the presence of shockwaves, where both, PJ′ and PS ′, or PF′ probabilities are all estimated to be high. That is due to the different shapes of the convolution kernels and according velocities that are considered for the determination of the phase. Since WMJs can propagate through other phases without interruption (Kerner, 2009), and probability PJ′ is supported by more distinctive criteria, it is reasonable to assume that these regions belong to the J phase rather than to one of the others. 651

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Consequently, in those cases dominance of the J phase over the other phases is modelled (comparable to the GASM where smoothed low velocities are given a higher weight). Final phase probabilities PP are set as:

PJ (t ,x ) = PJ′ (t ,x ),

(21)

PS (t ,x ) = PS′ (t ,x )·(1−PJ (t ,x )),

(22)

PF (t ,x ) = PF′ (t ,x )·(1−PJ (t ,x )).

(23)

3.3. Estimating phase velocities Up to this step, raw data has been used in order to reconstruct the traffic phases in space and time using characteristics of smoothing kernels and several other criteria that allow one to distinguish between them. For estimating traffic velocities it is very useful to know about the locations of these phases. They provide the information about the space and time for which velocity measurements are valid. E.g. assume a low velocity measurement that is estimated as part of a WMJ. In this case, the velocity measurement can be used to estimate the average traffic velocity of the WMJ. However, it does not provide information about the traffic velocity in an adjacent traffic free flow or synchronized flow phase. The phase-dependent velocities VpH are estimated as:

VpH (t ,x ) =

1 V V −1 (Pp,Φ Hp ,t ,x )

(24)

FCD

with −1 V FCD (t ,x ) =

1 VFCD (t ,x )

where V V −1 (Pp,ΦH p ,t,x) FCD volution kernel Φ Hp (A

(25)

−1 with the condenotes the convolution process (Eq. (5)) applied to the reciprocal trajectory velocities V VFCD

minimal velocity of 3 km/h is assumed). Instead of a normal weighting w 0 used previously (Eq. (11)), phase probabilities Pp are used as the weights for the input data. Effectively, measurements that are assigned to phase p have a greater influence on VpH than those that are assigned to another traffic phase. From smoothing the reciprocal velocities it follows that the smoothing process resembles a harmonic mean instead of an arithmetic one, which accounts for precision in travel time reconstruction (van Lint, 2010). Note that the convolution kernels Φ Hp can be parametrized differently from the convolution kernels Φp . The reason is that the shapes of Φp account for the characteristic propagation of phases, such as the propagation of WMJ phases upstream and the stationary character of S phases. Φ Hp , on the other hand, influences how data inside its respective phase is smoothed. A distinction between the parameters allows for a reconstruction of a stationary synchronized flow phase using a large kernel ΦS and a subsequent estimation of velocities inside the phase with a smaller kernel ΦSH which reconstructs minor shockwaves (narrow moving jams (Kerner, 2009)) occurring inside the phase. 3.4. Aggregating probabilities and velocities At this stage, all phase probabilities Pp , respective velocities VpH and probability PU are estimated. These quantities as well as the best-guess velocity VU are aggregated into the final velocity estimate VE as a weighted average:

VE (t ,x ) =

PU (t ,x )·VU (t ,x ) + ∑p ∈ Ω Pp (t ,x )·VpH (t ,x ) PU (t ,x ) + ∑p ∈ Ω Pp (t ,x )

. (26)

Other variants such as maximal values, harmonic averages, etc. are potential candidates for further studies. 4. Evaluation In the following, a study is conducted in order to assess the accuracy of the PSM estimating traffic speeds with FC data. The data was collected during a traffic jam on German freeway A99 on July 15th, 2014. First, the characteristics of the data and the discretization of the PSM are described briefly. In order provide an intuition about the method, resulting probabilities and velocities from the PSM are illustrated. In a subsequent comparison, the estimation quality of the PSM, the GASM, and an isotropic smoothing algorithm are compared. Next, a runtime analysis gives insights into the efficiency of the method. Finally, the sensitivity of the parameters on the estimation accuracy is presented. 4.1. Data Available data consists of timestamps and GPS positions sampled anonymously by individual vehicles with sampling times between 20 s and 30 s. Each vehicle is equipped with a filter mechanism that decides whether or not a sampled GPS positions is sent to the central server. Reason for this filter mechanism is privacy since individuals cannot be tracked along their journey. The filter mechanism continuously compares the vehicle’s velocity with the speed limit of the road or, if available, reported traffic speeds by a 652

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 5. Raw trajectory data reported for a congestion pattern on 15-July 2014 on A99 in west-bound direction. Gray triangles mark the positions of on- and off-ramps.

traffic information service. If over a period of time, up to a few minutes, only small deviations in the order of 10–30% occur, the recently sampled positions are retained and not reported. In case of greater deviations, positions and corresponding timestamps are transmitted. Effectively, the central server receives most trajectory data in case of traffic disturbances. Due to high dynamics in congested regions, reported traffic information by providers do not match the vehicles’ velocities exactly over time, such that only few trajectories are filtered out. This in effect considers the drivers’ privacy while providing detailed velocity information in case traffic is disturbed. In Fig. 5 raw trajectory data for the congestion on A99 in westbound direction are displayed. The pattern shows different characteristics often occurring in congested freeway traffic (compared to traffic patterns described in Kerner (2009) and Helbing et al. (2009)). At around 7:45 am a moving jam emerged that evolved into a WMJ and induced a traffic breakdown at the on-ramp at kilometer 10. The WMJ propagated further upstream and induced another traffic breakdown at the neighboring bottleneck. The pattern evolved into a General Pattern (GP) expanding over two bottlenecks where the downstream fronts of synchronized flow phases are fixed slightly downstream of the on-ramp positions. In the pinch region of the downstream synchronized flow phase a few WMJs originated and propagated upstream. 4.2. Discretization In order to apply estimation methods, time and space dimension are discretized into intervals of length ΔL = 50 m and ΔT = 10 s called grid cells. Average velocities between GPS positions and respective timestamps are computed and the vehicle occupied regions are determined (see ection 2). Each of the Nc trajectories is written as a set of tuples tr = {(t ,x ,v,ψ)1,…,(t ,x ,v,ψ)ntr } where each tuple represents the mean velocity of the vehicle at space-time (t ,x ) referring to the center of the grid cell in space and time, and cell

Fig. 6. Discrete region Ψ for the convolution process working as weights for the convolution processes. High values indicate a high coverage with data.

653

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 7. Resulting phase probabilities Pp for the Free Flow (left), Synchronized Flow (mid) and Wide Moving Jam phase (right).

occupation ψ ∈ [0,1]. If two velocities of different trajectories are assigned to the same cell, their mean value is determined and the occupation values are added. Fig. 6 depicts the discrete function Ψ(t ,x ) that stems from the individual occupations of the reported trajectories. According to average time headways measured in real traffic, headway T0 is set to one second and x 0 to 6 m (Krbalek et al., 2001; Brackstone et al., 2002). At this step, one benefit of the definition of occupation region Ψ(t,x) becomes visible. Discretizing this function allows to decouple the data from the chosen grid size. If the grid cells are larger in time and/or space, passing vehicles occupy smaller regions of a cell and the weighting decreases. On the other hand, if ΔL is smaller than the occupied region of a vehicle, the vehicle can occupy more than one cell at the same time. The traffic conditions are assumed to be similar on all lanes. As a result, all the collected data is mapped to one lane. This is a simplification since vehicles are able to overtake in free flow but in congested flow the assumption holds more or less. Besides the trajectory data, the continuous 2D convolution (Eqs. (3) and (4)) needs to be discretized. An efficient implementation of the discrete 2D convolution using the Fast Fourier Transform (FFT) is described thoroughly in Schreiter et al. (2010). 4.3. Applied PSM Applying the PSM to the data depicted in Fig. 5, the following intermediate and final results are obtained. Parameters of the method are given and discussed in the following section. After computing criteria probabilities Pp1 and Pp2 and their interactions, final phase probabilities Pp are available (Fig. 7). Raw data is weighted with these phase probabilities and smoothed using respective kernels, providing velocity estimates VpH (Fig. 8). The final velocity estimate VE (Fig. 9 left) results from aggregating phase probabilities, phase velocities, uncertainty PU , and a fallback velocity VU . The corresponding quality estimate Q (t ,x ) (Fig. 9 right) describes the trustworthiness of a computed velocity. It is close to the value one if data is nearby and converges to zero if few data is in the proximity of (t ,x ) . 4.4. Comparison algorithms For an objective evaluation, a comparison with other algorithms is conducted. Therefore, two common approaches that do not require density or flow data are taken into consideration. One is the GASM. Here, the GASM is applied as described in Treiber et al. (2010) and van Lint and Hoogendoorn (2009) with an adaption to sparse FCD as presented in Rempe et al. (2016). The adaption describes how sparse FC data and a best-guess velocity can be fused in order to provide a continuous velocity estimate if no data is nearby. The weighting ratio of FC data to the velocity fallback is 1000:1. Parametrization is chosen according to Treiber et al. (2010). As far as we know, isotropic smoothing is still widely used in many traffic related applications. Thus, besides the PSM and GASM, which incorporate findings from empirical traffic, also an isotropic smoothing approach is applied that smooths data mostly in time and slightly in space.

Fig. 8. Velocities VpH of Free Flow (left), Synchronized Flow (mid) and Wide Moving Jam phase (right). Regions where the respective phase probabilities are below 0.1% are colored gray.

654

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 9. Resulting velocity estimate VE (left) and quality estimate Q (right).

The applied PSM parameters are listed in Table 2. The velocity thresholds vpthres have already been motivated in ection 3.2. The parameters for the sigmoid functions λp are set equal. The propagation direction of WMJ phase fronts vJdir is set similar to Kerner (2009). Propagation direction vSdir and vFdir describe the expected propagation velocity of phase fronts of synchronized flow and free flow phases. According to the Three-Phase traffic theory, synchronized flow phases are rather stationary or, if not, their downstream fronts usually stick to adjacent upstream bottlenecks and remain there. Therefore, vSdir is set to 0 km/h. Since the downstream synchronized flow phase front is the upstream free flow phase front, which consequently has a stationary character as well, we set vFdir also to 0 km/h. The other parameters are mainly the parameters τ and σ , which influence the decay of the kernel function in time and space. The greater the value the slowlier the kernel decays, such that more distant measurements influence the estimation of a velocity estimate at space-time (t ,x ) . The kernels Φp have been parametrized with respect to typically observed properties of congestion: WMJs often have a high spatial extent but their width in time is limited (Kerner, 2009; Schoenhof and Helbing, 2007). Consequently, σJ is chosen relatively high and τJ relatively low. Due to the stationary character of the synchronized flow phase and the nearby free flow phase, τS,F is set significantly greater than τJ and σS,F significantly lower. The kernels Φ Hp describe the propagation of shockwave inside the phases. WMJs and synchronized flow phases are congested traffic phases. Since shockwaves in traffic are well-understood phenomena in traffic dynamics (Richards, 1956; Kerner, 2009; Newell, ,H is set according to these publications. vFdir is set along the same line of reasoning. In order to reconstruct 1993; Mika et al., 1969), vJdir ,S minor shockwaves also called narrow moving jams in synchronized flow with short temporal width and possibly shorter spatial extent than a fully developed WMJ, τJH,S is set similar to τJ and σJH,S is set smaller than σJ . The kernel parameters τFH and σFH are set to average values compared to all other parameters. Varying these two parameters did influence the estimation accuracy of this example significantly (see ection 5.4).

5. Results In the following, estimation accuracy of the PSM, the GASM and the isotropic smoothing approach are compared qualitatively and quantitatively. Subsequently, the results of a runtime analysis and a parameter study are presented.

Table 2 PSM parameters. Parameter⧹Phase

J

S

F

vpthres

30 km/h

65 km/h

55 km/h

λp

0.5 h/km

vpdir

−18 km/h

0 km/h

τp

30 s

250 s

σp

500 m

150 m

vpdir ,H

−18 km/h

70 km/h

τ pH

30 s

100 s

σpH

200 m

100 m

655

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 10. Ground Truth (GT) as an aggregation of trajectory data (left). Sparse vehicle trajectories for the comparison of estimation algorithms (right).

5.1. Qualitative results Fig. 10 left illustrates the Ground Truth (GT) VGT (t ,x ) as the aggregation of all trajectory velocities in a resolution of 60 s and 100 m. In order to compare several estimators, a subset of trajectories is drawn randomly from all reported trajectories. The accuracy of an algorithm is evaluated as its ability to reconstruct the GT given only the aforementioned subset. Fig. 10 right depicts one random subset of trajectories that serves as input VFCD (t ,x ) for the following traffic speed estimations. In Fig. 11 the velocity estimate VE (t ,x ) computed with the PSM and the absolute error |VE (t ,x )−VGT (t ,x )| to the GT are depicted. The error plot reveals that most regions of the estimation match well with the GT. Both, the stationary congestion at the bottlenecks at kilometers 8 and 11 as well as the WMJs are accurately reconstructed. Significant differences are marked in the plot: For instances, at (a) the moving jam is reconstructed inaccurately. At (b) the stationary congestion at kilometer 2 starts earlier than estimated and is therefore underestimated. At (c), the PSM overestimates the extent of the congestion at the upstream front for a short range in time. In Fig. 12 the velocity estimates computed with the naive isotropic smoothing (upper left) and the GASM (upper right) are illustrated. In order to compare the reconstruction accuracy of the PSM and these algorithms directly, the differences of the estimation accuracies are computed:

VΔ = |VE ,1−VGT |−|VE ,2−VGT |.

(27)

This error term returns positive values if the error of the first velocity estimate VE,1 dominates, and is negative if the second velocity estimate is more inaccurate. In Fig. 12 VE,2 denotes the PSM such that blue colored regions indicate that the PSM reconstructed the velocity more accurately than the comparison method. The plot corresponding to the isotropic smoothing method reveals regions around the moving jams (a), (b), (c) where isotropic smoothing results in large errors as it does not account for moving traffic patterns. At (d) this method estimates the stationary traffic slightly better than the PSM. The GASM manages to estimate the dissolving jam at (a) better than the PSM. On the other hand, it does not reconstruct the stationary traffic patterns well (b), (c), (d), but estimates it as free flow. Furthermore, the moving jams are extrapolated further than they actually occurred (e). Overall, the difference plots reveal that both algorithms result in larger erroneous regions than the PSM, which combines the strengths of both comparative algorithms: Reconstruct stationary as well as moving traffic patterns. The following quantitative comparison shall provide further results with respect to varying data input.

Fig. 11. Velocity estimate computed using the PSM and the sparse vehicle trajectories (left). Error of the estimate compared with the ground truth (right).

656

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 12. Velocity estimate using an isotropic smoothing method (upper left) and the GASM (upper right). Difference of estimation errors of both algorithms compared to the estimation error of the PSM (lower left and lower right).

5.2. Quantitative results In order to evaluate an algorithm quantitatively, the set of Nc trajectories is divided into a training set and a test set. The training set is used to estimate velocities VE , and the test set is used to evaluate the accuracy of that estimate. The size of test set STest = {tr1,tr2,…,trNTest } relates to the total number of trajectories as: (28)

NTest = αNc

with 0 < α < 1. The size of NE of the training (estimation) set SE is the (1−α ) part of all trajectories. Additionally, that part is varied with another factor 0 < β < 1 that is used to simulate different data densities:

NE = β (1−α ) Nc.

(29)

For the sake of interpretation, NE is normalized with the analyzed time interval T , resulting in the data density δE :

δE =

NE . T

(30)

Thus, δE denominates the number of trajectories per hour that are used for the estimation process. In Huber et al. (2014) a variety of quality measures are proposed and compared. Since a deviation between estimate and ground truth is more critical in congested regimes than in free flow, the Inverse Mean Absolute Error (IMAE) as metric is applied:

IMAEVE =

1 ∑tr ∈ STest |tr|

∑

∑

tr ∈ STest (t ,x ,v ) ∈ tr

1 1 − . VE (t ,x ) v

(31)

Fig. 13 depicts the IMAE with respect to data densities δE between 10 traces/h and 90 traces/h for the isotropic smoothing approach, the GASM, and the PSM. The presented error values are the average IMAEs after 100 iterations with randomly assigned training set and test set in order to ensure robustness of the results. The results show that all approaches achieve more accurate results with higher data densities. The traffic-motivated methods GASM and PSM achieve significantly more accurate results for all data densities than the isotropic smoothing. Compared to the other two algorithms, the PSM shows the lowest error values for all data densities and thus, can be entitled as the most accurate algorithm in this study. 5.3. Mega-jams The previous study focused on the reconstruction of a GP (Kerner, 2004) incuding a pinch region in which WMJs emerge. In this 657

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 13. Resulting IMAEs resulting from an isotropic smoothing approach, the GASM and the PSM for the estimation of velocities with respect to varying training data densities δE .

section the pattern mega-jam or, also-called HCT, is highlighted (Kerner, 2004; Helbing et al., 2009) (see ection 3.2.1). Although this pattern is rare (Schoenhof and Helbing, 2007), it is important for all types of ITS applications since travel-times increase significantly. In the Three-Phase theory, a mega-jam is described as ‘a wide moving jam with an extremely great width growing continuously over time’ (Kerner, 2004). Although this pattern is classified as a WMJ, it shows some irregularities: Instead of individual, distinct waves with constant downstream front velocity, a mega-jam consists of multiple WMJs which merge into one phase region. While a typical WMJ emerges as a narrow moving jam in a pinch region (in a synchronized flow phase) and develops slowly into a WMJ, there is no extended pinch region in a mega-jam, but all moving jams emerge directly at the bottleneck position. Resultantly, a fixed downstream phase front can be observed as long as the bottleneck is active. Fig. 14 (left) displays a congestion pattern that is a mixture of a GP until 6 pm and, subsequently, a mega-jam due to a laneclosure. Fig. 15 displays the reconstruction results using the PSM and the GASM of the training set (Fig. 14 right) with the same parametrization as used before. The first observation is that the downstream front of the mega-jam is more accurate using the PSM. This can be explained with two effects: (1) specifically for cases like these, the identification of the WMJ phase involves two smoothing kernels which prevent a spill-over in downstream direction (see ection 3.2.1.2). (2) Close to the downstream bottleneck there is a short transition region where vehicles start to accelerate. Since velocities exceed the WMJ speed this region is reconstructed as a synchronized flow phase. Data is smoothed in time-direction and the resulting velocity estimate ‘straightens’ the jagged estimate of the WMJ phases. However, the estimate of the front still lacks accuracy, which is noticeable at around 8 pm when part of the jam is under-estimated. Another observation is that both algorithms fail to provide accurate estimations of the upstream front. Note that the GASM always over-estimates the downstream congestion front in this case. The reason is the following: The original formulation of the GASM (Treiber and Helbing, 2003) smooths velocities which results in a jagged downstream front that varies around the real

Fig. 14. Raw trajectory data (left) and a subset of trajectories for reconstruction purposes (right) of a mega-jam congestion pattern on A99 on 21-July 2014 on A99 in westbound direction.

658

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 15. Reconstruction using the PSM (left) and the GASM (right).

congestion front. In this case, an extended version of the GASM is implemented (van Lint, 2010), which, in accordance with a harmonic mean, smooths the inverted velocities in order to improve travel-time accuracy. The consequence of this variation is that congested regions with very low velocities are overestimated in its size. Generally, the PSM and the GASM lack accuracy reconstructing this type of pattern, or respectively, more raw data is required in order to obtain accurate results. The reason is the substantial difference of the mega-jam or HCT pattern compared to the other congestion patterns: Other congestion patterns consist of regions with reduced vehicle speeds where perturbations cause the emergence of WMJs (or oscillating traffic), which propagate upstream. In contrast, mega-jams have queueing character in which vehicles move very slowly and the downstream front is fixed. Although the PSM involves a few measures to enhance the reconstruction accuracy compared to the GASM, its fundamental assumption that regions of low traffic speeds propagate upstream does not match the nature of this pattern. In spite of the rare occurrence of this pattern (Schoenhof and Helbing, 2007), it is quite relevant for ITS. Therefore, future work should focus on the development of extensions of the PSM (or GASM) that allow to reconstruct megajams more accurately. 5.4. Runtime analysis Besides accuracy, efficiency is an important aspect of a traffic state estimation method. Especially since real-time applications require an algorithm that processes sensor information of a large network as quickly as possible in order to broadcast up-to-date traffic information. But also for processing greater amounts of historic sensor data, efficient estimation algorithms are necessary. In order to provide insights into the efficiency of the PSM, this section presents its computational complexity and average runtimes required for scenarios with respect to varying network sizes. When applied to a space-time domain discretized into A intervals in time and B intervals in space ( A,B ∈  ), an implementation of the PSM consists of element-wise matrix operations and 2D convolution processes. With respect to the runtime, it is irrelevant whether grid cells contain data or not. Element-wise matrix operations have a complexity of O (AB ) . As shown in detail in Schreiter et al. (2010), the convolution process can be implemented efficiently using the Fast Fourier Transform (FFT) (Cooley and Tukey, 1965). The resulting complexity is O (AB logAB ) . For complexity considerations, the term with the highest order dominates. With C = AB denoting the total number of grid cells in the domain, the overall complexity of the PSM is:

O (C log(C )).

(32)

This class of complexity is close to linear complexity. Furthermore, when it is applied to increasing problem sizes (number of grid cells), the required number of computations increases slightly stronger than linearly. If hardware-specific effects such as memorymanagement are neglected, then also processing times are expected to grow approximately linearly. The computational complexity is important in order to estimate whether an algorithm is eligible to apply to huge problems, but does not give much further insight in the real computation times. Therefore, in the following, the computation times for a real-time scenario with increasing domain size are presented. We assume a domain with ΔT = 30 s and A = 60 such that data within a time interval of thirty minutes are available for traffic state estimation. ΔL is set to 50 m, and B , denoting the number of cells in space direction, is varied. We assume that the PSM is called iteratively (every time step) in the real-time case, such that the smoothing kernels can be preprocessed. The time that it takes to process a matrix with velocity data VFCD into the velocity estimate VE is measured. Each measurement is repeated for 50 times to ensure robustness of the results. The computations are run on a notebook with i7-4800 processor with 2.7 GHz and 8 GB memory. Fig. 16 illustrates the processing time with respect to an increasing problem size. As expected from the complexity considerations, runtime increases approximately linearly, with slight deviations at low runtimes. With a cell count of approximately 6·106 (5000 km) a significant jump of computation times occurs. The reason for the occurrence is because the required memory exceeds the available memory of the notebook. Despite considering the average runtime of 50 iterations, the curve is still noisy. It is caused by the FFT 659

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 16. Runtime of the PSM in a real-time scenario with respect to an increasing problem size.

implementation that internally expands matrices to 2 A and 2B cells, such that jumps in the computation time can be observed. All in all, we can already conclude that our un-optimized implementation of the PSM allows processing traffic data in network sizes of up to 1000 km in about 3 s. Further optimizations of the code as well as distributed computing are expected to handle far larger network sizes in real-time. 5.5. Sensitivity analysis The PSM involves several parameters (see Table 2) that need to be set. In order to understand the influence of certain parameters on the estimation accuracy better, in this section the sensitivity of the estimation result with respect to different parameters is presented. This is supposed to facilitate the adoption and parametrization of the PSM to other scenarios. Some of the parameters, such as the smoothing directions vpdir and vpdir ,H as well as the velocity thresholds between the phases are motivated by empirical traffic characteristics and reasonable values have been discussed in ection 4.4. However, the kernel parameters τp , τ pH , σp and σpH are not yet clear. The parameters time headway T0 and minimal distance x 0 have been subject of several studies publishing distributions of observed values in real traffic. Though, the microscopic car-following model we assume in order to determine Ψ(t ,x ) is a simplification. If it turns out that parameters T0 and x 0 influence the estimation result significantly, more sophisticated models might be necessary. A method that is commonly used to quantify the sensitivity of an input parameter on the cost function is the Variance-based Sensitivity Analysis (VSA) (Saltelli et al., 2007). Consider a function Y (X ) with X ∈ d as input vector. Then, the first-order sensitivity Si measures the effect of varying one input dimension Xi on the cost function Y :

Si =

Var (E (Y |Xi )) Var (Y )

(33)

where Var (·) denotes the variance and E (·) the expectation of a random variable. Note that ∑i Si < 1 in case there are interaction terms (Saltelli et al., 2007). In order to reduce the number of samples here the Fourier-Amplitude Sensitivity Test (FAST) is applied for the mentioned 10 parameters (Cukier et al., 1978; Saltelli et al., 1999). The sampling interval for time headway T0 is varied in the range [0.5 s,3.0 s], for x 0 in [2 m,20 m], for all τ values in [20 s,500 s] and for all σ in [100 m,1000 m]. Fig. 17 depicts the resulting first order sensitivities with 20,000 model evaluations of the scenario described in ection 5.1. The estimation accuracy appears to be less sensitive to variations of the parameters of the assumed car-following model. This supports our hypothesis that a simplified model is sufficient for the purpose of the PSM. Also some other parameters have very low first order sensitivities. Thus, varying these in the proposed realistic intervals does not impact estimation accuracy significantly. The highest impact on the estimation accuracy have parameters τS,F and τJH,S . When applying the PSM to other scenarios, an optimization of these two parameters is likely to produce the greatest gain in accuracy. 6. Discussion The presented results show the strengths of the PSM as an accurate and efficient algorithm for traffic speed estimation with sparse velocity data. In the qualitative comparison, it achieves to reconstruct different empirical traffic patterns, such as moving jams and stationary congested traffic as described in different traffic theories (Kerner, 2009; Schoenhof and Helbing, 2007) The qualitative 660

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Fig. 17. Global sensitivities of the parameters involved in the PSM.

results are supported by the accurate results obtained from a quantitative analysis where the PSM is applied to varying data densities. In these studies it outperforms a naive method and the GASM. An intuitive explanation for the better results is that the distinction between three traffic phases means incorporating more empirical traffic knowledge into the reconstruction method than it is currently done with the GASM. Effectively, the higher complexity of the algorithm enables higher estimation accuracies. 6.1. Critical review For now, the presented qualitative and quantitative evaluation are based on one congestion pattern. Since many congestion patterns on freeways consist of synchronized flow phases with fixed downstream fronts and occasionally occurring WMJs, it is expected that the PSM also provides good results in many other scenarios. However, a more comprehensive study is still required. Varying parametrization is another consideration which still has to be taken into account. Neither the PSM, nor the GASM or the naive method have a formally optimized set of parameters. For the GASM, Rempe et al. (2016) shows that parameters of the smoothing kernel have a significant impact on the estimation accuracy. Furthermore, for different input data densities, different parameter sets achieve the most accurate results. The same characteristic is expected from the PSM and the naive method with respect to their parametrizations. Therefore, in order to obtain quantitatively secured values describing the full potential of each method, for each scenario and data density an optimization would be necessary. A further point of discussion concerns the evaluation methodology where FCD is split into training and test sets and where they are used for estimation and evaluation respectively. One may argue that the data is biased, since vehicles do not send their positions continuously and the fleet of equipped vehicles are mainly private cars but no trucks. Consequently, the estimation and evaluation would not be an accurate representative for the macroscopic velocity. In fact, this is what actually happens in free flow phase. Because of the definition of the free flow state, the individual vehicles’ velocities differ significantly, especially among private cars and trucks. Therefore, there is always a high uncertainty of the macroscopic traffic velocity if not all vehicle velocities are reported. Additionally, the described filter mechanism often retains data, in case of free flow, with the result that less data is available. Effectively, which probably holds for all velocity estimators based on probe data the resulting free flow velocity has a lowered trustworthiness. This contrasts with the congested phases. By definition, due to high traffic density the velocities among all vehicles are synchronized in these phases. Velocities of individual vehicles can be seen as good approximations of the macroscopic speed. Consequently, the fact that only private vehicles report data does not impact the estimate of the macroscopic velocity in congested regimes. One can conclude that velocity estimation with probe data benefits mainly in congested states. Fortunately, that is the traffic condition which is mostly interesting for traffic information services and traffic control. The choice of the IMAE as error estimator, penalizing errors in free flow much less than errors corresponding to congested phases, emphasizes these considerations. Finally, the PSM exploits many properties of the three phases, but also neglects some. For instance, the downstream front of a synchronized flow phase is expected to maintain at the bottleneck. That does not fully agree with Kerner’s theory, since there are socalled Moving Synchronized Pattern (Kerner, 2009). Effectively, the reconstruction of these patterns is less accurate. Though, these patterns usually do not persist long. They often evolve into WMJs or attach to an upstream adjacent bottleneck. 6.2. A Comment on the relation to two-phase models There is an ongoing discussion whether two-phase or three-phase traffic models constitute more valid models of real traffic flow. The presented method assumes the validity of the Three-Phase traffic theory and builds a procedure that integrates empirically determined characteristics of the three phases. However, it is also possible to motivate the present approach assuming the validity of a modern two phase model: At an active bottleneck the downstream congestion front is usually pinned to the bottleneck position. Traffic is in the convectively unstable regime (Helbing et al., 1999) which is a congested state where traffic speeds are significantly higher than zero. In this state, perturbations originate that may grow into oscillating traffic further upstream of the bottleneck. The features of this state match with the properties of the smoothing kernel that is applied in order to identify the synchronized flow phase, i.e. a propagation velocity vpdir = 0 km/h that models the locality of the congestion, and vpdir ,H ≈ −18 km/h that models the emergence of perturbations. To conclude, the promising results of this method must not be interpreted as a proof that the Three-Phase 661

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

theory is superior over two-phase models. Rather, the comprehensive description of properties of congestion patterns and phases as provided by the Three-Phase theory forms a strong basis for the presented method. 7. Conclusion and outlook In this paper, a novel freeway traffic estimation algorithm called the Phase-based Smoothing Method (PSM) was presented. The PSM estimates macroscopic traffic speeds in time-space domain from sparse probe data. The approach is to identify the three traffic phases Wide Moving Jams, synchronized flow and free flow using a probabilistic model consisting of a set of fuzzy criteria. A subsequent phase-dependent smoothing of data allows incorporating several empirical characteristics of traffic into the estimation algorithm. In a qualitative and quantitative study using real FC data collected during a congestion on a German freeway, the PSM was compared to the Generalized Adaptive Smoothing Method and isotropic smoothing. Results show that due to its distinctive treatment of moving and stationary congested patterns the method achieves a more accurate reconstruction of real world congestion patterns. Besides the gain in accuracy, the PSM provides some features relevant for practical implementation. First, it computes a quality estimate alongside the velocity estimate. That generated quality is valuable information in order to quantify the trustworthiness of a given velocity estimate. Second, since it only requires gridded velocity information, it can be applied to basically any traffic velocity data. Similar to Treiber et al. (2010), a mechanism to fuse different data sources with different qualities is possible. Third, the exclusive use of smoothing operations make the approach efficient, which enables a large-scale application in real-time. Fourth, a parameter study reveals that the variation of only two parameters controls most of the estimation accuracy of the method, which facilitates adoption of the method to new scenarios. For now, the PSM is developed to be applied retrospectively where all reported data is available at once. In an online application where measurements are received with a temporal delay, slight modifications are required. It needs to be considered that current traffic speeds are predictions rather than reconstructions. Adaptation and evaluation of the method for online use remain for future work. References Aw, A., Rascle, M., 2000. Resurrection of second order models of traffic flow. SIAM J. Appl. Math. 60 (3), 916–938. Bekiaris-Liberis, N., Roncoli, C., Papageorgiou, M., 2016. Highway traffic state estimation with mixed connected and conventional vehicles. IEEE Trans. Intell. Transp. Syst. 1–14. Brackstone, M., Sultan, B., McDonald, M., 2002. Motorway driver behavior: studies on car following. Transp. Res. Part F: Traffic Psychol. Behav. 5 (1), 31–46. Cooley, J.W., Tukey, J.W., 1965. An algorithm for the machine calculation of complex Fourier series. Math. Comput. 19 (90), 297. Cukier, R.I., Levine, H.B., Shuler, K.E., 1978. Nonlinear sensitivity analysis of multiparameter model systems. J. Comput. Phys. 16, 1–42. Daganzo, C.F., 1996. The nature of freeway gridlock and how to prevent it. In: Proceedings of the 13th International Symposium on Transportation and Traffic Theory, pp. 629–646. Duret, A., Yuan, Y., 2017. Traffic state estimation based on Eulerian and Lagrangian observations in a mesoscopic modeling framework. Transp. Res. Part B: Methodol. 101, 51–71. Edie, L.C., 1963. Discussion of Traffic Stream Measurements and Definitions. Port of New York Authority. Hall, D., Agyemang-Duah, K., 1991. Freeway capacity drop and the definition of capacity. Transp. Res. Rec.: J. Transp. Res. Board 20–28. Helbing, D., Hennecke, A., Treiber, M., 1999. Phase diagram of traffic states in the presence of inhomogeneities. Phys. Rev. Lett. 82 (21), 4360. Helbing, D., Treiber, M., Kesting, A., Schoenhof, M., 2009. Theoretical vs. empirical classification and prediction of congested traffic states. Eur. Phys. J. 2009 (B 69), 583–598. Herrera, J.C., Bayen, A.M., 2010. Incorporation of Lagrangian measurements in freeway traffic state estimation. Transp. Res. Part B: Methodol. 44 (4), 460–481. Huber, G., Bogenberger, K., Bertini, R., 2014. New methods for quality assessment of real time traffic information. Transp. Res. Board 93. Jacobs, R.A., 1995. Methods for combining experts’ probability assessments. Neural Comput. 7 (5), 867–888. Jiang, Z., Chen, X.M., Ouyang, Y., 2017. Traffic state and emission estimation for urban expressways based on heterogeneous data. Transp. Res. Part D: Transp. Environ. 53, 440–453. Kerner, B.S., 2004. The Physics of Traffic. Springer. Kerner, B.S., 2009. Introduction to Modern Traffic Flow Theory and Control. Springer, Berlin, Heidelberg. Kerner, B.S., Rehborn, H., Aleksic, M., Haug, A., 2004. Recognition and tracking of spatial–temporal congested traffic patterns on freeways. Transp. Res. Part C: Emerg. Technol. 12 (5), 369–400. Krbalek, M., Seba, P., Wagner, P., 2001. Headways in traffic flow: remarks from a physical perspective. Phys. Rev. E: Stat., Nonlin, Soft Matter Phys. 64 (6), 066119. Laval, J.A., 2007. Linking synchronized flow and kinematic waves. In: Traffic and Granular Flow’05, pp. 521–526. Lighthill, M., Whitham, G., 1955. On kinematic waves II. A theory of traffic flow on long crowded roads. Proc. Roy. Soc. A 229, 317–345. Mika, H.S., Kreer, J.B., Yuan, L.S., 1969. Dual mode behavior of freeway traffic. Highway Res. Rec. 279, 1–12. Newell, G.F., 1993. A simplified theory of kinematic waves in highway traffic, part II, queueing at freeway bottlenecks. Transp. Res. Part B: Methodol. 27 (4), 289–303. Palmer, J., Rehborn, H., Gruttadauria, I., 2011. Reconstruction quality of congested freeway traffic patterns based on Kerner’s three-phase traffic theory. Int. J. Adv. Syst. Meas. 4, 168–181. Rempe, F., Franeck, P., Fastenrath, U., Bogenberger, K., 2016. Online freeway traffic estimation with real floating car data. In: Intelligent Transportation Systems, ITSC, pp. 1838–1843. Richards, P., 1956. Shock waves on the highway. Oper. Res. 4, 42–51. Saltelli, A., Ratto, M., Andres, T., Campolongo, F., Cariboni, J., Gatelli, D., Saisana, M., Tarantola, S., 2007. Global Sensitivity Analysis. The Primer. John Wiley, Sons Ltd. Saltelli, A., Tarantola, S., Chan, K.P.-S., 1999. A quantitative model- independent method for global sensitivity analysis of model output. Technometrics 41 (1), 39. Schoenhof, M., Helbing, D., 2007. Empirical features of congested traffic states and their implications for traffic modeling. Transp. Sci. 135–166. Schreiter, T., van Lint, H., Treiber, M., Hoogendoorn, S., 2010. Two fast implementations of the adaptive smoothing method used in highway traffic state estimation. In: 13th International IEEE Conference on Intelligent Transportation Systems. ITSC, pp. 1202–1208. Seo, T., Bayen, A.M., Kusakabe, T., Asakura, Y., 2017. Traffic state estimation on highway: a comprehensive survey. Ann. Rev. Control. Shafer, G., 1976. A Mathematical Theory of Evidence. Limited Paperback Editions, vol. 42 Princeton University Press, NY. Suzuki, H., Nakatsuji, T., Nanthawichit, C., 2003. Application of probe vehicle data for real-time traffic state estimation and short-term travel time prediction on a freeway. Transp. Res. Rec.: J. Transp. Res. Board 1855, 49–59.

662

Transportation Research Part C 85 (2017) 644–663

F. Rempe et al.

Treiber, M., Helbing, D., 2003. An adaptive smoothing method for traffic state identification from incomplete information. Interface Transp. Dyn. 32, 343–360. Treiber, M., Kesting, A., 2013. Traffic Flow Dynamics. Springer, Berlin, Heidelberg. Treiber, M., Kesting, A., Wilson, R., 2010. Reconstructing the traffic state by fusion of heterogeneous data. Comput.-Aided Civ. Infrastruct. Eng. 2011 (26), 408–419. van Hinsbergen, C.P.I.J., Schreiter, T., Zuurbier, F.S., van Lint, J.W.C., van Zuylen, H.J., 2012. Localized extended kalman filter for scalable real-time traffic state estimation. IEEE Trans. Intell. Transp. Syst. 13 (1), 385–394. van Lint, H., Djukic, T., 2014. Applications of kalman filtering in traffic management and control. New Direct. Informat. Optimiz. Logist. Prod. 59–91. van Lint, H., Hoogendoorn, S., 2009. A robust and efficient method for fusing heterogeneous data from traffic sensors on freeways. Comput. Aided Civ. Infrastruct. Eng. 24, 1–17. van Lint, J., 2010. Empirical evaluation of new robust travel time estimation algorithms. Transp. Res. Rec.: J. Transp. Res. Board 2160, 50–59. Wang, Y., Papageorgiou, M., 2005. Real-time freeway traffic state estimation based on extended Kalman filter: a general approach. Transp. Res. Part B Methodol. 39 (2), 141–167. Work, D.B., Tossavainen, O.-P., Blandin, S., Bayen, A.M., Iwuchukwu, T., Tracton, K., 2008. An ensemble kalman filtering approach to highway traffic estimation using GPS enabled mobile devices. In: 47th IEEE Conference on Decision and Control, pp. 5062–5068. Work, D.B., Tossavainen, O.-P., Jacobson, Q., Bayen, A.M., 2009. Lagrangian sensing: traffic estimation with mobile devices. In: American Control Conference, pp. 1536–1543. Yuan, Y., van Lint, J.W.C., Wilson, R.E., van Wageningen-Kessels, F., Hoogendoorn, S.P., 2012. Real-time Lagrangian traffic state estimator for freeways. IEEE Trans. Intell. Transp. Syst. 13 (1), 59–70.

663