A twodimensional datadriven model for traffic flow on highways
Abstract
Based on experimental traffic data obtained from German and US highways, we propose a novel twodimensional firstorder macroscopic traffic flow model. The goal is to reproduce a detailed description of traffic dynamics for the real road geometry. In our approach both the dynamics along the road and across the lanes is continuous. The closure relations, being necessary to complete the hydrodynamics equation, are obtained by regression on fundamental diagram data. Comparison with prediction of onedimensional models shows the improvement in performance of the novel model.
Msc
90B20; 35L65; 35Q91; 91B74
Keywords
Traffic flow, macroscopic model, two dimensional model, trajectory data, validation, data fitting
1 Introduction
The mathematical modeling of vehicular traffic flow uses different descriptions and we refer to [11, 30, 57] for some review papers. Besides microscopic and cellular models there has been intense research in continuum models where the temporal and spatial evolution of car densities is prescribed. Based on the level of detail there are gaskinetic or mesoscopic models (e.g., [32, 34, 39, 46, 56, 31]) and macroscopic models being fluiddynamics models (e.g., [8, 10, 13, 21, 22, 25, 26, 42, 43, 47, 50, 54, 55, 60, 66, 68]). Among the (inviscid) macroscopic models one typically distinguishes between firstorder models based on scalar hyperbolic equations and secondorder models comprised of systems of hyperbolic equations. The pioneering work of the first case is the Lighthill and Whitham [50] and Richards [60] model (LWR). While a specific example of the second case is the Aw and Rascle [8] and Zhang [68] model (ARZ). Depending on the detailed level of description of the underlying process different models have been employed and tested against data. In recent publications it has been argued that the macroscopic models provide a suitable framework for the incorporation of online traffic data and in particular of fundamental diagram data [4, 9, 23]. While microscopic models are nowadays widely used in traffic engineering, continuum models have been studied mathematically, but very little work has been conducted on their validation with traffic data [4, 6, 14].
So far, most of the proposed continuum models are for single lane vehicular traffic dynamics. However, the data for fundamental diagrams is taken from interstate roads with multiple lanes [2, 3] and can be used for deriving or testing models for real road geometry. Multilane models belong to this class. Typical modeling of multilane traffic uses a spatially onedimensional model (1D) of either LWR or ARZ type for each lane. The lanechanging of cars is then modeled by interaction terms (the sources on the righthand sides of the equations) using empirical interaction rates, see e.g. [44, 45, 46]. The interaction modeling is typically assumed to be proportional to local density on current and desired lane. A fluid dynamics model describing the cumulative density on all lanes is proposed in [16, 17, 64], where a twodimensional (2D) system of balance laws is obtained by analogy with the quasigasdynamics (QGD) theory. Here, the authors model 2D dynamics assuming that vehicles move to lanes with a faster speed or a lower density and the evolution equation for the lateral velocity is described by the sum of the three terms proportional local density and mean speed along the road.
A major problem of the approaches described above is to estimate from data the interaction rate or the great number of coefficients and parameters. Therefore, here, we propose a different approach: we treat also lanes as continuum and postulate a dynamics orthogonal to the driving direction. The precise form of the dynamics is established through comparisons with fundamental diagrams obtained from trajectory data recorded on a road section of the A3 German highway near Aschaffenburg. Thus, the experimental measurements allow us to derive a model being able to take into account the realistic dynamics on the real road geometry without prescribing heuristically the behavior of the flow of vehicles.
The contribution and the organization of this paper is summarized below.

Derivation and the presentation of historic fundamental diagrams data for the dynamics of traffic across the lanes (see Section 2). In fact, the German dataset provides the twodimensional timedependent positions of vehicles while crossing the road section. Therefore, in addition to the classical fundamental diagrams widely studied in the literature [5, 41, 48] and used for deriving onedimensional datafitted macroscopic models [23, 24], we can also generate diagrams for the dynamics across the lanes. Although twodimensional experimental traffic measurements are already available in the literature, this is, to our knowledge, the first time that they are used to study the dynamics orthogonal to the movement of vehicles;

Design of a new datafitted twodimensional firstorder model and the analysis of its mathematical properties (see Section 3). The historic data are therefore used to develop the novel macroscopic model defining the flux functions by means of a datafitting approach. The closures are necessary to complete the macroscopic equation and taking them using the experimental data allows to describe the real dynamics of the flow;

