Congestion Reduction via Personalized Incentives

Ali Ghafelebashi 0000-0001-8339-7960,
Meisam Razaviyayn 0000-0003-4342-6661,
and Maged Dessouky 0000-0002-9630-6201
This study was funded by a grant from the National Center for Sustainable Transportation (NCST), supported by the U.S. Department of Transportation’s University Transportation Centers Program. The contents of this project reflect the views of the authors, who are responsible for the facts and the accuracy of the information presented herein. This document is disseminated in the interest of information exchange. The U.S. Government and the State of California assumes no liability for the contents or use thereof. Nor does the content necessarily reflect the official views or policies of the U.S. Government and the State of California. This paper does not constitute a standard, specification, or regulation. This paper does not constitute an endorsement by the California Department of Transportation (Caltrans) of any product described herein.

With rapid population growth and urban development, traffic congestion has become an inescapable issue, especially in large cities. Many congestion reduction strategies have been proposed in the past, ranging from roadway extension to transportation demand management. In particular, congestion pricing schemes have been used as negative reinforcements for traffic control. In this project, we study an alternative approach of offering positive incentives to drivers to take different routes. More specifically, we propose an algorithm to reduce traffic congestion and improve routing efficiency via offering personalized incentives to drivers.  We exploit the wide-accessibility of smart devices to communicate with drivers and develop an incentive offering mechanism using individuals’ preferences and aggregate traffic information. The incentives are offered after solving a large-scale optimization problem in order to minimize the total travel time (or minimize any cost function of the network such as total Carbon emission). Since this massive size optimization problem needs to be solved continually in the network, we developed a distributed computational approach. The proposed distributed algorithm is guaranteed to converge under a mild set of assumptions that are verified with real data. We evaluated the performance of our algorithm using traffic data from the Los Angeles area. Our experiments show congestion reduction of up to 11% in arterial roads and highways.

§ INTRODUCTION

Today, traffic congestion is one of the most prevalent issues in large metropolitan areas, resulting in lowered quality of life for residents and economic losses. According to INRIX the United States' economy suffered $88 billion in losses and the drivers lost 99 hours due to traffic congestion in 2019 [1]. In addition to direct economic losses, traffic congestion can worsen air quality and adversely affect health conditions. According to the Transportation Research Board, vehicle emission is the main cause for air pollution [National Research Council (US). Committee for the Evaluation of the Congestion Mitigation and Air Quality Improvement Program and National Research Council, 2002, Health Effects Institute. Panel on the Health Effects of Traffic-Related Air Pollution, 2010, Zhang and Batterman, 2013, Hennessy and Wiesenthal, 1999, Selzer and Vinokur, 1974].

In this paper, the focus is on changing drivers' behavior by offering incentives in order to reduce traffic congestion. Perhaps one of the closest strategies to our solution is pricing mechanisms in the literature. Road pricing policies, such as assigning a fee or tax for driving on a highway/road, have been widely studied in theory and practice [Pigou, 1920, Knight, 1924, Verhoef et al., 2008, Van Amelsfort, 2009, Bouchelaghem and Omar, 2018, Zhang et al., 2013, Kachroo et al., 2016, Farokhi and Johansson, 2014, Groot et al., 2014, Cao et al., 2020]. In this approach, congestion reduction is expected as a result of discouraging people to use congested roads. Such pricing mechanisms could be dependent on different factors such as time [Zheng et al., 2016], distance [Daganzo and Lehe, 2015], or vehicle characteristics [Zhang et al., 2019].
While pricing is a promising approach from a market point of view, issues such as equity barriers complicate the implementation of congestion pricing/taxation schemes [Knockaert et al., 2012, Levinson, 2010, Martens et al., 2012, Raux and Souche, 2004, Ieromonachou et al., 2006, Hensher and Bliemer, 2014]. In addition, complexities and uncertainties in designing pricing mechanisms have prevented policymakers from implementing advanced congestion pricing schemes [Gu et al., 2018].
Tradable credits (TCs) or tradable mobility permits (TPMs) are another token-based pricing mechanism [Verhoef et al., 1997, Viegas, 2001, Raux, 2004, Fan and Jiang, 2013, Thøgersen and Møller, 2008, Wang et al., 2014]. The theoretical advantages of such tradable credits have been shown in [Tsekeris and Voß, 2009, Nie, 2015, Wu et al., 2012, Grant-Muller and Xu, 2014]. While such cap-and-trade programs have been implemented in some economic sectors, such as airport slot allocation [Fukui, 2010], it has not been implemented for individual-level personal travels and daily commutes [Dogterom et al., 2017] due to the design complexities of such token markets [Azevedo et al., 2018].

Lately, researchers have paid more attention to positive incentive policies. Based on the psychological theory of reactance, rewarding desirable behavior could work better than penalizing undesirable behavior [Brehm, 1966]. Moreover, rewarding is a more popular policy than a taxation approach [Knockaert et al., 2007]. While the effectiveness of rewarding in changing the individual's behavior have been shown in [Kreps, 1997] and [Berridge, 2001], there is a limited number of studies on the effectiveness of rewarding policies in the transportation area. Among these studies, the INSTANT project [Merugu et al., 2009] and the CAPRI project [Yue et al., 2015] have provided positive incentives to motivate commuters to avoid peak time for their travel. They have shown the effectiveness of a rewarding policy in congestion reduction. Another series of studies has been done in the Netherlands on the effectiveness of monetary rewards to avoid rush hour driving [Bliemer et al., 2009, Ben-Elia et al., 2011, Ben-Elia and Ettema, 2011, Ben-Elia and Ettema, 2011]. Xianbao et al. [Hu et al., 2015] provided information about the travel time for each offered departure time, and alternative time, so the driver does not depend on their experience in choosing the alternatives.

Recently, [Azevedo et al., 2018] offered token form incentives for different travel choices such as route, travel modes, and ride-sharing. The proposed model learns individuals' decisions and adapts to their preferences based on their travel history. While these policies were successful in short-term experiments, they do not necessarily result in permanent changes in the travelers' behavior [Kumar et al., 2016].
[Papadopoulos et al., 2018] and [Kordonis et al., 2019] provide route and an incentive or fee to volunteer truck drivers to improve the overall traffic condition in a budget balanced mechanism. [Papadopoulos et al., 2021] considers VOT (Value of Time) in the mechanism to make it personalized.

There have been different choices used as the incentive in transportation studies. [Fuji and Kitamura, 2003] and [Bamberg et al., 2003] used free bus tickets as an incentive to study the changes in the frequency of traveling by bus. In Australia, an early bird ticket program was offered to relieve the problem of rail overcrowding during peak hours [Currie, 2011]. Free WiFi and discounted fares have been effective for Beijing commuters to motivate them to avoid the rush hours of the morning [Zhang et al., 2014]. In the Tripod project [Azevedo et al., 2018], they offered token form incentives that can be redeemed for goods and services at local businesses. In the CAPRI project [Yue et al., 2015], users of the app collected points, and they could trade 100 points for $1 or use the points to play a game. Knockaert et al. [Knockaert et al., 2012] provided smartphones to users who could receive money for their credits or keep the smartphone by collecting enough credits.

The preferences of the drivers can be considered in the incentive offering platform. Mohan et al. [Mohan et al., 2019] divided effective factors on drivers' decision into two categories of static factors and dynamic factors. The static factors, which are fixed for each person, include the number of available transportation options and the distance of nearby transit. On the other hand, the dynamic factors include weather and travel purpose. They identified the factors from interviews and surveys, and concluded that personalizing can be advantageous for travel assistance.
Also, the drivers' preferences can be learned through interaction with the individual [Zhu et al., 2019, Azevedo et al., 2018]. Personalized incentives can be used to offer a unique alternative for each individual. The goal is to make the offer close to the individual's preferences and maximize the contribution to the network [Azevedo et al., 2018] or to minimize the cost of incentives by maximizing the probability of acceptance of the alternative route [Zhu et al., 2019].

In this paper, we study the problem of offering personalized incentives to maximize the global utility of the network. Although, previous studies such as [Bliemer et al., 2009] consider static rewards for static options like teleworking, biking, and walking, our model assigns different rewards for different alternative routes for different drivers based on the traffic condition and personalization factors. Consequently, we have a much larger system (and a resulting optimization problem).
The implementation of the proposed model could be through a smartphone app where the traffic data can be used to offer incentives to drivers. In addition, smartphones will help the central planner to distribute the computational load for finding the optimal incentive offering strategy. The rest of this paper is organized as follows. In Section <ref>, we present our models. Then to solve our optimization problem efficiently, we propose an efficient distributed algorithm. Results of our numerical experiments are presented in Section <ref> using data from the Los Angeles area.

§ INCENTIVE OFFERING MECHANISMS

Let us model the structure of the traffic network with a  directed graph $\cG = (\cV,\cE)$.
Here $\cV$ is the collection of all major intersections and ramps, which form the set of nodes in the graph. We use the set of edges $\cE$ to capture the connectivity of the nodes in the graph. Two different nodes are adjacent in the graph if it is possible to directly go from one to another without passing over any other node. The direction of an edge between two nodes is based on the direction of the road from which we can go from one point to another. We also use the notation $|\cE|$ to denote the total number of road segments/edges in our network (i.e. the cardinality of the set $\cE$). A route is a collection of adjacent edges that starts from one node and ends in another. We use the one-hot encoding scheme to denote the routes. In other words, a given route is represented by a vector $\mathbf{r} \in \{0,1\}^{|\cE|}$. Here, the $k$-th entry of vector $\mathbf{r}$ is one if the $k$-th edge is a part of route $\mathbf{r}$ and it is zero, otherwise.

Let  $\bT = \{1,\ldots,T\}$ denote the time horizon of interest assuming the system is currently at time $t = 1$. For any $t\in \bT$, we use the random vector $\bv_t \in \mathbb{R}^{|\cE|}$ to represent the traffic volume on different road segments at time $t$. The $k$-th entry of $\bv_t$ shows the total number of vehicles of road segment $k$ at time $t$. Notice that the offered incentives can change the drivers' behavior who are using the platform in the future and thus affecting the vector $\bv_t$.

We use $\cN$ to denote the set of drivers that we can influence their behavior through offering incentives. For any driver $n \in \cN$, let $\cR_n \subseteq \{0,1\}^{|\cE|}$ denote the set of possible route options for going from its origin to its destination. Let $\cI_n$ be the set of possible incentives we can offer to driver $n \in \cN$.
We also use the binary variable $s_i^{\br,n}\in\{0,1\}$ to represent the offered incentives. For any driver $n\in \cN$ and incentive $i \in \cI_n$, the variable $s_i^{\br,n} = 1$  if incentive $i$ is offered to driver $n$ to take route $r$; and $s_i^{\br,n} = 0$ otherwise. We assume that we incentivize each driver with only one offer, i.e., $\sum_{\br \in \cR_n} \sum_{i \in \cI_n}s^{\br,n}_i = 1$. Given any incentive offered to the drivers, we model the decision of the drivers stochastically. In particular, we assume after offering incentives, each driver $n$ chooses route $\br$ with a certain probability which depends on the amount of incentive, the route, and the driver's preferences, as described below.

The route preferences of drivers depend on different factors such as route travel time, gender, age, and particularly the (monetary) incentive provided to the drivers in our context. Such dependence can be learned using standard machine learning approaches in the presence of data [Xiong et al., 2019]. In this project, we rely on the model developed in [Xiong et al., 2019] for our preference modeling. We simplify their model by ignoring the less predictive features and only consider two major features: value of incentive and travel time. In particular, we assume that,  given incentive $i\in \cI_n$ to driver $n$, the driver chooses route $\br$ with probability

\begin{equation} \label{eq:DriverPreference}
p^{\br,n}_{i} = P(\widehat{T}_{\br},i),
\end{equation}

where $\widehat{T}_{\br}$ is the estimate of the travel time for route $\br$ provided by the incentive offering platform. Notice that when drivers make their routing decisions, they do not know the exact travel time $T_\br$ for route $\br$, but instead they rely on the estimate $\widehat{T}_{\br}$ provided by the system. Here, we make an implicit assumption that the drivers do not consider their own judgement about the travel time in their decision. However, if such individual biases for drivers exist, the system can learn it over time  using standard preference learning techniques. Modeling the drivers' behavior in a probabilistic fashion has its own benefits. The decision of a driver for a given incentive amount depends on many factors such as age, gender, income as also studied in [Xiong et al., 2019]. It is even likely that the driver's decision may depend on the driver's “state of mind” at the time that the incentive is offered. Thus, the features that influence the driver's decision are not completely known to the central planner. In such a setting, probabilistic models can be a better fit for modeling the system. For this reason, in the general area of “recommendation systems” in machine learning and statistics, probabilistic models have been widely used to model the behavior of individual users (drivers in our setting) [Zheng et al., 2016].
In addition, we do not assume any traffic control by modeling the probability of drivers' acceptance. Traffic control might be more effective but it needs an authority with the power of changing traffic which is not required in our framework. In our model, drivers can disregard the offers at any time but offers change the probability of accepting drivers' routing choices.

