# Sustainable optimal transport in multilayer networks

Abdullahi Adinoyi Ibrahim<EMAIL_ADDRESS>Max Planck
Institute for Intelligent Systems, Cyber Valley, Tübingen 72076, Germany
Daniela Leite<EMAIL_ADDRESS>Max Planck Institute for
Intelligent Systems, Cyber Valley, Tübingen 72076, Germany  Caterina De Bacco
<EMAIL_ADDRESS>Max Planck Institute for Intelligent
Systems, Cyber Valley, Tübingen 72076, Germany

###### Abstract

Traffic congestion is one of the major challenges faced by the transportation
industry. While this problem carries a high economical and environmental cost,
the need for an efficient design of optimal paths for passengers in multilayer
network infrastructures is imperative. We consider an approach based on
optimal transport theory to route passengers preferably along layers that are
more carbon-efficient than the road, e.g., rails. By analyzing the impact of
this choice on performance, we find that this approach reduces carbon
emissions considerably, compared to shortest-path minimization. Similarly, we
find that this approach distributes traffic more homogeneously, thus
alleviating the risk of traffic congestions. Our results shed light on the
impact of distributing traffic flexibly across layers guided by optimal
transport theory.

##  I Introduction

Traffic congestion is a major problem in the transportation industry, with
significant economic and environmental repercussions. The impacts of the
environmental cost such as carbon emissions and other air pollutants, on
public health can be sizable and need to be properly studied [1]. Combining
different transportation modalities, as in multilayer networks, can mitigate
congestion and thus improve urban sustainability [2]. Modeling traffic
congestion on multilayer networks is crucial to investigate the efficiency and
cost of operating such infrastructures [3]. Addressing this problem requires
extracting what paths passengers take from source to destination, information
that can then be used to analyze traffic patterns. Many route extraction
methods are based on shortest-path minimization [4, 5, 6, 7] or assignment
strategy [8]. However, the shortest paths (i.e. selfish routing) might not
always be the optimal path in a congested network [9, 10, 11], hence the need
for coordinated traffic congestion. In addition, empirical results have shown
that passengers may not always consider the shortest route [12, 13]. While
efforts have been made to go beyond shortest-path minimization using the
cavity method or message-passing algorithms [14, 10, 15, 16, 17], these
approaches are only valid in single-layer networks. In multilayer networks,
several works focus more on analyzing the properties of passenger flows rather
than proposing models to extract trajectories. They consider random walks [4]
or shortest-path optimization [6, 5, 7] to extract flows, thus necessarily
influencing the results of subsequent analysis based on these strategies.
Fewer models have been targeting transport optimization in multilayer networks
[3]. For instance, Ref. [8] proposed a flow-assignment strategy on multilayer
networks, while Ref. [18] developed a recurrent algorithm for communication
networks.

A principled and efficient approach for extracting optimal paths of passengers
in networks is optimal transport (OT) theory [19, 20, 21, 22]. This approach
has been applied recently to multilayer networks [23], where the key idea is
to flexibly tune between different cost functions in each of the different
layers, thus capturing the specificity of each type of infrastructure. For
instance, a road network is more sensitive to traffic congestion than a rail
one, while the infrastructure of a rail network may be more costly to build.
Our work builds from these ideas by adapting this model to study and evaluate
optimal paths on multilayer networks under different scenarios. The goal of
our work is to study the trajectories of optimal paths and compare them with
those extracted from standard approaches relying on shortest-path minimization
to identify key properties that are better optimized if one considers the
multilayer character of the network. Our main contribution is threefold:
First, we consider an optimal transport-based approach to extract optimal
paths for passengers in multilayer networks, contrarily to standard approaches
based on shortest-path minimization. Second, we propose a variant of this OT-
based method that interpolates between OT and shortest-path minimization.
While the extracted paths of the two OT-based models are longer than those
obtained by shortest-path minimization, the rail layer is used by more
passengers. Finally, we show that by using the optimal routes extracted by OT-
based algorithms, passengers are more likely to encounter little or no traffic
while less CO2, leading to a reduced environmental cost. Our empirical results
on synthetic and real data show the need for approaches that exploit the
multilayer nature of multimodal transportation networks.

##  II Optimal transport for traffic distribution in multilayer networks

We denote a multilayer network as a graph denoted as
$G(\\{\mathcal{V}_{\alpha}\\}_{\alpha},\\{\mathcal{E}_{\alpha}\\}_{\alpha},\\{\mathcal{E}_{\alpha}\gamma\\}_{\alpha\gamma})$,
where $\mathcal{V}_{\alpha}$, $\mathcal{E}_{\alpha}$, and
$\mathcal{E}_{\alpha\gamma}$ denote the set of nodes, edges in layer $\alpha$,
and inter-layer edges between layers $\alpha$ and $\gamma$;
$\alpha=1,\dots,L$, where $L$ is the number of layers. We denote the number of
nodes and edges as $N$ and $E$, and we assume that edges have length
$l_{e}>0$, which determines the cost of traveling through them. We consider
the case of a two-layer network, but all results are valid for a higher number
of layers. We denote the two layers as $\alpha,\gamma$ and consider a road
network for $\alpha$ and a rail network for $\gamma$, as explained in detail
in Sec. III. We show an example of this structure in Fig. 1.

