I Introduction
Road networks suffer from various types of equilibrium inefficiency due to selfish routing. Of particular interest is the Braess paradox, a counterintuitive phenomenon describing scenarios in which building new road links results in higher traffic delays at equilibrium. Since its introduction in 1968 in [braess1968paradoxon], the Braess paradox has been studied extensively to find efficient ways to predict, detect and prevent its existence.
Early results, such as [murchland1970braess], [leblanc1975algorithm], [fisk1979more], [stewart1980equilibrium], [frank1981braess], focus on the classic diamondshaped, fournode network with a single ODpair. A later study [steinberg] extends the analysis to general traffic networks to predict the occurrence of the paradox; however, the applicability of the results is limited by a restrictive assumption that all routes with nonzero flows in the original network are also utilized after the addition of a road link. Moreover, [steinberg] as well as related theoretical papers [dafermos1984some] and [taguchi1982braess] consider the special case when exactly one road link is built or removed from the network. Restricting the problem to a single link removal simplifies it to a polynomial time search, but limits the potential improvements.
Another approach to anticipating the Braess paradox is a network topology analysis. The study [Milchtaich] shows that a twoterminal network is immune to Braess paradox if and only if it is seriesparallel. In addition, the paper extends the result to account for any Pareto inefficiency in a twoterminal network. References [cenciarelli2016graph] and [chen2015excluding] extend the characterization from undirected graphs to directed graphs and allow for multiple commodities. Another theoretical study, [fujishige], explores the concept of matroids to identify networks that are immune to the Braess paradox. Unfortunately, few transportation networks exhibit a matroid structure; therefore, the applicability is limited in practice. A further shortcoming of the theoretical studies mentioned above is that they present structures that are immune to the Braess paradox, but do not provide tools to modify existing networks to eliminate this paradox and improve the efficiency of the equilibrium.
The computational study [ji2014detecting] proposes a mathematical programming method to detect the Braess paradox in a given network. The problem is formulated as a bilevel structure and then transformed into a singlelevel mixed integer program. By setting tolls to links and analysing the resulting network latency, the algorithm detects links that cause the Braess paradox and penalizes them to imitate road closure.
Instead of links, in our paper we search for ‘Braess routes’ (routes that cause the Braess paradox) with a greedy optimization algorithm and remove them sequentially from the navigation system. This leaves the drivers with a subset of routes that are immune to the paradox, which they are free to choose from. We believe that removing routes is advantageous over removing links, as removing a link adversely affects all other routes going through this link. Among other scenarios, our approach enables removing throughtraffic from residential roads, while allowing them to continue to serve the residents. One shortcoming of the route removal approach, however, is that customers might simply switch to another navigation system if their freedom of choice is limited. Further research is needed to address this issue, such as splitting populations into “selfish” and “altruistic” to model different levels of cooperation or creating a system of benefits to encourage drivers’ participation.
Unlike other methods for detecting and eliminating the Braess paradox, our model accounts for intersection structures and queues. A more accurate graph representation achieved by fitting functions for link delays and queue delays from data allows us to make theoretical methods applicable to realworld networks. We validate our approach on an extended graph model of a road network from North Bethesda, Montgomery County, Maryland. We slightly modified the geometry of this network with additional links to be able to test a bigger number of ODpairs and routes. We were able to demonstrate up to 12% delay reduction on this extended network. Despite this improvement, we do not claim the proposed route removal strategy is optimal. Indeed, as shown in [roughgarden], the problem of finding the optimal subnetwork is NPhard. Instead we trade optimality with computational tractability and applicability to realworld networks.
Ii Network Specification
We consider a network with several routes available to each ODpair. Every route is represented by a set of consecutive links connecting the origin to the destination.
Iia Link delay function
We estimate the time delay that a vehicle experiences on each link with a Bureau of Public Roads (BPR) function
[united1964traffic]:(1) 
where
is the freeflow time (),
and are parameters that depend on the link’s properties,
is the link throughput capacity (in ).
To compute the values of parameters and for a particular link, we followed the method described in [kucharski2017estimating]
, which requires speeddensity data points for a curvefitting algorithm that outputs appropriate parameters. We used Simulation of Urban Mobility (SUMO) opensource software to build a test case based on the modified North Bethesda, Montgomery County, Maryland network around the intersections of Montrose Rd and Montrose Pkwy (Fig.
1). In the absence of historical data, we generated traffic data from the simulation environment and estimated the delay functions. The derivation is detailed in Appendix A.IiB Queue delay estimation
Link delay functions do not account for the queues, which accumulate when the link throughput capacity is insufficient for the incoming flow. To address this issue and model queue delay we introduce “phantom” links, i.e. links that have no physical analogue in the realworld or simulation, but exist solely to account for additional delay related to queues. We insert “phantom” links into routes between consecutive edges incoming to intersections and edges leaving the same intersections (Fig. 2).
The analysis in Appendix B yields queue delay function:
(2) 
where
is expected constant delay due to intersection structure,
and are parameters that depend on the link’s properties,
is the link’s saturation rate.
The proposed function is a continuous nondecreasing nonnegative piecewise linear function. Therefore, it can be readily used in our route detection algorithm.
Iii Algorithm
Iiia Wardrop Equilibrium computation
To compute the Wardrop Equilibrium state we use a convex optimization problem that utilizes Beckmann potentials [beckmann1956studies] for delay functions. The objective function of the optimization problem to minimize is the sum of Beckmann potentials, i.e. integrals of delay functions, across all links:
(3) 
where
is the set of links in the network,
is the Beckmann potential of the link (),
is the delay function of the link ,
is the flow on the link (),
is the vector of route flows,
is the routing matrix defined as
(4) 
The first constraint ensures nonnegative flows on routes:
(5) 
where
is the set of routes corresponding to the ODpair ,
is the set of ODpairs.
The second constraint guarantees that flows on routes corresponding to the same ODpair sum up to the demand for that pair:
(6) 
where is the demand on the ODpair .
Therefore, the problem is:
(7) 
The solution of the problem (7) is the vector of route flows at equilibrium, which results in the total network delay:
(8) 
IiiB Elimination procedure
The objective of our algorithm is to find subsets of links/routes such that the remaining network has reduced latency. To achieve this we remove links and routes that cause the Braess paradox, which we refer to as ‘Braess’ links/routes, from the original system. Removing a route implies removing it from the set of options suggested by the navigation system (e.g., Google Maps); removing a link implies either physically closing down the road segment or reducing its capacity with tolls or signaling. As explained in the Introduction, route removal is preferable in practice; however, for the completeness of the study we discuss several approaches for route and link elimination, as they are easily obtained from the equilibrium computation method of the previous section. An important constraint in these approaches is to keep the connectivity of the network unchanged, i.e. every ODpair from the original network must retain at least one route in the reduced network.
In each approach we attribute values to the links/routes with the help of the optimization problem (7). First we obtain the equilibrium network delay, , for the original system. Then, we tentatively remove a link/route from the network (remove corresponding rows and columns from the routing matrix (4)) and solve (7) again, deriving the new equilibrium network delay, . Having both the original and the new network delays, the value of the removed link/route is:
(9) 
Links/routes with negative values are associated with the Braess paradox, because removing them reduces latency.
Greedy Single Link Removal