In the next subsections, we present our model and formulation in more detail. For the convenience of the reader, the list of notations defined here and later in the manuscript is presented in Appendix <ref>. We present our framework under two different scenarios: First, for simplicity of presenting the ideas, we study the case where it is possible to bring traffic flow below the network capacity. Then, we study the high demand scenario where there is no feasible strategy to bring the demand below the network capacity.

§.§ Scenario I: Operating Below Network Capacity

Let us first for simplicity assume that there exists a solution that all road segments operate below the capacity. Hence, for that solution, we can assume that the travel time will be based on the free flow traffic. As this section shows, this assumption will result in a mixed integer linear programming optimization which can be solved efficiently using standard solvers.

Given (<ref>), the expected value of the volume vector $\mathbf{v}_t$ can be computed as:

\begin{equation} \label{eq:Averagev_t}
\mathbb{E} \left[\mathbf{v}_t \right]= \sum_{n\in \mathcal{N}} \sum_{i \in \mathcal{I}_n} \sum_{\mathbf{r} \in \mathcal{R}_n} s_i^{\br,n} p^{\mathbf{r},n}_i \boldsymbol{\beta}_{\mathbf{r},t}
\end{equation}

where the vector $\boldsymbol{\beta}_{\br, t}, \in \mathbb{R}^{|\mathcal{E}|}$ shows the probability of being at different links of the network at time $t\in \mathbf{T}$, conditioned on the fact that driver $n$ is on route $\br$. For more details about the vector $\boldsymbol{\beta}_{\br, t}$, please refer to [Ma and Qian, 2018].

In order to minimize the drivers' total travel time while keeping the volume below the road segment capacity vector $\mathbf{v}_0$, we need to solve the following optimization problem

\begin{equation} \label{eq:LP2}
\small
\begin{split}
\min_{\{s_i^{\br,n}\}} \quad & \sum_{t\in \bT}\sum_{n\in \mathcal{N}} \sum_{r\in \cR_n} \sum_{i \in \mathcal{I}_n} s_i^{\br,n} p^{\mathbf{r},n}_i \boldsymbol{\beta}_{\mathbf{r},t}^{\intercal} \boldsymbol{\omega} \\
\st \quad & \sum_{n\in \mathcal{N}} \sum_{i \in \mathcal{I}_n} \sum_{\mathbf{r} \in \mathcal{R}_n} s_i^{\br,n} p^{\mathbf{r},n}_i \boldsymbol{\beta}_{\mathbf{r},t}   \leq \bv_0, \quad\forall t  \in \mathbf{T} \\
& \sum_{n\in \mathcal{N}} \sum_{i \in \mathcal{I}_n} \sum_{\br \in \cR_n} s_i^{\br,n} \eta_i \leq \Omega \\
& \sum_{\br \in \cR_n} \sum_{i \in \mathcal{I}_n}s_i^{\br,n} = 1\quad\forall n\in \mathcal{N} \\
& s_i^{\br,n} \in \{0,1\} \quad \forall n\in\mathcal{N}, \forall i\in \mathcal{I}_n, \forall \br \in \cR_n
\end{split}
\normalsize
\end{equation}

where $\boldsymbol{\omega} \in \mathbb{R}^{|\mathcal{E}|}$ is the vector of free flow travel time of the links, $\Omega$ is the total available budget, and $\eta_i$ is the cost of offering incentive $i$. To keep this optimization problem tractable, we rely on the assumption of a large number of vehicles in each road segment and approximate the random quantity $\bv_t$ with its average $\mathbb{E}[\bv_t]$ provided in equation (<ref>). Notice that the objective function is equal to

\begin{equation}
\small
\min_{\{s_i^{\br,n}\}} \sum_{n\in \mathcal{N}} \sum_{r\in \cR_n} \sum_{i \in \mathcal{I}_n} s_i^{\br,n} p^{\mathbf{r},n}_i \sum_{t\in \bT} \boldsymbol{\beta}_{\mathbf{r},t}^{\intercal} \boldsymbol{\omega} \nonumber
\normalsize
\end{equation}

in which $\sum_{t\in \bT} \boldsymbol{\beta}_{\mathbf{r},t}^{\intercal} \boldsymbol{\omega}$ is the expected travel time of driver $n$ driving on route $\br$.

Problem (<ref>) is a mixed integer linear program which can be solved via standard  solvers such as Gurobi, AMPL, GAMS, and CPLEX. We use Gurobi in our experiments because of its powerful LP solver.

§.§ Scenario II:  Operating Above Network Capacity

In this subsection, we assume that the demand is elevated; thus, there is no incentive offering strategy that can bring the traffic flow below the network capacity. In such a scenario, we still can “improve” the congestion via incentivizing individual drivers. Our goal is to optimize a utility of the system as a criterion to compare the traffic condition after incentivizing. To make the formulation more specific, we  use total travel time as the utility function. It is worth noting that while we use this utility, following our steps, one can use any other utility functions such as carbon emissions or energy consumption.

To compute the total travel time of the system, we sum the travel time of the drivers of all the links over all time periods:

\begin{equation} \label{ObjFunct1}
\small
\begin{split}
F_{tt}(\hat{\bv}) = & \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|} \hat{v}_{\ell, t} \delta_{\ell, t}(\hat{v}_{\ell, t})  \\
\end{split}
\normalsize
\end{equation}

where $\delta_{\ell, t}$ is the travel time of link $\ell$ at time $t$ (which itself is a function of the volume). Here, $\hat{\bv}$ is the vector of volume of links in which $\hat{v}_{\ell, t}$ is the $(|\cE|\times t + \ell)^{th}$ element of vector $\hat{\bv}$ representing the volume of the $\ell^{th}$ link at time $t$.

To understand the impact of our offered incentives, we estimate the drivers' decision based on the provided incentives, which in turn results in estimating the volume of the links in the horizon of interest. Given these estimated volume values, we  estimate the travel time in the links as described below.

Travel time value $\delta$: There are different functions that capture the relation between travel time and volume. For example, the link congestion function developed by the Bureau of Public Roads (BPR) [United States. Bureau of Public Roads, 1964] defines a nonlinear relation between the volume and travel time of the road segments:

\begin{equation}
\small
\begin{aligned}
f_{\text{BPR}}(v) = t_0 \left(1+0.15\left(\frac{v}{w}\right)^4\right) \nonumber
\end{aligned}
\normalsize
\end{equation}

where $f_{\text{BPR}}(v)$ is the travel time of the drivers on the link given the assigned traffic volume $v$; the parameter $t_0$ is the free flow travel time of the link; $v$ is the assigned traffic volume of the link; and $w$ is the practical capacity of the link. We learn $t_0$ and $w$ using historical traffic data. In order to estimate the total travel time of the system, we need to estimate the volume vector $\hat{\bv}$, which we discuss next.

Volume vector $\hat{\bv}$: To compute the volume vector, we need to know the routing  decision of the drivers to be able to (approximately) estimate their location at different times. Clearly, the drivers' decision is a function of the offered incentives. In other words, the location of a driver is dependent on the incentive that we assign to them, because the likelihood of various decisions changes with different incentives.
Let us first explain our notations for the offered incentives: For each driver, we have a one hot encoded vector describing which route has been incentivized and how much reward has been assigned to it. Thus, for each driver we have a binary vector $\bs_n \in \{0, 1\}^{|\cR|\cdot |\cI|}$ in which only one element has a value of one and it corresponds to the route and the incentive amount that we offer. As we need one vector for each driver, we can aggregate all our incentivization strategy in a matrix $\bS \in \{0, 1\}^{(|\cR|\cdot|\cI|) \times |\cN|}$. Naturally, routes that are not relevant to that OD pair of a driver will get a value of zero in the corresponding incentive vector (since we cannot offer those routes to the driver).

To understand the drivers' responses to our offered incentives, we need to estimate the probability of acceptance of incentivized routes under different incentives including zero incentive (i.e, no incentive). To model this probability, we use the utility function developed in [Xiong et al., 2019] and compute the probability of acceptance of each offered incentive (by using a Softmax function on top of the utility). While the model in [Xiong et al., 2019] takes many parameters (such as gender, age, and education of the driver) as input, in our model and numerical experiments we only consider static parameters of the travel time and the reward value  to generate the probability of acceptance of a given incentive/reward.
However, our framework is modular and we can use any prediction model that can estimate drivers' behavior given an incentive amount. We can use any personalized routing model that can learn drivers’ behavior such as a neural network.
Let $\bP\in [0, 1]^{|\cR| \times (|\cR|\cdot|\cI|)}$ be a matrix encoding the information of probability of picking different routes given the offered (route, incentive) pairs. Thus, the vector $\bP \bS \mathbf{1} \in \mathbb{R}^{|\cR|\times 1}$ shows the expected number of vehicles in each route.

Given the number of vehicles in each route, the location of each driver for the next time horizon can be modeled in a probabilistic fashion. For this purpose, we rely on the model developed in [Ma and Qian, 2018] where a specific matrix $\bR  \in [0, 1]^{(|\cE|\cdot|\bT|) \times |\cR|} $ is proposed to estimate the probability of the presence of a driver in a given road segment at a specific time in the future (assuming that the driver is picking a specific route). We can compute matrix $\bR$ by running a simulation model if we have enough computation power. In our experiments in subsection <ref>, we rely on the historical data in computing matrix $\bR$. Thus, the vector

\begin{equation} \label{xx}
\small
\begin{split}
\hat{\bv} =  & \bR \bP \bS \mathbf{1} \in \bbR^{(|\cE|\cdot|\bT|)\times1} \nonumber\\
\end{split}
\normalsize
\end{equation}

represents the expected number of vehicles in all the links at each time slot. Substituting the expression of $\hat{\bv}$ in (<ref>), we get

\begin{equation} \label{eq:F_tt}
\small
\begin{split}
F_{tt}(\hat{\bv}) =  & \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|} (\bA\bS\mathbf{1})_{\ell, t} \delta((\bA\bS\mathbf{1})_{\ell, t}) \\
= & \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|} (\ba_{\ell, t}\bS\mathbf{1}) \delta(\ba_{\ell, t}\bS\mathbf{1}) \\
\end{split}
\normalsize
\end{equation}

where $\ba_{\ell, t}$ is the row of matrix $\bA = \bR \bP$ which corresponds to link $\ell$  at time $t$.
Thus in order to minimize the total travel time in the system via providing incentives to drivers, we need to solve the following optimization problem:

\begin{equation}
\label{eq:ObjFunctFinal}
\small
\begin{split}
\min_{\bS} \quad & \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|} (\ba_{\ell, t}\bS\mathbf{1}) \delta(\ba_{\ell, t}\bS\mathbf{1}) \\
\st \quad & \bS^{\intercal} \mathbf{1} = \mathbf{1},\quad
\bc^{\intercal} \bS \mathbf{1} \leq \Omega\\
& \bD \bS \mathbf{1} = \bq,\quad
\bS \in \{0, 1\}^{(|\cR||\cI|) \times |\cN|}
\end{split}
\normalsize
\end{equation}

where $\bc \in \bbR_{+}^{|\cR|\cdot|\cI|}$ is the vector of cost of incentives assigned to each route, $\bD \in \{0, 1\}^{K \times (|\cR|\cdot|\cI|)}$ is the matrix of incentive assignment to the OD pairs, and $\bq \in {\bbR}^{K \times |\cI|}$ is the vector of the number of drivers for each OD pair. Here, $K$ is the number of OD pairs. We explain the constraints in more details below:

Constraint 1 ($\bS^{\intercal} \mathbf{1} = \mathbf{1}$): This constraint simply states that we only  assign  one incentive to each driver.

Constraint 2 ($\bc^{\intercal} \bS \mathbf{1} \leq \Omega$): This is our budget constraint. The vector $\bc \in \bbR^{|\cR|\cdot|\cI|}$ represents the cost of the different rewards assigned to each driver. $\Omega$ is the total budget.