$\ \alpha$$\ \gamma$ Figure 1: Multilayer structure with $N_{\alpha}=15$ and
$N_{\gamma}=4$. The network edges are represented by continuous lines (magenta
and brown) and the two-edge path by dashed lines. The thicker magenta nodes
represent stations belonging to both layers.

We consider passengers traveling through the networks and distinguish them by
their origin and destination (traffic demands) stations $(o_{i},t_{i})$, where
$o_{i},t_{i}\in\mathcal{V}=\cup_{\alpha}V_{\alpha}$. We denote as
$\mathcal{S}=\left\\{(o_{i},t_{i})\right\\}$ the set of all origin-destination
pairs, and $|\mathcal{S}|=M$ denotes their number.

We briefly describe the model of [23] to find optimal paths in multilayer
networks using optimal transport theory. It considers two main quantities on
network edges: fluxes $F_{e}$ of passengers traveling through edge $e$ and
conductivities $\mu_{e}$ determining $F_{e}$ passing through an edge $e$. To
keep track of the various routes that passengers have, a multi-commodity
approach is considered [22, 24] in which passengers are distinguished based on
their entry station $i\in\mathcal{S}$. With this approach, the flux $F_{e}$ is
an $M$-dimensional vector, where entries $F_{e}^{i}$ denote the flux of
passengers of type $i$ traveling on edge $e$. We assume the fluxes are
determined by pressure potentials $p_{u}^{i}$ and $p_{v}^{i}$ defined on nodes
as follows:

$F_{e}^{i}:=\frac{\mu_{e}}{l_{e}}\left({p_{u}^{i}-p_{v}^{i}}\right),\quad
e=(u,v)\quad,$ (1)

where $l_{e}$ is the length of edge $e$. Kirchhoff’s law is imposed on network
nodes to properly enforce mass conservation. Finally, the dynamics assumes
that the conductivity $\mu_{e}$ depends on flux $F_{e}$ as follows:

$\dot{\mu}_{e}=\mu_{e}^{\beta_{q_{e}}}\,\frac{\sum_{i\in\mathcal{S}}(p_{u}^{i}-p_{v}^{i})^{2}}{l_{e}^{2}}-\mu_{e},\quad\forall
e\in\mathcal{E}\quad,$ (2)

where $q_{e}$ encodes the layer to which the edge $e$ belongs. The parameter
$0<\beta_{q_{e}}<2$ determines the type of optimal transport problem one aims
to solve: $0<\beta_{q_{e}}<1$ discourage traffic congestion,
$1<\beta_{q_{e}}<2$ encourages path consolidation into few highways, while
$\beta_{q_{e}}=1$ is shortest-path-like. Interpreting the conductivities as
quantities proportional to the size of an edge, this dynamics enforces a
feedback mechanism such that the edge size increases if the flux through that
edge increases, and it decreases otherwise.

It can be shown [22, 23] that the stationary solutions of Eq. 2 minimize the
multilayer transport cost function:

$J_{\beta}=\sum_{\alpha=1}^{L}\sum_{e\in\mathcal{E}_{\alpha}}l_{e}||F_{e}||_{2}^{\Gamma(\beta_{\alpha})}\quad,$
(3)

where $\Gamma(\beta_{\alpha})=2(2-\beta_{\alpha})/(3-\beta_{\alpha})$ for all
$\alpha$, and the 2-norm is calculated over the $M$ entries of each $F_{e}$.
Intuitively, solving the system of Eqs. 1 and 2 and Kirchhoff’s law is
equivalent to finding the optimal trajectories of passengers in a multilayer
network, where optimality is given with respect to the transport cost in Eq.
3. We refer to this OT-based algorithm as MultiOT.

###  II.1 MultiOTsp: interpolating between OT and shortest-path minimization

The paths extracted by MultiOT will encourage path consolidation along with
the rail network and traffic minimization on the road network. Empirically, we
observe that this model tends to distribute passengers of the same type (i.e.
same origin and destination) along various routes, as shown in Fig. 2. While
most of these passengers take the shortest among these routes, some distribute
on longer ones to prevent traffic congestion. This suggests an alternative
algorithm that interpolates between MultiOT and shortest-path minimization to
select only the main relevant routes for each origin-destination pair among
those extracted by MultiOT. This can be done by inputting the solution of
MultiOT for each passenger type $i$ into a weighted shortest-path algorithm
with edge weights defined as:

$w_{e}=\frac{l_{e}}{|F_{e}^{i}|}\quad,$ (4)

where the fluxes $F_{e}$ are those extracted from MultiOT. All the passengers
of type $i$ are then routed along the output path. We call this algorithm
MultiOTsp and show a pseudo-code in Algorithm 1. The advantage of this weight
function $w_{e}$ is that weakly used edges are consolidated on others that are
on the optimal path (according to OT) of many passengers. These edges thus
become more desirable when designing and individual “consensus” OT-based path
that takes into account both path length and optimal fluxes. The paths
selected by MultiOTsp rely strongly on how the fluxes are selected in the
first place to determine the weights $w_{e}$. As the $F_{e}$ are calculated by
considering all the passengers simultaneously (using MultiOT), the final
optimal trajectories of MultiOTsp are significantly distinct from those
obtained by shortest-path minimization, which are independent from the
surrounding environment. We show an example of this in Fig. 2.

Algorithm 1 MultiOTsp

1:Graph $G(V,E,W)$, set $\mathcal{S}$ of origin-destination pairs,
$\beta=(\beta_{1},\dots,\beta_{L})$

2:Fluxes $\left\\{F_{e}\right\\}_{e}$

3:function MultiOTsp($G,\mathcal{S},\beta$)

4: $\left\\{F_{e}\right\\}_{e}\leftarrow$ MultiOT($G,\mathcal{S},\beta$)

5: for $i=1,\dots,M$ do

6: $\left\\{F_{e}^{i}\right\\}_{e}\leftarrow$ weighted
Dijkstra$(G,\mathcal{S},w)$ with $w_{e}=l_{e}/|F_{e}^{i}|$

7: end for

8: $F_{e}=(F_{e}^{1},\dots,F_{e}^{M})$, $\forall\,e$

9:end function

In the following, we study the trajectories of optimal paths extracted by the
three approaches: MultiOT, MultiOTsp, and shortest-path minimization (SP). We
use the implementation in [25] for MultiOT, while for SP we use the Dijkstra
algorithm [26].

Figure 2:  Example trajectories. We show the trajectory of one type of
passenger (black edges) whose origin and destination stations are the green
and magenta nodes, respectively. We also highlight the total fluxes on edges,
solutions of the OT problem including all other passengers, for a total of
$M=300$. Blue and red edges denote road ($\alpha$) and rail ($\gamma$) layers,
respectively. Edge widths are proportional to the amount of passengers
travelling through an edge. The exact width of the black edge has been either
increased (for SP) or reduced (for OT-based methods) in order to distinguish
the flux of this type of passengers from the overall trajectories. Origin-
destination pairs have been selected so that $80\%$ of the passengers are
directed towards a central node; $\beta_{\alpha}=0.5$ and
$\beta_{\gamma}=1.9$.

##  III Empirical results

To investigate the relevant properties of the optimal paths extracted by the
different algorithms, we simulate a variety of realistic traffic scenarios.
Specifically, we generate a dataset of synthetic 2-layer planar networks,
where $\alpha$ simulates a road network and $\gamma$ simulates a rail network
(e.g., a tram). The layer $\alpha$ is constructed first, by randomly placing
$N_{\alpha}$ nodes in $[0,1]\times[0,1]$, and extracting its Delaunay
triangulation [27]. We then select among them a subset of $N_{\gamma}$ nodes
to build the layer $\gamma$ with an analogous procedure, thus ensuring that
the two layers are connected. In total, in this construction the multilayer
network has $N=N_{\alpha}$ nodes and resembles the situation in which all the
stations in the second layer also have access to the road network. Notice that
other constructions are possible, but this choice does not impact the validity
of our model. In our simulations, we set $N_{\alpha}=300$ and $N_{\gamma}=60$.
We extract 20 different networks and $100$ random samples of origin-
destination pairs for each of them, for a total of 2000 realizations for each
parameters’ configuration. With this, we aim at capturing different
transportation scenarios in the two layers, as rail networks are less subject
to traffic congestion but more costly to build, while we can state the
opposite for road networks. MultiOT (and thus MultiOTsp) can capture these
differences by suitably tuning $\beta$ in each layer: to discourage traffic
congestion in the road layer, we set $\beta_{\alpha}=0.5$ and vary
$\beta_{\gamma}$ in $0<\beta_{\gamma}<2$ to study various scenarios. In
realistic scenarios, passengers have different origins and destinations; see
Fig. 2 for an example. As we may expect in many urban scenarios that the most
frequent destination is located in the city center, we assign to each
passenger type its destination by default to be a central node. Then, to
explore alternative scenarios where destinations are more heterogeneous, we
consider a rewiring probability $p=[0,1]$ to rewire its destination at random.
Specifically, for each passenger type, we rewire its destination to a random
node with probability $p$. Hence, $p=0.0$ corresponds to having a monocentric
destination where all passengers move towards a central node and $p=1.0$
corresponds to selecting all passengers’ destinations at random. We consider
$p=\\{0.2,0.5,0.8\\}$, but we show results for $p=0.5$, as the qualitative
behavior is similar to that for the others, see Appendix A. These settings
exhibit three important properties of the OT-based algorithms.