Find the link with the smallest value .

If , remove the link and revert to step (1).

Otherwise, terminate.
As stated earlier, removing a link affects the entire set of routes going through this link, which compromises the effectiveness of the algorithm. Moreover, this approach is slower than some algorithms discussed further.
Link Combination Removal

Compute the values of all possible combinations of links and find the one with the smallest value .

If , remove the corresponding combination.

Otherwise, the network is Braess paradoxfree.
Since the number of link combinations grows exponentially with the number of links, this method is computationally expensive and does not scale well to large networks.
LinkRoute Combination Removal This approach modifies the first method and tries to address the issue with subset of routes elimination. Instead of removing links completely,

For each link find the subset of routes utilizing this link that results in the minimal network delay.

Remove the routes that were not present in at least one optimal configuration.
Unlike previous algorithms, this method removes routes, but does it indirectly by working with linkroute combinations.
Greedy Single Route Removal

Find the route with the smallest value .

If , remove the route and revert to step (1).

Otherwise, terminate.
This approach, unlike the first method, deals with routes directly, which is faster since in practice the number of commonly used routes is smaller than the number of links. An ODpair can potentially have a large number of possible routes; however the dominant part of the drivers uses the a very limited subset of those routes. It also allows to identify the occurrence of the Braess paradox after the first iteration: if the route with the minimal value has zero initial flow at equilibrium, then the network is Braess paradoxfree.
Route Combination Removal

Compute the values of all possible combinations of routes and find the one with the smallest value .

If , remove the corresponding combination.

Otherwise, the network is paradoxfree.
This approach is similar to the second approach and has the same limitation due to computational complexity.
All five approaches were tested to identify the fastest and most accurate method for further implementation and verification. The improvements achieved on the sample network were comparable. Therefore, we report below the results for the fourth method (Greedy Single Route Removal) due to its computational speed and advantages in implementation.
Iv Simulation
Series of simulations were conducted to demonstrate the improvement, evaluate the prediction accuracy and estimate the computational precision of our algorithm when applied to physical networks. The testing was designed to model the behavior of a realworld traffic system, therefore a complex structure of ODpairs and corresponding routes was built upon (Fig. 1). The parameters of the simulation model are presented in Table I.
Parameters  Values 