Constraint 3 ($\bD \bS \mathbf{1} = \bq$): This constraint makes sure that we offer the correct number of rewards for the routes between OD pairs. Recall that $\bS \mathbf{1}$ represents the (expected) number of drivers that have been offered different routes given different rewards. We use matrix $\bD$ to sum the number of drivers that received different reward offers for routes between the same OD pair. $\bq$ is the vector of the actual number of drivers that are travelling between OD pairs and $\bD \bS \mathbf{1}$ must be equal to $\bq$.

Constraint 4 ($\bS \in \{0, 1\}^{(|\cR||\cI|) \times |\cN|}$):  This  constraint imposes binary  structure on our decision parameters. In other words, 0 is not choosing an incentive and 1 is selecting the incentive.

To illustrate our model and the above constraints, we provide an example in Appendix <ref>.

§.§ Algorithm for Offering Incentives and A Distributed Implementation

The optimization problem (<ref>) is of large size while it needs to be solved in almost real time (or  hourly if the drivers send their travel information to the central planner every hour before their trip) in the network. However, due to the existence of binary variable $\bS$, solving this problem efficiently is difficult[We conjecture that  problem (<ref>) is NP-hard to solve  since it is a special instance of polynomial optimization with discrete variables and there does not appear to be any special structure in function $f$ to reduce its complexity.].
In order to develop an efficient “approximate” solver for (<ref>), we first relax the binary constraint in (<ref>) and replace it with the relaxed convex constraint $ \bS \in [0, 1]^{(|\cR||\cI|) \times |\cN|}$, leading to the relaxed formulation

\begin{equation}
\small
\label{eq:ObjFunctRelaxed}
\begin{split}
\min_{\bS} \quad & \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|} (\ba_{\ell, t}\bS\mathbf{1}) \delta(\ba_{\ell, t}\bS\mathbf{1}) \\
\st \quad & \bS^{\intercal} \mathbf{1} = \mathbf{1},\quad
\bc^{\intercal} \bS \mathbf{1} \leq \Omega\\
& \bD \bS \mathbf{1} = \bq,\quad
\bS \in [0, 1]^{(|\cR||\cI|) \times |\cN|}.
\end{split}
\end{equation}

The constraints in the above optimization problem are clearly convex. By substituting $\ba_{\ell, t}\bS\mathbf{1}$ by $\bga_{\ell, t}$, the objective function becomes a summation of monomial functions with positive coefficients. Moreover, $\bga_{\ell, t}$ is an affine mapping of the optimization variable $\bS$. Since our domain is the nonnegative orthant and monomials are convex in this domain, the objective function is convex. This convexity will allow us to explore the use of standard solvers such as CVX [Grant and Boyd, 2014]. However, these solvers rely on methods such as interior point methods [Grant and Boyd, 2008] which requires $O(n^3)$ number of iterations with $n$ being the number of variables.
This heavy computational complexity prevents us from applying standard solvers.
In our context, each driver is equipped with a smartphone and; thus, we can distribute the computational burden of solving (<ref>) among the drivers. In what follows,  we propose a simple reformulation of the problem leading to a distributed algorithm for solving (<ref>). To present our algorithm, let us start by reformulating (<ref>) as

\begin{equation}
\label{eq:ObjFunctRelaxed2}
\small
\begin{split}
\min_{\bga,\bu,\bS,\bW,\bH,\bz,\beta} \quad & \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|} \bga_{\ell, t} \delta(\bga_{\ell, t}) \\ &- \frac{\tilde{\lambda}}{2} \sum_{r=1}^{|\cR|} \sum_{i=1}^{|\cI|}  \sum_{n=1}^{|\cN|}  \bH_{r, i, n}(\bH_{r, i, n}-1)\\
\st \quad & \bS \mathbf{1} = \bu,\quad
\bW^{\intercal} \mathbf{1} = \mathbf{1} \\
& \bD \bu = \bq,  \quad
\bA \bu = \bga  \quad\quad\quad \\
& \bH = \bS,  \quad
\bW = \bS \quad\quad\quad \\
& \bc^{\intercal} \bu + \beta = \Omega,\quad
\beta \geq 0  \\
&  \bH \in [0, 1]^{(|\cR|\cdot|\cI|) \times |\cN|}.\quad\quad\quad\\
\end{split}
\normalsize
\end{equation}

This formulation is amenable to the ADMM method [Boyd et al., 2011, Gabay and Mercier, 1976, Glowinski and Marroco, 1975, Hong et al., 2016, Huang et al., 2021, Barazandeh et al., 2021], which has a natural distributed implementation.
Our ADMM formulation (<ref>) shows that this computation burden can be distributed among drivers' cell phones. This distributed optimization/federated learning framework can have other standard advantages of federated learning/distributed systems [Li et al., 2020, Lowy et al., 2022]. For example, when proper privacy preserving mechanisms (such as differential privacy [Dwork, 2006]) are utilized, we can guarantee the privacy of drivers since they can participate in the optimization procedure without completely sharing their data.
The steps of this algorithm is summarized in Algorithm <ref> and the details of the derivation of its steps is provided in Appendix <ref>.
Distributed Incentivization via ADMM
[1]
Input: Initial values: $\bga^0$, $\bS^0$, $\bH^0$, $\bW^0$, $\bu^0$, $\beta^{0}$, $\bla^0_1 \in \mathbb{R}^{|\cR|\cdot|\cI| \times 1}$,  $\bla^0_2 \in \mathbb{R}^{|\cN| \times 1}$,  $\bla^0_3 \in \mathbb{R}^{K \times 1}$,  $\bla^0_4 \in \mathbb{R}^{|\cE| \cdot \bT \times 1}$,  $\bLa^0_5 \in \mathbb{R}^{|\cR|\cdot|\cI| \times |\cN|}$,  $\lambda^0_6 \in \mathbb{R}$,  $\bLa^0_7 \in \mathbb{R}^{|\cR|\cdot|\cI| \times |\cN|}$, Dual update step: $\rho$, Number of iterations: $\tilde{T}$.
$t = 0, 1, \ldots, \tilde{T}$

$\bu^{t+1} = (\rho\bI +\rho \bD^{\intercal}\bD + \rho\bA^{\intercal}\bA + \rho\bc\bc^{\intercal})^{-1}(\bla_{1}^{t} + \rho\bS^{t}\mathbf{1} - \bD^{\intercal}\bla_3^{t} + \rho\bD^{\intercal}\bq - \bA^{\intercal}\bla_4^{t} + \rho\bA^{\intercal}\bga^{t} - \bc(\lambda_6^{t} + \beta^{t} - \Omega))$

$\bW^{t+1} = (\rho\mathbf{1}\mathbf{1}^{\intercal} + \rho\bI)^{-1}(\rho\mathbf{1}\mathbf{1}^{\intercal} + \rho\bS^{t} - \bLa_7^{t} - \mathbf{1}\bla_2^{t \intercal}) $

$\bH^{t+1} = \mathbbm{1}({\rho>\tilde{\lambda}}) \Pi {\left(\left(\frac{1}{\rho-\tilde{\lambda}}\right)(\rho\bS^{t} - \bLa_5^{t} - \frac{\tilde{\lambda}}{2})\right)}_{[0, 1]} + \mathbbm{1}({\rho<\tilde{\lambda}}) \Pi {\left(\left(\frac{1}{\rho-\tilde{\lambda}}\right)(\rho\bS^{t} - \bLa_5^{t} - \frac{\tilde{\lambda}}{2})\right)}_{\{0, 1\}}$

$\bS^{t+1} = (\rho\bu^{t+1}\mathbf{1}^{\intercal} + \bLa_5^{t} + \rho\bH^{t+1} + \bLa_7^{t} + \rho\bW^{t+1} - \bla_1^{t}\mathbf{1}^{\intercal})(\rho\mathbf{1}\mathbf{1}^{\intercal} + 2\rho\bI)^{-1}$

$\ell = 0, 1, \ldots, |\cE|$
$\hat{t} = 1, \ldots, |\bT|$
$\gamma^{t+1}_{\ell, \hat{t}} = \argmin\limits_{\gamma_{\ell, \hat{t}}}   \gamma_{\ell, \hat{t}} \delta(\gamma_{\ell, \hat{t}})  + \bla^{t}_{4, (\ell, \hat{t})}(\ba_{\ell, \hat{t}} \bu^{t} - \gamma_{\ell, \hat{t}}) + \frac{\rho}{2} (\ba_{\ell, \hat{t}} \bu^{t} - \gamma_{\ell, \hat{t}})^{2}$
$\beta^{t+1} = \Pi\left(\Omega - \bc^{\intercal}\bu^{t+1} -\frac{1}{\rho}\lambda_6^{t} \right)_{\bbR_{+}}$

$\bla^{t+1}_{1} = \bla^{t}_{1} + \rho (\bS^{t+1}\mathbf{1} - \bu^{t+1})$

$\bla^{t+1}_{2} = \bla^{t}_{2} + \rho (\bW^{t+1 \intercal} \mathbf{1} -\mathbf{1})$

$\bla^{t+1}_{3} = \bla^{t}_{3} + \rho (\bD \bu^{t+1} - \bq)$

$\bla^{t+1}_{4} = \bla^{t}_{4} + \rho (\bA \bu^{t+1} - \bga^{t+1})$

$\bLa^{t+1}_{5} = \bLa^{t}_{5} + \rho (\bH^{t+1}  - \bS^{t+1})$

$\lambda^{t+1}_{6} = \lambda^{t}_{6} + \rho (\bc^{\intercal} \bu^{t+1} + \beta^{t+1} - \Omega)$

$\bLa^{t+1}_{7} = \bLa^{t}_{7} + \rho (\bW^{t+1} - \bS^{t+1})$

Return: $\bS^{\tilde{T}}$

In Algorithm <ref>, $\Pi(\cdot)_{[0,1]}$ is the projection operator that projects each entry of the input matrix to the interval $[0,1]$ and $\Pi(\cdot)_{\bbR_{+}}$ is the projection operator that projects each entry of the input matrix on to $\bbR_{+}$. For a better convergence of Algorithm <ref>, we randomly permute the order of update of two blocks of ADMM at each iteration. The first block consists of variables $\bu, \bW$, and $\bH$ and the second block includes $\bS, \gamma_{\ell, \hat{t}}$, and $\beta$. In Algorithm <ref> the steps <ref>, <ref>, and <ref> are computationally cumbersome due to the size of the matrices $\bW,\bH,$ and $\bS$. However, notice that each column of the matrices $\bW,\bH$, and $\bS$ corresponds to single driver and hence the computation corresponding to each column can be performed in parallel on smartphone devices of the drivers. Moreover, since the steps are not coupled, they can be solved in parallel on the driver's smart devices. Theorem <ref> guarantees the convergence of our ADMM algorithm.

Algorithm <ref> finds an $\epsilon$-optimal solution of problem (<ref>) in $O(1/\epsilon)$ iterations [He and Yuan, 2012].

In optimization problem (<ref>) (and consequently (<ref>)), all solutions $\bS^{\ast}$ with a fixed value of $\bS^{\ast} \mathbf{1}=\bu^{\ast}$ lead to the same objective as long as $\bS^{\ast \intercal} \mathbf{1}= \mathbf{1}$. Hence, this convex problem can have infinite number of solutions (with many of them not even close to binary). Therefore, in order to find (approximately) binary solutions, we add the following regularizer to the objective function in (<ref>):

\begin{equation}
\Re(\bH_{r, i , n}) = - \frac{\tilde{\lambda}}{2} \bH_{r, i , n}(\bH_{r, i , n}-1)
\end{equation}

where $\tilde{\lambda} \in \bbR_{+}$ is the regularization parameter and $\bH_{r, i, n} \in [0, 1]$. This regularizer forces the elements of matrix $\bH$ to be as close as possible to the binary domain $\{0, 1\}$.

While Algorithm <ref> returns the solution of the optimization problem (<ref>), this problem (<ref>) is a relaxation of the original problem (<ref>). Hence, the obtained solution in Algorithm <ref> must be utilized to obtain a feasible point in $\eqref{eq:ObjFunctFinal}$. For this step, we solve the following mixed integer (linear) problem

\begin{equation}
\label{eq:linearModelS}
\small
\begin{split}
\min_{\bS} \quad & \| \bS\mathbf{1} - \bu^{\ast}\|_{1} \\
\st \quad & \bS^{\intercal} \mathbf{1} = \mathbf{1},\quad
\bc^{\intercal} \bS \mathbf{1} \leq \Omega\\
& \bD \bS \mathbf{1} = \bq, \quad
\bS \in \{0, 1\}^{(|\cR||\cI|) \times |\cN|}
\end{split}
\normalsize
\end{equation}