##  IV Longer lengths but higher rail network usage

Shortest-path optimization is utilized to minimize the total path length taken
by passengers, hence we expect MultiOT and MultiOTsp to underperform SP on
this task. In fact, the performance of OT-based algorithms is expected to
decrease as $\beta_{\gamma}$ increases, as shown in Fig. 3 by the average path
length
$\langle\,\l\,\rangle=\frac{1}{M}\sum_{e\in\mathcal{E}}l_{e}||F_{e}||_{1}$
over the one obtained from a shortest-path algorithm.

Figure 3: Average total path length ratio. We show the ratio of the average
total path length to the one extracted from SP. We set $p=0.5$,
$\beta_{\alpha}=0.5$ and vary $0<\beta_{\gamma}<2$. The results are averaged
over 20 different network realizations with 100 randomly selected origin-
destination pairs for each network realization. The markers and error bars are
averages and standard deviations.

This is expected given that higher $\beta_{\gamma}$ encourages more traffic to
be routed towards the rail network at the cost of increased distance to cover,
as the rail network has fewer and more distant nodes to reach than a road
network. We then measure how passengers are distributed in the two layers by
defining a coupling coefficient, a known concept to describe how well two
layers are linked [6]. We define:

$\lambda=\frac{1}{M}\sum_{i\in\mathcal{S}}\left(\frac{\sum_{e\in\mathcal{E}_{\gamma}}|F_{e}^{i}|}{\sum_{e\in\mathcal{E}_{\alpha}\cup\mathcal{E}_{\gamma}}|F_{e}^{i}|}\right)\quad,$
(5)

where the numerator inside the parentheses contains only the flux in the rail
layer so we can distinguish how many passenger types effectively use that
layer in their trajectories. This definition is valid for two-layer networks,
such as the empirical networks studied here. However, one can appropriately
generalize it for networks with more than two layers. The usage of the rail
layer increases monotonically for both OT-based algorithms, as shown in Fig.
4, with MultiOTsp reaching higher usage values. This suggests that the
shortest-path routes selected from the possible paths output by MultiOT are
composed of a significant amount of rail edges. This also shows that the raw
solutions output of MultiOT consider paths more distributed across the road
layer, as qualitatively observed in Fig. 2.

Figure 4: Coupling between layers. We show the coupling coefficient as defined
in Eq. 5. All other settings remain the same as in Fig. 3.

##  V Multilayer optimization can decrease carbon consumption

As more passengers take longer paths while being encouraged to use the rail
network, they also consume less carbon per unit of length. The question is
whether the increased length can be properly compensated by the decrease in
carbon consumption. We tested this on the same paths extracted to plot Fig. 3
by measuring the average $CO_{2}$ consumption per passenger as:

$\langle\,CO_{2}\,\rangle=\frac{1}{M}\sum_{e\in\mathcal{E}}r_{q_{e}}l_{e}||F_{e}||_{1}\quad,$
(6)

where $r_{\alpha}$ is the carbon emission rate in layer $\alpha$. This has a
dimension of unit of mass (e.g. $g$) per passenger per unit of length (e.g.
pkm). For instance, a bus on average generates $101.87g/pkm$ [28] while a
train generates $28.39g/pkm$ [29]. Hence, defining $r_{\alpha}$ as the rate of
the road layer and considering buses traveling on it, we can set
$r_{\gamma}=28.39r_{\alpha}/101.87=0.28\,r_{\alpha}$. These values can be
changed accordingly with more specific values if a traffic manager has precise
statistics of vehicles’ types traveling on the network. By leveraging optimal
transport with a bias towards shortest paths, MultiOTsp is able to decrease
carbon consumption the most compared to SP, measured by the ratio of its
$\langle\,CO_{2}\,\rangle$ over that produced by SP. A minimum is reached for
$1.1\leq\beta_{\gamma}\leq 1.5$ where MultiOTsp produces $25\%$ fewer
emissions than a shortest-path routing algorithm, as shown in Fig. 5. This
important result is a consequence of flexibly tuning the cost to be optimized
in each layer, as allowed by $\beta$ in Eq. 3. In particular,
$\beta_{\gamma}>1$ encourages paths to consolidate on the rail layer, while
$\beta_{\alpha}=0.5$ controls for traffic congestion on the road layer. The
fact that the minimum consumption of MultiOTsp has not been realized at the
highest value $\beta_{\gamma}=1.9$, where the paths are consolidated into the
fewest rail routes, further suggests that there is a trade-off between keeping
the path lengths short while directing more passengers towards the rail layer.
In fact, at $\beta_{\gamma}=1.9$, as the number of passengers redirected
towards the second layer increases, they also have to take longer routes, thus
emitting more carbon. A value of $\beta=1.3$ results in a nice tradeoff
between these two competing behaviors in terms of carbon emission. On the
contrary, MultiOT shows a monotonic decreasing behavior with a minimum reached
at $\beta=1.9$, but still higher than that emitted by SP. This is a
consequence of the higher number of possible paths that passengers can take as
routed by MultiOT, which are by default longer than those obtained by
MultiOTsp and use more edges of the road layer. As a consequence, the longer
length does not seem to justify the higher usage of the rail layer.