Number of ODpairs  12 
Number of routes per ODpair  15 
Total number of routes  38 
Number of simulated vehicles  500  3200 
A side benefit of the route removal approach is that it may help small sideroads in residential areas that usually suffer from congestion due to drivers attempting to avoid busy freeways. In the era of Google maps, Waze and other navigation systems, which try to discover shortcuts and guide vehicles through neighborhoods to free up main roads, residents of these regions face busy traffic consequences and spend significantly more time than usual on short trips [Mac19].
The testing procedure consists of the following steps:

Solve the problem (7) for the original network to obtain the equilibrium flow distribution, , and equilibrium network delay, .

Derive the route configuration resulting in minimal delay via route elimination algorithm and reduce the network.

Solve the problem (7) for the reduced network to obtain the new equilibrium flow, , and network delay, . The theoretical delay reduction is: .

Feed into the SUMO to obtain the simulation network delay, , as a sum of travel times of all vehicles.

Feed into the SUMO to obtain the new simulated network delay, . The simulation delay reduction is .

Compare the theoretical and simulated delay reductions and , and estimated total delays and .
Every simulation set was run with different initial demands corresponding to ODpairs (last row of Table I). Congestion scenarios were of particular interest, because they allowed us to test the developed queue delay estimation model.
It is important to mention that our model does not account for spillbacks (full occupation of a link that causes the queue propagation to the upstream edge). Therefore, upper bounds on demands were applied to avoid spillbacks. Further research is required for an appropriate spillback representation.
V Results
In this section we present the results of our algorithm for several scenarios. The first two rows in Table II show the size of the simulated traffic. We addressed conditions associated with different times of day.
Setup 1  Setup 2  Setup 3  Setup 4  

Demand  Low  Medium  High  High 
Number of vehicles  500  1600  2600  3300 
Improvement ()  0%  3.2%  11.8%  8.1% 
Improvement Diff.  0%  5.3%  15.3%  10.7% 
Network Delay Diff.  1.7%  8.1%  9.5%  9.2% 
We are interested in three main metrics when analyzing the efficiency of our approach, which are presented in Table II. The first metric (Row 3) is the theoretical improvement in the network delay, i.e. total travel time saved by all vehicles after implementing our algorithm. In the low demand scenario, the network was Braess paradoxfree, and no travel time reduction was achieved. For moderate and congested traffic the improvement ranges from 3.2% to 11.8%, which is reasonable considering the insignificant effort required to modify the network configuration.
The second metric (Row 4) shows the difference between the theoretically predicted delay improvement and the corresponding simulation result. For the freeflow setup there was no Braess paradox and no comparison was needed. For medium and high demands, the difference varies between 5.3% and 15.3%. Taking into account the improvement value itself, we can conclude that the algorithm makes a relatively accurate prediction. Additionally, the theoretical solution always resulted in actual simulation improvement, i.e., .
The third metric (Row 5) is the difference between the predicted network delay and the simulated one. According to the results, we managed to keep the deviation under 10% for all scenarios. Furthermore, some cases demonstrated almost identical values for the network delays. The accuracy of the prediction suggests that the proposed graph representation accounting for queue delays yields a reasonable model to reflect real traffic conditions.
Vi Conclusion
We presented a greedy optimization algorithm to detect and sequentially remove the routes that contribute to the Braess paradox. The algorithm assigns values to existing routes in terms of how much the network delay would increase if they were removed, so that routes with negative values indicate a Braess paradox. The algorithm removes the route with the most negative value at each iteration. The resulting reduced network configuration demonstrates up to 12% travel time reduction. The proposed extended graph representation to account for queue delays provides a major improvement in how well the theoretical predictions match the simulated ones. We are interested in further studying the possible methods of modeling spillbacks and other congestionrelated phenomena.
Appendix A Link Delay Function Derivation
In this section we present a detailed link delay function derivation from empirical data points. We treated flow data for each link individually, simulating traffic on one edge at a time. To capture a wide range of flow values, we generated a new random flow from the interval [0,3000] every 200 seconds. The provided data were sufficient to construct BPRfunctions for 90% of links (Fig. 3). Tuning simulation parameters covered additional 70% of the remaining links (Fig. 4). For the rest of the failed edges, we used parameter values corresponding to delay functions of upstream successful links with minor modifications. We assumed these parameters are likely to have close values and, therefore, can be almost interchangeable. However, if this assumption is false the discrepancy between ground truth delay and estimated delay produced by inaccurate parameter substitution of one link is insignificant in a network scale.
The negligibility of this difference follows from the fact that road links constituting a network are relatively short in general. One edge in the graph network representation is usually displayed by a series of short connected links in the SUMO simulation network (Fig. 1). The travel time contribution of one link is in the order of few seconds, which is insignificantly small to deviate from the ground truth. Moreover, since the number of poorlyfitted links makes up less than 2% of the number of all links, we can conclude that the proposed design is sufficiently accurate.
Additionally, the link capacities can be extracted from fundamental diagrams built upon the collected data (in ) (Fig. 2(b)). The first method is to set the capacity estimate to the maximal recorded flow value. Another approach is to use a piecewiselinear approximation to derive the capacity as a ycoordinate of the intersection of two approximation lines.
Appendix B Queue Delay Function Derivation
In this section we present a detailed queue delay function derivation based on the queue formation analysis. Queues occur when the flow () on the link exceeds the saturation rate () of this link. The saturation rate is the upper bound on the number of vehicles able to leave the link within a period of time. Therefore, the queue formation depends on the difference between the inflow and the outflow on a particular link. Depending on the link type and its relative position with a specific intersection, we can distinguish several possible options for saturation rate estimation:

If the link incomes to a signalized intersection, its corresponding saturation rate equals to the maximum number of vehicles () able to pass the intersection on green light within one cycle normalized to one hour:
(10) where, is the cycle duration (in seconds).

If there is a STOP sign at the end of the link, its corresponding saturation rate is:
(11) where, is the delay (in seconds) a vehicle experiences when forced to stop at the STOP sign on an empty road. This value depends primarily on the speed limit on the link and has a small fluctuation from edge to edge.

If the link is free from any traffic regulation causing vehicles to stop or slow down, the saturation rate equals to the link physical throughput capacity.

The link is the secondary link, having no priority on an unsignalized intersection. Up to this point, the queue delay is given as which is a function of the flow on the corresponding link. In this case, however, the delay incurred by the intersection structure on the secondary link is also a function of the flow of the primary link. Although this violates the fundamental assumption of delay being a function of the corresponding link flow only, if the delay was symmetric for both the primary and secondary links, we could have still implemented this in our optimization problem and computed the equilibrium points. However, primary link sees no delay from the intersection hence the delay function is asymmetric which makes it hard if not impossible to compute equilibrium points the same way. Therefore, this scenario is not represented in our model and is subject to further research.
Similar to the link delay function derivation, we need to use realworld or simulated empirical data to estimate the parameters for queue delay function. We present the detailed function computation procedure for only the first scenario, which features a link incoming to a signalized intersection. Other scenarios utilize a slightly modified approach, which we will mention at the end of this section.
Knowing the traffic light cycle length () and keeping in mind the queue size dependency on the inflowoutflow difference, we estimate the one cycle queue growth rate as:
(12) 
where is the throughput of the link.
Based on the flow value, we can distinguish two scenarios:

The flow is smaller than the saturation rate (). In this case, , all vehicles are able to pass the intersection and no queue is forming. To account for a potential red phase arrival, instead of setting the queue delay to zero, we choose a specific constant , which is the expected value of a delay due to phase change:
(13) where is the duration of the red phase.

The flow is greater than the saturation rate (). In this case, and the queue is growing. To find the growth rate, the flowthroughput dependency and, thus, the saturation rate for a particular link are required.
Based on the queue formation model presented earlier, queue delay function depends on the saturation rate of the link. To find the onecycle throughput capacity of the link we simulate flows of various values in ascending order for a short period of time (200 seconds) each and record the number of vehicles leaving the link (entering the intersection) within one cycle. The flow point at which the throughput linear growth stops (Fig. 4(a)) corresponds to the onecycle saturation rate of the link. To obtain the value for onehour period, scale the result by . The throughput of the link (in ) is a piecewise function of the flow which has the form:
(14) 
For the flows smaller than the saturation rate, the throughput is a simple linear function, because all vehicles are able to leave the link. Otherwise, the throughput equals to the sum of the saturation rate and an additional term, which never exceeds 2. This term accounts for rare scenarios with one or multiple abnormally fastmoving vehicles that manage to exceed the usual saturation rate.
Feeding the derived throughput function back into equation (12), we obtain a onecycle gain in numbers of vehicles to the queue due to excessive flow. It is important to note that the imposed queue delay varies among vehicles within the same flow, and drivers in the head of a heavy traffic would spend significantly less time in the queue, than the ones in the tail. To avoid complications in the optimization problem formulation, an average delay can be assigned to every vehicle on the link. Plotting the computed delay results in a piecewise linear function (Fig. 5(a)).
The last step is to fit a linear function to the nonconstant area of the graph (Fig. 5(b)) and learn the parameters and resulting in equation (2).
The presented approach can also be used to determine queue delay functions for links ending up with the STOP sign. Both the expected delay and the cycle length would be artificially set to from equation (11), since every vehicle without exception is required to stop at the STOP sign, which takes it exactly seconds to do.