where $\bu^{\ast}$ is the optimal solution obtained by Algorithm <ref>. We can use off-the-shelf solvers such as Gurobi to solve (<ref>).

The BPR function of Algorithm <ref> in Section <ref> can capture both Scenario I in Section <ref> and Scenario II in Section <ref>. However, the computational requirements for the free-flow case in Scenario I in Section <ref> is less expensive compared to that of the congested case in Scenario II in Section <ref>. Thus, the model (<ref>) in Scenario I in Section <ref> is an alternative when computational resources are limited.

§ NUMERICAL EXPERIMENTS

We evaluate the performance of our algorithms using data from the Los Angeles area. The Los Angeles region is ideally suited for being the validation area because there are multiple routes connecting most OD pairs. Additionally, researchers at the University of Southern California have developed the Archived Data Management System (ADMS) that collects, archives, and integrates a variety of transportation datasets from Los Angeles, Orange, San Bernardino, Riverside, and Ventura Counties. ADMS includes access to real-time traffic data from 9500 highway and arterial loop detectors with measurements every 30 seconds and 1 minute respectively.

Due to the lack of access to the drivers' routing information, we need to estimate the origin-destination (OD) matrix from the network flow information. Rows and columns of the OD matrix correspond to the origin and destination points respectively. For OD matrix $A$, the element $A_{(i,j)}$ is the number of drivers going from point $i$ to point $j$. The OD matrix estimation problem is severely under-determined [Van Zuylen and Willumsen, 1980, Cascetta, 1984, Bell, 1991]. There are two categories of OD matrices: static and dynamic [Bera and Rao, 2011]. Due to the high resolution of our data, most of the existing dynamic OD estimation (DODE) methods become computationally inefficient. In addition, we do not have prior data of the OD matrix which many studies consider as given data [Krishnakumari et al., 2019, Carrese et al., 2017, Nigro et al., 2018, Kim et al., 2014] and we do not have access to prior observations of the OD matrix. Given these barriers, we relied on the algorithm proposed by [Ma and Qian, 2018]. This algorithm performs without employing any prior OD matrix information.

§.§ Simulation Model

In our numerical experiments, we integrate different datasets and models to evaluate the performance of our algorithms. First, we extract the speed data, volume data, and sensor information including the location of sensors from the Archived Data Management System (ADMS). Then, we use the distances of sensors, extracted from the location of sensors using Google Maps API, to create the graph of the network. We created three sets of graph networks corresponding to the regions depicted in Fig. <ref>, Fig. <ref>, and Fig. <ref>. In the next step, the speed data, volume data, and the network graph are used for estimating OD pairs by the algorithm provided in [Ma and Qian, 2018]. The total number of estimated incoming drivers for all three experiments is presented in Fig. <ref>. For each OD pair, we find up to 4 different routing options. In particular, we start by the shortest path for each OD pair. Then, we remove the edges in this path and go with the second shortest path, and we continue this process until we find 4 different routes between the origin and destination (or no other routes exist). We use the model in [Xiong et al., 2019] to compute the acceptance probability for the different offers on the different routes for each individual driver. The parameters used for computing the probability are static values provided by [Xiong et al., 2019]
and we only calibrate some of the parameters because we are only using historical data and we do not have access to the drivers' features such as age or gender.
In this paper, we do not learn the route choice model of the drivers so the parameters of the probability model are fixed but it is possible to adapt our routing model to the drivers’ preferences by observing the drivers’ behavior.
We run three different experiments that model the road network at different scales.  In Experiment I, we model an arterial region (Fig. <ref>) but includes surface streets. For Experiment II, we model a large  network of highways (Fig. <ref>). For Experiment III, we model a moderate region (Fig. <ref>) which is a subset of the region in Experiment II.
We use travel time savings as our metric for performance evaluation in all experiments. Beside travel time savings, we also include the monetary value of traffic reduction based on the Value of Time (VOT) as an alternative metric. Our base Value of Time (VOT) is derived from the estimation of [Castillo, 2020] which is $2.63 per minute or $157.8 per hour.
As we discussed in section <ref>, model (<ref>) assumes that there exists a traffic flow solution operating below the network capacity. When this assumption is not satisfied, our model results in an “approximate” solution. To evaluate the validity of this approximation in heavily congested networks, we run model (<ref>) for heavily congested networks (Experiments I and II) with many alternative routes so that we can reasonably reduce the congestion level. As we see in Experiments I and II, this model can provide a reasonable approximation in both arterial (Experiment I) and highways (Experiment II) and leads to congestion reduction even when the final result is above the system capacity. Experiment III, on the other hand, shows a limitation of model (<ref>) where incentives offered does not necessarily translate to reduced congestion. This is because of the heavy congestion and the lack of availability of enough alternative routes to reduce congestion (so that the final solution is far away from the free flow traffic and a linear approximation of travel time is no longer accurate enough). In this scenario, our numerical experiments demonstrate the superiority of model (<ref>) in reducing congestion even in such a heavily congested network. We were not able to run model (<ref>) for Experiments I and II due to the large number of nodes in the network. However, relying on edge computation, this model could be solved efficiently in practice as we discussed in subsection <ref>.
To solve model (<ref>), we use the Gurobi solver in all experiments. Also, we solve model (<ref>) in Experiment III utilizing Gurobi and MOSEK to compare their results with Algorithm <ref>.
Gurobi and MOSEK are the state of the art off-the-shelf commercial solvers of linear and mixed integer optimization problems.
To have a better balance between accuracy and required time for solving the problem, we set the relative mixed integer programming optimality gap at 0.01 for both Gurobi and MOSEK in the experiments.
The comparison between ADMM, Gurobi, and MOSEK is shown in Experiment III.
While we only provide incentives to the drivers that enter the system in the first time interval,  our incentive offering mechanism considers estimations of the traffic flows in the next time intervals.
The selected drivers for incentivization are from the same cohort. We randomly select a group of drivers between 7 AM and 7:15 AM. Then, we use the selected drivers to compare the performance of the model with different budget values on the total travel time for 7 AM to 8 AM.
While our formulation is static, it
can be applied in a dynamic environment if
solved frequently in the network in order to offer incentives to the drivers.

§.§ Experiment I

In Experiment I, we check the performance of model (<ref>) using the ADMS data for May 5th, 2018 with the incentive set $\cI=\{\$0,\$2,\$10\}$. The studied region, which is depicted in Fig. <ref>, includes the data of 301 sensors. Based on the ADMS data, we created a graph with 41 nodes, 139 links, and $105.5$ miles of road. OD points are located on intersections and close to the ramp of the highways. The number of OD pairs is 1681 and there are 4278 paths between them in total. We assume 7494 drivers enter the system between 7 AM and 8 AM and we consider 1805 drivers entering the system in the first 15 minutes for incentivization.
Results of model (<ref>) at 100% penetration rate are presented in TABLE <ref>. In this table, “Total travel time” shows the travel time computed via the BPR function after offering incentives. When the budget is increased from $1000 to $10,000, the percentage of total travel time decrease is improved by almost a factor of 2 from $4.07\%$ to $7.37\%$. The row “Cost” in the table shows the amount of the budget that was used. In all cases, almost all the budget is used. The value of saved time is much larger than budget except in the budget of $10000. Note that a budget of zero is the case of no incentive.

Studied region in Experiment I.

4c|Budget ($)
$0$         $100$             $1000$             $10000$
1|c|Cost ($)     $0$         $100$             $1000$             $9996$
1|c|Value of saved
time ($)     $0$         $775$             $4930$             $8927$
1|c|Total travel
time (hour)     $767$        $762$             $736$             $710$

Experiment I: Linear model (<ref>).

Notice that model (<ref>) may result in an infeasible optimization problem (particularly in a heavily congested network). Hence, we included a parameter $\alpha$ in our formulation as the multiplier of the allowed capacity. In other words, we use $\alpha\times \bv_0$ instead of $\bv_0$ in model (<ref>). We only consider this multiplier during the computation of incentives; however, during the computation of total travel time, we use the original capacity. We assumed zero dollar incentive in our probabilistic model for drivers that are not receiving an incentive.

of drivers
the system     Budget ($)     % of
drivers     Average
in total
travel time
2*7402     1000     6.67%     $2.00     4.07%
10000     13.64%     $9.78     7.37%

Comparison of $1000 and $10000 budget in Experiment I.

TABLE <ref> shows that increasing the budget results in higher percentage of drivers to whom we offered the incentive and a higher average amount of offered incentives. In addition, we observe that even offering incentives to 6.67% of the drivers (with an average of $2.00  monetary incentive per driver) can reduce the total travel time by 4.07%. If approximately 13.64% of the drivers are incentivized with an average of $9.78 per driver, the total travel time can be reduced by almost 7.37%. For more details about the distribution of offered incentives to drivers in Experiment I, please see TABLE <ref> in the Appendix.

Fig. <ref> shows the effect of the penetration rate (percentage of drivers who are considered in incentivization) on travel time decrease. By reducing the penetration rate, we experience smaller travel time decrease because of the flexibility of the model in selecting drivers decreases. Although reducing the penetration rate adversely affects the incentivization, the model focuses on available drivers for reducing the travel time. For more details on the numbers provided in Fig. <ref>, please see TABLE <ref> and TABLE <ref> in the Appendix.

Effect of the penetration rate on percentage of travel time decrease in Experiment I.

§.§ Experiment II

In Experiment II, we evaluate the performance of our methods for the region depicted in Fig. <ref> with 753 sensors under two different possible sets of incentives:

$\bullet$ $\cI_{1}=\{\$0, \$2,\$10\}$

$\bullet$ $\cI_2=\{\$0,\$1,\$2,\$3,\$5,\$10\}$

This region only includes data of highway sensors with 25 OD points and 32 links which includes $707.6$ miles of road. The number of OD pairs is 625, and there are 1331 paths between them in total. We assume 15093 drivers enter the system between 7 AM and 8 AM. Our incentivization model considers 4126 drivers entering the system in the first 15 minutes.
The results of our experiment at 100% penetration rate are presented in TABLE <ref> for incentive set $\cI_1$, and in TABLE <ref> for incentive set $\cI_2$. The value of saved time is much larger than budget for all budget values for both incentive set and it can go up to 97 times of the cost.

Studied region in Experiment II.

4c|Budget ($)
$0$         $100$             $1000$             $10000$
1|c|Cost ($)     $0$         $100$             $1000$             $9998$
1|c|Value of saved
time ($)     $0$         $6378$             $12919$             $15712$
1|c|Total travel
time (hour)     $5023$        $4983$       $4941$           $4924$

Experiment II: Linear model (<ref>) for incentive set $\cI_1$.

4c|Budget ($)
$0$         $100$             $1000$             $10000$
1|c|Cost ($)     $0$         $100$             $1000$             $10000$
1|c|Value of saved
time ($)     $0$         $9791$             $18499$            $28535$
1|c|Total travel
time (hour)     $5023$        $4961$             $4906$             $4843$

Experiment II: Linear model (<ref>) for incentive set $\cI_2$.

of drivers
the system               Budget
($)      % of
drivers     Average
amount     Reduction
in total
travel time
1|c|2*Exp. II-1
1000     2.97%     $2.23     1.63%
1|c|              10000     16.42%     $4.03     1.98%

1|c|2*Exp. II-2     2*15093     1000     3.87%     $1.71     2.33%
1|c|             10000     15.97%     $4.15     3.60%

Comparison of $1000 and $10000 budget in Experiment II.

In addition to confirming our previous observations in Experiment I, Experiment II shows the  diversity gain related to the incentive set $\cI_1$  (see the “Reduction in total travel time” column in TABLE <ref>). In other words, more choices in the incentive set provides more flexibility for the algorithm, resulting in the total travel time reduction.
For more details about the distribution of offered incentives to drivers in Experiment II, please see TABLE <ref> and TABLE <ref> in the Appendix. By examining Experiments I and II, we observe that more alternative routes leads to more gain in travel time reduction.

§.§ Experiment III

In Experiment III, we compare the performance of linear model (<ref>) and ADMM model (<ref>) using the incentive set $\cI=\{\$0,\$2,\$10\}$. The region considered in our analysis is depicted in Fig. <ref>. This region includes the data of 293 sensors. Based on the ADMS data, we created a graph with 12 nodes, 32 links, and $288.1$ miles of road. The number of OD pairs is 144 and there are 270 paths between them in total. The estimated total number of drivers incoming to the system between 5 AM to 9 AM by the OD estimation algorithm is depicted in Fig. <ref> (c).
In our simulations, we assume 8220 drivers enter the system between 7 AM and 8 AM. Our incentivization model considers 2248 drivers entering the system in the first 15 minutes.
Results of model (<ref>) and model (<ref>)  at 100% penetration rate are presented in TABLE <ref>.