Figure 5: Carbon emission ratio. We show the ratio of the average carbon
emissions as defined in Eq. 6 to the one obtained by SP. All other settings
remain the same as in Fig. 3.

##  VI Traffic congestion

All the results of the previous section were interpreted with the assumption
that the flow of passengers is regular, even on high-traffic edges. Instead,
if we account for traffic slowing down the flow on edges with a high density
of travelers, those vehicles emit more carbon while they keep their engines on
longer. The routes suggested by MultiOT are less sensitive to this, hence we
also expect a lower carbon emission than that shown in Fig. 5 when accounting
for traffic. We thus measure traffic load on edges as

$\displaystyle T_{e}=\frac{1}{n}||F_{e}||_{1}\quad,$ (7)

where $n$ is the total number of passengers and measures the Gini coefficient
Gini$(T_{e})\in[0,1]$ as a global network metric of inequality of how traffic
is distributed on the network [30], with a Gini close to 1 meaning high
inequality in flow assignment along edges. As the road layer is the one more
sensitive to potential traffic bottlenecks, we consider only the traffic on
road edges and denote with Gini$(T_{e}^{(\alpha)})$ the Gini coefficient
calculated using only $e\in\mathcal{E}_{\alpha}$. As expected, MultiOT has
more balanced traffic than the other two algorithms, as shown in Fig. 6. While
congestion increases with $\beta_{\gamma}$, even at the maximum $\beta=1.9$
the Gini coefficient is lower than that of SP. The reason for this increase is
that paths consolidate more on those fewer road edges that allow a connection
to the rail layer, as can be seen on the example optimal routes in Fig. 6, a
behaviour also observed in previous studies [31, 32]. This is exacerbated in
MultiOTsp, as one can notice that central road edges are overly trafficked
when many passengers exit the rail to reach the final destination in the
center. This also causes the Gini coefficient of MultiOTsp to be higher than
that of SP. In other words, few central edges cause most of the traffic for
MultiOTsp. This can be partially alleviated by increasing $p$ towards 1 as
fewer destinations are directed towards the network center, although this may
become an unrealistic assumption in urban scenarios. Alternatively, one can
simply add rail stations in the center, so that passengers do not have to
commute one extra mile to reach their final destination, a scenario that we
explore below in the case of a real network.

Figure 6: Traffic distribution. Top) Gini coefficient of the traffic on the
road layer $\alpha$. Bottom) optimal trajectories. Edge widths are
proportional to $||F_{e}||_{1}$ averaged over 10 samples of origin-destination
configurations; $p=0.5$, $\beta_{\alpha}=0.5$. When $\beta_{\gamma}>1$, OT-
based methods are consolidated into fewer edges in layer $\gamma$. The SP on
the other hand is not affected by this parameter, showing flows of passengers
on more edges compared to OT-based which consolidate into fewer edges. Blue
and red edges correspond to road ($\alpha$) and rail ($\gamma$) layers,
respectively.

##  VII Real multilayer network

Next, we examine these properties on a real 2-layer network of the city of
Bordeaux [33], where $\alpha$ and $\gamma$ represent the bus and tram
networks, respectively. Similar to the synthetic network, we compare the
performances of OT-based algorithms with SP on this network. We set $p=0.2$ to
consider the situation where the majority of passengers are directed towards
the city center, a central node coinciding with a tram station, and we extract
100 realizations of origin-destination pairs.