validation of the novel 2D macroscopic model via timedependent trajectory data and the definition of a systematic methodology to study and to compare the predictive accuracy with respect to the 1D LWR model (see Section 4).
2 Dataset description and fundamental diagrams
We use a set of experimental data recorded on a German highway. Precisely, we have twodimensional trajectory data collected on a meter stretch of the westbound direction of A3 highway near Aschaffenburg. Laser scanners detect the twodimensional positions of each vehicle at time on the road segment with a temporal resolution of seconds for a total time of approximately minutes. Here the position is in driving direction, the position is across lanes. During the time observation, the laser scanners record the trajectories of vehicles.
The road section consists of three lanes and an outgoing ramp. However, we only consider the stretch as if there is no ramp. In fact, the data show that the flow on the ramp does not influence the traffic conditions, namely the amount of traffic on the ramp is not significant. Taking into account only the three main lanes, the road width is meter. The stretch we are considering is not straight but it has a turn with a small radius of curvature. We have taken into account this feature when cleaning the experimental data since the curvature could mainly falsify the data in direction. Precisely, the cleaning procedure is based on the knowledge of the curvature of the roadsection. We clean the timedependent positions of the vehicles on the road by adjusting them by a factor which allows to not consider movement in if it is only due to the curvature of the road section.
As pointed out in the Introduction, in this paper we are interested in the study of macroscopic traffic models. In other words, instead of looking at the motion of each single vehicle, we wish to “zoom out” to a more aggregate level by treating traffic as a fluid. Therefore, with the aim of proposing a novel datafitted 2D macroscopic model, from the microscopic experimental data we need to recover the macroscopic quantities, namely the density (measured as number of vehicles per kilometer), the flux (measured as number of vehicles per hour) and the mean speed (measured as kilometer per hour) of the flow.
To this end, firstly, we observe that the microscopic positions of vehicles at each time are sufficient to recover the microscopic velocities of vehicles. In particular, since the road section is relatively short, we compute the velocity, both in  and direction, of each vehicle by using a linear approximation in the least squares sense of its positions, and respectively, on the road during the time interval. In other words we assume that the vehicle velocity is constant during the crossing of the road section and is exactly the slope of the linear fit. Thus, we associate at each vehicle the vector of the microscopic velocities . The maximum detected speed in direction is about kilometer per hour which means about seconds to travel the meters of the road section.
The timedependent microscopic positions and the microscopic velocities of vehicles are used to compute the macroscopic data as we describe in the following. Clearly, since we are aimed to develop a twodimensional datafitted firstorder macroscopic model, the derivation of macroscopic quantities, such as the flux and the mean speed, should be done for each direction, along the road (direction) and across the lanes (direction), separately.
The macroscopic density gives information on the congestion level of the road section. It is usually expressed in number of vehicles per unit length (here kilometers) and therefore it ignores the concept of traffic composition. This is not restrictive for our purpose of deriving a twodimensional firstorder macroscopic model for traffic. The modeling of the heterogeneous composition of vehicles is studied in multipopulation models, e.g., in [12] at the macroscopic level and in [59, 58] at the kinetic level, where the concept of density is replaced by the rate of occupancy. In order to compute the macroscopic density, we first fix a sequence of equally spaced discrete times such that , and , where is the final observation time in the dataset (here 20 minutes). Then we count the number of vehicles on the road at each discrete time defining
where is the length of the road section expressed in the unit length. Finally, we consider a moving mean by aggregating with respect a certain time period , with and including consecutive observations. This temporal average leads to values of the density
In particular, in this paper we take second and then we aggregate the data over the time period seconds.
Observe that, clearly, the density does not depend on the direction we are looking at. The computation of the flux and of the mean speed is, instead, a little bit more complex.
In our approach, we first compute the mean speeds of the flow. Consider the same sequence of equally spaced discrete times introduced above. Then we average the microscopic velocities and of all vehicles being on the road at a fixed time with respect to the number of vehicles , so that we define
Using and we then compute the fluxes at each discrete time by means of the hydrodynamics relation
and, as done for the density, we consider a temporal average by aggregating with respect the same time period , leading to the following expressions of the macroscopic fluxes
for . Finally, using again the hydrodynamics relation, we get the mean speeds of the flow as
For a more detailed discussion on the computation of macroscopic quantities from microscopic data, we refer to [38, 51].
The diagrams showing the relations between the vehicle density and the fluxes , or the mean speeds , are called fundamental diagrams and speeddensity diagrams, respectively. They represent the basic tools for the analysis of traffic problems operating in a homogeneous steady state or equilibrium conditions.
In Figure 1 we show the diagrams resulting from the German dataset: the top row shows the relations and , while the bottom row shows the relations and . Observe that the dataset provides during the time period several levels of congestion but we never observe bumpertobumper conditions. In fact, the maximum density is about vehicles per kilometer.
In addition to the classical fundamental relations and we also show the diagrams and . We highlight that datasets providing 2D trajectories are already available, e.g., see [2, 3]. But the attempt of taking into account the study of the dynamics across the lanes, and thus considering also the data in direction, is, to our knowledge, a novelty in the mathematical literature on traffic flow.
The qualitative structure of such diagrams is defined by the properties of different regimes, or phases, of traffic. For a description of the diagrams in direction we refer, e.g., to [5, 41, 48]. Clearly, the diagrams in direction show a different quantitative and qualitative behavior with respect to the classical ones. Firstly, we observe that the values of the flux and of the mean speed are about smaller than the values in direction. This is obvious since the velocity of vehicles along the road is higher then the lateral velocity and thus is not dominant with respect . In other words, we are looking at two behaviors occurring at different scales. But this does not mean that the behavior in direction can be neglected. In fact, an analysis of the trajectories shows that about the of the total vehicles crosses a lane while traveling the road section. See Figure 2. In [36] we also show, using the same German dataset described above, that the level of potential conflicts is strongly affected by the lane changing. To this end we computed the average risk on the two directions of the flow by means of the TimeToCollision metric and we observe that the risk due to lane changing is higher.
Moreover, notice that and have positive and negative values since across the lanes vehicles are free to travel in the two directions, towards right and left. Precisely, we assume that positive speeds represent the motion towards the leftmost lane, instead negative speeds represent the motion towards the rightmost lane.
Remark 1.
As pointedout above, the fundamental diagrams in Figure 1 are obtained by averaging over a time period of seconds macroscopic data computed each second. If we take more frequent time observations (e.g. every seconds) the structure of the diagrams does not change. We only observe that the data become thicker. Compare the topleft panel in Figure 1 with the left panel in Figure 3. This behavior is due to the fact that we consider a linear approximation of the vehicle trajectories, i.e. we assume that the vehicle velocities are constant. Thus, considering time observations every second does not keep out information.
In contrast, the time for the data aggregation slightly influence the structure of the fundamental diagrams. Compare the topleft panel in Figure 1 with the right panel in Figure 3. In the following we will consider the diagrams obtained with an aggregation time period of seconds which is widely used in the engineering literature on traffic flow.
3 Twodimensional LWRtype model
One dimensional firstorder macroscopic traffic models are based on the continuity equation
(1) 
which gives the conservation of vehicles on the road segment . In (1), the vehicle density is , and the vehicle velocity field is , where is the position along the road, and is time.
The simplest macroscopic traffic model, the LWR [50, 60] model, is obtained by assuming a functional relationship between and , i.e., . This turns equation (1) into a scalar hyperbolic conservation law
(2) 
where the flux is given by the flow rate function . Because the LWR model (2) is a closed model consisting of a single equation, it is denoted a firstorder model. The velocity function is commonly assumed to be decreasing in with for some maximal vehicle density . Here, is assumed to be the density in bumpertobumper conditions and its value is given in Section 3.1 below.
The strict functional relationship between and is called closure law and is loosened in the socalled secondorder models, which augment (1) by an evolution equation for the velocity field, see [8, 23].
The onedimensional model (2) describes the flow of vehicles in the simple case of a singlelane road or, if the road has multiple lanes (in a given direction), it considers these aggregated into the scalar field quantities and . Nevertheless, the dynamics of traffic on a multilane highway could be more complex and is strongly influenced by the motion of vehicles across the lanes. For this reason, we take into account the intrinsic multidimensional characteristic of traffic flow by extending model (2) to the twodimensional firstorder macroscopic model
(3) 
where and are the fluxes in the two possible directions of the flow and , are the speed along the road and the lateral speed, respectively. The quantities and are the length and the width of the road, respectively. Clearly, one expects that . As in the onedimensional model, we have the following two closures
The velocity function is the same speed introduced in the onedimensional LWR model (2) and thus it is obvious to assume that it has the same properties discussed previously. A heuristic description of and as function of the density is not immediate since it depends strongly on the preference on the drivers as well as general imposed traffic rules. It is natural to assume that the lateral speed is for and . In fact, in the first case vehicles would travel towards the rightmost lane since the road is free (according to the traffic rules), while in the second case they cannot change lanes and thus travel in direction since the lanes are almost congested. Note that contrary to the direction the speed in direction can be negative. In the following section we propose a functional relation obtained from data.
We finally stress the fact that the twodimensional model (3) is able to take into account the dynamics of traffic on a multilane highway but actually it is not a multilane model. In fact, notice that equation (3) is continuous in . Instead a multilane model requires to treat lanes as discrete object and, thus, to develop a system of balance laws in which the source terms describe the mass exchange between the lanes.
3.1 Macroscopic closures and datafitting
For the dynamics along the road, namely in direction, several laws have been considered in the literature: popular examples of flow rate functions are the Greenshields’ flux [28], in which is a quadratic function, and the NewellDaganzo flux [20, 52], in which is a piecewise linear function. These different choices of functions lead to wellposed firstorder models. Many closure laws were proposed in the literature, for further discussions we refer, e.g., to the book [62].
A natural way to derive closure laws is to construct a fitting of the experimental data. Although this approach ignores the scattered behavior of data, we expect to characterize key properties of the traffic flow (as the critical density, the maximum flow, ). For comparisons between models using classical closure laws and datafitted models, see [23, 24] based on the NGSIM dataset [2] and on the RTMC dataset [3].
We are considering the German dataset described in Section 2 and we are proposing a twodimensional firstorder macroscopic model. Then, in order to get the closure laws to complete equation (3) we proceed by constructing the best fitting via a least squares fit to the data computed in Section 2 and showed in Figure 1. The closures for the firstorder macroscopic model (3) must represent these data via singlevalued functions and .
Since the stagnation density is not represented well via data, we prescribe it as a fixed constant, given by the ratio between the number of lanes and the typical vehicle length of meters, plus of additional safety distance, so that
where this value is obtained by considering the unit distance of 3000 meters for the three lanes (1 kilometer per lane).
As visible in the fluxdensity diagram in the top left panel of Figure 1, the data tend to exhibit a relatively linear increasing relationship between and for low densities. In turn, for higher densities, a significant spread is visible, i.e., a single value corresponds to many different flow rates . For the data in direction then we employ the same approach presented in [23] and in [24] by selecting a threeparameter family of smooth and strictly concave flow rate curves as
(4) 
where
Each flow rate function in this family vanishes for and . The three free parameters allow for controlling three important features of the fluxdensity diagram in direction: the value of maximum flow rate (mainly determined by ), the critical density (mainly controlled by ), and the curvature (dominated by ).
In case of the direction a significant spread of the data is visible also at lower densities and they show both positive and negative values since the motion is allowed in two directions: towards the leftmost and the rightmost lane. Moreover, we observe only for , proving the fact that vehicles tend to travel towards the rightmost lane in the freeflow regime. To take into account these features we have to propose a different flow rate curve with respect to (4). Precisely, in this case we choose a simple twoparameter family of smooth functions as follows
(5) 
Each flow rate function in this family vanishes for . The two free parameters allow for controlling the speed in the freeflow regime (determined by ) and the shape of the of the curve due to the data (mainly controlled by ).
Remark 2.
Clearly, a more complex flow rate curve (5) may be postulated. The simple choice (5) is justified by the behavior of the diagrams provided by the A3 German highway which do not give information on how the data behave for higher density values. In fact, the only realistic apriori assumption for the congested regime is that .
From the three and the twoparameter family of flow rate curves (equations (4) and (5), respectively), the closures and are selected in such a way they are the closest, in a leastsquares sense, to the experimental data points , and , respectively. Thus we solve
(6) 
The minimization problems (6) are solved numerically by using the Matlab solver fmincon which finds the minimum of constrained nonlinear functions. For the German dataset the solver provides the following values for the free parameters