Studied region in Experiment III.

4c|Budget ($)
$0$        $100$        $1000$        $10000$
(Linear)     Cost ($)     $0$      $100$        $1000$         $10000$
1|c|                       Value of
time ($)     $0$         $-7265$             $1498$             $-25118$
1|c|                      Total travel time (hour)      3005.22        3051        2995       3164
(ADMM)     Cost ($)      $0$       $100$          $1000$         $9954$
1|c|                       Value of
time ($)     $0$         $4105$             $37016$             $53646$
1|c|                       Total travel time (hour)      3005        2979        2771        2665

Experiment III: Linear model (<ref>) and model (<ref>).

Total estimated number of drivers entering the system (in 15-minute intervals). (a) Experiment I, (b) Experiment II, and (c) Experiment III.

of drivers
the system               Budget
($)      % of
drivers     Average
amount     Reduction
in total
travel time
1|c|2*Exp. III
Model (<ref>)     2*8220
1000     6.08%     $2.00     0.32%
1|c|                                            10000     15.86%     $7.71     -5.30%

1|c|2*Exp. III
Alg. <ref>     2*8220     1000     5.26%     $2.31     7.81%
1|c|                                            10000     16.17%     $7.49     11.31%

Comparison of $1000 and $10000 budget in Experiment III.

TABLE <ref> shows a major shortcoming of model (<ref>). Although the objective function in model (<ref>) reduces the total free flow travel time, the actual travel time is not based on the free flow travel time. Hence, model (<ref>) does not necessarily decrease the total travel time when the available budget increases as can be seen in TABLE <ref>. However, the decrease in total travel time is guaranteed in model (<ref>) for larger budget using Gurobi, MOSEK, and Algorithm <ref>. Model (<ref>) directly minimizes the travel time based on the BPR function so it captures the nonlinear relation between travel time and volume. Moreover, when the volume is greater than the capacity, model (<ref>) which is a more accurate model representing the traffic network produces better results.
These phenomena can be observed in Fig. <ref> and Fig. <ref>. We have not reported the travel time decrease by MOSEK because its results are less than $0.6\%$ different with the results of Gurobi. Although both Gurobi and MOSEK find slightly better solutions for model (<ref>), it can take take up to 35 hours for Gurobi and up to 80 hours for MOSEK to solve the problem. However, Algorithm <ref> takes at most 1.12 hours. Algorithm <ref> by offering incentives to 16.17% of the vehicles at 100% penetration rate in Experiment III (with an average of $7.49 monetary incentive per driver)  can reduce the total travel time by 11.31% using model (<ref>). For a budget of $\$1000$, model (<ref>) has 7.92% larger reduction in percentage of travel time compared to model (<ref>) although both offer almost the same amount of incentive on average to the almost same percentage of drivers. For a $10000 budget, not only model (<ref>) does not compete with model (<ref>) but also it increases the total travel time. Despite this shortcoming of model (<ref>), it is still useful in low congested areas and when the computational resources are limited. In most cases, model (<ref>) found a reasonable approximate solution except when the budget was $10000.
The computation time of model (<ref>) is 2.6 minutes, but model (<ref>) requires up to 14.12 hours to run if we employ  serial computation. Utilizing parallel computation as described in section <ref>, we can reduce the computational time to at most 1.12 hours.
If we decompose matrices A and D beforehand among the drivers in our parallel computation, we can further reduce the computational time to 0.11 hours (or 6 minutes). The value of saved time using Algorithm <ref> is much larger than budget for all budget values and it can go up to 41 times of the cost.
For more details about the distribution of the offered incentives to the drivers in Experiment III, please see TABLE <ref>, TABLE <ref>, and TABLE <ref> in the Appendix.
The effect of the penetration rate on travel time decrease in Experiment III for model (<ref>) is depicted in Fig. <ref> and Fig. <ref>. The behavior is similar to our observation in Fig. <ref> which was for Model (<ref>). For more details of the numbers provided in Fig. <ref> and Fig. <ref>, please see TABLE <ref>, TABLE <ref>, TABLE <ref>, and TABLE <ref> in the Appendix.

Effect of the penetration rate on percentage of travel time decrease in Experiment III, model (<ref>), Algorithm <ref>.

Effect of the penetration rate on percentage of travel time decrease in Experiment III, model (<ref>), Gurobi solver.

§ CONCLUSION

In this paper, we developed mathematical models and proposed algorithms for offering personalized incentives to drivers to reduce congestion in the network. In this framework, drivers share their origin-destination and routing information with a central planner. Based on this information, the central planner then offers incentives to drivers to incentivize/enforce a socially optimal routing strategy. The incentives are offered based on solving large-scale optimization problems in our framework. In our framework, we bring together prior works to model the behavior of drivers in response to the offered incentives as well as the resulting congestion reduction in the network where no traffic control is required.
We paid special attention to minimizing the total travel time of the network. In addition, we showed that this problem can be solved in a distributed fashion where some of the computations are performed on individual drivers’ smart devices. Finally, we evaluated the performance of our models and algorithms using Archived Data Management System (ADMS) data. Our experiments showed that the proposed framework can lead up to a 11% decrease in the total travel time of the system during rush hour times.

In this work, the incentives are only offered to alter the routing decision of the drivers. In future work, it is crucial to look at the effect of offering incentives to change the mode or time of the drivers’ trips. These options will bring additional flexibility to the model, which in turn will result in further congestion reduction. To compute the drivers' acceptance probability, we can include more aspects of drivers' characteristics and features into account such as gender, age. and salary. In addition, we can utilize preference learning in computing drivers' acceptance probability if we have access to the data of preferences of the drivers. All the codes for this project can be found in [86].

[1]
INRIX 2019 Global Traffic Scorecard.
Accessed: 2021-02-16.

[National Research Council (US). Committee for the Evaluation of the Congestion Mitigation and Air Quality Improvement Program and National Research Council, 2002]
National Research Council (US). Committee for the
Evaluation of the Congestion Mitigation and Air Quality
Improvement Program and National Research Council.
The congestion mitigation and air quality improvement program:
assessing 10 years of experience, volume 264.
Transportation Research Board, 2002.

[Health Effects Institute. Panel on the Health Effects of Traffic-Related Air Pollution, 2010]
Health Effects Institute. Panel on the Health Effects of
Traffic-Related Air Pollution.
Traffic-related air pollution: a critical review of the
literature on emissions, exposure, and health effects.
Number 17. Health Effects Institute, 2010.

[Zhang and Batterman, 2013]
Kai Zhang and Stuart Batterman.
Air pollution and health risks due to vehicle traffic.
Science of the total Environment, 450:0 307–316,

[Hennessy and Wiesenthal, 1999]
Dwight A Hennessy and David L Wiesenthal.
Traffic congestion, driver stress, and driver aggression.
Aggressive Behavior: Official Journal of the International
Society for Research on Aggression, 250 (6):0 409–423,

[Selzer and Vinokur, 1974]
Melvin L Selzer and Amiram Vinokur.
Life events, subjective stress, and traffic accidents.
American Journal of Psychiatry, 1310 (8):0
903–906, 1974.

[Pigou, 1920]
Arthur C Pigou.
The economics of welfare Macmillan and Co.
London, United Kingdom, 1920.

[Knight, 1924]
Frank H Knight.
Some fallacies in the interpretation of social cost.
The Quarterly Journal of Economics, 380 (4):0
582–606, 1924.

[Verhoef et al., 2008]
Erik Verhoef, Michiel C.J. Bliemer, Linda Steg, and Bert van Wee, editors.
Pricing in Road Transport.
Number 4192 in Books. Edward Elgar Publishing, July 2008.
ISBN ARRAY(0x4e494af0).
URL <https://ideas.repec.org/b/elg/eebook/4192.html>.

[Van Amelsfort, 2009]
Dirk Hendrik Van Amelsfort.
Behavioural responses and network effects of time-varying road

[Bouchelaghem and Omar, 2018]
Siham Bouchelaghem and Mawloud Omar.
Reliable and secure distributed smart road pricing system for smart
IEEE Transactions on Intelligent Transportation Systems,
200 (5):0 1592–1603, 2018.

[Zhang et al., 2013]
Guohui Zhang, Xiaolei Ma, and Yinhai Wang.
Self-adaptive tolling strategy for enhanced high-occupancy toll lane
IEEE Transactions on Intelligent Transportation Systems,
150 (1):0 306–317, 2013.

[Kachroo et al., 2016]
Pushkin Kachroo, Saumya Gupta, Shaurya Agarwal, and Kaan Ozbay.
Optimal control for congestion pricing: Theory, simulation, and
IEEE Transactions on Intelligent Transportation Systems,
180 (5):0 1234–1240, 2016.

[Farokhi and Johansson, 2014]
Farhad Farokhi and Karl H Johansson.
A study of truck platooning incentives using a congestion game.
IEEE Transactions on intelligent transportation systems,
160 (2):0 581–595, 2014.

[Groot et al., 2014]
Noortje Groot, Bart De Schutter, and Hans Hellendoorn.
Toward system-optimal routing in traffic networks: A reverse
stackelberg game approach.
IEEE Transactions on Intelligent Transportation Systems,
160 (1):0 29–40, 2014.

[Cao et al., 2020]
Zhiguang Cao, Hongliang Guo, Wen Song, Kaizhou Gao, Liujiang Kang, Xuexi Zhang,
and Qilun Wu.
Improving the performance of transportation networks: A
semi-centralized pricing approach.
IEEE Transactions on Intelligent Transportation Systems,
220 (10):0 6353–6364, 2020.

[Zheng et al., 2016]
Nan Zheng, Guillaume Rérat, and Nikolas Geroliminis.
Time-dependent area-based pricing for multimodal systems with
heterogeneous users in an agent-based environment.
Transportation Research Part C: Emerging Technologies,
62:0 133–148, 2016a.

[Daganzo and Lehe, 2015]
Carlos F Daganzo and Lewis J Lehe.
Distance-dependent congestion pricing for downtown zones.
Transportation Research Part B: Methodological, 75:0
89–99, 2015.

[Zhang et al., 2019]
Shu Zhang, Ann M Campbell, and Jan F Ehmke.
Impact of congestion pricing schemes on costs and emissions of
commercial fleets in urban areas.
Networks, 730 (4):0 466–489, 2019.

[Knockaert et al., 2012]
Jasper Knockaert, Yin-Yen Tseng, Erik T Verhoef, and Jan Rouwendal.
The Spitsmijden experiment: A reward to battle congestion.
Transport Policy, 24:0 260–272, 2012.

[Levinson, 2010]
David Levinson.
Equity effects of road pricing: A review.
Transport Reviews, 300 (1):0 33–57, 2010.

[Martens et al., 2012]
Karel Martens, Aaron Golub, and Glenn Robinson.
A justice-theoretic approach to the distribution of transportation
benefits: Implications for transportation planning practice in the United
Transportation Research Part A: Policy and Practice,
460 (4):0 684–695, 2012.

[Raux and Souche, 2004]
Charles Raux and Stéphanie Souche.
The acceptability of urban road pricing: A theoretical analysis
applied to experience in Lyon.
Journal of Transport Economics and Policy (JTEP), 380
(2):0 191–215, 2004.

[Ieromonachou et al., 2006]
Petros Ieromonachou, Stephen Potter, and James P Warren.
Evaluation of the implementation process of urban road pricing
schemes in the United Kingdom and Italy.
European Transport, 32:0 49–68, 2006.

[Hensher and Bliemer, 2014]
David A Hensher and Michiel CJ Bliemer.
What type of road pricing scheme might appeal to politicians?
viewpoints on the challenge in gaining the citizen and public servant vote by
staging reform.
Transportation Research Part A: Policy and Practice,
61:0 227–237, 2014.

[Gu et al., 2018]
Ziyuan Gu, Zhiyuan Liu, Qixiu Cheng, and Meead Saberi.
Congestion pricing practices and public acceptance: A review of
Case Studies on Transport Policy, 60 (1):0
94–101, 2018.

[Verhoef et al., 1997]
Erik Verhoef, Peter Nijkamp, and Piet Rietveld.
Tradeable permits: their potential in the regulation of road
transport externalities.
Environment and Planning B: Planning and Design, 240
(4):0 527–548, 1997.