The tram in this city travels through its own reserved lanes, independently
from other vehicles. Hence, although the two layers are physically located
next to each other, there is no mixing of fluxes from the two layers on edges.
This may differ in other real situations. While in principle our model is best
suited for independent usage of the space by the various layers (e.g. road and
subway or the case studied here), the results shown here may still apply if we
assume that the physical presence of the tram only marginally impacts the
traffic in the road layer, compared to other vehicles. Specifically, the
combination of high enough capacity (number of passengers that can fit into a
tram) and lower frequency of trams than other vehicles on the road may allow
us to assume independence between the tram and road layer. In fact, while the
tram may have many passengers traveling at any given time along an edge, these
are all entering inside the same wagons. Thus, the space occupied by the tram
is limited by its physical shape. Instead, the same amount of passengers would
need to distribute in many different cars, thus occupying much more space,
potentially creating congestion. In general, in situations in which trams have
reserved lanes that cars cannot enter, our treatment applies without any
further assumption.

We find that MultiOTsp produces $25\%$ fewer carbon emissions than SP for
$1.1\leq\beta_{\gamma}\leq 1.5$, as shown in Fig. 7, similar to what observed
on synthetic networks. MultiOT has a minimum at $\beta_{\gamma}=1.9$, but the
emissions are higher than SP. We argue that also in this case this is due to
the assumption that the flow of vehicles is smoothly moving, with no traffic
congestion causing velocity to decrease, and thus causing emissions to
increase nearby traffic bottlenecks. To assess this hypothesis, we investigate
the distribution of fluxes on the road layer by measuring the traffic $T_{e}$
on edges in layer $\alpha$. We find that indeed MultiOT has path trajectories
more homogeneously distributed across the road layer as measured by the Gini
coefficient plotted in Fig. 8 along with an example solution, potentially
lowering the number of traffic jams. As we can see from an example solution in
the same figure, the two OT-based variants distribute passengers in higher
amounts along the tram network, thus lowering the road’s usage, while SP makes
use of the tram mainly in the vicinity of the central node. We can further
notice how MultiOT uses the road with higher intensity than MultiOTsp, but the
road edges have less traffic than those obtained by SP. As for MultiOTsp, the
road edges with the most traffic are those nearby tram stations, like those in
the upper left corner in the figure.

Figure 7: Carbon emissions on the Bordeaux network. Ratio of the average
carbon emission over that of SP. Here, we set $p=0.2$, thus favoring
monocentric destinations. The grey-dashed line shows the minimum value
obtained by MultiOTsp, corresponding to $\beta_{\gamma}=1.3$. Inset: zoom-in
of MultiOTsp values for $1.1\leq\beta_{\gamma}\leq 1.9$.

To better quantify the potential impact of traffic congestion as a proxy for
the potential increase in CO2 emission, we consider a measure of transport
cost used before in similar problems [14, 15] and defined as

$J_{\alpha}=\sum_{e\in\mathcal{E_{\alpha}}}l_{e}||F_{e}||_{1}^{2}\quad,$ (8)

where the exponent $2$ discourage traffic. This should not be confused with
the definition of Eq. 3, which is the one used to extract the optimal paths in
our model, i.e. to solve the OT problem. Specifically, Eq. 3 uses different
$\beta$ for edges in different layers. In particular, it allows us to
encourage both path consolidation in one layer and path distribution in
another. Instead, Eq. 8 only discourage traffic, as the exponent is greater
than 1. In addition, Eq. 3 considers the norm-2 of the passengers’ flows,
while Eq. 8 considers norm-1. While the latter is more intuitive, as it is the
total number of passengers traveling along an edge, the former admits rigorous
theoretical guarantees for OT to converge to an optimal solution. This does
not apply to a cost function using norm-1; see [24] for a detailed discussion.

As seen in Fig. 9, both OT-based algorithms outperform SP as $\beta_{\gamma}$
increases, meaning that passengers traveling on paths generated by the OT-
based algorithms will generally record less road traffic congestions compared
with the paths extracted by SP. Assuming that velocity decreases along
congested edges, we conjecture that this would result in MultiOT having lower
carbon emissions than SP.

Figure 8:  Traffic distribution on the Bordeaux network. We set $p=0.2$ and
$\beta_{\alpha}=0.5$. Left) Gini coefficient calculated using the traffic on
the road network. Error bar (standard deviation) are smaller than marker size.
Center-Right) Traffic distribution for each of the algorithms with
$\beta_{\gamma}=1.9$. Red and blue edges denote tram and road layers,
respectively. Node and edge sizes are proportional to the amount of passengers
traveling through them. The results are averaged over 100 samples of origin-
destination pairs. Figure 9: Transport cost on the road layer $\alpha$ of
Bordeaux network. The cost is defined as in Eq. 8; here $p=0.2$.

##  VIII Discussion & Conclusion