, and ;

and .
For the parameters in directions we do not prescribe a range, i.e. . While we require that and . The constrained minimization problems are quickly solved by the Matlab fmincon function. The CPU time needed is about seconds in each direction. The relative errors obtained with the above optimal parameters are
In Figure 4, the red curves represent the leastsquares fits to the given data points computed in Section 2. These functions are used to close the twodimensional firstorder macroscopic model (3) and to validate the model in the next section.
The mathematical properties of the proposed flux in are similar to classical LWR type models. Therefore, a detailed discussion is skipped. For the numerical scheme below we remark that the conservation law (3) with choice (5) is strictly hyperbolic. Moreover, the optimal parameters and lead to a single inflection point function of the flux function and therefore conservation law (3) still give rise to only simple waves (either shock or rarefaction waves) also in direction.
4 Numerical simulations
In this section we study the predictive accuracy of the 2D LWRtype model (3) with respect to measurement data. In particular, we show that model (3) is more accurate than its 1D version (2) in which we choose as closure the flow rate function (4).
To this end, we firstly present the scheme used to numerically solve model (3). Then, we specify how continuous field quantities are constructed from the trajectory data, following the approach described in [23] and [24] for 1D datafitted models.
4.1 Numerical scheme
In the following, we simply describe the numerical scheme for solving the twodimensional model (3). For the onedimensional one (2), the same numerical scheme is used, clearly neglecting the computation of transport term in .
In order to approximate the solution of (3), we use the dimensional splitting method or method of fractional steps, [49, 65]. We split (3) into
(7a)  
(7b) 
and for each problem we apply a finite volume approximation. To this end we divide the spatial domain into cells , , , such that and , . Thus .
We consider a semidiscrete finite volume scheme and denote by
the cell average of the exact solution in the cell at time and its numerical approximation. By integrating each equation (7) over , dividing by , using the midpoint rule and finally a stage Strong Stability Preserving RungeKutta method (SSPRK) with Butcher’s tableau and time step , we get the fully discrete scheme
(8a)  
(8b) 
giving the approximation of the solutions at time , where is the initial time. Notice that in the sweeps we start with the cell averages at time and solve onedimensional problems with fixed updating to . In the sweeps we then use the values as data for solving the problems with fixed, which results in . Here,
for , are the numerical fluxes approximating and , respectively. We consider and as local LaxFriedrichs fluxes. We could also use the Godunov scheme which is less diffusive. However, we choose LaxFriedrichs for its ease. Instead, and are the reconstructions at the cell interfaces, at left and right sides, from the stage values and , respectively. The stage values of the cell averages are evolved by
where the th stage value is assumed to be at time .
Without using additional tools, the scheme described above is firstorder accurate. In order to get a secondorder scheme the following ingredients are necessary. The reconstruction at the interfaces from the stage values is performed using a piecewise linear reconstruction in each direction. To guarantee the nonoscillatory nature of the reconstruction, we apply a nonlinear limiter for the computation of the slopes and here we use the minmod slope limiter. Thus, e.g.,
where
and the minmod function is defined as
The dimensional splitting (8a)(8b) is only firstorder accurate. See [27]. For a secondorder scheme the Strang splitting technique [1] has to be employed. This method consists in a slight different application of the equations (8a)(8b). More precisely, equation (8a) is used to obtain the update up to time , i.e. with time step . This datum is then used in equation (8b). Finally, the datum resulting from (8b) is used to compute the approximation of the solution at time by means of equation (8a) starting from the time level , thus with time step . For further details we refer to [49, 1].
A timestepping of (at least) secondorder is mandatory for all subproblems described in the Strang splitting. Here, as RungeKutta scheme we take the Heun’s method [37] whose coefficients , and are defined in the following Butcher tableau
The time step is chosen in such a way it satisfies the CFL condition [18]. In particular, if not explicitly specified, in the following we will consider as CFL and the time step is
where the maximum of the derivative of the flux functions is computed on the density profile at initial time.
For the numerical solution of the onedimensional LWR model (2) we consider the natural onedimensional version of the secondorder finite volume scheme presented above.
The scheme described above is a secondorder scheme and for our purposes is sufficient. We choose to employ a dimensional splitting technique since it is conceptually easy to understand, allowing to take advantage of using classical onedimensional methods for conservation laws. There are also methods for multidimensional conservation laws that are intrinsically multidimensional, see e.g. [7]. These methods should be used to get more accurate numerical schemes and in this case highorder spatial reconstructions [19, 63] combined with highorder RungeKutta schemes have to be considered.
Although the scheme presented here is already wellknown, for the sake of completeness in Figure 5 we show the order of convergence for the case of the twodimensional linear scalar conservation law
with a Gaussian initial datum
up to time (left panel of Figure 5) and with a double sine wave initial datum
up to time (right panel of Figure 5). In both cases we use periodic boundary conditions. It is wellknown that the scheme still remains secondorder accurate also for nonlinear conservation laws on smooth solutions.
4.2 Treatment of experimental data
The numerical implementation of the macroscopic traffic models (2) and (3), require the knowledge of continuously in space initial data. Since we are aimed to compare the predictive accuracy of the two models against the dataset described in Section 2, here we specify how continuous field quantities can be constructed from trajectory data. In particular, we follow the same approach used in [23] and [24]. The same idea is applied for computing the reference data in order to compare the model predictions.
From the German dataset we have the twodimensional trajectories of vehicles , with a temporal resolution 0.2 seconds, that is essentially continuous in time. However, at each time, the vehicle positions are discrete. In order to incorporate this data as initial condition into the continuous models (2) and (3), we must generate a function , for , that is defined everywhere on the road segment. This approach also allows to compare the model accuracy against the experimental data, constructing at a certain final time the continuous function from the discrete positions of vehicles with .
The construction of density functions from discrete samples is a statistic problem. We employ a kernel density estimation (KDE) approach, with a fixed Gaussian kernel. The specific KDE approach used here is described in [24] for the case of traffic models and it is called the ParzenRosenblatt window method [53, 61]: Assume that at time we have the positions of vehicles on the road. This data are interpreted as a finite sample of some (unknown) density function. The goal is to reconstruct a kernel density estimator from the discrete information that is close to. At each time instant , KDE starts with a twodimensional comb function
where is the twodimensional Dirac delta function and the number of vehicles on the road section at time . Thus the function accounts for the positions of vehicles on the road at time . Clearly, cannot be used as initial condition of numerical simulations but we need to define its smoothed version. To this end, we consider a twodimensional Gaussian kernel
(9) 
and we define the density function at time as
(10) 
Here and are the bandwidths in  and direction respectively. There are several works dealing with the optimal choice of the bandwidth in the KDE approach, e.g., see [15, 40]. Here we take the same values already used in [23, 24] in which the bandwidths are chosen in such a way equally distant vehicles lead to an almost constant density profile. It is clear that this technique for choosing the bandwidths does not depend on the type of data but only on the road section dimensions. Therefore, we recompute the values given in [23, 24] for the case of a meter road and then we get the kernel widths meter and meter.
Finally, we note that the road section is meters and therefore vehicles travel from the initial point to the exiting one in about seconds, at the maximum speed.
4.3 Validation of the model
In the following, we validate the presented twodimensional firstorder macroscopic model (3) by comparing the evolution of the model with the corresponding measured trajectories. Also, we compare the predictive accuracy of the model with respect to its onedimensional version (2).
The deviation between predicted and real traffic states quantifies the model error. Thus we choose the spatial discretization sufficiently fine, namely meter.
In order to quantify the deviation of the model predictions from the real data, we proceed as follows. Firstly, we compute the continuous density that defines the starting condition at a fixed initial time as in (10). Then, we numerically evolve the density profile up to a final time using the numerical scheme defined in Section 4.1 and applied to the model (3). The numerical simulation gives the data output . The continuous reference solution at time is constructed from the real data by means of the density estimation defined in (10). From this function we obtain discrete values , and we finally compute the prediction error as
(11) 
4.3.1 Predictive accuracy against trajectory data
We study the predictive accuracy of the 2D model (3) with respect to the trajectory data provided by the data. In the first test we simply study the accuracy of the model without possible spurious errors included by the treatment of boundary data. We choose an initial time and using the kernel density estimation approach we compute the density profile. Then we evolve it up to a final time , such that seconds, in order to guarantee that the simulation is not influenced by outgoing boundary conditions. Moreover, at the same time we wish to reduce errors due to the numerical scheme and therefore we choose a very small CFL condition. Finally, we compute the difference between the simulated profile and the real density profile at final time, normalizing with respect to the maximum value, c.f. Figure 6 and in Figure 7.
Clearly, the boundary data are important for computing long time simulations. To this end, we extrapolate the incoming and the outgoing boundary data by artificially extending the trajectory data in computational cells outside the domain. In fact, recall that the trajectories of vehicles are approximated by means of a least squares linear approximation of their positions on the road section (see Section 2), thus we are able to detect the cars in the ghost part of the computational domain at a fixed time. Then, using the KDE technique, we compute the twodimensional density in the ghost cells which is used as boundary condition. More precisely, the extrapolation and the computation of the density in the ghost cells is based on the following procedure:

as described in Section 2, starting from the knowledge of the time dependent positions of each vehicle, we have considered the linear approximation of the data and , where and would represent the minimum and the maximum time, respectively, such that for each . Thus, for each vehicle , we have the linear trajectories
(12) minimizing

once the linear trajectories, and thus the coefficients and are known , we can use equations (12) to extrapolate the position along the road of a vehicle at time such that , computing

using step 2, we count and identify the vehicles that at time are in a position such that the KDE approach (10) applied to these data produces a nonzero density profile in the ghost cells. We use this density as boundary data.
Notice that, since we choose a very fine space discretization, the above extrapolation is supposed to be not too far beyond the known data.
In Figure 8 we study the predictive accuracy of the 2D model (3) for seconds taking into account the boundary data. We choose different time periods for the simulations in Figure 6 and in Figure 7. The error is computed every seconds on the whole domain using the norm error, see equation (11).
4.3.2 Comparison between the 1D and the 2D model.
We compare now the 2D model (3) with respect its 1D version (2) in order to estimate the benefit of a refined model compared with a commonly used averaged onedimensional model.
We select different initial conditions, characterized by different densities on the road. Then, we evolve the initial density profiles up to different final times seconds, with . In Figure 9 we compare the errors (the norm error as in (11)) on the density produced by the two macroscopic models (2) (red data) and (3) (blue data). The results show that the 2D model produces smaller errors and therefore it can results in more realistic evolution of traffic conditions. This is mainly due to the fact that, for the 1D model (2), the kernel density estimation as well as the definition of the prediction error is done following the same approach of [23]. Hence, we project the positions on the same coordinate and this may result in an overestimation of the density. In the projected point of view, we could have two or more vehicles being near each other leading to, in the kernel density estimation approach, higher values of density traveling thus with a lower speed. But in the realistic data vehicles could be distant due to different positions.
In order to better evaluate the performances between the 2D and the 1D model, in Figure 10 we also compare them in the prediction of traveling time, only for the case , thus for the same initial time used in the middle left panel of Figure 9. The traveling time for both models is computed up to the final times seconds, with , as
where and are the speed functions computed using the flux functions and in equation (4) and equation (5), respectively, with the optimal parameters determined in Section 3.1. The exact traveling time is instead computed by applying the kernel density estimation approach to find a continuous estimation for the fluxes at time as
where is defined in (9) and , are the microscopic velocities in  and direction of vehicle being on the road section at time . Finally, a continuous estimation for the speeds is computed as
where is given by (10) and the exact traveling time, for the 2D case, is
The traveling time for the 1D model is instead computed by applying the kernel density estimation approach for the flux following the same approach of [23] and thus, as done above, by projecting the positions on the same coordinate. In this way, we get the estimation for the flux in direction and then we can compute the estimation for the velocity which is used to find the exact traveling time.
4.3.3 Dependence on the fitting parameters.
We study the dependence of the presented results on the choice of the fitting parameters, those being crucial in the derivation. In particular, we are interested in the magnitude of the error changes due to the variation of the parameters defining the closure in direction, see equation (5).
We consider the same initial conditions as studied in Figure 6 and in Figure 7. In both cases we consider values of the fitting parameters and sampled from the intervals and , respectively. Then, we evolve the density profile with the 2D model (3) for seconds and for each of the pairs . Finally, we compute the errors at final time as in (11). We recall that the error for the optimal pair is for the time period and for the time period . Notice that the different parameters do not modify the errors strongly and therefore the presented procedure is robust against those variations. In order to quantify this consideration, in Figure 11 we show the relative difference between and . We observe that the maximum differences are of order and and therefore of the order of the numerical scheme.
5 Conclusions and Outlook
In this paper we proposed a twodimensional scalar macroscopic model to describe traffic flow on multilane roads. Therefore, the equation generalizes the onedimensional LWR model. We prescribed the closure laws describing the two flux functions by using a datafitting technique with respect to experimental measurements on a German highway.
Since laser sensors provide the twodimensional trajectory of vehicles, we recovered also the fundamental diagram for traffic behavior across lanes. To our knowledge, this the first time that the resulting behavior of the flow across lanes is taken into account. This is possible thanks to the particular traffic rules on European highways which lead to a nonnaive dynamics in the orthogonal direction to the movement of vehicles. In fact, if we consider experimental data on US highway, where there is no obligation to overtake on left lanes, the resulting behavior across lanes is naive. For instance, see Figure 12 in which we show the fundamental and speeddensity diagrams in direction computed from NGSIM data on US101 [2]. The mean dynamics of the flow seems to suggest and therefore the 2D model (3) reduces to the 1D LWRtype model (2).
On the German dataset, numerical examples show the validity of the macroscopic modeling when comparing with experiments. In particular, the numerical comparison with trajectory data shows that the twodimensional scalar model already outperforms a corresponding laneaveraged onedimensional model. From an application point of view, in future works we plan to investigate now the effect of regulations on lane reduction using the twodimensional setting as well as higherorder models. In fact, it is expected that as in [23] the additional degree of freedom in the modeling allows for a better adjustment of the models to data. Further, in [33] we study a secondorder twodimensional macroscopic model as limit of a twodimensional microscopic followtheleader model which is validated using the German dataset. Instead, in [35] we study a twodimensional hybrid kinetic model which allows to take into account the twodimensional phenomena of traffic flow in the kinetic theory showing, as application, that it is able to reproduce the US datasets.
Acknowledgments
This work has been supported by HE5386/1315 and DAAD MIUR project. We also thank the ISAC institute at RWTH Aachen, Prof. M. Oeser and MSc. F. Hennecke for kindly providing the trajectory data.
References
 [1]
 [2] Federal Highway Administration US Department of Transportation. Next Generation Simulation (NGSIM). http://ops.fhwa.dot.gov/trafficanalysistools/ngsim.htm.
 [3] Minnesota Department of Transportation. mn/DOT Traffic Data. http://data.dot.state.mn.us/datatools.
 [4] Mobile millennium project: http://traffic.berkeley.edu.
 [5] 75 Years of the Fundamental Diagram for Traffic Flow Theory: Greenshields Symposium, 2011. Transportation Research Board, Circular EC149.
 [6] S. Amin et al. Mobile century — Using GPS mobile phones as traffic sensors: A field experiment. In 15th World Congress on Intelligent Transportation Systems, New York, Nov. 2008.
 [7] D. AregbaDriollet and R. Natalini. Discrete kinetic schemes for multidimensional conservation laws. SIAM J. Numer. Anal., 37(6):1973–2004, 2000.
 [8] A. Aw and M. Rascle. Resurrection of “second order” models of traffic flow. SIAM J. Appl. Math., 60(3):916–938 (electronic), 2000.
 [9] A. M. Bayen and C. G. Claudel. LaxHopf based incorporation of internal boundary conditions into HamiltonJacobi equation. Part I: Theory. IEEE Trans. Automat. Contr., 55(5):1142–1157, 2010.
 [10] A. M. Bayen and C. G. Claudel. Convex formulations of data assimilation problems for a class of HamiltonJacobi equations. SIAM J. Control Optim., 49(2):383–402, 2011.
 [11] N. Bellomo and C. Dogbé. On the modeling of traffic and crowds: A survey of models, speculations, and perspectives. SIAM Rev., 53:409–463, 2011.
 [12] S. BenzoniGavage and R. M. Colombo. An populations model for traffic flow. European J. Appl. Math., 14(5):587–612, 2003.
 [13] F. Berthelin, P. Degond, M. Delitala, and M. Rascle. A model for the formation and evolution of traffic jams. Arch. Ration. Mech. Anal., 187:185–220, 2008.
 [14] S. Blandin, G. Bretti, A. Cutolo, and B. Piccoli. Numerical simulations of traffic data via fluid dynamic approach. Appl. Math. Comput., 210(2):441–454, 2009.
 [15] R. Cao, A. Cuevas, and W. G. Manteiga. A comparative study of several smoothing methods in density estimation. Comput. Stat. Data. An., 17(2):153–176, 1994.
 [16] B. N. Chetverushkin, N. G. Churbanova, I. R. Furmanov, and M. A. Trapeznikova. 2D micro and macroscopic models for simulation of heterogeneous traffic flows. In J. C. F. Pereira and Adélia Sequeira, editors, Proceedings of the ECCOMAS CFD 2010, V European Conference on Computational Fluid Dynamics, Lisbon, Portugal, 2010.
 [17] B. N. Chetverushkin, N. G. Churbanova, A.B. Sukhinova, and M.A. Trapeznikova. Congested traffic simulation based on a 2D hydrodynamical model. In B.A. Schrefler and U. Perego, editors, Proceedings of 8th World Congress on Computational Mechanics and 5th European Congress on Computational Methods in Applied Science and Engineering, WCCM8 & ECCOMAS 2008, Barcelona, Spain, 2008.
 [18] R. Courant, K O. Friedrichs, and H. Lewy. Über die partiellen differenzengleichungen der mathematischen physik. Mathematische Annalen, 100(1):32–74, 1928.
 [19] I. Cravero, G. Puppo, M. Semplice, and G. Visconti. CWENO: uniformly accurate reconstructions for balance laws. Math. Comp. In press. arXiv:1607.07319.
 [20] C. F. Daganzo. The cell transmission model: A dynamic representation of highway traffic consistent with the hydrodynamic theory. Transportation Research Part B, 28(4):269–287, 1994.
 [21] C. F. Daganzo. Requiem for second order fluid approximations of traffic flow. Transp. Res. B, 29B:277–286, 1995.
 [22] C. F. Daganzo. In traffic flow, cellular automata = kinematic waves. Transp. Res. B, 40:396–403, 2006.
 [23] S. Fan, M. Herty, and B. Seibold. Comparative model accuracy of a datafitted generalized AwRascleZhang model. Netw. Heterog. Media., 9(2):239–268, 2014.
 [24] S. Fan and B. Seibold. A comparison of datafitted first order traffic models and their second order generalizations via trajectory and sensor data. ArXiv eprints, August 2012.
 [25] M. Garavello and B. Piccoli. Traffic flow on networks, volume 1 of AIMS Series on Applied Mathematics. American Institute of Mathematical Sciences (AIMS), Springfield, MO, 2006. Conservation laws models.
 [26] P. Goatin. The Aw–Rascle vehicular traffic flow model with phase transitions. Math. Comput. Modeling, 44(3):287–303, 2006.
 [27] S. K. Godunov. Finite difference methods for numerical computation of discontinuous solutions of the equations of fluid dynamics. Mat. Sbornik, 47:271–295, 1959.
 [28] B. D. Greenshields. A study of traffic capacity. Proc. Highway Res., 14:448–477, 1935.
 [29] A. Harten. High resolution schemes for hyperbolic conservation laws. J. Comp. Phys., 49:357–393, 1983.
 [30] D. Helbing. Traffic and related selfdriven manyparticle systems. Reviews of Modern Physics, 73:1067–1141, 2001.
 [31] R. Herman and I. Prigogine. Kinetic theory of vehicular traffic. Elsevier, New York, 1971.
 [32] M. Herty and R. Illner. Analytical and numerical investigations of refined macroscopic traffic flow models. Kinet. Relat. Models, 3(2):311–333, 2010.
 [33] M. Herty, S. Moutari, and G. Visconti. Macroscopic modeling of multilane motorways using a twodimensional secondorder model of traffic flow. 2017.
 [34] M. Herty and L. Pareschi. FokkerPlanck asymptotics for traffic flow models. Kinet. Relat. Models, 3(1):165–179, 2010.
 [35] M. Herty, A. Tosin, G. Visconti, and M. Zanella. Hybrid stochastic kinetic description of twodimensional traffic dynamics. 2017.
 [36] M. Herty and G. Visconti. Analysis of risk levels for traffic on a multilane highway. 2017.
 [37] K. Heun. Neue Methoden zur approximativen Integration der Differentialgleichungen einer unabhängigen veränderlichen. Z. Math. Phys, 45:23–38, 1900.
 [38] S. P. Hoogendoorn, Faculty of Civil Engineering Delft University of Technology, and Geosciences. Traffic Flow Theory and Simulation: CT4821. TU Delft, 2007.
 [39] R. Illner, A. Klar, and T. Materne. VlasovFokkerPlanck models for multilane traffic flow. Commun. Math. Sci., 1(1):1–12, 2003.
 [40] M. C. Jones, J. S. Marron, and S. J. Sheather. A brief survey of bandwidth selection for density estimation. J. Am. Statist. Assoc., 91(433):401–407, 1996.
 [41] B. S. Kerner. The Physics of Traffic. Understanding Complex Systems. Springer, Berlin, 2004.
 [42] B. S. Kerner and P. Konhäuser. Cluster effect in initially homogeneous traffic flow. Phys. Rev. E, 48:R2335–R2338, 1993.
 [43] B. S. Kerner and P. Konhäuser. Structure and parameters of clusters in traffic flow. Phys. Rev. E, 50:54–83, 1994.
 [44] A. Klar and R. Wegener. A hierarchy of models for multilane vehicular traffic. I. Modeling. SIAM J. Appl. Math., 59(3):983–1001 (electronic), 1999.
 [45] A. Klar and R. Wegener. A hierarchy of models for multilane vehicular traffic. II. Numerical investigations. SIAM J. Appl. Math., 59(3):1002–1011 (electronic), 1999.
 [46] A. Klar and R Wegener. Kinetic derivation of macroscopic anticipation models for vehicular traffic. SIAM J. Appl. Math., 60(5):1749–1766 (electronic), 2000.
 [47] J. P. Lebacque. Les modeles macroscopiques du traffic. Annales des Ponts., 67:24–45, 1993.
 [48] W. Leutzbach. Introduction to the Theory of Traffic Flow. Springer, New York, 1988.
 [49] R. J. LeVeque. FiniteVolume Methods for Hyperbolic Problems. Cambridge University Press, 2002.
 [50] M. J. Lighthill and G. B. Whitham. On kinematic waves. II. A theory of traffic flow on long crowded roads. Proc. Roy. Soc. London. Ser. A., 229:317–345, 1955.
 [51] S. Maerivoet and B. De Moor. Traffic flow theory. Technical Report 05154, Katholieke Universiteit Leuven, 2005.
 [52] G. F. Newell. A simplified theory of kinematic waves in highway traffic II: Queueing at freeway bottlenecks. Transp. Res. B, 27:289–303, 1993.
 [53] E. Parzen. On estimation of a probability density function and mode. Ann. Math. Statist., 33:1065–1076, 1962.
 [54] H. J. Payne. Models of freeway traffic and control. Simulation Council, 1971.
 [55] H. J. Payne. FREFLO: A macroscopic simulation model for freeway traffic. Transportation Research Record, 722:68–77, 1979.
 [56] W. F. Phillips. A kinetic model for traffic flow with continuum implications. Transportation Planning and Technology, 5:31–138, 1979.
 [57] B. Piccoli and A. Tosin. Vehicular traffic: A review of continuum mathematical models. In R. A. Meyers, editor, Encyclopedia of Complexity and Systems Science, volume 22, pages 9727–9749. Springer, New York, 2009.
 [58] G. Puppo, M. Semplice, A. Tosin, and G. Visconti. Fundamental diagrams in traffic flow: the case of heterogeneous kinetic models. Commun. Math. Sci., 14(3):643–669, 2016.
 [59] G. Puppo, M. Semplice, A. Tosin, and G. Visconti. Analysis of a multipopulation kinetic model for traffic flow. Commun. Math. Sci., 15(2):379–412, 2017.
 [60] P. I. Richards. Shock waves on the highway. Operations Res., 4:42–51, 1956.
 [61] M. Rosenblatt. Remarks on some nonparametric estimates of a density function. Ann. Math. Statist., 27:832–837, 1956.
 [62] M. Rosini. Macroscopic models for vehicular flows and crowd dynamics: theory and applications. Springer, Basel, Switzerland, 2013.
 [63] C. W. Shu. High Order Weighted Essentially Nonoscillatory Schemes for Convection Dominated Problems. SIAM REVIEW, 51(1):82–126, 2009.
 [64] A.B. Sukhinova, M.A. Trapeznikova, B. N. Chetverushkin, and N. G. Churbanova. TwoDimensional Macroscopic Model of Traffic Flows. Mathematical Models and Computer Simulations, 1(6):669–676, 2009.
 [65] E. F. Toro. Riemann Solvers and Numerical Methods for Fluid Dynamics. Springer, Berlin, 2009.
 [66] R. T. Underwood. Speed, Volume, and Density Relationships: Quality and Theory of Traffic Flow. Technical report, Yale Bureau of Highway Traffic, 1961.
 [67] B. van Leer. Towards the ultimate conservative difference scheme III. Upstreamcentered finitedifference schemes for ideal compressible flow. J. Comp. Phys., 23(3):263–275, 1977.
 [68] H. M. Zhang. A nonequilibrium traffic model devoid of gaslike behavior. Transport. Res. BMeth., 36:275–290, 2002.