[Viegas, 2001]
Jose M Viegas.
Making urban road pricing acceptable and effective: searching for
quality and equity in urban mobility.
Transport Policy, 80 (4):0 289–294, 2001.

[Raux, 2004]
Charles Raux.
The use of transferable permits in transport policy.
Transportation Research Part D: Transport and Environment,
90 (3):0 185–197, 2004.

[Fan and Jiang, 2013]
Wenbo Fan and Xinguo Jiang.
Tradable mobility permits in roadway capacity allocation: Review and
Transport Policy, 30:0 132–142, 2013.

[Thøgersen and Møller, 2008]
John Thøgersen and Berit Møller.
Breaking car use habits: The effectiveness of a free one-month
Transportation, 350 (3):0 329–345, 2008.

[Wang et al., 2014]
Guangmin Wang, Ziyou Gao, Meng Xu, and Huijun Sun.
Models and a relaxation algorithm for continuous network design
problem with a tradable credit scheme and equity constraints.
Computers & Operations Research, 41:0 252–261, 2014.

[Tsekeris and Voß, 2009]
Theodore Tsekeris and Stefan Voß.
Design and evaluation of road pricing: state-of-the-art and
methodological advances.
NETNOMICS: Economic Research and Electronic Networking,
100 (1):0 5–52, 2009.

[Nie, 2015]
Yu Marco Nie.
A new tradable credit scheme for the morning commute problem.
Networks and Spatial Economics, 150 (3):0
719–741, 2015.

[Wu et al., 2012]
Di Wu, Yafeng Yin, Siriphong Lawphongpanich, and Hai Yang.
Design of more equitable congestion pricing and tradable credit
schemes for multimodal transportation networks.
Transportation Research Part B: Methodological, 460
(9):0 1273–1287, 2012.

[Grant-Muller and Xu, 2014]
Susan Grant-Muller and Meng Xu.
The role of tradable credit schemes in road traffic congestion
Transport Reviews, 340 (2):0 128–149, 2014.

[Fukui, 2010]
Hideki Fukui.
An empirical analysis of airport slot trading in the United
Transportation Research Part B: Methodological, 440
(3):0 330–357, 2010.

[Dogterom et al., 2017]
Nico Dogterom, Dick Ettema, and Martin Dijst.
Tradable credits for managing car travel: a review of empirical
research and relevant behavioural approaches.
Transport Reviews, 370 (3):0 322–343, 2017.

[Azevedo et al., 2018]
Carlos Lima Azevedo, Ravi Seshadri, Song Gao, Bilge Atasoy, Arun Prakash
Akkinepally, Eleni Christofa, Fang Zhao, Jessika Trancik, and Moshe
Tripod: sustainable travel incentives with prediction, optimization,
and personalization.
In Transportation Research Board 97th Annual Meeting, 2018.

[Brehm, 1966]
Jack Williams Brehm.
A Theory of Psychological Reactance.
New York: Academic Press, 1966.

[Knockaert et al., 2007]
J Knockaert, M Bliemer, D Ettema, D Joksimovic, A Mulder, J Rouwendal, and
D Van Amelsfort.
Experimental design and modelling Spitsmijden.
Utrecht, Consortium Spitsmijden, 2007.

[Kreps, 1997]
David M Kreps.
Intrinsic motivation and extrinsic incentives.
The American Economic Review, 870 (2):0
359–364, 1997.

[Berridge, 2001]
Kent C Berridge.
Reward learning: reinforcement, incentives, and expectations.

[Merugu et al., 2009]
Deepak Merugu, Balaji S Prabhakar, and N Rama.
An incentive mechanism for decongesting the roads: A pilot program in
In Proc. of ACM NetEcon Workshop. Citeseer, 2009.

[Yue et al., 2015]
Jia Shuo Yue, Chinmoy V Mandayam, Deepak Merugu, Hossein Karkeh Abadi, and
Balaji Prabhakar.
Reducing road congestion through incentives: a case study.

[Bliemer et al., 2009]
Michiel Bliemer, Matthijs Dicke-Ogenia, and Dick Ettema.
Rewarding for avoiding the peak period: A synthesis of three studies
in the Netherlands.
In European Transport Conference 2009. Citeseer, 2009.

[Ben-Elia et al., 2011]
Eran Ben-Elia, Dick Ettema, and Hennie Boeije.
Behaviour change dynamics in response to rewarding rush-hour
avoidance: A qualitative research approach.

[Ben-Elia and Ettema, 2011]
Eran Ben-Elia and Dick Ettema.
Changing commuters’ behavior using rewards: A study of rush-hour
Transportation Research Part F: Traffic Psychology and
Behaviour, 140 (5):0 354–368, 2011a.

[Ben-Elia and Ettema, 2011]
Eran Ben-Elia and Dick Ettema.
Rewarding rush-hour avoidance: A study of commuters’ travel
Transportation Research Part A: Policy and Practice,
450 (7):0 567–582, 2011b.

[Hu et al., 2015]
Xianbiao Hu, Yi-Chang Chiu, and Lei Zhu.
Behavior insights for an incentive-based active demand management
International Journal of Transportation Science and
Technology, 40 (2):0 119–133, 2015.

[Kumar et al., 2016]
Vivek Kumar, Chandra R Bhat, Ram M Pendyala, Daehyun You, Eran Ben-Elia, and
Dick Ettema.
Impacts of incentive-based intervention on peak period traffic:
experience from the Netherlands.
Transportation Research Record, 25430 (1):0
166–175, 2016.

[Papadopoulos et al., 2018]
Aristotelis-Angelos Papadopoulos, Ioannis Kordonis, Maged Dessouky, and Petros
Coordinated freight routing with individual incentives for
IEEE Transactions on Intelligent Transportation Systems,
200 (9):0 3397–3408, 2018.

[Kordonis et al., 2019]
Ioannis Kordonis, Maged M Dessouky, and Petros A Ioannou.
Mechanisms for cooperative freight routing: incentivizing individual
IEEE Transactions on Intelligent Transportation Systems,
210 (5):0 2155–2166, 2019.

[Papadopoulos et al., 2021]
Aristotelis-Angelos Papadopoulos, Ioannis Kordonis, Maged M Dessouky, and
Petros A Ioannou.
Personalized pareto-improving pricing-and-routing schemes for
near-optimum freight routing: An alternative approach to congestion pricing.
Transportation Research Part C: Emerging Technologies,
125:0 103004, 2021.

[Fuji and Kitamura, 2003]
Satoshi Fuji and Ryuichi Kitamura.
What does a one-month free bus ticket do to habitual drivers.
Transportation, 30:0 81–95, 2003.

[Bamberg et al., 2003]
Sebastian Bamberg, Icek Ajzen, and Peter Schmidt.
Choice of travel mode in the theory of planned behavior: The roles of
past behavior, habit, and reasoned action.
Basic and Applied Social Psychology, 250 (3):0
175–187, 2003.

[Currie, 2011]
Graham Currie.
Free fare incentives to shift rail demand peaks–medium-term impacts.
Technical report, 2011.

[Zhang et al., 2014]
Zheng Zhang, Hidemichi Fujii, and Shunsuke Managi.
How does commuting behavior change due to incentives? an empirical
study of the Beijing Subway System.
Transportation Research Part F: Traffic Psychology and
Behaviour, 24:0 17–26, 2014.

[Mohan et al., 2019]
Shiwali Mohan, Matthew Klenk, and Victoria Bellotti.
Exploring how to personalize travel mode recommendations for urban
In IUI Workshops, 2019.

[Zhu et al., 2019]
Xi Zhu, Feilong Wang, Cynthia Chen, and Derek D Reed.
Personalized incentives for promoting sustainable travel behaviors.
Transportation Research Part C: Emerging Technologies, 2019.

[Xiong et al., 2019]
Chenfeng Xiong, Mehrdad Shahabi, Jun Zhao, Yafeng Yin, Xuesong Zhou, and Lei
An integrated and personalized traveler information and incentive
scheme for energy efficient mobility systems.
Transportation Research Part C: Emerging Technologies, 2019.

[Zheng et al., 2016]
Yin Zheng, Bangsheng Tang, Wenkui Ding, and Hanning Zhou.
A neural autoregressive approach to collaborative filtering.
In International Conference on Machine Learning, pages
764–773. PMLR, 2016b.

[Ma and Qian, 2018]
Wei Ma and Zhen Sean Qian.
Estimating multi-year 24/7 origin-destination demand using
high-granular multi-source traffic data.
Transportation Research Part C: Emerging Technologies,
96:0 96–121, 2018.

[United States. Bureau of Public Roads, 1964]
United States. Bureau of Public Roads.
Traffic assignment manual for application with a large, high
speed computer, volume 37.
US Department of Commerce, Bureau of Public Roads, Office of
Planning, Urban Planning Division, 1964.

[Grant and Boyd, 2014]
Michael Grant and Stephen Boyd.
CVX: Matlab software for disciplined convex programming, version
<http://cvxr.com/cvx>, March 2014.

[Grant and Boyd, 2008]
Michael Grant and Stephen Boyd.
Graph implementations for nonsmooth convex programs.
In V. Blondel, S. Boyd, and H. Kimura, editors, Recent Advances
in Learning and Control, Lecture Notes in Control and Information Sciences,
pages 95–110. Springer-Verlag Limited, 2008.
<http://stanford.edu/ boyd/graph_dcp.html>.

[Boyd et al., 2011]
Stephen Boyd, Neal Parikh, and Eric Chu.
Distributed optimization and statistical learning via the
alternating direction method of multipliers.
Now Publishers Inc, 2011.

[Gabay and Mercier, 1976]
Daniel Gabay and Bertrand Mercier.
A dual algorithm for the solution of nonlinear variational problems
via finite element approximation.
Computers & Mathematics with Applications, 20
(1):0 17–40, 1976.

[Glowinski and Marroco, 1975]
Roland Glowinski and A Marroco.
Sur l'approximation, par éléments finis d'ordre un, et la
résolution, par pénalisation-dualité d'une classe de
problèmes de dirichlet non linéaires.
ESAIM: Mathematical Modelling and Numerical
Analysis-Modélisation Mathématique et Analyse Numérique,
90 (R2):0 41–76, 1975.

[Hong et al., 2016]
Mingyi Hong, Zhi-Quan Luo, and Meisam Razaviyayn.
Convergence analysis of alternating direction method of multipliers
for a family of nonconvex problems.
SIAM Journal on Optimization, 260 (1):0
337–364, 2016.

[Huang et al., 2021]
Tianjian Huang, Prajwal Singhania, Maziar Sanjabi, Pabitra Mitra, and Meisam
Alternating direction method of multipliers for quantization.
In International Conference on Artificial Intelligence and
Statistics, pages 208–216. PMLR, 2021.

[Barazandeh et al., 2021]
Babak Barazandeh, Ali Ghafelebashi, Meisam Razaviyayn, and Ram Sriharsha.
Efficient algorithms for estimating the parameters of mixed linear
regression models.
arXiv preprint arXiv:2105.05953, 2021.

[Li et al., 2020]
Tian Li, Anit Kumar Sahu, Ameet Talwalkar, and Virginia Smith.
Federated learning: Challenges, methods, and future directions.
IEEE Signal Processing Magazine, 370 (3):0
50–60, 2020.

[Lowy et al., 2022]
Andrew Lowy, Ali Ghafelebashi, and Meisam Razaviyayn.
Private non-convex federated learning without a trusted server.
arXiv preprint arXiv:2203.06735, 2022.

[Dwork, 2006]
Cynthia Dwork.
Differential privacy.
In International Colloquium on Automata, Languages, and
Programming, pages 1–12. Springer, 2006.

[He and Yuan, 2012]
Bingsheng He and Xiaoming Yuan.
On the O(1/n) convergence rate of the Douglas–Rachford
alternating direction method.
SIAM Journal on Numerical Analysis, 500 (2):0
700–709, 2012.

[Van Zuylen and Willumsen, 1980]
Henk J Van Zuylen and Luis G Willumsen.
The most likely trip matrix estimated from traffic counts.
Transportation Research Part B: Methodological, 140
(3):0 281–293, 1980.

[Cascetta, 1984]
Ennio Cascetta.
Estimation of trip matrices from traffic counts and survey data: a
generalized least squares estimator.
Transportation Research Part B: Methodological, 180
(4-5):0 289–299, 1984.

[Bell, 1991]
Michael GH Bell.
The real time estimation of origin-destination flows in the presence
of platoon dispersion.
Transportation Research Part B: Methodological, 250
(2-3):0 115–125, 1991.