Designing and extracting optimal passenger flows in a transportation network
is crucial for reducing traffic congestion and environmental costs. Methods
based on shortest-path optimization are optimal in terms of reducing the
average shortest path length to reach destination, but they may fail in terms
of other relevant transportation metrics. In addition, passengers do not
always follow the shortest route [12], hence the need for alternative
approaches to extract path trajectories and investigate their properties in
multilayer networks. We present two models based on optimal transport theory
that can flexibly tune the amount of traffic routed in the different layers to
encourage usage of rail networks while reducing traffic on the road. As a
result, optimal trajectories extracted with these methods significantly
decrease the amount of carbon emissions compared to shortest-path
minimization, while also being more robust to traffic congestions. In
particular, we found that MultiOTsp, by interpolating between optimal
transport and shortest-path minimization, can achieve the lowest amount of
carbon emissions under the hypothesis of smooth flow of passengers in a
network. Instead, MultiOT, based purely on optimal transport, distributes
paths more homogeneously, thus being potentially more robust against increased
carbon emissions when accounting for passengers’ flow slowing down along
traffic bottlenecks. This can be tested quantitatively in real scenarios by
having access to empirical data of different velocities during traffic
congestion, along with detailed velocity limits imposed by regulations in
different parts of the network. One could potentially compare the theoretical
results with the empirical ones observed from real data as in [34].
In general, we show that models based on optimal transport can be used to
design optimal routes for passengers in a multilayer network, and we
investigate scenarios beyond those obtained by using standard shortest-path
algorithms. In this work, we assumed fixed origin-destination pairs, but one
can further generalize this analysis by considering dynamical traffic demands
that change in time. This would require suitably adapting the models studied
in this work to account for this, for instance borrowing ideas from [35, 36,
37, 38, 39]. Similarly, we did not explore here the possibility of traffic
diversions due to road blockages or changing conditions in the network
structure [40, 41, 42, 43]. Studying the robustness of the methods
investigated in this work to these scenarios would be an interesting subject
for future work. Finally, it would be interesting to investigate more complex
scenarios with more than two layers, possibly on a larger scale than that of a
unique urban scenario. To facilitate future analysis, we provide an open
source implementation of our code at [25].

### Acknowledgments

The authors thank the International Max Planck Research School for Intelligent
Systems (IMPRS-IS) for supporting A.A.I and D.L.

##  Appendix A Additional results varying $p$

We set $p=\\{0.2,0.5,0.8\\}$ to capture different traffic demand scenarios,
where $p=0.2$ and $p=0.8$ correspond to having the majority and minority of
the passengers with a monocentric destination. We show in Fig. 10 the
performance of the algorithms in terms of the same metrics investigated in the
main manuscript. All displayed results have the same settings described in
Sec. III.

Figure 10: Additional results on synthetic data for varying $p$. The results
are averages and standard deviations over $20$ different network realizations
with $100$ independent samples of origin-destination pairs on each network
realization. Other parameters used here are $\beta_{\alpha}=0.5$,
$N_{\alpha}=300$, $N_{\gamma}=60$.