[Bera and Rao, 2011]
Sharminda Bera and KV Rao.
Estimation of origin-destination matrix from traffic counts: the
state of the art.
European Transport\Trasporti Europei, 2011.

[Krishnakumari et al., 2019]
Panchamy Krishnakumari, Hans van Lint, Tamara Djukic, and Oded Cats.
A data driven method for OD matrix estimation.
Transportation Research Part C: Emerging Technologies, 2019.

[Carrese et al., 2017]
Stefano Carrese, Ernesto Cipriani, Livia Mannini, and Marialisa Nigro.
Dynamic demand estimation and prediction for traffic urban networks
adopting new data sources.
Transportation Research Part C: Emerging Technologies,
81:0 83–98, 2017.

[Nigro et al., 2018]
Marialisa Nigro, Ernesto Cipriani, and Andrea del Giudice.
Exploiting floating car data for time-dependent origin–destination
matrices estimation.
Journal of Intelligent Transportation Systems, 220
(2):0 159–174, 2018.

[Kim et al., 2014]
JinYoung Kim, Fumitaka Kurauchi, Nobuhiro Uno, Takeshi Hagihara, and Takehiko
Using electronic toll collection data to understand traffic demand.
Journal of Intelligent Transportation Systems, 180
(2):0 190–203, 2014.

[Castillo, 2020]
Juan Camilo Castillo.
Who benefits from surge pricing?
Available at SSRN 3245533, 2020.

[86]
Github codes.
Accessed: 2021-04-18.

§ LIST OF NOTATIONS

The following symbols are used in this paper.

* $\cG$: Directed graph of the traffic network

* $\cV$: Set of nodes of graph $\cG$ which correspond to major intersections and ramps

* $\cE$: Set of edges of graph $\cG$ which correspond to the set of road segments

* $|\mathcal{E}|$: Total number of road segments/edges in the network $\cG$ (i.e. the cardinality of the set $\mathcal{E}$)

* $\mathbf{r}$: Route vector

* $\mathbf{T}$: Time horizon

* $|\mathbf{T}|$: Number of time units (i.e. the cardinality of $\mathbf{T}$)

* $\mathbf{v}_0$: Capacity vector of road segments

* $\mathbf{v}_t$: Volume vector of road segments at time $t$

* $\mathcal{N}$: Set of drivers

* $|\mathcal{N}|$: Number of drivers (i.e. the cardinality of the set $\mathcal{N}$)

* $\mathcal{R}_n$: Set of possible route options for driver $n$

* $\mathcal{R}$: Total set of possible route options for all drivers

* $|\mathcal{R}|$: Number of possible route options (i.e. the cardinality of the set $\mathcal{R}$)

* $\mathcal{I}_n$: Set of possible incentives to offer to driver $n$

* $\mathcal{I}$: Total set of possible incentives to all drivers

* $|\mathcal{I}|$: Number of possible incentives (i.e. the cardinality of the set $\mathcal{I}$)

* $s_i^{\mathbf{r},n}$: Decision parameter indicates whether incentive $i$ is offered to driver $n$ for route $\br$

* $p_{i}^{\mathbf{r},n}$: The probability of acceptance of route $\mathbf{r}$ by driver $n$ given incentive $i$

* $\widehat{T}_{\br}$: The estimate of the travel time for route $\br$ provided by the incentive offering platform

* $T_{\br}$: The exact travel time for route $\br$

* $\boldsymbol{\beta}_{\mathbf{r},t}$: The vector of location of driver that is traveling a route $\br$ at time $t$

* $\eta_i$: The cost of incentive $i$

* $F_{tt}(.)$: Total travel time function

* $\delta_{\ell, t}$: Travel time of link $\ell$ at time $t$

* $\hat{\bv}$: The vector of volume of links at different times in the horizon

* $\hat{v}_{\ell, t}$: The $(|\cE|\times t + \ell)^{th}$ element of vector $\hat{\bv}$ representing the volume of $\ell^{th}$ link at time $t$

* $t_0$: The free flow travel time of the link

* $v$: The traffic volume of the link

* $w$: The practical capacity of the link

* $\bs_n$: The binary decision vector for one driver in which only one element has value of one and it corresponds to the route and the incentive amount that we offer

* $f_{BPR}(.)$: BPR function

* $\bS$: Decision matrix

* $\bR$: The matrix of location of a driver

* $\bP$: Route choice probability matrix

* $\bD$: The matrix of incentive assignment to OD pairs

* $\bq$: The vector of number of drivers for each OD pair

* $\bc$: The vector of cost of incentives assigned to each route

* $\Omega$: Budget

* $\boldsymbol{\omega}$: The vector of free flow travel time of links

* $\ba_{\ell, t}$: The row of matrix $\bA = \bR \bP$ which corresponds to link $\ell$  at time $t$

* $K$: The number of OD pairs

* $e$: An edge of graph $\cG$ which corresponds to a road segments in the traffic network

§ DETAILS OF ALTERNATING DIRECTION METHOD OF MULTIPLIERS (ADMM)

In this section, we explain the details of solving (<ref>) utilizing ADMM method [Boyd et al., 2011].

\begin{align}
& \mathcal{L}(\bga, \bS, \bH, \bW, \bu, \beta)  \nonumber\\
\triangleq & \; F_{tt}(\bga) + \mathbb{I}_{[0, 1]^{(|\cR|\cdot|\cI|) \times |\cN|}}(\bH)+\mathbb{I}_{\mathbb{R}^{+}}(\beta) \nonumber\\
&  + \langle\bla_1, \bS\mathbf{1}-\bu\rangle + \langle\bla_2, \bW^{\intercal}\mathbf{1} - \mathbf{1}\rangle + \langle\bla_3, \bD\bu-\bq\rangle \nonumber\\
&  + \langle\bla_4, \bA\bu - \bga\rangle +  \langle\bLa_5, \bH - \bS\rangle \nonumber\\
&   + \lambda_6(\bc^{\intercal}\bu  + \beta - \Omega) + \langle\bLa_7, \bW - \bS\rangle    \nonumber\\
& + \frac{\rho}{2}||\bS\mathbf{1} - \bu||^{2} + \frac{\rho}{2}||\bW^{\intercal}\mathbf{1}  -\mathbf{1}||^{2} \nonumber\\
&  + \frac{\rho}{2}||\bD\bu - \bq||^{2}  + \frac{\rho}{2}||\bA\bu - \bga||^{2}    \nonumber\\
& + \frac{\rho}{2}||\bH - \bS||^{2} +  \frac{\rho}{2}(\bc^{\intercal}\bu + \beta - \Omega)^{2} \nonumber\\
&    + \frac{\rho}{2}||\bW - \bS||^{2} - \frac{\tilde{\lambda}}{2} \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|}  \bH_{\ell, t}(\bH_{\ell, t}-1) \nonumber\\
\end{align}

be the augmented Lagrangian function of (<ref>) with the set of Lagrange multipliers $\{\bla_1, \bla_2, \dots, \bLa_7\}$ and $\rho > 0$ be the primal penalty parameter. Then, ADMM solves (<ref>) by the following iterative scheme

\begin{align*}
\bu^{t+1} = & \argmin\limits_{\bu} \quad  \langle\bla_1^{t}, \bS^{t+1}\mathbf{1} - \bu\rangle +  \langle\bla_3^{t}, \bD\bu - \bq\rangle   \\
&  + \langle\bla_4^{t}, \bA\bu-\bga\rangle + \lambda_6(\bc^{\intercal}\bu  + \beta - \Omega)  \\
& + \frac{\rho}{2}||\bS^{t+1}\mathbf{1} - \bu||^{2} + \frac{\rho}{2}||\bD\bu - \bq||^2 \\
& + \frac{\rho}{2}||\bA\bu - \bga||^2  +  \frac{\rho}{2}(\bc^{\intercal}\bu + \beta - \Omega)^{2}\\
\bW^{t+1} = & \argmin\limits_{\bW} \quad  \langle\bla_2^{t}, \bW^{\intercal}\mathbf{1} - \mathbf{1}\rangle + \langle\bLa_7^{t}, \bW - \bS^{t+1}\rangle \\
& + \frac{\rho}{2}||\bW^{\intercal}\mathbf{1} - \mathbf{1}||^2  +  \frac{\rho}{2}||\bW - \bS^{t+1}||^2  \\
\bH^{t+1} = & \argmin\limits_{\bH}  \quad \mathbbm{1}({\rho>\tilde{\lambda}})\mathbb{I}_{[0, 1]^{(|\cR|\cdot|\cI|) \times |\cN|}}(\bH)\\
& + \mathbbm{1}({\rho<\tilde{\lambda}})\mathbb{I}_{\{0, 1\}^{(|\cR|\cdot|\cI|) \times |\cN|}}(\bH) \nonumber\\
& + \langle\bLa_5^{t}, \bH- \bS^{t+1}\rangle + \frac{\rho}{2}||\bH-\bS^{t+1}||^{2} \nonumber\\
& - \frac{\tilde{\lambda}}{2} \sum_{\ell=1}^{|\cE|} \sum_{t=1}^{|\mathbf{T}|}  \bH_{\ell, t}(\bH_{\ell, t}-1) \\
\bS^{t+1} = & \argmin\limits_{\bS}  \quad \langle\bla_1^{t}, \bS\mathbf{1}-\bu^{t}\rangle + \langle\bLa_5^{t}, \bH^{t} - \bS\rangle      \\
&  + \langle\bLa_7^{t}, \bW^{t} - \bS\rangle +  \frac{\rho}{2}||\bS\mathbf{1} - \bu^{t}||^{2}  \\
& + \frac{\rho}{2}||\bH^{t}- \bS||^{2} + \frac{\rho}{2}||\bW^{t} - \bS||^{2}  \\
\beta^{t+1} = & \argmin\limits_{\beta} \quad  \mathbb{I}_{+}(\beta) + \lambda_6^{t}\left(\bc^{\intercal}\bu^{t+1} + \beta - \Omega\right) \\
& + \frac{\rho}{2}\left(\bc^{\intercal}\bu^{t+1} + \beta - \Omega\right)^2  \\
\bla^{t+1}_{1} = & \bla^{t}_{1} + \rho \left(\bS^{t+1}\mathbf{1} - \bu^{t+1}\right)   \\
\bla^{t+1}_{2} = & \bla^{t}_{2} + \rho \left(\bW^{t+1 \intercal} \mathbf{1} -\mathbf{1}\right)  \\
\bla^{t+1}_{3} = & \bla^{t}_{3} + \rho \left(\bD \bu^{t+1} - \bq\right)  \\
\bla^{t+1}_{4} = & \bla^{t}_{4} + \rho \left(\bA \bu^{t+1} - \bga^{t+1}\right)  \\
\bLa^{t+1}_{5} = & \bLa^{t}_{5} + \rho \left(\bH^{t+1}  - \bS^{t+1}\right)  \\
\lambda^{t+1}_{6} = & \lambda^{t}_{6} + \rho \left(\bc^{\intercal} \bu^{t+1} + \beta^{t+1} - \Omega\right)  \\
\bLa^{t+1}_{7} = & \bLa^{t}_{7} + \rho \left(\bW^{t+1} - \bS^{t+1}\right)
\end{align*}

The primal update rules can be simplified as

\begin{align*}
\gamma^{t+1}_{\ell, \hat{t}} = & \argmin\limits_{\gamma_{\ell, \hat{t}}}   \gamma_{\ell, \hat{t}} \delta(\gamma_{\ell, \hat{t}}) + \bla^{t}_{4, (\ell, \hat{t})}(\ba_{\ell, \hat{t}} \bu^{t} - \gamma_{\ell, \hat{t}})\\
& + \frac{\rho}{2} (\ba_{\ell, \hat{t}} \bu^{t} - \bga_{\ell, \hat{t}})^{2}, \quad \forall \ell, \forall \hat{t}  \\
\bS^{t+1} = &\; \frac{1}{\rho} (-\bla_1^{t} \mathbf{1}^{\intercal} + \bLa_{5}^{t} + \bLa_{7}^{t} + \rho \bu^{t} \mathbf{1}^{\intercal} + \rho \bH^{t}\\
& + \rho \bW^{t})(\mathbf{1}\mathbf{1}^{\intercal} + 2\bI)^{-1}  \\
\bH^{t+1} =  &\; \mathbbm{1}(\rho>\tilde{\lambda})\Pi {\left(\left(\frac{1}{\rho-\tilde{\lambda}}\right)(\rho\bS^{t} - \bLa_5^{t} - \frac{\tilde{\lambda}}{2})\right)}_{[0, 1]}\\ & +\mathbbm{1}(\rho<\tilde{\lambda})\Pi {\left(\left(\frac{1}{\rho-\tilde{\lambda}}\right)(\rho\bS^{t} - \bLa_5^{t} - \frac{\tilde{\lambda}}{2})\right)}_{\{0, 1\}}  \\
\bW^{t+1} = &\; \frac{1}{\rho} (\bI + \mathbf{1}\mathbf{1}^{\intercal})^{-1} (-\mathbf{1}\bla_2^{t \intercal} - \bLa_7^{t} + \rho \mathbf{1}\mathbf{1}^{\intercal} + \rho \bS^{t+1})  \\
\bu^{t+1} = &\; \frac{1}{\rho} (\bI + \bD^{\intercal} \bD + \bA^{\intercal}\bA + \bc \bc^{\intercal})^{-1} (\bla_1^{t} - \bD^{\intercal} \bla_3^{t} - \bA^{\intercal} \bla_4^{t}\\
& + \rho\bS^{t+1}\mathbf{1} + \rho \bD^{\intercal}\bq + \rho \bA^{\intercal}\bga^{t+1} - \lambda_6 \bc - \beta \rho\bc + \Omega\rho \bc)\\
\beta^{t+1} = & \; \Pi\left( \frac{1}{\rho} (-\lambda_6^{t} -\rho \bc^{t \intercal} \bu^{t} + \rho \Omega)\right)_{\bbR_{+}}
\end{align*}

§ AN EXAMPLE OF THE MODEL AND NOTATIONS

In this section, we present a small example of a network to illustrate our model and notations. Consider the  network

Network example $\cG_1$.

where $\cV = \{\nu_1, \nu_2, \nu_3\}$ is the set of nodes and $\cE = \{e_1, e_2, e_3\}$ is the set edges (roads). Details of the links and attributes are represented in TABLE <ref>. The (origin, destination) pair is ($\nu_1$, $\nu_3$). There are two routes going from origin to destination as illustrated in TABLE <ref>. The time horizon set is $\bT = \{1, 2, 3\}$ and each time is $0.2$ hour. To estimate the location of drivers at each time, we need matrix $\bR \in [0, 1]^{9\times6}$ as follows

\begin{align*}
\bR =
\begin{blockarray}{ccccccc}
\small
& \makecell{t_1\\=1\\ \br_1} & \makecell{t_1\\=1\\ \br_2} & \makecell{t_1\\=2\\ \br_1} & \makecell{t_1\\=2\\ \br_2} & \makecell{t_1\\=3\\ \br_1} & \makecell{t_1\\=3\\ \br_2} \\
\begin{block}{c(cccccc)}
\makecell{t_2=1, e_1} & 1 & 1 & 0 & 0 & 0 & 0 \\
\makecell{t_2=1, e_2} & 0 & 0 & 0 & 0 & 0 & 0 \\
\makecell{t_2=1, e_3} & 0.5 & 0 & 0 & 0 & 0 & 0 \\
\makecell{t_2=2, e_1} & 0 & 0 & 1 & 1 & 0 & 0 \\
\makecell{t_2=2, e_2} & 0 & 1 & 0 & 0 & 0 & 0 \\
\makecell{t_2=2, e_3} & 0.5 & 0 & 0.5 & 0 & 0 & 0 \\
\makecell{t_2=3, e_1} & 0 & 0 & 0 & 0 & 1 & 1 \\
\makecell{t_2=3, e_2} & 0 & 0 & 0 & 1 & 0 & 0 \\
\makecell{t_2=3, e_3} & 0 & 0 & 0.5 & 0 & 0.5 & 0 \\
\end{block}
\end{blockarray}
\end{align*}

where $t_1$ is the entrance time of the driver and $t_2$ is the driver's arrival time at the road. In model (<ref>), the column vector $\boldsymbol{\beta_{\br, t}}$ corresponds to the columns of matrix $\bR.$

(Mile)     Speed
(mph)     Travel time
1|c|$e_1$      5     50     0.1
1|c|$e_2$      10     50     0.2
1|c|$e_3$      5     50     0.1

Set of edges.

$\br$     Graph
1|c|Route 1
$e_1 \rightarrow e_3$      $\br_1=\begin{bmatrix} 1 \\ 0 \\ 1\end{bmatrix}$
1|c|Route 2
$e_2 \rightarrow e_3$     $\br_2 =\begin{bmatrix} 0 \\ 1 \\ 1\end{bmatrix}$

Set of routes.

Assume there are two  drivers in the system and $\cN = \{d_1, d_2\}$. We want to offer rewards from the set $\cI = \{\$0, \$5\}$ to control the traffic. To estimate the probability of choosing routes given an offered incentive at a time, we use matrix $\bP \in [0, 1]^{6\times12}$ when incentive $i$ is offered:

\begin{align*}
\bP_{t_i} =
\begin{blockarray}{ccccc}
& \multicolumn{2}{c}{\;\;\;\text{No incentive}} & \$5 \rightarrow \br_1 & \$5 \rightarrow \br_2  \\
\begin{block}{c(cccc)}
t=1, \br_1 & \;\;\; 0.50 & 0.50 & 0.97 & 0.03 \\
t=1, \br_2 & \;\;\; 0.50 & 0.50 & 0.03 & 0.97 \\
t=2, \br_1 & \;\;\; 0.50 & 0.50 & 0.97 & 0.03 \\
t=2, \br_2 & \;\;\; 0.50 & 0.50 & 0.03  & 0.97 \\
t=3, \br_1 & \;\;\; 0.50 & 0.50 & 0.97 & 0.03 \\
t=3, \br_2 & \;\;\; 0.50 & 0.50 & 0.03  & 0.97 \\
\end{block}
\end{blockarray}
\\,\forall i \in \{1, 2, 3\}
\end{align*}

\begin{align*}
\bP = \begin{bmatrix}
\bP_{t_1} &  \bP_{t_2} &  \bP_{t_3} \\
\end{bmatrix}
\end{align*}

Probability matrices for all three times are equal because the speed is the same in all three times. We compute the probability of choosing route $k$ given that $\$i'$ is offered for route $j'$ by

\begin{equation}
\small
\begin{split}
& \mathbb{P}(\br=k, i=(\$i' \rightarrow \text{route} \; j')) \\
= & \frac{\exp{( -0.086 tt_k + 0.7 i' \mathbb{I}_{k=j'})}}{\exp{(-0.086tt_{j'} + 0.7 i')} + \sum_{j\neq j'} \exp{(-0.086tt_j)}}
\end{split}
\normalsize
\end{equation}

where $tt_j$ is the travel time of route $j$. We use [Xiong et al., 2019] to extract these coefficients.

§ DETAILS OF THE NUMERICAL EXPERIMENTS

2c|2*                                                                                2*Budget     3c|Incentive
2c|                                                                                1c|$0      1c|$2     $10
1|c|8*c]@c@[origin=c]90Penetration Rate     2*25%                           $1000                      1c|7242     1c|191     61
1|c|                                                                                $10000                   1c|7198     1c|14      282
1|c|                                                                                2*50%                           $1000                      1c|7463     1c|414     17
1|c|                                                                                $10000                     1c|6975     1c|0       519
1|c|                                                                                2*75%                           $1000                     1c|6994     1c|500     0
1|c|                                                                                $10000                     1c|6717     1c|0       777
1|c|                                                                                1l|2*100%     $1000                      1c|6994     1c|500     0
1|c|                                                                                1l|                           $10000                     1c|6472     1c|28      994

Distribution of the offered incentives in Experiment I with different penetration rates.

2*Budget     3c|Incentive
$0         $2          $10
$1000           14645        435           13
$10000           12614         1849          630

Distribution of the offered incentives in Experiment II for incentive set $\cI_1$ with penetration rate of 100%.

2*Budget     6c|Incentive
$0      $1      $2      $3      $5      $10
$1000       14509      351      152      30      51      0
$10000       12682     184      305      832      838      252

Distribution of the offered incentives in Experiment II for incentive set $\cI_2$ with penetration rate of 100%.

2*Budget     3c|Incentive
$0         $2          $10
$1000           7720        500           0
$10000           6916         380          924

Distribution of the offered incentives in Experiment III for model (<ref>) with penetration rate of 100%.

2c|2*                                                                                2*Budget     3c|Incentive
2c|                                                                                1c|$0      1c|$2     $10
1|c|8*c]@c@[origin=c]90Penetration Rate     2*25%                           $1000                      1c|8040     1c|100     80
1|c|                                                                                $10000                     1c|7866     1c|115     239
1|c|                                                                                2*50%                           $1000                      1c|7894     1c|285     41
1|c|                                                                                $10000                     1c|7153     1c|94      973
1|c|                                                                                2*75%                           $1000                      1c|7821     1c|374     25
1|c|                                                                                $10000                     1c|6960     1c|335     925
1|c|                                                                                1l|2*100%     $1000                      1c|7780     1c|415     17
1|c|                                                                                1l|                           $10000                     1c|6891     1c|417     912

Distribution of the offered incentives in Experiment III with different penetration rates for model (<ref>), Algorithm <ref>.

2c|2*                                                                                2*Budget     3c|Incentive
2c|                                                                                1c|$0      1c|$2     $10
1|c|8*c]@c@[origin=c]90Penetration Rate     2*25%                           $1000                      1c|7468     1c|110     78
1|c|                                                                                $10000                     1c|7329     1c|9       320
1|c|                                                                                2*50%                           $1000                      1c|7980     1c|175     65
1|c|                                                                                $10000                     1c|7576     1c|59      585
1|c|                                                                                2*75%                           $1000                      1c|7972     1c|185     63
1|c|                                                                                $10000                     1c|7246     1c|78      896
1|c|                                                                                1l|2*100%     $1000                      1c|7896     1c|280     44
1|c|                                                                                1l|                           $10000                     1c|7022     1c|248     950

Distribution of the offered incentives in Experiment III with different penetration rates for model (<ref>), Gurobi.

2*Budget     4c|Penetration Rate
1c|25%      1c|50%      1c|75%      100%
$100                       1c|3.61      1c|3.90      1c|3.93      4.91
$1000                      1c|16.27     1c|24.77     1c|31.52     31.24
$10000                     1c|20.41          1c|41.05     1c|50.09     56.57

Effect of the penetration rate on travel time decrease (hour) in Experiment I.

2*Budget     4c|Penetration Rate
1c|25%       1c|50%       1c|75%       100%
$100                       1c|0.47%     1c|0.51%     1c|0.51%     0.64%
$1000                      1c|2.12%     1c|3.23%     1c|4.11%     4.07%
$10000                     1c|2.66%     1c|5.35%     1c|6.63%     7.37%

Effect of the penetration rate on percentage of travel time decrease in Experiment I.

2*Budget     4c|Penetration Rate
1c|25%       1c|50%       1c|75%       100%
$100                       1c|39.41      1c|43.78      1c|28.91      26.01
$1000                      1c|156.32     1c|215.23     1c|220.24     234.57
$10000                     1c|164.80     1c|255.77     1c|313.13     339.96

Effect of the penetration rate on travel time decrease (hour) in Experiment III, model (<ref>), Algorithm <ref>.

2*Budget     4c|Penetration Rate
1c|25%       1c|50%       1c|75%        100%
$100                       1c|1.31%     1c|1.46%     1c|0.96%      0.87%
$1000                      1c|5.20%     1c|7.16%     1c|7.33%      7.81%
$10000                     1c|5.48%     1c|8.51%     1c|10.42%     11.31%

Effect of the penetration rate on percentage of travel time decrease in Experiment III, model (<ref>), Algorithm <ref>.

2*Budget     4c|Penetration Rate
1c|25%       1c|50%       1c|75%       100%
$100                       1c|70.75      1c|77.67      1c|86.99      88.51
$1000                      1c|201.04     1c|288.96     1c|362.90     402.51
$10000                     1c|215.95     1c|361.06     1c|473.70     536.59

Effect of the penetration rate on travel time decrease (hour) in Experiment III, model (<ref>), Gurobi.

2*Budget     4c|Penetration Rate
1c|25%       1c|50%        1c|75%        100%
$100                       1c|2.35%     1c|2.58%      1c|2.89%      2.95%
$1000                      1c|6.69%     1c|9.62%      1c|12.08%     13.39%
$10000                     1c|7.19%     1c|12.01%     1c|15.76%     17.86%

Effect of the penetration rate on percentage of travel time decrease in Experiment III, model (<ref>), Gurobi.