## References

  * Landrigan  _et al._ [2018] P. J. Landrigan, R. Fuller, N. J. Acosta, O. Adeyi, R. Arnold, A. B. Baldé, R. Bertollini, S. Bose-O’Reilly, J. I. Boufford, P. N. Breysse, _et al._ , The lancet 391, 462 (2018).
  * Orozco  _et al._ [2021] L. G. N. Orozco, L. Alessandretti, M. Saberi, M. Szell, and F. Battiston, arXiv preprint arXiv:2111.02152 (2021).
  * Wu  _et al._ [2020] J. Wu, C. Pu, L. Li, and G. Cao, Digital Communications and Networks 6, 58 (2020).
  * De Domenico  _et al._ [2014] M. De Domenico, A. Solé-Ribalta, S. Gómez, and A. Arenas, Proceedings of the National Academy of Sciences 111, 8351 (2014).
  * Lampo  _et al._ [2021] A. Lampo, J. Borge-Holthoefer, S. Gómez, and A. Solé-Ribalta, Physical Review Research 3, 013267 (2021).
  * Morris and Barthelemy [2012] R. G. Morris and M. Barthelemy, Physical review letters 109, 128703 (2012).
  * Solé-Ribalta  _et al._ [2016] A. Solé-Ribalta, S. Gómez, and A. Arenas, Physical review letters 116, 108701 (2016).
  * Gao  _et al._ [2019] L. Gao, P. Shu, M. Tang, W. Wang, and H. Gao, Physical Review E 100, 012310 (2019).
  * Lazar  _et al._ [2020] D. A. Lazar, S. Coogan, and R. Pedarsani, IEEE Transactions on Automatic Control 66, 2664 (2020).
  * Po  _et al._ [2021] H. F. Po, C. H. Yeung, and D. Saad, Physical Review E 103, 022306 (2021).
  * Mehr and Horowitz [2018] N. Mehr and R. Horowitz, in  _2018 IEEE Conference on Decision and Control (CDC)_ (IEEE, 2018) pp. 1788–1793.
  * Quercia  _et al._ [2014] D. Quercia, R. Schifanella, and L. M. Aiello, in  _Proceedings of the 25th ACM conference on Hypertext and social media_ (2014) pp. 116–125.
  * Zheng  _et al._ [2018] Z. Zheng, Z. Huang, F. Zhang, and P. Wang, EPJ Data Science 7, 1 (2018).
  * Yeung and Saad [2012] C. H. Yeung and D. Saad, Physical Review Letters 108, 208701 (2012).
  * Yeung  _et al._ [2013] C. H. Yeung, D. Saad, and K. M. Wong, Proceedings of the National Academy of Sciences 110, 13717 (2013).
  * Altarelli  _et al._ [2015] F. Altarelli, A. Braunstein, L. Dall’Asta, C. De Bacco, and S. Franz, PloS one 10, e0145222 (2015).
  * De Bacco  _et al._ [2014] C. De Bacco, S. Franz, D. Saad, and C. H. Yeung, Journal of Statistical Mechanics: Theory and Experiment 2014, P07009 (2014).
  * Zhou  _et al._ [2013] J. Zhou, G. Yan, and C.-H. Lai, EPL (Europhysics Letters) 102, 28002 (2013).
  * Bonifaci  _et al._ [2012] V. Bonifaci, K. Mehlhorn, and G. Varma, Journal of theoretical biology 309, 121 (2012).
  * Baptista  _et al._ [2020] D. Baptista, D. Leite, E. Facca, M. Putti, and C. De Bacco, Scientific reports 10, 1 (2020).
  * Baptista and De Bacco [2020] D. Baptista and C. De Bacco, arXiv preprint arXiv:2012.12758 (2020).
  * Lonardi  _et al._ [2021a] A. Lonardi, E. Facca, M. Putti, and C. De Bacco, Physical Review Research 3, 043010 (2021a).
  * Ibrahim  _et al._ [2021] A. A. Ibrahim, A. Lonardi, and C. D. Bacco, Algorithms 14, 189 (2021).
  * Lonardi  _et al._ [2022] A. Lonardi, M. Putti, and C. De Bacco, Scientific reports 12, 1 (2022).
  * [25] https://github.com/cdebacco/multiot.
  * Dijkstra  _et al._ [1959] E. W. Dijkstra _et al._ , Numerische mathematik 1, 269 (1959).
  * Guibas and Stolfi [1985] L. Guibas and J. Stolfi, ACM transactions on graphics (TOG) 4, 74 (1985).
  * DEFRA [2017] DEFRA, Department for environment, food and rural affairs. greenhouse gas reporting: Conversion factors 2017 (2017).
  * EEA [2017] EEA, Energy efficiency and specific co2 emissions (2017).
  * Dixon  _et al._ [1987] P. M. Dixon, J. Weiner, T. Mitchell-Olds, and R. Woodley, Ecology 68, 1548 (1987).
  * Chodrow  _et al._ [2016] P. S. Chodrow, Z. Al-Awwad, S. Jiang, and M. C. González, PloS one 11, e0161738 (2016).
  * Strano  _et al._ [2015] E. Strano, S. Shai, S. Dobson, and M. Barthelemy, Journal of The Royal Society Interface 12, 20150651 (2015).
  * Kujala  _et al._ [2018] R. Kujala, C. Weckström, R. K. Darst, M. N. Mladenović, and J. Saramäki, Scientific data 5, 1 (2018).
  * Taillanter and Barthelemy [2021] E. Taillanter and M. Barthelemy, Journal of the Royal Society Interface 18, 20210391 (2021).
  * Yeung [2019] C. H. Yeung, Physical Review E 99, 042123 (2019).
  * Lonardi  _et al._ [2021b] A. Lonardi, E. Facca, M. Putti, and C. De Bacco, arXiv preprint arXiv:2112.10620 (2021b).
  * Hu and Cai [2013] D. Hu and D. Cai, Phys. Rev. Lett. 111, 138701 (2013).
  * Corson [2010] F. Corson, Phys. Rev. Lett. 104, 048703 (2010).
  * Hu  _et al._ [2012] D. Hu, D. Cai, and A. V. Rangan, PLoS ONE 7, 1 (2012).
  * Tai and Yeung [2021] T. S. Tai and C. H. Yeung, Physical Review E 104, 024311 (2021).
  * Kleineberg  _et al._ [2017] K.-K. Kleineberg, L. Buzna, F. Papadopoulos, M. Boguñá, and M. Á. Serrano, Physical review letters 118, 218301 (2017).
  * Wandelt  _et al._ [2018] S. Wandelt, X. Sun, D. Feng, M. Zanin, and S. Havlin, Scientific reports 8, 1 (2018).
  * Wang  _et al._ [2018] Z. Wang, D. Zhou, and Y. Hu, Physical Review E 97, 032306 (2018).