an electric blackout should be promptly signalled to warn vehicles entering the affected area;
the same may be said for emergencies due to floods.

According to the circumstances, city requirements may: (i) frequently change on a temporal basis, and (ii) be relaxed or tightened to meet changing Service Level Agreements (SLAs). In technical terms, approaching the complex needs of a city means fast supporting the elicitation of new service requirements and satisfying the ever-changing non-functional ones. This section discusses these requirement classes and other helpful general features a city platform should be equipped with.

Moreover, city complexity requires factorising the infrastructures and related behaviours
in a unifying meta-model to reduce the platform complexity by focusing on a specific set of
functional and architectural choices and quickly extending it to address new and changing
needs. By harnessing the power of abstraction, these tasks can be eased by adopting an
underlying unified meta-model. The considerations reported in Chapter 4, motivated us to
adopt a graph-based meta-model, as will be further detailed below.

5.1.1 Complex network modelling of city infrastructures

Natural and human-engineered systems can be represented as a set of *complex time-changing*
entities that interact with one another through *multiple complex time-changing* relationships
(both entities and relationships can appear, disappear, and change their attributes over
time), that originate *complex networks*.

We consider as complex network a graph that has numerous non-trivial structural features (e.g., multiple attributes on edges and nodes, multiple edges between the same pair of nodes, attributes that change value over time, etc.) and can represent different real-world systems such as city, social, biological, computer, and technological networks [215, 216]. Moreover, complex networks can also include other types of features such as (i) multiple layers of semantically-different interactions between the modelled entities, and (ii) multiple intra-layer sub-networks, each corresponding to a group of entities with specific properties (e.g., entities within a delimited geographical area, entities belonging to a specific organisation, etc.).

PROMENADE has been designed to perceive city networks as *dynamic complex networks* in which *nodes* (vertices) represent different homogeneous or heterogeneous city entities, and *links* (edges) correspond to the diverse relationships existing between these entities. Although a platform centred around the concept of complex network can be exploited to monitor and study numerous *complex dynamic systems* in disparate contexts, we mostly focused on urban technical networks (and specifically transportation), which can be organised in three levels of abstraction (see Fig. 5.1):

* At level 0, the bare physical space, which is represented by different (interdependent) city complex systems, like, for instance, transportation infrastructures (road networks, buses, trams, and metros), telecommunication/phone networks, electric grids, person networks, etc.
* At level 1, the deployed IoT, Edge, Fog, and Cloud nodes, which are responsible for providing the necessary technological resources to collect, transfer, store, and process data coming from multiple real-time sources integrated with the complex systems from level 0, at high-frequency, fine-resolution, and large spatio-temporal scales [217, 218].
* At level 2, complex network models and the different application services for knowledge extraction, data visualisation, machine learning, etc.

* At level 1, the deployed IoT, Edge, Fog, and Cloud nodes, which are responsible for providing the necessary technological resources to collect, transfer, store, and process data coming from multiple real-time sources integrated with the complex systems from level 0, at high-frequency, fine-resolution, and large spatio-temporal scales [217, 218].

* At level 2, complex network models and the different application services for knowledge extraction, data visualisation, machine learning, etc.

At the topmost level (Fig. 5.1), the city-system under analysis is abstracted as a complex multi-modal, multi-layer attributed time-varying network where nodes can have a hierarchical organisation (forming multiple physical or logical layers). These different layers can be utilised to model specific perceptions of the same system or a combination of different systems. For example, if we consider transportation systems, the obtained network can be
<figure><img src="image_63.png" /><figcaption>**Figure 5.1:** Collecting and leveraging big smart city data as complex dynamic networks.</figcaption></figure>

a multi-layer combination of road and urban transit networks, including buses, trams, and
metros. The considered city-system (each of its layers) can also be divided into numerous
geographical *areas* or logical *groups* that identify specific sets of nodes and links (Fig. 5.1:
level 2).

In addition, both nodes and links are characterised by multiple attributes belonging to two types: (i) *static* attributes representing the monitored system (e.g., road network topology and its basic features), and (ii) *dynamic* attributes obtained through different continuous operations such as data sampling, aggregation, and computation (e.g., average travel speed observations).

*   **Static attributes:** to ensure the proper functioning of the platform, nodes and links must have a non-empty set of static descriptive characteristics that depend on the considered application. Along with these attributes, both nodes and links can have any additional number of attributes, depending on the entities they represent. For instance, nodes and links must have a unique identity (*nodeID* and *linkID*). As a second example, when entities of level 0 (see Fig. 5.1) are non-mobile (have a fixed location) in the physical space, their representative nodes can be characterised by two or three static parameters, depending on whether a 2D or 3D space is considered: *longitude, latitude, altitude*. In some cases (e.g., road networks), links can also be located using the nodes that delimit them: *fromNode*, *toNode* (which are the *IDs* of the *adjacent* nodes of each link<fnref n="3" />).
*   **Dynamic attributes:** each operation performed on nodes and links (e.g., data sampling, graph metric computation, etc.) to measure a dynamic attribute must include at least two elements: (i) the *identity* of the entity to which the operation is associated and (ii) the *timestamp* at which this operation was executed. This means that nodes and links must be uniquely referenced through an *identifier* (ID). An example of attributes falling into this category is the dynamic location of mobile entities of level 0.

*   **Dynamic attributes:** each operation performed on nodes and links (e.g., data sampling, graph metric computation, etc.) to measure a dynamic attribute must include at least two elements: (i) the *identity* of the entity to which the operation is associated and (ii) the *timestamp* at which this operation was executed. This means that nodes and links must be uniquely referenced through an *identifier* (ID). An example of attributes falling into this category is the dynamic location of mobile entities of level 0.

<footnote n="3">³We recall that any two nodes connected by an edge or any two edges connected by a node are said to be adjacent.</footnote>
5.1.2 Functional requirements

The graph-based meta-model can be exploited to implement the diverse inter-connected
systems composing a city in its entirety. By running specific standard or custom graph
algorithms, the platform enables to (i) analyse a city partially or as a whole, (ii) analyse a
single subsystem or a combination of them, and (iii) rapidly tackle emerging problems by
setting properties on the graphs.

According to this main pillar, we identified a (minimal) set of functional requirements a
graph-based smart city platform should satisfy:

*   *Static data upload and update:* the platform must be able to upload and update the essential static data (initial graph structure), which represent the main entities of the targeted city system and their relationships.

*   *Dynamic data acquisition and storage:* the platform should be able to continuously ingest node and link dynamic attributes from different sources. The collected dynamic data must be efficiently stored to be easily accessed via multiple frequent queries.

*   (Graph-based) data processing support: the gathered big data may be useless if not transformed into meaningful and valuable information. Therefore, the platform must provide the necessary general-purpose (stream and batch) processing solutions to close the semantic gap between raw data collected from the environment and the information needed by the targeted services. At the same time, the processing middleware should be flexible enough to easily integrate additional processing services. Moreover, the platform must also provide graph-based data processing support. This set of functions or graph algorithms (such as depth-first search, breadth-first search, Dijkstra single-source shortest path, betweenness centrality, degree centrality, etc.) can be offered by specific graph-based software or implemented using the general-purpose processing engine.

*   **Data visualisation:** the platform must be able to provide a user-friendly, map-based Graphical User Interface (GUI) to visualise city infrastructures, their digital representations, as well as the collected dynamic data, and to enable interaction with the platform itself.

5.1.3 Non-functional requirements

The adoption of a graph meta-model for a large-scale system requires satisfying some non-functional requirements as the model grows in complexity, such as:

*   *Time efficiency:* the objective is to offer smooth, timely online and offline services to final users. To do so, in addition to the possible pre-processing of Edge and Fog nodes, the different implemented cloud software components must be highly time-efficient in terms of data ingestion, storage, processing, and querying (both live-streaming and historical);

*   Availability: the proposed architecture with all its components and services must be fault-tolerant. That is, regardless of software and hardware failures, the system must be able to continue delivering its requested functionalities, and recovery must be instantaneous and transparent (without data loss or service interruption);

*   Scalability: the proposed components must be highly scalable so that they can maintain the required low-latency performance in front of both an increasing data volume (large-size or large-amount of data attributes associated with graph entities) and increasing velocity (high-throughput data);

*   Flexibility: the system should be easily changed by adding, updating, or removing software components at any time without interrupting the currently provided services. This point is fundamental as one cannot foresee and provide in advance all possible services that can be useful for a city.

By pairing DevOps, a methodology as well as a set of tools and practices that *merge* and *automate* development and IT operations, with the adoption of the cloud-native microservices architectural style [219, 220], the effort to implement these requirements can be
lessened.

### 5.1.4 Additional features

Other requirements depend on the applications and can be satisfied by orthogonal layers. In general, a plug-in mechanism should be provided to decouple the platform core from every end-user service that may be needed. The following are examples of possible additional features for some applications.

*   **Data privacy**: the collected IoT data can be personal and sensitive. They cannot be utilised in their arrival raw state and must be pre-processed before use (anonymised, filtered, etc.). Across Europe, as already mentioned in Sec. 3.4, to solve this issue and harmonise data protection and privacy, the EU has specified a set of strict regulations and rules, called GDPR<fnref n="4" />, which all digital businesses and services (using personal sensitive data) must comply with. To anonymise and filter the received data according to the consent of customers, the latter can be offered an online contract that lists all their preferences (e.g., who can access their data, for how long, etc.), and allows them to access and modify these preferences at any time. The solutions proposed in Sec. 3.4 based on k-anonymity are an example in this direction.

*   **Security**: whether in the case of collecting personal-sensitive information (such as personal identification and localisation) or not, security mechanisms (i.e., confidentiality, integrity, and availability of information) should be provided.

*   **Transaction integrity**: some applications, such as the ones that deal with citizens data for city services fruition, should ensure the collection of data as shared and immutable transactions in order to avoid repudiation of the originated data.

Even though we acknowledge the importance of the aforementioned additional requirements, these features and those not mentioned here were not specifically addressed in the prototypical implementation of the platform, being left for future work and platform extensions.

## 5.2 PROMENADE architecture

PROMENADE is a set of interconnected microservices that can be operationally divided into four *layers* and one *plane* (Fig. 5.2). The platform Initialization and Monitoring Plane (IMP) is a set of software components dedicated to initialise and monitor the platform and its different layers. It is responsible for uploading and updating the main static data of the monitored system. It is also responsible for monitoring the health of the cloud nodes, the hosted software components, their CPU, memory, and network usage, etc. All these features can be offered to end-users/administrators via a Web-GUI, through which they can initialise, monitor, and control the platform, its static and dynamic data, and its overall status.

The first three layers, which have been divided into stream and batch components, are in charge of big data ingestion, storage, and processing: (i) Stream data Ingestion (SI) and Historical data Ingestion (HI), (ii) Distributed Asynchronous producer-consumer Queue (DAQ) and Permanent Storage (PS), (iii) Stream Processing (SP) and Batch Processing (BP). The fourth Service offering layer (SVC) can query either or both the stream and permanent storage (DAQ and PS), which can contain raw and processed data, to provide the desired services to end-users. According to the Lambda architecture [221], PROMENADE handles both live streams and historical data. However, due to the nature of the applications, it mainly works on time windows of different sizes, often relying on Kappa architecture principles [222] at the implementation level.

<footnote n="4">⁴General Data Protection Regulation; applicable as of May 25, 2018, in all member states of EU (https://gdpr-info.eu/).</footnote>
<figure><img src="image_64.png" /><figcaption>**Figure 5.2:** Detailed architecture of PROMENADE.</figcaption></figure>

In the following, we detail the layers of PROMENADE in this order: (1) SI and HI, (2) DAQ and PS, (3) SP and BP, and finally (4) SVC.

### 5.2.1 Live-stream and historical-data ingestion

Once the main software components (e.g., microservices) have been deployed and the platform properly initialised with the underlying data (static attributes of nodes and links), it is ready to start ingesting dynamic data and performing its required tasks. To achieve this goal, the main role of this ingestion (SI and HI) layer is (i) gathering (stream and historical) data from the defined sources (using predefined or custom software components) and (ii) preparing this input by making it usable by the upper layers. The HI is a specific software component responsible for ingesting pre-recorded dynamic data from external historical sources, while the SI provides live stream data ingestion from IoT devices. For instance, if the adopted communication protocol is MQTT<fnref n="5" />, the SI component can be a high-performance MQTT broker responsible for receiving, decoding, formatting, and handing the output to the upper DAQ component [223].

<footnote n="5">⁵MQTT: The Standard for IoT Messaging, https://iot.mqqt.org</footnote>
5.2.2 Temporary and permanent storage

Raw and processed data are stored in two modes: *temporary stream queuing* and *permanent historical storage*.

* *Temporary stream queuing*: implemented through a **DAQ**, this component can mediate the interaction between every other platform component, thus enabling an asynchronous communication model. For example, the **DAQ** can receive both raw data coming from the bottom ingestion (**SI**) components and processed data generated by the top processing (**SP** and **BP**) components, making them consumable by any other interested component.

* **Permanent storage**: the **PS** component can be provisioned from three data sources: (i) initially by the IMP plane (with the structure and characteristics of the monitored system, i.e., static node and link attributes), (ii) continuously from the previous **DAQ** component located to its left in Fig. 5.2 (dynamic data continually flows from the temporary space to the dedicated permanent storage), or (iii) at any time from the historical data ingestion (**HI**) component located below it. For the received *static* and *dynamic* data to be efficiently accessed and analysed at any subsequent time, we designed them to be separately stored in two distinct non-ephemeral **PS** components (Fig. 5.2), and therefore undergo different operations. Concerning static data (i.e., the nodes and the links, as well as their static attributes), they can be created, updated, and deleted. This is not the case for dynamic data (i.e., collected or computed attributes); they are recorded as immutable append-only.

5.2.3 Stream and batch processing

Data are sent to the platform to be processed for immediate (stream) or subsequent (his-
torical) use. These two processing modes are based on the previously presented temporary
and permanent storage (Sec. 5.2.2).

*   **Stream processing**: the **SP** is responsible for ensuring low-latency, high-throughput retrieval and processing of continuous **IoT** data streams temporarily stored in the lower **DAQ** component. The **SP** is also responsible for data pre-processing. According to the considered application, the continuously ingested streams can be sent to this component to undergo all necessary beforehand cleansing, filtering, or aggregation. This stream processing component can retrieve data from both the stream and permanent (**DAQ** and **PS**) storage, but can only write (the obtained results) to the stream queue (**DAQ**). For instance, after pre-processing stream data, the modified data can be returned to the **DAQ**, and from there stored in the dedicated permanent storage (Fig. 5.2). This architectural choice provides *time* and *synchronisation decoupling*, which among other things help prevent backpressure issues.

*   **Batch processing**: the **BP** is responsible for processing historical data. Similar to **SP**, this component can also read data from both the permanent and temporary storage (**DAQ** and **PS**), and can only write to the stream **DAQ**. As examples of possible tasks assigned to this component, we mention searching, querying, and applying different machine learning algorithms (classification, clustering, anomaly detection, etc.) to the large available historical data.

5.2.4 Service offering

Services can perform data *processing* or data *reporting*, where the latter may depend on the former. In particular, *processing* services, deployed at the level of the stream and historical pipelines, further process the incoming data and either expose their insight to the upper layer or persist them back into the storage components. On the other side, *reporting* services, which are either consuming raw data or exploiting the *processing* services output, offer fully-fledged dashboards.
## 5.3 Platform implementation and its specialisation in road network monitoring

This section discusses the implementation of the platform core framework, called PROMENADE-v2.0 (Sec. 5.3.1), and then tackles the specialisation of this core for road network monitoring (Sec. 5.3.2). The empirical evaluation of the platform performance is not reported in this document for brevity: readers can refer to the full paper [214] for a detailed presentation of the conducted experiments. This study allowed us to assess the scalability of PROMENADE-v2.0, and to set the main architectural parameters that significantly impact the performances of the proposed solution.

### 5.3.1 Platform core implementation

PROMENADE-v2.0 is the evolution of two previous versions proposed in [224] and [225]. In the following, we first introduce the underlying software and hardware technologies on top of which this core has been implemented and then present its main components and their technological choices. It is worth underlining that the development of the platform is an ongoing process and that only a part of the components has been fully implemented in the prototype, as graphically detailed in Fig. 5.4.

**Underlying software and hardware technologies**

The deployment environment of PROMENADE-v2.0 is a software/hardware infrastructure composed of:

* An open-source middleware built around Docker<sup>6</sup>, Kubernetes, and Red Hat OpenShift (Fig. 5.3), the currently most popular DevOps tools to create and orchestrate containers.

<figure><img src="image_65.png" /><figcaption>Figure 5.3: Underlying software and hardware technologies of PROMENADE-v2.0.</figcaption></figure>

<footnote n="6">⁶Docker: https://www.docker.com/</footnote>
* A virtualization environment built atop OpenStack<sup>7</sup> that provides services for provisioning virtual machines (Nova), providing networking capabilities (Neutron) and managing disks (Cinder) to offer a block storage system.

* Five physical machines to host the OpenStack services and virtual machines. One of the physical machines operates as an OpenStack controller (with 1 CPU at 3.20GHz with 4 cores/4 threads), while the other four are configured to work as either Compute or BlockStorage nodes (each with 2 CPUs Intel Xeon, running at 2.20 or 2.40 GHz, more than a hundred of GB of RAM, TBs of hard disk, and two 1 Gb/s Ethernet adapters). The characteristics of the controller, physical hosts (Host 0, 1, 2, and 3) as well as their corresponding guest VMs (Master, Infra, Compute 1, 2, 3, 4, 5, and 6) are summarized in Fig. 5.3.

The overall infrastructure provides dynamically provisioned computing units, called Pods, that can contain one or multiple containers. To satisfy the desired non-functional requirements, the infrastructure is specialised to govern the way such Pods are built, deployed, and kept in their desired state. This is performed by different layers of controllers and mechanisms such as: ReplicaSet (self-healing, scalability, ...), Deployment (updates and rollbacks), StatefulSet, Jobs, Services, Routes, Custom Resource Definition (CRD), Operators, etc. We redirect the reader to the specific technical documentation to have more details about the concepts introduced above.

The main role played by each component of PROMENADE has been discussed in Sections 5.2.1 to 5.2.4. In the following, we discuss the implementation of these components using open-source software technologies, which have been chosen not only because of their popularity, stability, and efficiency, but also because they comply with the design objective of PROMENADE. Tab. 5.1 compiles and maps the components of PROMENADE-v2.0 to their corresponding Kubernetes objects.

## Live-stream ingestion

As Fig. 5.4 shows, Apache ActiveMQ Artemis<sup>8</sup> is the stream ingestion tool used to collect dynamic data from the lower IoT/Edge networks. To face any increase/decrease in the size of the monitored network, the deployed Artemis broker instances can be scaled up and down respectively (manually or automatically as needed). Also, to maintain the persistence and consistency of the data in case of possible crashes of the broker instances, a StatefulSet controller has been used (Tab. 5.1). These two points apply to any other component of the architecture that needs such control.

<table><caption>**Table 5.1:** Mapping of PROMENADE-v2.0 components to their corresponding K8s objects.</caption><thead><tr><td>Component</td><td>K8s object</td></tr></thead><tbody><tr><td>ActiveMQ Artemis</td><td>Stateful Set (STS) with 2 Persistent Volume Claims (PVCs) of 8 GiB each.</td></tr><tr><td>Neo4J</td><td>Three read-write cores (STS) and one read replica (Deployment). The three read-write cores consist of 1 leader and 2 followers with 3 PVCs of 10 GiB each.</td></tr><tr><td>Kafka</td><td>STS with 3 brokers/pods (0, 1, and 2) that have each one PVC of 10 GiB.</td></tr><tr><td>ZooKeeper</td><td>STS (0, 1, and 2) with 6 PVCs of 10 GiB each.</td></tr><tr><td>Kafka Center</td><td>Deployment.</td></tr><tr><td>kafka-Artemis source connector</td><td>Deployment.</td></tr><tr><td>Kafka-Mongodb sink ingestion/processing connector</td><td>Deployment.</td></tr><tr><td>MongoDB</td><td>Deployment with 1 PVC of 20 GiB.</td></tr><tr><td>Spark</td><td>The considered cluster is composed of one master. (CPU: 1 core, RAM: 1 GiB) and one worker (CPU: 32 cores, RAM: 200 GiB).</td></tr></tbody></table>

<footnote n="7"><sup>7</sup>OpenStack: https://www.openstack.org/</footnote>

<footnote n="8"><sup>8</sup>Artemis: https://activemq.apache.org/components/artemis/</footnote>
<figure><img src="image_66.png" /><figcaption>**Figure 5.4:** Detailed architecture of PROMENADE-v2.0. Gray boxes indicate that the corresponding components are not implemented in the current prototype of the platform.</figcaption></figure>

As discussed in Sec. 5.1.1, the considered city-system (at each layer) can be split into logical groups (e.g., geographical areas) of nodes and links. To implement this design, we considered a one-to-one mapping between the groups and the internal data structures of SI, DAQ, and PS components. For instance, a group corresponding to a geographical area in a city can be mapped to an Artemis queue (SI component). It is worth remembering that these groups can be defined via the initialisation component from the IMP plane. The implementation details about the one-to-one mapping of the DAQ and PS components will be presented in the next subsections. This partitioning configuration allowed us to achieve better performances.

## Temporary queuing storage

The popular Apache Kafka<fnref n="9" /> and its discovery service ZooKeeper have been considered to comply with the design of the queuing component. To deal with data persistence while providing fault-tolerance, a StatefulSet controller has been used for both Kafka and ZooKeeper, as is the case for Artemis (Tab. 5.1). This way, if one of the stream brokers goes down, it will be recreated anew while retaining the old session.

<footnote n="9">⁹Kafka: https://kafka.apache.org/</footnote>
In accordance with the *temporary* nature of this component, data is not permanently
stored: depending on the chosen retention policy, data can be retained for one hour, one
day, one week, *etc*. However, an important point about data-flow concerns the *source* and
*sink* connectors, which are respectively in charge of getting data in and out of Kafka. We
used three connectors to implement PROMENADE-v2.0 (Fig. 5.4): one customized source
connector (which links Kafka to Artemis) and two out-of-the-box sinks (which both connect
Kafka to MongoDB for ingestion and processing, respectively). The role of MongoDB is
explained in the next subsection.

The adopted Kafka source connector was implemented through an adapter that retrieves packets from Artemis queues, uniforms their format, and forwards them to Kafka topics. The two sink connectors work as follows: when Kafka receives data from the Artemis-based SI component, the first connector writes the received dynamic data of links and nodes (e.g., *average travel time* data of links) from Kafka to MongoDB, while the second is responsible for writing the computed node and link attributes (e.g., *betweenness centrality* of nodes or any other information computed at the SVC level) from Kafka to MongoDB.

As done for Artemis (SI component), to implement the logical grouping of city entities,
also this DAQ component requires a specific topic to be associated to each group.

Permanent historical storage

As detailed in Sec. 5.2.2, PROMENADE-v2.0 handles static and dynamic data with two
distinct non-ephemeral components. To implement this design, we opted for two popular
technologies: Neo4J<sup>10</sup> and MongoDB<sup>11</sup>, for static and dynamic data, respectively.

The former is an open-source NoSQL graph DBMS implemented in Java that supports
efficient graph-based queries of static data (structure and characteristics of the monitored
system). To this end, Neo4J uses Cypher, a declarative query language similar to SQL but
optimized for graph databases.

MongoDB is an open-source NoSQL document DBMS, developed in C++ around two
main concepts: *documents* and *collections*. A *collection* is a group of documents, and a
*document* is a data structure composed of field-value pairs (JSON objects). We used Mon-
goDB for dynamic data because samples in the form of JSON documents can be immediately
stored, while embedded documents and arrays reduce the need for expensive joins. More-
over, the dynamic schema support, which means that collections do not enforce any schema,
ensures high flexibility for handling different kinds of data.

For each logical group of city entities, we created specific MongoDB collections. Each of
these collections is ultimately in charge of storing the node and link samples coming from
the stream ingestion component, and the node and link attributes computed by the data
processing component.

Data processing

To satisfy the functional requirements listed in Sec. 5.1.2, PROMENADE-v2.0 offers both *general-purpose* and *graph-based* data processing support. The former has been implemented by leveraging a general-purpose engine for large-scale data processing. Specifically, as shown in Fig. 5.4, we opted for the well-known Apache Spark<sup>12</sup> framework, which can work on both data streams and data at rest. This way, by injecting appropriate (general-purpose or graph-based) algorithms into this data processing component, the platform can be easily extended to provide end-users with fine-grained services catering to the different problems that should be solved.

<footnote n="10">¹⁰Neo4J: https://neo4j.com/</footnote>

<footnote n="11">¹¹MongoDB: https://www.mongodb.com/</footnote>

<footnote n="12">¹²Spark: https://spark.apache.org/</footnote>
The basic functions of the graph-based support are provided by Neo4J (which is used for storing and querying static data). For instance, the library of Neo4J contains implementations of path discovery and centrality algorithms<fnref n="13" />. The first category helps to identify shortest paths and perform search tasks (e.g., depth-first search, breadth-first search, Djekstra single-source shortest path, etc.), whereas the second determines the importance of distinct nodes in a network (e.g., betweenness, degree, eigenvector centrality, etc.). Other examples of provided algorithms include community detection, similarity functions, topological link prediction, node embeddings, node classification, and link prediction.

### 5.3.2 Specialisation for Road Network Resilience

To assess the flexibility of the proposed platform architecture and prototype, we report in this section how PROMENADE-v2.0 has been specialised for real-time road traffic monitoring and management applications.

#### Data model

Instantiating the general platform means (i) instantiating the graph-based meta-model (of Sec. 5.1.1) into a specific model related to the considered specialisation, and (ii) providing the necessary dynamic (node and link) values for this graph through appropriate data sources, whereas the desired non-functional requirements are satisfied as a consequence of the application-level choices and the available hardware resources. In particular, we specialised PROMENADE-v2.0 for road traffic monitoring and management by (i) modelling the monitored road infrastructure as a one-layer directed multi-attributed graph, and (ii) consuming the constantly changing traffic conditions (dynamic observations) retrieved from its links.

Formally, the road network graph is denoted as $G_{RN} = (V, E, A, B)$, where $V$ is the set of vertices or road intersections, $E$ is the set of edges or road segments, $A$ is the set of attributes associated to nodes, i.e., $\{A_v | v \in V\}$ and $B$ is the set of attributes associated to the edges, i.e., $\{B_e | e \in E\}$. The attributes associated to nodes and edges in such a graph could be static (not frequently changing) or dynamic. More specifically, as shown in Tab. 5.2, nodes static attributes include: a *unique identifier* (osmID<fnref n="14" />), a *latitude-longitude* value, and their *area name*. On the other hand, edges static attributes include: *linkID*, *fromNode*, *toNode*, *street name*, *street length*, *speed limit*, *estimated free-flow speed*, *area name*, and *coordinates*; *fromNode* and *toNode* are the osmID of the adjacent nodes of an edge and *coordinates* is its geometry. The choice of the attributes listed above is dictated by the available fields in the datasets used in our case study evaluation. However, nodes and links can have other static attributes depending on the specific application. For instance, our dataset allows further characterising a link via the *width* (in meters), the *number of lanes*, the presence of *bridge* or *tunnel* for the associated road segment. These attributes were ignored in this work as irrelevant to our case study.

As an additional feature, the single layer of our graph, i.e., $G_{RN}$, is composed of multiple sub-networks $G_{RN}^z$, each representing a different geographical area $z \in \mathcal{Z}$ (from a given segmentation of the modelled city) identifiable via the *area name* parameter. More formally, $G_{RN}^z$ is the sub-network induced by the set of nodes $V^z \subseteq V$ and edges $E^z \subseteq E$ that are geographically located in the given area $z$.

The mandatory and non-mandatory *static* characteristics of road networks (Tab. 5.2) are injected into the implemented platform during the initialisation phase. Once configured and ready, the platform can then be continuously fed with *dynamic* attributes (in our case, *average travel time* of links, as from Tab. 5.3). The values of the link attributes can be thus used to produce insightful vulnerability indicators such as the betweenness centrality

<footnote n="13">¹³Neo4j Graph Algorithms: https://neo4j.com/developer/graph-data-science/graph-algorithms/</footnote>

<footnote n="14">¹⁴OpenStreetMap ID (https://www.openstreetmap.org/).</footnote>
<table><caption>**Table 5.2:** Road segments (links) and intersections (nodes): static data provided to PROMENADE-v2.0.</caption><thead><tr><th>Node attrib.</th><th>Type</th><th>Description</th></tr></thead><tbody><tr><td>osmID</td><td>long</td><td>OpenStreetMap ID<sup>14</sup>.</td></tr><tr><td>longitude</td><td>double</td><td>Longitude (x-coordinate in decimal degrees).</td></tr><tr><td>latitude</td><td>double</td><td>Latitude (y-coordinate in decimal degrees).</td></tr><tr><td>area name</td><td>String</td><td>Node area name.</td></tr><tr><th>Link attrib.</th><th>Type</th><th>Description</th></tr><tr><td>linkID</td><td>long</td><td>Link ID.</td></tr><tr><td>fromNode</td><td>long</td><td>OsmID of first adjacent node u.</td></tr><tr><td>toNode</td><td>long</td><td>OsmID of second adjacent node v.</td></tr><tr><td>name</td><td>String</td><td>Link (street) name.</td></tr><tr><td>length</td><td>double</td><td>Link length (meters).</td></tr><tr><td>speedLimit</td><td>int</td><td>Maximum speed (km/h); 0 if not applicable.</td></tr><tr><td>ffs</td><td>double</td><td>Estimated free-flow speed (km/h).</td></tr><tr><td>areaName</td><td>String</td><td>Link area name.</td></tr><tr><td>coordinates</td><td>List(double, double)</td><td>Sequence of (lat, long) points representing a street between two intersections (going from u to v).</td></tr></tbody></table>

<table><caption>**Table 5.3:** Road segments (links) and intersections (nodes): dynamic data collected and produced by PROMENADE-v2.0.</caption><thead><tr><th>Link attrib.</th><th>Type</th><th>Description</th></tr></thead><tbody><tr><td>linkID</td><td>long</td><td>Link ID.</td></tr><tr><td>timestamp</td><td>dd/mm/yyyy hh:mm:ss</td><td>Data aggregation time.</td></tr><tr><td>average tt</td><td>double</td><td>Average travel time (s).</td></tr></tbody></table><table><thead><tr><th>Node attrib.</th><th>Type</th><th>Description</th></tr></thead><tbody><tr><td>osmID</td><td>long</td><td>Node ID.</td></tr><tr><td>timestamp</td><td>dd/mm/yyyy hh:mm:ss</td><td>Data computation time.</td></tr><tr><td>BC</td><td>double</td><td>Betweenness centrality.</td></tr></tbody></table>

of nodes (see Sec. 4.2.2), which represents, in turn, a time-varying attribute of nodes (as reported again in Tab. 5.3).

## Stream data generation

The IoT networks correspond to the underlying level on top of which the platform was built (Fig. 5.4). Although being crucial, in the implemented version of our platform specialisation, they have been emulated by playing back an available offline dataset. This choice was motivated by two main reasons. First, the strong dependency of the IoT/Edge systems on the physical infrastructure of cities, and second, the flexibility of generating data as needed. The dataset used in the emulation to reproduce realistic road traffic conditions was the same as the one used in Sec. 4.2 for our experimental studies on vulnerability assessment related to the city of Lyon. This dataset was organised according to the structure presented in Tab. 5.2 and Tab. 5.3.

We implemented the emulation software visualised in Fig. 5.5: it can be (i) exploited and parameterised to stress and test the developed platform while considering extreme conditions (something that is not possible or very difficult in a physical environment), (ii) easily reused to study and analyse other complex networks for which offline historical data is available, and (iii) rapidly and straightforwardly substituted whenever physical IoT networks are available. Moreover, it does not reduce the overall validity of the platform. It is worth noting that the emulation focuses on data generation and does not consider edge-to-cloud communication delays, which were left out of scope for our prototype implementation.

From a technical point of view, the main task of the emulator is to consume the provided historical dynamic (raw) data (described in Tab. 5.4) and behave similarly to a set of Edge nodes.

Specifically, we coded the emulator to continuously produce the aggregate (average)
<figure><img src="image_67.png" /><figcaption>**Figure 5.5:** Architecture of the IoT data-driven emulation software.</figcaption></figure>

travel time on links (Tab. 5.3) based on measures of speed observed from individual vehicles. To this end, the emulator was designed (see Fig. 5.5) to be a multi-threaded Java application that is composed of three main layers: *sensors layer*, *Edge processing layer*, and *Edge control layer*:

*   **Sensors layer**: its role boils down to emulating data readings at the Edge. As shown in Fig. 5.5, it retrieves all samples (raw non-aggregated data of Tab. 5.4) and sends them as messages to the different queues of the Edge processing layer, based on the *area name* attribute from the link static data file of Tab. 5.2. In fact, the Edge processing layer hosts a queue for each geographic area. By this approach, we kept the geographic division inside the digital model by implementing it at the level of software components.

*   **Edge processing layer**: the role of each Edge node of this layer is to consume raw samples (vehicle speed readings from sensors or on-board connected GPS devices) that are continuously coming to its respective queue (from the emulated sensors layer) and operates as follows. During the aggregation interval (*e.g.*, 3 minutes in our case), an emulated Edge node just handles a given vehicle sample by (*i*) increasing the number of observations for the link crossed by the vehicle (*link.numObservations++*), (*ii*) computing the vehicle travel time (*vtt*)<fnref n="15" />, and (*iii*) adding the computed (*vtt*) value to the link statistics (*link.stats*). When the aggregation interval expires (after 3 minutes), the node computes and writes the aggregation result of its links (*e.g.*, *avgTravelTime*, *stdDevTravelTime*, *aggTimestamp*, *etc.*) to a JSON file that refers to a specific area represented by the queue

<footnote n="15">¹⁵vtt = (link.length * msg vehicleCoverage) / msg vehicleSpeed</footnote>
from which it consumes the raw data (Fig. 5.5).

* Edge control layer: mimics the behaviour of multiple *parallel* Edge nodes. More specifically, once the previous layer has finished its task of aggregating all data, this layer can start using its output JSON files. As shown in Fig. 5.5, this layer is composed of one main thread and *n* Edge threads (corresponding to the *n* areas of the city of Lyon). The role of the main thread is to control the online Edge threads by creating, starting, and synchronising them. Upon their creation, Edge threads sleep indefinitely and wait to be awakened. Once awakened, each Edge thread reads its corresponding aggregation file, sends data to Artemis, and then goes back to sleep (waiting for the next aggregation interval to come).

<table><caption>**Table 5.4:** Road segments: dynamic raw (non-aggregated) data of vehicles.</caption><thead><tr><th>Link attribute</th><th>Type</th><th>Description</th></tr></thead><tbody><tr><td>linkID</td><td>long</td><td>Link ID.</td></tr><tr><td>timestamp</td><td>dd/mm/yyyy<br>hh:mm:ss</td><td>Data collection time.</td></tr><tr><td>v_ID</td><td>int</td><td>Vehicle or observation ID.</td></tr><tr><td>v_speed</td><td>double</td><td>Measured vehicle speed (km/h).</td></tr><tr><td>visen coverage</td><td>double</td><td>How much of the road segment has been covered by the vehicle when the sample was taken (0 ≤ coverage ≤ 1).</td></tr></tbody></table>

The emulation process has been designed to be configurable. It can be *realistic* (samples generated according to their original timestamps) or *virtual* (samples generated at a very high or low speed to enable several analysis scenarios). To be able to do so, a temporal factor is utilised to allow control of the frequency at which samples are generated while retaining their respective order. Moreover, each emulated Edge node can either compute the average travel time (of each link during a given time interval, e.g., every 3 minutes), or it can behave differently depending on the considered application. For instance, it can act as a collector by forwarding raw samples to the Cloud platform.

We recall that the one-to-one mapping between geographical/logical partitions and queues implemented in the emulated Edge processing layer is also adopted in the SI, DAQ, and PS components of PROMENADE-v2.0.

## Implemented services

The specialisation of PROMENADE-v2.0 for road network monitoring provides two main classes of services: *real-time monitoring* and *traffic management*. We also implemented visualisation support for real-time presentation of the observed traffic conditions.

Real-time monitoring exploits the betweenness centrality as a metric to *continuously* detect critical elements of the network, i.e., those characterised by higher values of central-ity [213]. Introduced in [193], BC measures the importance of nodes or links for the flow of information in a graph. It is widely used to identify opinion leaders or influential people in social network analysis [226], vulnerabilities in computer networks [227], and threats from terrorist networks [228]. As discussed in Chapter 4, BC can also help detect critical intersections in transportation networks [213, 229, 230, 231, 232] and be used as a vulnerability indicator. The knowledge of critical nodes or links may enact appropriate control actions, provide relevant recommendations, and guide drivers through less vulnerable itineraries not traversing high-BC nodes.

Based on Brandes’ algorithms [192, 233], the BC solutions provided by Neo4J are inefficient [213]. Thus, we used the general-purpose processing engine available in PROMENADE (Spark) to run our algorithmic approach proposed in [213] for approximate BC computation on directed and weighted graphs. The latter has been further enhanced in [234] for the exact computation of BC on undirected and unweighted graphs (see Sec. 5.4). These algorithms have been integrated as a core feature in PROMENADE-v2.0 to periodically compute the BC values of large-scale time-varying networks.
<figure><img src="image_68.png" /><figcaption>**Figure 5.6:** GUI examples - Map with vulnerability indicator (BC). The hottest spots of the heatmap represent critical nodes with high BC values.</figcaption></figure>

The collected dynamic data, together with the computed *BC* values (sent to Kafka and then stored in MongoDB), can be exploited by other additional data processing and visualisation services. We used two different connectors for collecting traffic conditions (i.e., the average travel times) and for producing the real-time vulnerability indicators (i.e., the *BC* values of nodes). The first connector is based on a RESTful API and operates in a pull mode, whereas the second exploits WebSockets and works in a push mode to provide continuous updates. The Kafka broker handles these operations: whenever an update/event occurs, it triggers the process of retrieving the necessary data and serving subscriptions. Fig. 5.6 reports an example of the output produced by the platform. *BC* values are represented in two different moments of a day using a heatmap layer, showing the evolution of the most central nodes over the day depending on varying traffic conditions.

The details of the adopted *BC* algorithm for real-time traffic monitoring are described in the next section, focusing on the solution proposed for rapid exact computation on large-scale undirected and unweighted networks [234], leaving the interested reader to our papers [235, 213] for a detailed description of the approximate versions for weighted and directed networks.

The solutions implemented for real-time traffic management are instead reported in
Sec. 5.5.

**5.4 Fast Cluster-based Computation of Betweenness Cen-**
**trality in Large Graphs**

Despite its great potential, *BC* requires a considerable computation time, which often rep-
resents a barrier to applying this metric in large-scale contexts, especially with dynamic
graphs, whose topology or weights change frequently.

In recent years, the Floyd method [236], which requires $O(n^3)$ computation time, has
been overcome by the well-known Brandes’ algorithm [192]. Given a graph $G(V, E)$, it exhibits $O(n + m)$ space complexity, $O(nm)$ time complexity for unweighted graphs, and $O(nm + n^2 \log(n))$ for weighted ones, where $n = |V|$ is the number of nodes and $m = |E|$ the number of edges. However, the polynomial complexity of Brandes’ algorithm is still an obstacle for analysing very large networks. Such a problem becomes even more evident and limiting if centrality is used for real-time (or quasi-real-time) analysis of dynamic networks.

By building on our previous work on approximated BC computation for unweighted-undirected graphs [237, 232] and weighted-directed ones [235, 213], we developed, as part of Cecile Daniel’s PhD thesis, a novel algorithm (widely discussed in [238, 234]) which makes it possible to perform the exact calculation of BC on large, undirected graphs with a high speedup when compared to Brandes’ algorithm and a significant improvement over recent variants of Brandes’ algorithm based on clustering [239].

Undirected graphs are very common in real-world systems: examples include social networks, communication networks, protein interaction graphs, people interaction graphs, finite element meshes, and the topology of transport networks. Among these graphs, scale-free and Barabási-Albert graphs [240] represent an ideal candidate for the analysis of the proposed algorithm, since they model many real-world systems, such as the World Wide Web, the Internet and other computer networks, citation networks, social networks, airline networks, financial networks, etc. However, considering the application of our approach for real-time monitoring of road traffic networks, we also examined small-world networks by including in our benchmark the real-world road network of Lyon, France.

In the following, we introduce the notation used throughout the paper and briefly describe Brandes’ algorithm. Then, we present our algorithm and the most relevant performance evaluation results. The presentation of the literature review, the details of the algorithm implementations (a Scala map-reduce parallel and a sequential one), the mathematical proofs and detailed performance analyses are left out of this document for brevity but can be found in the aforementioned publications [238, 234].

### 5.4.1 Background

Let $G(V, E)$ be an undirected unweighted graph with $V$ representing the set of $n$ vertices (or nodes) and $E$ the set of $m$ edges (or links). Let $s, t \in V$ be two generic nodes of $G$. We denote by $e_{s,t}$ the edge connecting $s$ and $t$. The neighbours of a vertex $s$ are all vertices $u$ such that $e_{s,u} \in E$. The distance between $s$ and $t$, denoted by $d_G(s, t)$, is the length of the shortest path(s) connecting them in $G$. The number of shortest paths between $s$ and $t$ is denoted by $\sigma_{s,t}$, whereas the number of shortest paths between $s$ and $t$ that cross a generic node $v \in V$ is denoted by $\sigma_{s,t}(v)$. It is worth noting that since the graph is undirected, $d_G$ and $\sigma$ are symmetric functions, thus $d_G(s, t) = d_G(t, s)$, $\sigma_{s,t} = \sigma_{t,s}$ and $\sigma_{s,t}(v) = \sigma_{t,s}(v)$. Given a generic node $w \in V$, $P_s(w) = \{u \in V : e_{u,w} \in E, d_G(s, w) = d_G(s, u) + 1\}$ is the set of direct predecessors of vertex $w$ on shortest paths from $s$.

The Betweenness Centrality (BC) of a vertex $v \in V$ is defined as follows:

$$BC(v) = \sum_{s \neq v \neq t \in V} \frac{\sigma_{s,t}(v)}{\sigma_{s,t}} \qquad (5.1)$$

$BC(v)$ thus represents the fraction of shortest paths containing $v$ among all the shortest paths in the graph between any generic pair of nodes $s$ and $t$, summed over all possible pairs $s$ and $t$ with $s \neq v, s \neq t$ and $v \neq t$.

We refer to Tab. 5.5 for a summary of the notation used in the paper.
<table><caption>Table 5.5: The notation used with the BC exact computation algorithm.</caption><thead><tr><th>Notation</th><th>Description</th></tr></thead><tbody><tr><td><b>G</b></td><td>undirected unweighted input graph</td></tr><tr><td><b>G&#770;</b></td><td>a connected sub-graph of <b>G</b></td></tr><tr><td><b>V</b></td><td>set of vertices of <b>G</b> (|<b>V</b>| = <i>n</i>)</td></tr><tr><td><b>V&#770;</b></td><td>set of vertices of <b>G</b> inducing <b>G&#770;</b> (set of vertices of <b>G&#770;</b>)</td></tr><tr><td><b>V&#770;/G&#770;</b></td><td>set of vertices in <b>V</b> \ <b>V&#770;</b></td></tr><tr><td><b>V</b><sub>HSN</sub></td><td>set of vertices of HSN</td></tr><tr><td><b>E</b></td><td>set of edges of <b>G</b> (|<b>E</b>| = <i>m</i>)</td></tr><tr><td><i>e</i><sub>s,t</sub></td><td>edge connecting vertices <i>s</i> and <i>t</i></td></tr><tr><td><i>d</i><sub><b>G</b></sub>(<i>s</i>, <i>t</i>)</td><td>distance between vertices <i>s</i> and <i>t</i> in <b>G</b></td></tr><tr><td><i>d&#770;</i><sub><b>G</b></sub>(<i>s</i>, <i>t</i>)</td><td>normalised distance between vertices <i>s</i> and <i>t</i> in <b>G</b></td></tr><tr><td><i>σ</i><sub>s,t</sub></td><td>number of shortest paths between vertices <i>s</i> and <i>t</i></td></tr><tr><td><i>σ</i><sub>s,t</sub>(<i>v</i>)</td><td>number of shortest paths between vertices <i>s</i> and <i>t</i> which cross vertex <i>v</i></td></tr><tr><td><i>σ&#770;</i><sub>s,t</sub></td><td>normalised number of shortest paths between vertices <i>s</i> and <i>t</i></td></tr><tr><td><b>P</b><sub>s</sub>(<i>v</i>)</td><td>set of direct predecessors of vertex <i>v</i> on shortest paths from vertex <i>s</i></td></tr><tr><td><b>P</b><sub>s</sub>(<b>V</b>)</td><td>set of direct predecessors of vertices in <b>V</b> on shortest paths from vertex <i>s</i></td></tr><tr><td><i>BC</i>(<i>v</i>)</td><td>betweenness centrality of vertex <i>v</i></td></tr><tr><td><i>δ</i><sub>s,t</sub>(<i>v</i>)</td><td>pair-dependency of the pair of vertices (<i>s</i>, <i>t</i>) on the intermediary vertex <i>v</i></td></tr><tr><td><i>δ</i><sub>s,•</sub>(<i>v</i>)</td><td>dependency score of vertex <i>s</i> on vertex <i>v</i> due to all destination vertices</td></tr><tr><td><i>δ</i><sub>s,V&#770;</sub>(<i>v</i>)</td><td>dependency score of vertex <i>s</i> on vertex <i>v</i> due to all destination vertices in <b>V&#770;</b></td></tr><tr><td><b>C</b></td><td>set of clusters of <b>G</b></td></tr><tr><td><b>C</b><sub><i>i</i></sub></td><td>a generic cluster in <b>C</b></td></tr><tr><td><b>C</b>(<i>v</i>)</td><td>the cluster vertex <i>v</i> belongs to</td></tr><tr><td><b>C</b><sup>*</sup></td><td>set of extended clusters in <b>G</b></td></tr><tr><td><b>C</b><sup>*</sup><sub><i>i</i></sub></td><td>a generic extended cluster in <b>C</b><sup>*</sup></td></tr><tr><td><b>K</b></td><td>set of all the equivalence classes</td></tr><tr><td><b>K</b><sub><i>i</i></sub></td><td>an equivalence class</td></tr><tr><td><b>K</b><sub><b>C</b><sub><i>i</i></sub></sub></td><td>set of equivalence classes of cluster <b>C</b><sub><i>i</i></sub></td></tr><tr><td><b>P</b></td><td>set of all the pivots</td></tr><tr><td><i>k</i><sub><i>i</i></sub></td><td>pivot node of the equivalence class <b>K</b><sub><i>i</i></sub></td></tr><tr><td><b>EN</b></td><td>set of all the external nodes</td></tr><tr><td><b>EN</b><sub><b>C</b><sub><i>i</i></sub></sub></td><td>set of external nodes of cluster <b>C</b><sub><i>i</i></sub></td></tr><tr><td><b>BN</b></td><td>set of all the border nodes</td></tr><tr><td><b>BN</b><sub><b>C</b><sub><i>i</i></sub></sub></td><td>set of border nodes of cluster <b>C</b><sub><i>i</i></sub></td></tr><tr><td><b>BN</b><sub><b>C</b><sub><i>i</i></sub>(<i>s</i>, <i>t</i>)</sub></td><td>set of border nodes of cluster <b>C</b><sub><i>i</i></sub> on shortest paths from <i>s</i> ∈ <b>V</b><sub><b>C</b><sub><i>i</i></sub></sub> to <i>t</i> ∈ <b>V&#770;</b><sub><b>C</b><sub><i>i</i></sub></sub></td></tr><tr><td><i>b</i><sub><i>i</i></sub></td><td>a generic border node in <b>BN</b></td></tr><tr><td><i>δ</i><sup>γ</sup><sub>s,•</sub>(<i>v</i>)</td><td>global dependency score of <i>s</i> on <i>v</i> due to all <i>t</i> ∈ <b>V&#770;</b><sub><b>C</b><sub><i>s</i></sub></sub> (same as <i>δ</i><sup>γ</sup><sub>s,V&#770;<sub>C<sub><i>s</i></sub></sub></sub>(<i>v</i>))</td></tr><tr><td><i>δ</i><sup>γ</sup><sub>s,V&#770;<sub>C<sub><i>v</i></sub></sub></sub>(<i>v</i>)</td><td>global dependency score of <i>s</i> on <i>v</i> due to all <i>t</i> ∈ (<b>V&#770;</b><sub><b>C</b><sub><i>s</i></sub></sub> ∩ <b>V&#770;</b><sub><b>C<sub><i>v</i></sub></b></sub>)</td></tr><tr><td><i>δ</i><sup>γ</sup>(<i>v</i>)</td><td>sum of all the global dependency scores (global BC) on <i>v</i></td></tr><tr><td><i>δ</i><sup>λ</sup><sub>s,•</sub>(<i>v</i>)</td><td>local dependency score of <i>s</i> on <i>v</i> due to all <i>t</i> ∈ <b>V&#770;</b><sub><b>C</b><sub><i>s</i></sub></sub> = <b>V&#770;</b><sub><b>C<sub><i>v</i></sub></b></sub></td></tr><tr><td><i>δ</i><sup>λ</sup>(<i>v</i>)</td><td>sum of all the local dependency scores (local BC) on <i>v</i></td></tr><tr><td><i>δ</i><sup>ε</sup><sub>s,•</sub>(<i>v</i>)</td><td>dependency score of <i>s</i> on <i>v</i>, as external node, due to all <i>t</i> ∈ <b>V&#770;</b><sub><b>C</b><sub><i>s</i></sub></sub></td></tr><tr><td><i>δ</i><sup>ε</sup>(<i>v</i>)</td><td>sum of all the dependency scores on <i>v</i> as external node</td></tr></tbody></table>

## 5.4.2 Brandes’ algorithm

Brandes’ algorithm is the fastest known general-purpose sequential algorithm for computing BC. It is based on the notions of *pair-dependency* and *dependency score*. Let us consider two generic nodes $s, t \in V$. Given shortest paths counts $\sigma_{s,t}(v)$ and $\sigma_{s,t}$, the pair-dependency $\delta_{s,t}(v)$ of a pair $s, t$ on an intermediary node $v \in V$ is defined as follows:

$$ \delta_{s,t}(v) = \frac{\sigma_{s,t}(v)}{\sigma_{s,t}} \qquad (5.2) $$
The pair-dependency represents the fraction of shortest paths between s and t crossing v. The dependency score $\delta_{s,\bullet}(v)$ of a vertex s on a vertex v $\in \mathbf{V}$ is then defined as follows:

$$ \delta_{s,\bullet}(v) = \sum_{t \in \mathbf{V} \setminus \{s,v\}} \delta_{s,t}(v) \qquad (5.3) $$

BC can thus be redefined in terms of dependency score:

$$ BC(v) = \sum_{s \neq v \neq t \in \mathbf{V}} \frac{\sigma_{s,t}(v)}{\sigma_{s,t}} = \sum_{s \neq v \neq t \in \mathbf{V}} \delta_{s,t}(v) = \sum_{s \in \mathbf{V}} \delta_{s,\bullet}(v) \qquad (5.4) $$

The key observation of Brandes’ algorithm is that the dependency score obeys a recursive formula that considers the dependency scores of the successors of node v, i.e., $\delta_{s,\bullet}(w)$. In particular, for each $s \in \mathbf{V}$ we have:

$$ \delta_{s,\bullet}(v) = \sum_{w:v \in \mathbf{P}_s(w)} \frac{\sigma_{s,v}}{\sigma_{s,w}} \cdot (1 + \delta_{s,\bullet}(w)) \qquad (5.5) $$

Brandes’ algorithm runs in two phases, exploiting equation 5.5. For each (source) node $s \in \mathbf{V}$, in the first phase, a Single-Source Shortest-Paths (SSSP) algorithm, based on Breadth-First Search (BFS), is executed on $\mathbf{G}$ to find all the shortest paths rooted in $s$. In the second phase, dependency scores are accumulated by backtracking along the discovered shortest paths using the recursive relation in Eq. 5.5. In backtracking, nodes are visited in descending order of distance from the source. During these two phases, for each node $v \in \mathbf{V}$, the algorithm builds and exploits the following data structures: the set of direct predecessors $\mathbf{P}_s(v)$ on shortest paths from the source, the distance $d_{\mathbf{G}}(s, v)$ from the source, the number of shortest paths $\sigma_{s,v}$ from the source and the dependency score $\delta_{s,\bullet}(v)$ that accumulates the contribution of the source on node v due to all destinations during the back-propagation step.

### 5.4.3 Equivalence class

To reduce the number of SSSP explorations in the Brandes’ algorithm, and thus lower the BC computation time, we exploited the concept of equivalence class. The algorithm leverages structural properties of graphs to find classes of equivalent nodes: by selecting one representative node for each class, it becomes possible to calculate BC by significantly reducing the number of SSSP explorations and thus accelerating its computation.

Formally, given a connected sub-graph $\hat{\mathbf{G}}$ of $\mathbf{G}$ induced by the set of nodes $\mathbf{V}_{\hat{\mathbf{G}}} \subset \mathbf{V}$, we define an equivalence class $\mathbf{K}_i$ as any subset of nodes in $\mathbf{V}_{\hat{\mathbf{G}}}$ that produce the same dependency score on all nodes - and for destinations - outside sub-graph $\hat{\mathbf{G}}$ when used as sources for SSSP explorations.

By choosing only one representative node (called pivot) for each class, the correct dependency scores of nodes can be computed by multiplying the scores computed via the SSSP rooted in the pivot by the cardinality of the class. Let $k_i$ be a pivot of $\mathbf{K}_i$ and $v \notin \mathbf{V}_{\hat{\mathbf{G}}}$, i.e., a node outside sub-graph $\hat{\mathbf{G}}$. According to our definition, we have:

$$ \sum_{s \in \mathbf{K}_i} \sum_{t \notin \mathbf{V}_{\hat{\mathbf{G}}}} \delta_{s,t}(v) = |\mathbf{K}_i| \cdot \sum_{t \notin \mathbf{V}_{\hat{\mathbf{G}}}} \delta_{k_i,t}(v) $$

which, according to our notation, can be re-written as:

$$ \sum_{s \in \mathbf{K}_i} \delta_{s, \overline{\mathbf{V}_{\hat{\mathbf{G}}}}}(v) = |\mathbf{K}_i| \cdot \delta_{k_i, \overline{\mathbf{V}_{\hat{\mathbf{G}}}}}(v) \qquad (5.6) $$

with $k_i$ representing the pivot node of equivalence class $\mathbf{K}_i$. Equation 5.6 clearly shows that a low number of classes significantly reduces the computation time, by allowing to skip a high number of SSSP explorations.
### 5.4.4 Clustering and BC computation

A possible technique to identify equivalence classes is to consider reference nodes. Given a generic sub-graph $\hat{G}$, the reference nodes in $V_{\hat{G}}$ are those that need to be traversed to reach, via shortest paths from nodes in $V_{\hat{G}}$, any other node in $\overline{V_{\hat{G}}}$.

To easily identify reference nodes, we used clustering, and to increase the chances of identifying a low number of equivalence classes, we considered a clustering technique based on modularity (i.e., the Louvain method [241]), which allows reducing the amount of connections among groups of nodes belonging to different clusters, and, consequently, lowers the number of reference nodes to be considered for discovering equivalence classes.

#### Equivalence class with clustering

Let us assume a given graph $G$ is split into a set of clusters $C$, where a single cluster $C_i$ is a connected sub-graph of $G$ induced by a set of nodes $V_{C_i} \subset V$.

For each cluster $C_i \in C$, it is possible to identify a set of border nodes $BN_{C_i}$. A border node $b_i \in BN_{C_i}$ is a node belonging to $C_i$ and having at least one neighbor belonging to another cluster, as graphically presented in Fig. 5.7 (circled nodes are border nodes).

To discover equivalence classes, for each cluster $C_i$, it is possible to group nodes based on their distance and number of shortest paths to the border nodes. To this end, we can leverage the following theorem (see [234] for the proof).

**Theorem 5.4.1** Let $k \in \mathbb{R}^+$ and $l \in \mathbb{R}$, let $C_i$ be a generic cluster of graph $G$ with border nodes $BN_{C_i}$, and $s, p \in V_{C_i}$. If $\forall b_j \in BN_{C_i}, \sigma_{s,b_j} = k \cdot \sigma_{p,b_j}$ and $d_G(s, b_j) = d_G(p, b_j) + l$, then $\delta_{s, \overline{V_{C_i}}}(v) = \delta_{p, \overline{V_{C_i}}}(v), \forall v \in \overline{V_{C_i}}$.

In other words, any given pair of nodes $s, p$ belonging to the sub-graph induced by nodes in cluster $C_i$ (i.e., $s, p \in V_{C_i}$), produces the same dependency score on all nodes $v \in \overline{V_{C_i}}$ for destinations $t \in \overline{V_{C_i}}$ if the distances and the number of shortest paths from $s$ and $p$ to every border node of $C_i$ are the same, except for an additive or multiplicative factor, respectively.

From the previous theorem, we can derive the following corollary:

**Corollary 5.4.1** If $\forall b_j \in BN_{C_i}: \hat{\sigma}_{s,b_j} = \hat{\sigma}_{p,b_j}$ and $\hat{d}_G(s, b_j) = \hat{d}_G(p, b_j)$, then $\delta_{s, \overline{V_{C_i}}}(v) = \delta_{p, \overline{V_{C_i}}}(v), \forall v \in \overline{V_{C_i}}$.

where $\hat{d}_G(s, b_j)$ represents the normalised distance of the generic node $s$ to the generic border node $b_j$, defined as follows:

$$ \hat{d}_G(s, b_j) = d_G(s, b_j) - \min_{b_k \in BN_{C_i}} d_G(s, b_k) $$

and $\hat{\sigma}_{s,b_j}$ represents the normalised number of shortest paths from the generic node $s$ to the generic border node $b_j$, and is defined as:

$$ \hat{\sigma}_{s,b_j} = \sigma_{s,b_j} / \min_{b_k \in BN_{C_i}} \sigma_{s,b_k} $$

Normalised distance and normalised shortest paths simplify the identification of classes since they are characterised by the same vector of normalised distances and shortest paths as explained below with an example and the support of a graphical representation.

Let $G$ be the simple graph reported in Fig. 5.7(a), decomposed in three clusters, each separately shown in Fig. 5.7(b). We focus on the blue cluster, denoted as $C_1$, in order to illustrate the concept of equivalence class (see Tab. 5.6 and Fig. 5.8).
<figure><img src="image_69.png" /><figcaption>**Figure 5.7:** Example of clustering.</figcaption></figure>

<figure><img src="image_70.png" /><figcaption>**Figure 5.8:** Classes of equivalent nodes in the blue cluster $C_1$</figcaption></figure>

In $C_1$, nodes **1** and **2** are border nodes (also denoted as $b_1$ and $b_2$ in Tab. 5.6). The nodes of $C_1$ are related to $b_1$ and $b_2$ according to the properties detailed in Tab. 5.6: for each node the normalised distances and normalised number of shortest paths to the border nodes are reported. According to our previous definitions, nodes **3**, **4**, **6**, **14** and **5**, **1** can be grouped in two classes respectively, whereas node **2** is assigned to a singleton class. Nodes **1**, **2** and **14** are the pivots<fnref n="16" />.

<footnote n="16">¹⁶In our previous version of the algorithm, the pivots were chosen to minimise the error. Here, as we will explain later, any node in a class can be a pivot.</footnote>
<table><caption>**Table 5.6:** Normalised distances and normalised number of shortest paths for the blue cluster $C_1$. Note that the values of $\hat{d}$ and $\hat{\sigma}$ are computed considering external node 17, as from Fig. 5.9(b).</caption><thead><tr><th>node v</th><th>⟨d<sub>C<sub>1</sub></sub>(v, b<sub>1</sub>), v⟩</th><th>⟨d<sub>C<sub>1</sub></sub>(v, b<sub>2</sub>), v⟩</th><th>⟨σ<sub>vb<sub>1</sub></sub>, v⟩</th><th>⟨σ<sub>vb<sub>2</sub></sub>, v⟩</th></tr></thead><tbody><tr><td><b>1</b></td><td>0</td><td>2</td><td>1</td><td>2</td></tr><tr><td><b>2</b></td><td>2</td><td>0</td><td>2</td><td>1</td></tr><tr><td><b>3</b></td><td>0</td><td>0</td><td>1</td><td>1</td></tr><tr><td><b>4</b></td><td>0</td><td>0</td><td>1</td><td>1</td></tr><tr><td><b>5</b></td><td>0</td><td>2</td><td>1</td><td>2</td></tr><tr><td><b>6</b></td><td>0</td><td>0</td><td>1</td><td>1</td></tr><tr><td><b>14</b></td><td>0</td><td>0</td><td>1</td><td>1</td></tr></tbody></table>

## Cluster-based exact BC computation

The equivalence classes allow us to compute the dependency score on nodes - and for destinations - that do not belong to the same cluster of the source node, which means that the contributions computed via this approach are only *partial*. To obtain the total BC, we rewrite Eq. 5.4 as follows:

$$
\begin{align}
BC(v) &= \sum_{s \in \mathbf{V}} \delta_{s, \bullet}(v) \nonumber \\
&= \sum_{s \in \mathbf{V}} \sum_{t \in \mathbf{V}_{C(s)}} \delta_{s,t}(v) + \sum_{s \in \mathbf{V}} \sum_{t \notin \mathbf{V}_{C(s)}} \delta_{s,t}(v) \nonumber \\
&= \underbrace{\sum_{s \in \mathbf{V}_{C(v)}} \sum_{t \in \mathbf{V}_{C(v)}} \delta_{s,t}(v)}_{\text{sum of local dependency scores} = \delta^{\lambda}(v)} + \underbrace{\sum_{s \in \mathbf{V}} \sum_{t \notin \mathbf{V}_{C(s)}} \delta_{s,t}(v)}_{\text{sum of global dependency scores} = \delta^{\gamma}(v)} \tag{5.7} \\
&\quad + \underbrace{\sum_{s \notin \mathbf{V}_{C(v)}} \sum_{t \in \mathbf{V}_{C(s)}} \delta_{s,t}(v)}_{\text{sum of dependency scores on external nodes} = \delta^{\epsilon}(v)} \nonumber
\end{align}
$$

In Eq. 5.7, the term $\sum_{s \in \mathbf{V}} \sum_{t \in \mathbf{V}_{C(s)}} \delta_{s,t}(v)$ is decomposed in two terms, by splitting the external sum operation (i.e., $\sum_{s \in \mathbf{V}}$) into two parts, (i.e., $s \in \mathbf{V}_{C(v)}$ and $s \notin \mathbf{V}_{C(v)}$), given that $t \in \mathbf{V}_{C(v)}$ is equivalent to $t \in \mathbf{V}_{C(s)}$, when $s \in \mathbf{V}_{C(v)}$.

As a result, we can distinguish two main components in Eq. 5.7, *local* and *global* dependency scores. The third additional term is necessary to properly take into account the possible existence of shortest paths connecting nodes of the same cluster via nodes belonging to one or more different clusters, *i.e.*, *external nodes*.

We define the *local* dependency score of a node $s$ on a node $v$, $\delta_{s,\bullet}^{\lambda}(v)$, as the sum of pair dependency scores for which source $s$, the destinations and node $v$ belong all to the same cluster. We define the *local BC* of a node $v$, $\delta^{\lambda}(v)$, as the *BC* of $v$ computed on the sub-graph $C(v)$.

Local BC is computed using Brandes' algorithm inside each cluster<fnref n="17" />, which generates, as a by-product, additional information (i.e., the number of shortest paths and distances to border nodes). This information is later used to group nodes into equivalence classes and to fasten the computation of global dependency scores.

The *global* dependency score of a node $s$ on a node $v$, $\delta_{s,\bullet}^{\gamma}(v)$, is the sum of all the pair dependency scores for which destinations do not belong to the same cluster of source node $s$. The *global BC* of the generic node $v$, $\delta^{\gamma}(v)$, is thus the sum of the global dependency scores for source node $s$ ranging over the whole set of nodes $\mathbf{V}$.

The dependency score of a node $s$ on an external node $v$, i.e. $C(v) \neq C(s)$, noted as $\delta_{s,\bullet}^{\epsilon}(v)$, is the sum of all the pair dependency scores for which destinations belong to the

<footnote n="17">¹⁷As explained later, in the special case where there are external shortest paths in the cluster, the local BC is actually computed inside the *extended* cluster.</footnote>
same cluster of the source node $s$. We denote by $\delta^\epsilon(v)$ the sum of all the dependency scores on $v$, when $v$ is an external node and the sources and destinations are in the same cluster, different from $C(v)$.

This last term $\delta^\epsilon(v)$ is equal to zero when the clustering is ideal, i.e. when all the shortest paths between any pair of nodes of a cluster only contain nodes from that same cluster. When this condition is not fulfilled, multiple side effects due to the presence of external nodes have to be taken into account, as discussed below.

## External nodes/shortest paths

Given a cluster $C_i$, two nodes $s, t \in C_i$ and two border nodes $b_1, b_2 \in C_i$, there may exist shortest paths between $s$ and $t$ that exit $C_i$ through $b_1$, cross a certain number of nodes belonging to other clusters and then re-enter $C_i$ through $b_2$. We call these shortest paths external shortest paths and the nodes lying on them that do not belong to $C_i$, $\text{EN}_{C_i}$, external nodes of $C_i$. If the existence of such external shortest paths is neglected, BC computation will be affected by an error due to the incorrect values of the lengths and the counts of shortest paths between pairs of nodes inside the same cluster. Consequently, an error in the computation of the local BC, $\delta^\lambda$, and in the identification of equivalence classes will be introduced. This was one of the approximation errors that affected the previous versions of our algorithm [237, 232]. To remove this intra-cluster error, we borrowed the idea proposed by the authors in [239]. After clustering, we build a Hierarchical Sub-Network (HSN), i.e., a sub-graph of $G$ induced by the border nodes of all the clusters and nodes lying on the intra-cluster shortest paths between pairs of border nodes of the same cluster.

By retrieving all the shortest paths between pairs of border nodes of the same cluster via the HSN, we can identify possible external nodes for that cluster. Afterwards, we can extend each cluster with the related external nodes and use the extended clusters as sub-graphs to identify equivalence classes and pivots. Thus, local BC $\delta^\lambda$ can be correctly computed inside these extended clusters instead of the initial ones. The computation of local BC on the extended clusters also allows, as a by-product, computing the pair dependency scores on the external nodes. By accumulating the contributions on each external node from different clusters (in case a node is external to multiple clusters), it becomes trivial to compute the $\delta^\epsilon(v)$ term in Eq. 5.7.

Formally, an extended cluster $C_i^*$ of a cluster $C_i \in C$ is defined as a connected sub-graph induced by nodes $V_{C_i^*} = V_{C_i} \cup \text{EN}_{C_i}$.

<figure><img src="image_71.png" /><figcaption>**Figure 5.9:** Example of external node found through the HSN.</figcaption></figure>

We provide an illustrative example to better understand how the HSN is built and how it is used to form the extended clusters. Let us consider again the clustered graph from Fig. 5.7. In cluster $C_1$, nodes 1 and 2 are border nodes, while node 4 lies on the only intra-cluster shortest path between them. In cluster $C_2$, nodes 17 and 20 are border nodes and nodes 15, 21, 19 and 16 lie on the intra-cluster shortest paths between them. Finally,
in cluster $C_3$, there is only border node 8. All the aforementioned nodes build up the HSN (see Fig. 5.9(a)). If we now consider the shortest paths between border nodes 1 and 2 via the HSN, we notice that node 17 lies on a shortest path connecting the two former nodes. Consequently, it represents an external node of $C_1$ (see Fig. 5.9(b)).

Dependency score of pivots

From the equivalence class relationship described in Sec. 5.4.3, a pivot of such a class only represents the dependency scores on nodes $v$ - and for destinations $t$ - which do not belong to its own cluster. In fact, given a cluster $C_i \in C$ and all its equivalence classes $K_{C_i}$, from Eq. 5.6, we have:

$$ \sum_{s \in K_i} \delta_{s, \overline{V_{C_i}}}(v) = |K_i| \cdot \delta_{k_i, \overline{V_{C_i}}}(v) \quad \forall v \in \overline{V_{C_i}}, K_i \in K_{C_i}. \qquad (5.8) $$

This equation can be exploited to speed up the computation of BC, by building on Brandes’ algorithm and SSSP explorations, but only holds if $v \in \overline{V_{C_i}}$. Thus, it cannot be directly applied to correctly compute values of global BC when $v$ is in the same cluster of the source. Therefore, the algorithm requires a more elaborated approach to properly and efficiently calculate the contribution from the pivot of $K_{C_i}$ to the BC of nodes $v \in \overline{V_{C_i}}$<fnref n="18" />.

First of all, let us decompose the global dependency scores from Eq. 5.7 based on the cluster of node $v$ as follows:

$$ \delta^{\gamma}(v) = \sum_{s \notin \overline{V_{C(v)}}} \sum_{t \notin (\overline{V_{C(v)}} \cup \overline{V_{C(s)}})} \delta_{s,t}(v) + \sum_{s \notin \overline{V_{C(v)}}} \sum_{t \in \overline{V_{C(v)}}} \delta_{s,t}(v) + \sum_{s \in \overline{V_{C(v)}}} \sum_{t \notin \overline{V_{C(v)}}} \delta_{s,t}(v) \quad (5.9) $$

The previous equation can be further simplified by considering the following claim (see [234] for the proof):

**Claim 5.4.1** In undirected graphs:

$$ \sum_{s \in \overline{V_{C(v)}}} \sum_{t \notin \overline{V_{C(v)}}} \delta_{s,t}(v) = \sum_{s \notin \overline{V_{C(v)}}} \sum_{t \in \overline{V_{C(v)}}} \delta_{s,t}(v) \qquad (5.10) $$

By relying on Eq. 5.10, it becomes possible, for any given node $v$, to replace with zero the sum of the pair-dependencies $\delta_{s,t}(v)$ for which $s \in \overline{V_{C(v)}}$ and $t \in \overline{V_{C(v)}}$ (the third term in Eq. 5.9) and compensate later the lack of this term by doubling the sum of the pair-dependencies $\delta_{s,t}(v)$ for which $s \in \overline{V_{C(v)}}$ and $t \in \overline{V_{C(v)}}$ (the second term in Eq. 5.9). In other words, Eq. 5.9 can be rewritten for undirected graphs as follows:

$$ \delta^{\gamma}(v) = \sum_{s \notin \overline{V_{C(v)}}} \sum_{t \notin (\overline{V_{C(v)}} \cup \overline{V_{C(s)}})} \delta_{s,t}(v) + 2 \cdot \sum_{s \notin \overline{V_{C(v)}}} \sum_{t \in \overline{V_{C(v)}}} \delta_{s,t}(v) \quad (5.11) $$

With this formulation, the computation of global BC for any given node $v$ can be performed by only considering the contributions that source nodes have on nodes $v$ and for destinations $t$, where both $v$ and $t$ are outside the source’s cluster. This satisfies the requirement for applying equivalence classes to compute dependency scores. Therefore, we can now use pivots as sources to efficiently compute the exact global BC.

In particular, let $\delta_{s, \overline{V_{C(v)}}}(v)$ and $\delta_{s, \overline{V_{C(v)}}}(v)$ be the global dependency scores from node $s$ on node $v$ for destinations not belonging to $C(s)$, but belonging to $C(v)$, and the global

<footnote n="18">¹⁸In our previous version of the algorithm, we used Eq. 5.8 without taking into account the cluster of $v$ and the pivots were chosen to minimize the error. Here, we avoid such an error during the computation of global dependency scores by exploiting the properties of undirected graphs, as explained later.</footnote>
dependency score from node s on node v for destinations not belonging to C(s) and C(v). Eq. 5.11 can be rewritten as follows:

$$ \delta^{\gamma}(v) = \sum_{s \notin \mathbf{V}_{C(v)}} \left[ 2 \cdot \delta_{s, \mathbf{V}_{C(v)}}^{\gamma}(v) + \delta_{s, \overline{\mathbf{V}_{C(v)}}}^{\gamma}(v) \right] \quad (5.12) $$

Therefore, given a cluster $C_i \in C$ and all its equivalence classes $K_{C_i}$, we have:

$$ \forall v \notin \mathbf{V}_{C_i}, \mathbf{K}_i \in \mathbf{K}_{C_i}, \quad \sum_{s \in \mathbf{K}_i} \left( 2 \cdot \delta_{s, \mathbf{V}_{C(v)}}^{\gamma}(v) + \delta_{s, \overline{\mathbf{V}_{C(v)}}}^{\gamma}(v) \right) = |\mathbf{K}_i| \cdot \left( 2 \cdot \delta_{k_i, \mathbf{V}_{C(v)}}^{\gamma}(v) + \delta_{k_i, \overline{\mathbf{V}_{C(v)}}}^{\gamma}(v) \right) \quad (5.13) $$

Eq. 5.13 means that, during the back-propagation phase, we should distinguish between contributions due to destinations inside the same cluster of v and contributions due to destinations outside the cluster of v (and the cluster of s).

<figure><img src="image_72.png" /></figure>

(a) Problem with nodes belonging to clusters of pivots (e.g., 1, 2). Dependency scores from node 14 are set to null.

(b) Solution for undirected graphs. The contributions on nodes of the blue cluster (e.g., 1, 2) are corrected by doubling dependency scores from sources of other clusters for destinations inside the blue cluster).

<caption>**Figure 5.10:** Global SSSP explorations from pivots</caption>

For a better understanding of the formulas above, let us consider an illustrative example by leveraging again the clustered graph from Fig. 5.7 and the equivalence classes of cluster $C_1$ from Fig. 5.8.

The pivot node of the equivalence class composed of nodes {3, 4, 6, 14} is node 14. According to the proposed approach, we calculate the dependency scores from node 14 on all nodes of clusters $C_2$ (orange) and $C_3$ (green) and multiply them by 4, avoiding to calculate the dependencies scores from nodes 3, 4 and 6. This way, the computation time is divided by 4. However, while it is correct to multiply by 4 the dependency scores for nodes in $C_2$ and $C_3$, it is not for nodes belonging to the same cluster of the pivot (see Fig. 5.10(a)) since nodes 14, 3, 4, 6 of the class are indeed equivalent, but only with respect to border nodes of cluster $C_1$ (i.e., nodes 1, 2). Therefore, we cannot multiply by 4 the dependency scores on the nodes internal to cluster $C_1$ (i.e., 1, 2, 3, 4, 5, 6) since these scores are not the same when computed, for instance, from node 14 or node 4. To avoid the problem, we put these dependency scores to 0 during the SSSP explorations performed from pivot node 14. These neglected scores will be later compensated (via multiplication by two) during the SSSP explorations performed from the pivot nodes of $C_2$ and $C_3$ (see Fig. 5.10(b)).
**Back-propagation.** Differently from Brandes' algorithm, it is not possible to directly express the global dependency score<fnref n="19" /> of a node $v$, $\delta_{s,\bullet}^{\gamma}(v)$, in terms of the global dependency scores of $w$, $\delta_{s,\bullet}^{\gamma}(w)$, where $v \in \mathbf{P}_s(w)$. Indeed, when $\mathbf{C}(v) \neq \mathbf{C}(w)$ (i.e., when crossing a cluster), the set of destinations of $w$ which do not belong to $\mathbf{C}(w)$ can be composed of both destinations belonging to $\mathbf{C}(v)$ and destinations not belonging to $\mathbf{C}(v)$: for the former, the pair-dependencies have to be multiplied by 2, whereas for the latter no further operation is needed (see Eq. 5.13).

To overcome this problem, we decided to apply the classic recursive formula of Brandes' algorithm (Eq. 5.5) on a *vector* of contributions while propagating the global dependency scores $\delta_{s,\bullet}^{\gamma}(v)$. The dimensions of this vector of contributions correspond to the number of clusters, so that the contribution due to a destination $t$ is assigned to $\delta_{s,\mathbf{V}_{C(t)}}^{\gamma}(v)$. Formally, we have the following recursive formula:

$$ \forall \mathbf{C}_i \in \mathbf{C} \setminus \mathbf{C}(s) : \delta_{s, \mathbf{V}_{C_i}}^{\gamma}(v) = \sum_{w:v \in \mathbf{P}_s(w)} \frac{\sigma_{s,v}}{\sigma_{s,w}} * (1_{w \in \mathbf{C}_i} + \delta_{s, \mathbf{V}_{C_i}}^{\gamma}(w)), \quad (5.14) $$

where $1_{w \in \mathbf{C}_i}$ represents a Boolean variable equal to 1 if $w \in \mathbf{C}_i$, 0 otherwise<fnref n="20" />. At the end of the back-propagation phase, we put the dependency scores of nodes $v$ belonging to the same cluster of the (pivot) source node to 0, whereas the dependency scores of nodes belonging to the other clusters are computed using the following formula:

$$ \delta_{s, \overline{\mathbf{V}_{C(s)}}}(v) = 2 \cdot \delta_{s, \mathbf{V}_{C(v)}}^{\gamma}(v) + \sum_{\mathbf{C}_i \neq \mathbf{C}(v)} \delta_{s, \mathbf{V}_{C_i}}^{\gamma}(v) \quad (5.15) $$

Finally, according to Eq. 5.13, $\delta_{s,\bullet}^{\gamma}(v)$ is multiplied by the cardinality of the equivalence class $s$ belongs to.

### 5.4.5 Empirical Evaluation with Synthetic and Realistic Networks

We compared the execution times obtained with the implementation of our algorithm (named E1C-FastBC) to those obtained with other algorithms by using the Algorithmic Speedup (AS). Given two algorithms, $a_1$ and $a_2$, the AS of $a_1$ over $a_2$ with $p$ cores, noted as $AS_p^{a_1/a_2}$, is defined as $T_p^{a_2}/T_p^{a_1}$, where $T_p^{a_2}$ and $T_p^{a_1}$ are the computation times obtained with $p$ cores using algorithms $a_2$ and $a_1$, respectively. Hence, the larger the value of AS, the faster $a_1$ is compared to $a_2$ with the same computing resources. For example, $AS_p^{a_1/a_2} = 2$ means that the time taken by $a_1$ is half the time taken by $a_2$ and therefore that $a_1$ is two times faster than $a_2$.

In particular, we compared the E1C-FastBC algorithm, labelled with $\mathcal{E}$, with Brandes' algorithm, labelled as $\mathcal{B}$ and with the solution proposed in [239], labelled with $\mathcal{H}$. We chose this algorithm for comparison because it belongs to the same category as ours (cluster-based computation) and it addresses the problem of exact BC computation.

However, due to the unavailability of source/executable code for $\mathcal{H}$, we only considered the AS metric in sequential mode, by relying on the indications provided by the authors in the paper for its computation (see Eq. 7 in [239]).

To further explore the performance of our solution, we also analysed the efficiency of the E1C-FastBC algorithm, based on the canonical definition of speedup. Specifically, the speedup obtained with $p$ cores is defined as $S_p = T_s/T_p$, where $T_s$ is the computation time in sequential mode and $T_p$ is the computation time with $p$ cores. The efficiency with $p$ cores, noted as $E_p$, is then defined as $S_p/p^{21}$.

<footnote n="19">¹⁹$\delta_{s,\bullet}^{\gamma}(v)$ is equivalent to $\delta_{s,\overline{\mathbf{V}_{C(s)}}}(v)$</footnote>

<footnote n="20">²⁰This is the part of contribution due to $w$ as a destination</footnote>

<footnote n="21">²¹Efficiency $E_p$ is strictly, but non linearly, related to the AS metric when used to compare to different algorithms. From the definition of speedup: $T_p^{a_2} = T_s^{a_2}/S_p^{a_2}$, where $T_s^{a_2}$ and $S_p^{a_2}$ are the execution time</footnote>
In all reported tests, we checked the accuracy of our solution by always observing zero error on BC values.

Datasets and Experimentation Testbed

In our experimental tests, we considered both synthetic and real graphs.

For the first category, we focused on scale-free graphs generated using the implementation of the Barabási-Albert model provided by the Python library NetworkX. According to that model, a graph of *n* nodes is grown by attaching new nodes, one at a time, each with *m'* edges that are preferentially attached to existing nodes with high degree. In our case, *m'*, which is called the preferential attachment coefficient, was chosen equal to 1. This led us to consider graphs with *m* = *n* − 1 edges and an average degree approximately equal to 2<sup>22</sup>. This choice was motivated by the features of the current implementation of our algorithm that benefits of high modularity. In other words, this class of dataset was considered as the best-case scenario. However, as mentioned in the introduction, this does not limit the applicability of our solution because many real-world systems can be represented with the Barabási-Albert model. In particular, to analyse the algorithm in terms of performance and scalability, we generated graphs with different sizes (see Table 5.7).

<table><caption>**Table 5.7:** Topological information of synthetic & real graphs. The names of the graphs are given in the first column, whereas the number of nodes and edges are given in the second and third columns. *d*<sub>avg</sub> and *d*<sub>max</sub> are the average and max degree, respectively. *cc*<sub>avg</sub> is the average clustering coefficient.</caption><thead><tr><th></th><th>Graph</th><th>n</th><th>m</th><th>d<sub>avg</sub></th><th>d<sub>max</sub></th><th>cc<sub>avg</sub></th></tr></thead><tbody><tr><td rowspan="8">Synthetic</td><td>barabási-albert</td><td>6,250</td><td>6,249</td><td>1.999</td><td>126</td><td>0.000</td></tr><tr><td>barabási-albert</td><td>12,500</td><td>12,499</td><td>"</td><td>225</td><td>"</td></tr><tr><td>barabási-albert</td><td>25,000</td><td>24,999</td><td>"</td><td>344</td><td>"</td></tr><tr><td>barabási-albert</td><td>50,000</td><td>49,999</td><td>"</td><td>463</td><td>"</td></tr><tr><td>barabási-albert</td><td>100,000</td><td>99,999</td><td>"</td><td>1,138</td><td>"</td></tr><tr><td>barabási-albert</td><td>200,000</td><td>199,999</td><td>"</td><td>676</td><td>"</td></tr><tr><td>barabási-albert</td><td>400,000</td><td>399,999</td><td>"</td><td>1,142</td><td>"</td></tr><tr><td>barabási-albert</td><td>800,000</td><td>799,999</td><td>"</td><td>1,587</td><td>"</td></tr><tr><td rowspan="5">Real</td><td>web-webbase-2001 [242]</td><td>16,062</td><td>25,593</td><td>3.187</td><td>1,679</td><td>0.224</td></tr><tr><td>ego-twitter [243]</td><td>22,322</td><td>31,823</td><td>2.851</td><td>238</td><td>0.072</td></tr><tr><td>internet [242]</td><td>124,651</td><td>193,620</td><td>3.107</td><td>151</td><td>0.062</td></tr><tr><td>lyon-road-network</td><td>156,102</td><td>178,845</td><td>2.291</td><td>8</td><td>0.017</td></tr><tr><td>email-euAll [244]</td><td>224,832</td><td>339,925</td><td>3.024</td><td>7,636</td><td>0.079</td></tr></tbody></table>

For the second category, we focused on real graphs<sup>23</sup> available in public datasets. Table 5.7 reports all the graphs used, together with some relevant properties. In particular, for each graph, the table reports the average degree (*d*<sub>avg</sub>), the max degree (*d*<sub>max</sub>) and the average clustering coefficient (*cc*<sub>avg</sub>).

All the datasets, except the one related to the Lyon road network<sup>24</sup>, are scale-free graphs. For the sake of brevity, we only report results related to the synthetic graphs and the l yon-road-network dataset in this document.

in sequential mode and the speedup obtained with *p* cores of algorithm *a*<sub>2</sub>. Similarly, $T_p^{a_1} = T_s^{a_1}/S_p^{a_1}$. By these equations, the AS metric can be rewritten as: $AS_p^{a_1/a_2} = (T_s^{a_2}/T_s^{a_1}) \cdot (S_p^{a_1}/S_p^{a_2})$. Since for a given number of cores *p*, the ratio of the speedups equals the ratios of efficiency of the two algorithms, i.e., $S_p^{a_1}/S_p^{a_2} = E_p^{a_1}/E_p^{a_2}$, we thus have $AS_p^{a_1/a_2} = (T_s^{a_2}/T_s^{a_1}) \cdot (E_p^{a_1}/E_p^{a_2})$. The relationship between AS and $E_p$ suggests that if $a_1$ and $a_2$ have comparable efficiency with *p* cores, then the AS only depends on the ratio of the execution times of the two algorithms in sequential mode, thus providing interesting insights to analyse the two different solutions comparatively.

<footnote n="22"><sup>22</sup>Each link contributes 2 to the sum of the degrees. Thus, the sum of all degrees in the network is 2 · (N − 1). For large values of N, the average degree approaches 2.</footnote>

<footnote n="23"><sup>23</sup>For each graph, we extracted the largest connected component. Then, the latter was converted in an unweighted undirected graph.</footnote>

<footnote n="24"><sup>24</sup>The l yon-road-network dataset was supplied by the French National Institute of Geographic Information (IGN). https://www.ign.fr</footnote>
The tests were conducted on the PROMENADE-v2.0 platform prototype (Sec. 5.3.1), using the Spark BP component for our experiments. It was deployed on a server machine with 128 GB of RAM and 2 sockets Intel(R) Xeon(R) CPU E5-2680 v4 @ 2.40GHz, with 14 physical cores and 2 threads per core for a total of 28 logical cores per socket and 56 virtual cores in hyper threading, running Linux Debian as operating system. Both Brandes' algorithm and E1C-FastBC were implemented in Scala.

## Synthetic Graphs Analysis

Fig. 5.11 shows the algorithmic speedup of E1C-FastBC over Brandes’ algorithm, $AS_p^{\epsilon/\beta}$, obtained on the synthetic graphs in both sequential and parallel modes. In particular, we doubled the number of nodes from 25,000 to 800,000, and we considered a number of cores $p$ equal to 1, 5, 10, 15, 20 and 25. We estimated by log-log regression the computation times with Brandes’ algorithm for the graphs with 400,000 and 800,000 nodes, since executions would have required weeks to complete, whereas our algorithm completed in maximum 31.5 minutes and 1.64 hours, respectively (sequential mode).

<figure><img src="image_73.png" /><figcaption>**Figure 5.11:** Comparison with Brandes’ algorithm. Algorithmic speedup analysis - $AS_p^{\epsilon/\beta}$=[1,5,10,15,20,25]</figcaption></figure>

As shown in Fig. 5.11, $AS_p^{\epsilon/\beta}$ increases with the size of the graph, meaning that E1C-FastBC is not only faster than Brandes’ algorithm but its speedup grows with larger graphs. This is due to the fact that the computation of our algorithm is strongly dependent on the number of border nodes ($|BN|$), pivots ($|K|$) and external nodes ($|EN|$), in addition to the number of nodes ($n$) and edges ($m$). The first two variables increase slowly compared to the number of nodes and edges, while the third is almost always zero (only in one case it was equal to 2). The drawback is that $AS_p^{\epsilon/\beta}$ decreases as the number of cores increases. This behaviour is due to the fact that the Brandes’ algorithm is more efficient than E1C-FastBC (see Fig. 5.12).

This means that the ratio $E_p^{a_1}/E_p^{a_2}$ in the relationship between the AS metric and the efficiency is lower than 1. Consequently, the AS value in the sequential case is not preserved as the number of cores increases. However, as the following efficiency analysis will further clarify, this does not mean that E1C-FastBC is less scalable than Brandes’ algorithm, but rather that it needs very large graphs to better exploit the available computing resources. This statement is also confirmed by Fig. 5.11, which clearly shows that when the graph size is 400,000, a higher number of cores performs even better than a smaller one: in particular, we have that the $AS_p^{\epsilon/\beta}$ is better with 5 cores than with 1 core. To have a similar behaviour even for a number of cores greater than 5, we should consider larger graphs.

To better understand the performance of E1C-FastBC, we investigated its efficiency with respect to that of Brandes’ algorithm. Fig. 5.12 reports the results of the efficiency analysis
<figure><img src="image_74.png" /><figcaption>**Figure 5.12:** Comparison with Brandes' algorithm. Efficiency analysis - $E_p=[1,5,10,15,20,25]$</figcaption></figure>

performed for the two algorithms. In both cases, it is possible to observe that: (i) the efficiency decreases as the number of cores increases and (ii) for a given number of cores, it increases as the number of nodes increases.

However, it is worth highlighting that in the efficiency analysis, we use different but overlapping ranges of values for the number of nodes. In particular, for our solution, we select larger graphs since we aim to show that our algorithm scales well, especially with very large graphs. In fact, the efficiency trend is almost the same in the two cases reported in Fig. 5.12. Moreover, given the maximum values of the number of nodes for the two algorithms (800,000 for ours, 200,000 for Brandes'), efficiency values are approximately the same with 5 cores (i.e., the first considered parallel configuration) but significantly diverge as the number of cores increases. In particular, efficiency of *E1C-FastBC* decreases with a higher rate than that of Brandes.

The reason for this behaviour lies in the reduced amount of computation required by our solution, which can be satisfied with limited parallelism (i.e., number of cores). Indeed, pivots allow to significantly decrease the number of (modified) Brandes' SSSP explorations performed on the whole graph, which represent the heaviest part of the whole computation, thus reducing the workload of each core.

Our solution also introduces another benefit: it allows to mitigate the variability of the computation times due to the different topological characteristics of the graphs and to the partitioning of data performed by Spark during executions. Indeed, there may exist some partitions of the RDDs characterized by a high concentration of nodes that generate the most complex shortest path trees. The time required to process these partitions directly impacts the time required to process the whole RDD, since partitions are processed in parallel. However, Spark tasks process each single partition sequentially. This aspect, combined with the fact that the number of partitions of an RDD is always equal to the number of cores and the default partitioning scheme of Spark distributes data evenly across the partitions, explains the punctual efficiency drops that can be observed in the plot related to Brandes' algorithm, when using graphs with 50,000 and 100,000 nodes and a low number of cores (i.e., the cavity that can be observed for 50,000 and 100,000 nodes with 5-10 cores in the 3d plot of Fig. 5.12(b)).

Fig. 5.13(a) reports the algorithmic speedup of *E1C-FastBC* over Brandes' algorithm, $AS_{p=1}^{\mathcal{E}/B}$, alongside with the algorithmic speedup of the approach in [239] over Brandes, $AS_{p=1}^{\mathcal{H}/B}$, on synthetic graphs and in sequential settings. $AS_{p=1}^{\mathcal{H}/B}$ is analytically computed based on Eq. 7 provided in [239].

Using such equation, it is possible to observe that: (i) $AS_{p=1}^{\mathcal{H}/B}$ depends on the number
of clusters (|C|) and the average degree ($d_{avg}$), and (ii) when $|C| + 2 \gg d_{avg}/2$, it can be approximated with $d_{avg}/2$. Therefore, since for synthetic graphs the average degree is constant and the number of clusters increases with the number of nodes, $AS_{p=1}^{H/B}$ is always approximately equal to 1 (the average degree is 2). In particular, the higher the number of clusters, the closer to 1 the $AS_{p=1}^{H/B}$. This means that the algorithm proposed in [239] cannot improve that of Brandes. Conversely, ours improves Brandes by a large multiplicative factor. We can thus conclude that our solution always outperforms the one in [239] with synthetic graphs.

## Real Graphs Analysis

Fig. 5.13(b) reports the results of the analysis of $AS_{p=1}^{\mathcal{E}/B}$ and $AS_{p=1}^{H/B}$ carried out on real graphs. $AS_{p=1}^{H/B}$ is computed again using Eq. 7 provided in [239]. In all cases, our solution outperforms the one in [239].

To further confirm the considerations on the scalability of our solution, reported in the previous section, we analyse in the following both the algorithmic speedup and efficiency values of E1C-FastBC on the *lyon-road-network* graph, for which we observed a very high number of pivots (about 60% of the number of nodes) and the lowest algorithmic speedup factor. As shown in Fig. 5.14(a), the AS is always greater than 1, thus confirming the usefulness of our solution, although the reported values are not comparable to those obtained on synthetic graphs (see Fig. 5.11) with a similar number of nodes (100,000 and 200,000). In spite of this, the algorithm becomes more scalable and efficient than in the case of synthetic graphs with 100,000 and 200,000 nodes due to the increased amount of computation resulting from the higher number of border nodes, pivots and external nodes (see Fig. 5.14(b)).

In conclusion, our performance analysis confirmed the possibility to use E1C-FastBC for faster computation of BC and motivated its integration in the PROMENADE-v2.0 platform core. Specifically, we observed computation times in the order of few minutes to few seconds for analysing a large road network (more than 150,000 nodes and 170,000 links) such as that of Lyon, France, depending on the available hardware, and with limited parallelism. These figures indicate that our solution can effectively support (quasi) real-time monitoring of large networks, thereby representing a useful tool in the context of resilience enhancement of rapidly varying networks.

<figure><img src="image_75.png" /><figcaption>(a) Algorithmic speedup analysis with synthetic graphs in sequential mode - $AS_{p=1}^{\mathcal{E}/B}$ and $AS_{p=1}^{H/B}$. $AS_{p=1}^{H/B}$ is always equal to 1 for synthetic graphs (see Subsection 5.4.5). We therefore add only one single zoomed bar in the figure just to make it visible to the reader.</figcaption></figure>

<caption>(b) Algorithmic speedup analysis with real graphs in sequential mode - $AS_{p=1}^{\mathcal{E}/B}$ and $AS_{p=1}^{H/B}$</caption>

**Figure 5.13:** Comparison with the HSN solution for synthetic (a) and real graphs (b).
<figure><img src="image_76.png" /><figcaption>**Figure 5.14:** Performance analysis of the E1C-FastBC algorithm versus Brandes on the real graph
lyon-road-network.</figcaption></figure>

Moreover, our solution could benefit from customised clustering methods for more effec-
tive identification of border nodes in (synthetic and real) graphs with different topologies
and further improve its scalability. The efficiency of E1C-FastBC could be also improved by
means of a better mapping of the algorithm on distributed resources, when data-parallelism
is exploited, and by improving locality especially when different Spark executors are used.

5.5 Multi-agent-based Large-Scale Traffic Control

As a second case study of the PROMENADE platform, we proposed a multi-agent macro-
scopic solution for large-scale traffic control in urban environments [245].

Traffic control strategies are fundamental to reducing congestion and pollution, but they
require to address multiple challenges to be effective, especially if the goal is to address
large traffic networks encompassing multiple and heterogeneous road axis and multi-modal
mobility options. This is typical in urban areas of large metropolitan modern cities. More-
over, besides efficiency and scalability, the control strategy should respect user privacy by
relying on a minimum to null amount of personal information, preserving accessibility by
not penalising specific travellers or city areas, and introducing minimal burden to individu-
als with respect to the travel experience. Furthermore, large-scale control strategies should
help achieve resilience objectives, i.e., by reducing circulation in highly-vulnerable areas and
enforcing effective rerouting in the presence of accidents or for emergency evacuation.

Traditional approaches to traffic control include perimeter control [246] and optimal
route guidance [247]. Despite their proven effectiveness in controlled environments and
some successful small-scale on-the-field application, such approaches usually suffer in real-
world settings, either in terms of application to large metropolitan networks (i.e., they tend
to simply move congestion outside the controlled areas) or in terms of user acceptability as
they depend completely upon the users following the rerouting suggestions. Moreover, such
solutions typically rely on centralised approaches, which significantly reduces their scalability
[248].

Based on these considerations, we designed and experimented, in the last part of Ce-
cile Daniel's thesis, a novel decentralised control strategy, named COoperative Multi-agent
FramewOrk for large-scale Routing-based Traffic control (COMFORT). The solution aims
at homogeneously distributing road traffic between zones of a large urban network. The
idea is to periodically assess, in a real-time like fashion, the traffic conditions locally to each zone, by relying on traffic sensors (such as loop detectors) and the MFD theory [249]. Local information about traffic congestion is shared between zones via gossiping to guarantee rapid reaction to congestion formation and indirectly transmitted to vehicles through rerouting suggestions. The choice of being rerouted and the possible new route computation are performed at the vehicle level, ensuring privacy by design<fnref n="25" /> and leaving the final choice of rerouting to the driver<fnref n="26" />.

Differently from [248] and by fully adhering to the principles of the PROMENADE architecture in terms of scalability requirements, we adopted a fully decentralised approach. The latter is based on a cooperative Multi-Agent System (MAS) paradigm partially inspired by the work of Lequay et al. [250] related to reducing energy consumption in an electricity distribution network. We designed and implemented a simulation-based prototype of COMFORT as part of the PROMENADE platform.

According to the MAS paradigm, vehicles were considered in COMFORT as one specific kind of agent. Geographical zones were also considered as agents. They process information related to the accumulation of vehicle agents within the corresponding spatial region and accomplish specific interactions (mainly control actions, but also information sharing) with vehicles and neighbouring zone agents. This inherently defines a multi-layer framework. The main interactions between the two considered kinds of agents are reported in the activity diagram of Fig. 5.15.

<figure><img src="image_77.png" /><figcaption>**Figure 5.15:** Activity diagram of the interactions between agents in COMFORT.</figcaption></figure>

The zone agents are associated with attributes that represent the main time-dependent variables involved in the control strategy: (i) total flexibility, $f_i^{tot}(t)$, a time-dependent function inversely related to the congestion state of the zone, i.e., the difference between the

<footnote n="25">²⁵In our approach, vehicles do not communicate their final destination or the route they are following. Instead, they only signal the next area they are willing to traverse to the zone they are currently in.</footnote>

<footnote n="26">²⁶It is worth noting that COMFORT assumes a share of travelling vehicles to be connected via the urban infrastructure or via a software mobile application available to the driver. Connected vehicles are the only vehicles that can receive recommendations for rerouting from the area they are travelling through. On the contrary, traffic conditions (estimated via sensors) include both connected and non-connected vehicles. To model different scenarios of connectivity and user acceptability, we tested different values for the probability of a vehicle accepting a rerouting suggestion.</footnote>
critical accumulation (derived from the MFD of the zone) and the current observation of traffic conditions in the zone<fnref n="27" />; (ii) actual flexibility, $f_i(t)$, representing the evolution of traffic accumulation in the zone, computed over consecutive discrete time steps and corresponding to the number of vehicles actually accepted by the given zone during the last time step; (iii) grade, $g_i(t)$, expressing the ability of the zone to cooperate with its neighbours by reliably engaging its flexibility over time, i.e., by accepting a sufficiently high number of vehicles with respect to its total flexibility. The core concept of the control strategies is to generate rerouting advice for vehicles before they enter a congested area. This is achieved by computing a rerouting timeout, which depends on a combination of the following variables: the level of congestion (i.e., lower total flexibility) and the reliability (i.e., the grade) of the neighbouring zone. Essentially, the more congested and less reliable a given neighbouring zone is, the sooner a vehicle expected to enter it will be advised to compute an alternative route, avoiding that zone. A detailed description of the aforementioned variables and of the rerouting strategy can be found by the interested reader in [245].

The preliminary prototype of COMFORT has been evaluated through simulation on a grid Manhattan network with 3,712 nodes and 10,324 edges, partitioned into 25 (5X5) zones, enabling a comparison of our control strategy performance depending on the zone size and density. We considered the Total Travel Time Reduction (%TR), expressed in percentage, to evaluate the impact of our strategy and compare the results to a no-control baseline:

$$ \%TR = \frac{TTT_{nocontrol} - TTT_{control}}{TTT_{nocontrol}} \qquad (5.16) $$

with $TTT_i$ representing the total travel time, i.e., the sum of the duration of all trips generated during the simulation, when using a specific control strategy $i$. A positive %TR indicates the control strategy has been capable of improving the reference no-control scenario, whereas a negative %TR points out a pejorative impact of the control strategy with respect to the reference situation. An additional element to be considered to evaluate the effectiveness of the control strategy is the number of generated reroutings. We hypothesised that rerouting represents a cost for the driver, and thus the strategy should be tuned to maximising %TR with a low number of vehicle rerouting.

In the Manhattan scenario, we measured a consistent improvement in terms of %TR over the reference no-control situation, with respect to different combinations of hyper-parameters. The best configuration allowed achieving a %TR up to 16%, for a total number of approximately 2500 reroutings (with respect to 16,743 vehicles considered in the simulation).

Similar analyses were performed also in the realistic case of the road network of Lyon, France (see Fig. 5.16). We observed an improvement in terms of %TR up to 18% in the best hyper-parameter configuration. Very few reroutings were necessary to improve the reference situation. The maximum percentage of rerouted vehicles was 6.7%, corresponding to 4,700 rerouted vehicles.

The results reported in this section are encouraging and represent a first step towards an efficient and privacy-aware real-time traffic control solution for large-scale urban networks. Additional efforts are needed to refine the control strategy and prove the applicability and efficiency of the proposed solution. In particular, the prototype requires an extensive evaluation in terms of the actual acceptability rate of the rerouting suggestions by travellers as well as the critical amount of connected vehicles necessary to significantly reduce the overall urban congestion. In that respect, we are considering a real-world deployment of the COMFORT prototype on the PROMENADE platform, which would require the real-time collection of traffic flow data, jointly with the development of a mobile application to provide a real-time rerouting service to vehicles. We are also currently exploring the potential integration of the proposed approach with our existing solutions for efficiently computing

<footnote n="27">²⁷A zone has a high flexibility when it can accept a high number of vehicles with respect to its critical capacity.</footnote>
<figure><img src="image_78.png" /><figcaption>**Figure 5.16:** The northern region network of the conurbation of Lyon (3rd, 6th district of Lyon and Villeurbanne). The network contains 1,883 nodes, 3,383 links and is divided in 17 zones. Each zone is represented by a distinct colour.</figcaption></figure>

betweenness centrality. The primary idea involves incorporating an additional term, representing a time-varying zone-criticality indicator, into the rerouting trigger. This indicator is computed by aggregating the BC values of nodes within neighbouring zones. The aim is to mitigate the pressure on areas that could become more prone to grid-locks at any given moment due to a high traffic concentration and limited alternative routes available to bypass emerging traffic bottlenecks or capacity drops (e.g., due to accidents or other disruptive events). Lastly, to further minimise reaction time and introduce proactive measures in managing traffic congestion and disruptions, we are working on developing traffic forecasting solutions (see Chapter 6) that could be integrated into the control strategy.
# Chapter 6

# Conclusion and Perspectives

In this document, I have attempted to summarise the main results I have attained in the past seven years of research within the interdisciplinary and challenging domain of human mobility resilience. This contribution has been presented in terms of three main methodolog-ical axes: (i) mining of large and heterogeneous data, (ii) modelling and analysis of complex and large networks, and (iii) development of real-time big data platforms and algorithms.

These axes have been explored in the context of four different aspects of resilience: (i) learning from the past, (ii) anticipating what to expect, (iii) monitoring the critical properties of the system, and (iv) responding with the help of proper decision-making tools.

More specifically, the contribution that I consider the most relevant among the different results achieved are the following:

* We demonstrated the possibilities offered by passively collected mobile phone data, and, particularly, Network Signalling Data, to accurately estimate travel demand and traffic performance indicators, as well as precisely reconstruct mobility trajectories in urban environments.

* We have developed state-of-the-art approaches to k-anonymise trips reconstructed from mobile phone data, with minimal loss of information and reduced computation time, paving the way to the possibility of making the resulting datasets available to the community in an open, privacy-aware aggregate format (dynamic Origin-Destination matrices).

* We have combined solutions from two different but increasingly close fields - data science and complex network theory - to identify metrics and approaches that can help perform real-time monitoring of transportation networks and anticipate the occurrence of disruptive events.

* We have combined historical data analysis with numerical programming to identify resilient and sustainable on-demand flexible mobility services (i.e., Park-and-Ride systems), which can complement the existing transport offer with a solution that can withstand a large set of likely-occurring disruptive scenarios.

* We have leveraged distributed computing, cloud and big data technologies to build the foundation of a novel extensible platform for real-time monitoring and control of large-scale networks, and, more in general, complex and interdependent urban systems.

* We have developed a state-of-the-art algorithm for the fast computation of Betweenness Centrality, a fundamental metric for discovering opinion leaders and vulnerable components of different kinds of networks.
Once more, I am grateful to all the researchers I have had the opportunity to collaborate
with over the past years, as their invaluable contributions deserve a substantial portion of
credit for these achievements.

In the rest of this chapter, I will highlight some ongoing works, presented as short-term
perspectives (Sec. 6.1), that will further develop the main findings described above. The
planned short-term research has fairly well-identified objectives and methodological bases,
and I, therefore, expect conclusive results in the next few months and years. In section 6.2,
I will instead outline some broader scientific directions that are likely to drive my research
in the longer term, and for which the boundaries and precise methodologies are inevitably
less clear. Nevertheless, I will try to identify more specific research themes or questions to
help implement my research plan and facilitate the achievement of concrete results in the
long term.

## 6.1 Short-term Perspectives

In the short term, my primary objective is to finalise and extend the work conducted within
various research projects I have recently been or am still involved in. This section provides
an overview of this ongoing work, capturing its current state and outlining my immediate
future directions. These activities revolve around three main topics: (i) travel demand
analysis and modelling with a special focus on multi-modal demand variations in the presence
of disruptive events; (ii) accurate short-term traffic forecasting, in multi-modal settings
and with respect to atypical events; and (iii) synthetic generation of presence and multi-
modal mobility data. It is worth underlining that these themes all relate to some extent
to enhancing the resilience of multi-modal transport networks, as will be further detailed
below.

The first and second research subjects are strongly interrelated, as they globally aim to
develop accurate models of multi-modal transport network dynamics by advancing the un-
derstanding of travel demand changes during non-recurrent events with respect to recurrent
ones. These advancements will ultimately contribute to the state-of-the-art in the field of
mobility management during atypical events. For the first topic, we are leveraging unsu-
pervised machine learning solutions to identify and characterise travel demand variations in
the presence of disrupted mobility episodes by mining heterogeneous data sources covering
large periods of observations. Concerning the second one, we are exploring supervised deep
learning frameworks to enhance the quality of short-term traffic prediction. Specifically, we
are focusing on graph neural networks and transformer-based architectures, which represent
the state-of-the-art in the field. The initial results of these studies are reported in the two
subsections of Sec. 6.1.1.

The third research subject focuses on the generation of synthetic data by employing sta-
tistical inference techniques to fuse multiple available sources of real data, such as mobile
telephony data, population statistics, survey data, and OD matrices. The goal of this re-
search is to obtain rich and reliable synthetic population, travel demand and transport-mode
choice data that closely mirror real-world conditions, as reported in the input data. These
synthetic dataset can have useful applications in terms of resilience, enabling researchers to
conduct in-depth and realistic vulnerability analyses and experiments without compromising
user privacy. Ultimately, this research shall contribute to developing more robust and sus-
tainable transportation design and management strategies. Our current approach, further
detailed in Sec. 6.1.2, not only has the potential to enhance the accuracy of simulation tools
but also ensures the privacy and anonymity of individuals represented in the generated data.

Parts of the following text are drawn from papers that are currently being written or
under review and might be published as conference proceedings or journal papers between
the writing of this chapter and the public defence of my HDR. Whenever possible, references
to pre-print or preliminary versions are provided within the text.
6.1.1 Short-term Multi-modal Traffic Forecasting: looking for the
atypical

Forecasting traffic-related variables, such as flow, speed or travel demand, plays a critical role in traffic monitoring and control as it allows for anticipating changes and proactively developing mitigation strategies to handle potential disruptions and ensure resilient mobility. In our ongoing research, we aim to tackle the challenge of understanding, modelling and predicting fluctuations in traffic variables during highly dynamic and non-recurrent situations, such as peak hours, holidays, and special or disruptive events. We believe that this topic represents an under-explored area of research that should occupy a central position at the intersection of traffic modelling and transport resilience engineering.

To achieve this objective, we are currently investigating the dynamics of multi-modal
demand during atypical events as well as the integration into deep-learning forecasting ar-
chitectures of multi-modal travel demand data and context-specific information, such as
weather conditions, planned road closures, accidents, and public events. Our aim is to de-
velop innovative machine learning solutions capable of accurately forecasting<fnref n="1" /> multi-modal
traffic variations in both regular and atypical traffic scenarios. Ultimately, such prediction
tools can be integrated into scenario analyses and decision-making solutions aimed to op-
timise the response of the transport network to stress situations. Existing literature on
this topic is limited [251]. Traditional forecasting methods are typically designed for stable
time frames and specific geographical areas, neglecting rarer and more uncertain situations.
Extreme examples of such dynamic conditions that we would like to address include the
COVID-19 pandemic, population protests, and exposure to climate-change-related events
like heavy rainfall, floods and heatwaves.

The preliminary outcomes of our research on the topic are reported below and have
been mostly obtained in the context of the ANR PRCE MOBITIC project and the PhD
theses of Ali Benam and Romain Rochas. These results represent our first steps towards
answering the following research questions: (i) Can we automatically identify and model
categories of irregular travel demand behaviours and use them to further improve short-
term multi-modal traffic prediction? (ii) Which deep-learning frameworks are most suitable
for short-term traffic prediction? (iii) Are graph neural network solutions appropriate for
predicting multi-modal travel demand? (iv) What level of performance can be achieved in
perturbed and unstable mobility scenarios?

Moreover, the study of travel volume variations and modal shifts due to abnormal events is at the core of the research that will be carried out within the scope of Paul Denailly's postdoctoral fellowship. This research is funded by the ANR MOBITIC project and supervised in collaboration with the GRETTIA lab at the University Gustave Eiffel, focusing on the identification and understanding, via data fusion techniques, of the impacts on mobility induced by the COVID-19 pandemic.

Multi-source data-driven identification and classification of multi-modal travel demand variations during atypical events

We have developed a preliminary methodology to automatically detect and analyse multi-
modal travel demand variations following disruptive events. Currently, the methodology can
be applied to perform anomaly detection at a network-wide scale and is partially inspired
by the mobile-phone signature-based approach described in Sec. 2.1.2.

Previous research [252, 253] has examined the passenger responses to disruptions in
transport systems using surveys and mobile or public transport ticketing data to analyse the

<footnote n="1">¹We consider traffic forecasting as related to predicting, e.g., OD pairs travel demand or average speed at the level of road sections, zones or the entire transport network, depending on the specific case study. We refer to short-term prediction with respect to prediction horizons in the order of minutes or hours, by relying on historical time series of several hours, days or a week.</footnote>
demand dynamics and shifts between two specific modes during certain events. However, current literature lacks an integrated, automated approach that goes beyond single-kind disruption events and accounts for multiple urban transport modes. Our work provides a general data-driven methodology for the identification and categorisation of complex multi-modal demand dynamics during disruptive events.

The developed methodological framework comprises four main steps, as detailed in Figure 6.1, and has been evaluated in a case study related to the multi-modal transport demand for the city of Lyon, France. A brief summary is provided below, while more details can be found in [254].

<figure><img src="image_79.png" /><figcaption>Figure 6.1: Methodology flowchart</figcaption></figure>

In step 1, a compressed representation of the typical demand for each transport mode is constructed by unifying multi-source datasets. Traffic observations are aggregated hourly, per each mode. A unified dataset is thus produced with demand counts, denoted as $D_{\text{demand}} = \{N_t^q\}$. $N_t^q$ represents the observed demand count for mode q on the generic 1-hour time slot t from the data availability period. In our case study, $q \in Q = \{\text{bus, tram, metro, bike, car}\}$.

In step 2, the methodology defines a rolling signature $S_{(m,d,h)}^q$ that represents the expected demand for the generic mode q over a given reference period (e.g., a given week m, on the generic hour h of day of the week d). The definition of signature is based on a compressed representation combining the hourly mean ($\mu$) and standard deviation ($\lambda$) counts computed from past and future demand observations. The temporal support ($W_m^K, d, h$) used for the aggregation includes a subset of weeks ($W_m^K \subset M$), with K indicating the number of past and future weeks (with respect to week m) used for the aggregation, selected weekdays ($d \in \Delta$), and selected hours of the day ($h \in H$).

$$S_{(m,d,h)}^q = \mu_{(W_m^K, d, h)}^q \pm \lambda_{(W_m^K, d, h)}^q \quad (6.1)$$

An example of a signature for each transport mode is reported in Fig. 6.2.

In step 3, the continuous variable $\delta_t^q$ is computed as formulated in Equation 6.2. $\delta_t^q$ represents the standardised deviation of the observed demand from the expected demand for a particular mode at a given hour. By pivoting the data on the deviance variable $\delta_t^q$, we produce a dataset $D_{\text{deviance}}$ with transport modes for columns, hours for rows, and standardised deviance of expected demand for values.

$$\delta_t^q = \frac{N_t^q - \mu_{(W_m^K, d, h)}^q}{\sigma_{(W_m^K, d, h)}^q} \quad (6.2)$$
<figure><img src="image_80.png" /><figcaption>**Figure 6.2:** The weekly signature of counts and observations for a week in March 2019</figcaption></figure>

In step 4, as we are solely looking through anomalous hours for their categorisation, we omit hours in which demand in all modes falls inside the expected signature in a new dataset denoted as $D_{anomaly}$. This dataset is produced by filtering out regular hours; thus, it only includes hours where at least one mode's absolute value of deviance $|\delta_t^q|$ is larger than the range amplitude of the signature $\lambda_{(W_m^K, d, h)}^q$. We then apply clustering to this dataset to identify groups of time slots corresponding to categories of anomalies regarding multi-modal demand dynamics.

For the evaluation of the methodology, we focused on road traffic GPS observations, public transport ticket validations (tap-in for trams, metro, and bus lines) and shared-bike trips in the greater Lyon area, France, between 2019 and 2020. These sources of data have been simplistically leveraged as proxies of travel demand for the corresponding mode (i.e., motorised vehicles, transit and cycling).

For each identified cluster, we selected all hours belonging to the cluster to study how that cluster manifested a particular multi-modal dynamic. We averaged the deviance $\delta_t^q$ for each mode in the given cluster and produced a multi-modal demand profile for each cluster. We then used radar plots to illustrate these multi-modal profiles. Additionally, we labelled the dates and the hours of anomalies belonging to each cluster to explore their temporal correspondence to traceable events.

Fig. 6.3 shows three clusters of anomalous data using radar plots to illustrate the average deviance of each mode in each cluster. Cluster 1 indicates lower public transport demand on national holidays, Cluster 2 shows decreased shared-bike use during periods of higher precipitation, and Cluster 3 represents metro service disruptions with passengers shifting to other modes. The findings were validated using news and weather data, supporting the explainability of the multi-modal dynamics within the clusters.

The results presented in this section instil confidence in the potential to develop a comprehensive framework for detecting and analysing the spatio-temporal dynamics of the multi-modal travel demand. Ultimately, our goal is to inspect the retrieved anomalous dynamics and infer a general model of modal shifts in the presence of the associated non-recurrent
<figure><img src="image_81.png" /><figcaption>**Figure 6.3:** Average deviance of each mode per cluster (three largest identified clusters).</figcaption></figure>

events, at various spatio-temporal granularities. While the current approach focuses on network-wide resolution, we are exploring the feasibility of defining reliable signatures at the level of public transport lines and geographical areas (e.g., TAZ, IRIS zones, etc.) with promising initial results. By defining spatially-anchored multi-modal signatures, it should become possible to explore smaller-scale anomalies in travel demand, such as those due to localised events and disruptions. This approach is expected to facilitate the study of demand variation propagation in both space and time, thus offering the opportunity to identify areas, transport modes, or lines that are most vulnerable, i.e., susceptible to disruptions, and with less accessibility, i.e., with reduced possibility of modal shifts. These areas should be the target of transport policy aimed at reinforcing the mobility offer of the city by preserving equal access to transportation means.

Our short-term plans also involve integrating the signature-based approach with deep-learning forecasting models, which will be presented in the next subsection. By training such models with historical data related to multi-modal travel demand variations and contextual details of the corresponding anomalies, an ambitious objective would be to create forecasting tools capable of anticipating the increase or decrease in travel demand for one or multiple modes when specific disruptions are likely to occur or known to have recently occurred in an area of the city or network-wide. This research shall extend the state-of-the-art in deep learning-based traffic forecasting solutions for transport resilience, by enabling accurate predictions in both recurrent and non-recurrent scenarios and proactive decision-making in response to disruptive events.

## Travel demand prediction during atypical events

The importance of traffic forecasting applies to all modes of transportation, including public transport and soft modes such as private and shared cycling. By accurately knowing such demand and anticipating its variations over time, transport operators can plan, develop and adapt the city mobility infrastructure to accommodate varying needs, promoting active and sustainable transportation, and reducing dependence on cars. Furthermore, incorporating context information into traffic prediction models is crucial for improving accuracy and enhancing the usefulness of forecasts. Factors like weather conditions, such as rain, snow, or extreme temperatures, can significantly impact travel demand and traffic patterns. For example, during inclement weather, travellers may choose to abandon active modes and switch to public transportation or private cars, leading to shifts in demand and traffic congestion. By integrating meteorological data into forecasting models, transportation authorities can anticipate these changes and implement appropriate measures, such as adjusting public transportation schedules or providing real-time information to connected vehicles. Special events, such as concerts, sports games, or festivals, can also have a significant impact on
traffic patterns and demand. Predicting the influence of these events allows transportation authorities to prepare for increased demand, provide additional transportation services, or implement temporary traffic management measures to minimise disruptions. In a similar way, incorporating real-time data on accidents, road closures, or construction activities can improve forecasting accuracy, allowing authorities to inform travellers about alternative routes to mitigate congestion.

With these objectives in mind, we are currently investigating the performance of state-of-the-art Graph Neural Networks (GNN) and Graph Convolutional Neural Networks (GCNN)-based solutions for travel demand prediction. As a preliminary step, we have focused on single-mode OD travel demand prediction for shared bikes, with the goal to highlight the limitations of such solutions when it comes to anticipate variations of demand in the presence of disruptive events. Bike-sharing systems serve as connection on the multi-modal transport network [255], allowing for reduced travel times, cost-effectiveness, and a smaller spatial and environmental footprint compared to traditional modes of transportation. While demand forecasting in this domain has historically relied on statistical and machine learning techniques, the advent of deep learning, and, more recently, GNNs, has ushered in a new era of forecasting algorithms, promising improved accuracy. Various models exist for bike-sharing demand prediction, ranging from station level forecast [256], to cluster-based prediction [257], and many other variations deeply detailed in surveys such as [258]. However, only a limited number of studies have focused on OD travel demand prediction [259]. OD prediction poses unique challenges, as the final destinations are usually unknown, and OD matrices generated are sparse [260]. Moreover, forecasting demand in the case of atypical phenomena, such as non-recurrent events, remains a challenge [261]. Bike-sharing demands exhibit temporal and calendar dependencies tied to factors such as weekdays, public holidays, workdays, and school holidays [262]. Furthermore, they are influenced by weather conditions such as wind, humidity, temperature, and, particularly, rainfall [262].

We have chosen to investigate a state-of-the-art prediction algorithm, the Spatio-Temporal Encoder-Decoder Residual Multi-Graph Convolutional network (ST-ED-RMGC), recently proposed by Ke et al. [260] to perform bike-sharing OD demand forecast. ST-ED-RMGC consists in an encoder-decoder based model, where the encoder is composed by a temporal module which takes into account the temporal dependencies of all OD pairs, and a spatial module with several Residual Multi-graph Convolutional networks (RMGCs), which take several adjacency matrix and the graph OD demand as inputs. The RMGC combines a residual module with a multi-graph convolution to capture the spatial correlation between OD pairs. The residual module is introduced to tackle the issue of gradient explosion/gradient vanishing in complex deep networks. The multi-graph convolution applies graph convolution on stacked weighted adjacency matrices.

By focusing on the prediction of bike-sharing trip count between IRIS zones of the city of Lyon, we have explored the possibility of incorporating contextual information such as time, weather conditions, and multi-modal data (e.g., road car flow data) in order to improve the forecasting accuracy of the original ST-ED-RMGC algorithm. The details of the proposed enhanced ST-ED-RMGC architecture can be found in [263, 264].

Both the extended and the original ST-ED-RMGC algorithms have been trained on a dataset including hourly trip counts from January 8th, 2019 to November 9th, 2019. We considered 2500 ODs of the selected territory covering the cities of Lyon and Villeurbanne, France. For performance evaluation, we have considered a test dataset covering the 7 a.m. and 9 p.m. period for all days between December 1<sup>th</sup>, 2019, and March 10<sup>th</sup>, 2020. Alongside this test dataset, named *global scenario* in the following, we also identified specific *weather-related scenarios* to assess the accuracy of bike-sharing demand prediction under degraded weather conditions. All these scenarios were employed to make a comparative evaluation between the original ST-ED-RMGC, used as a baseline, and our enhanced model, which incorporates contextual data.

To construct the weather-related scenarios, we selected subsets from the global scenario
based on two factors computed from the available weather-related historical time series: the
observed daily cumulative rainfall level (dcr) and the hourly amount of rain (hr). Since the
choice of bike sharing as a mode of travel can be greatly influenced by ongoing, expected or
recently concluded rainfall episodes, we considered in our scenarios: peak hours with rain
(denoted as hr > 0), peak hours without rain (denoted as hr = 0), and peak hours with a
significant level of daily cumulative rainfall. The metrics used to evaluate the models were
the Mean Squared Error (MSE) and the MAPE. With no surprise, in scenarios without
rain (hr = 0, dcr = 0) and in the global scenario, the variations of the ST-ED-RMGC model
with weather features (both historical and forecast) exhibit only marginal improvement
compared to the baseline ST-ED-RMGC model. However, in rainfall scenarios (hr > 0), all
models incorporating weather information outperform the baseline, showing up to 20% and
27% improvements in MSE and MAPE, respectively. Notably, rainfall duration emerged
as a relevant predictor of bike-sharing demand during degraded weather conditions. The
initial findings of this study emphasise the significance of the considered contextual data,
particularly weather and calendar-related information, in forecasting bike-sharing demand,
especially during adverse weather conditions.

Despite such initial promising results, further investigation is needed, particularly regard-
ing non-recurrent events beyond weather-related occurrences. First of all, this study serves
as a stepping stone towards extending the proposed approach to encompass multi-modal
travel demand prediction. We are currently working on modifying the original architecture
for incorporating data on public transportation usage, car flow, and speed data. This inte-
gration aims to develop a versatile framework that could capture inter-modal travel demand
relationships jointly with contextual data to perform accurate multi-modal travel forecast-
ing. Although the literature on this subject is relatively limited, the work from Liang et
al. [265] presents an interesting research direction. The authors of this paper combine multi-
relational spatio-temporal GNN with attention-based mechanisms to capture the complex
spatio-temporal patterns across multiple modes and accurately predict multi-modal demand.

Another notable challenge that we aim to address in the short-term is the scarcity of
atypical events in the available datasets. To address this issue, an intriguing research di-
rection involves exploring transfer-learning-based approaches to test the trained model on
diverse traffic-related datasets from various cities and transport modes, thereby expanding
its applicability.

As another research direction, we are attempting to leverage additional data sources,
particularly mobile phone app consumption data, to develop innovative traffic forecasting
solutions capable to better capture dynamic contextual information on both human presence
and human mobility to anticipate traffic variations. Mobile phone app consumption data
offer valuable insights into the activities individuals perform at the specific locations where
they use their mobile phones, such as browsing the web, utilising navigation apps, streaming
audio/video content, or engaging with social media platforms. Sudden fluctuations in the
overall usage of a particular service could serve as an indicator of the initiation of an event,
which, in turn, may result in subsequent traffic state or demand changes. For instance, a
sudden surge in social network usage might be indicative of individuals witnessing an attack
or an extreme event, which could potentially lead to a mass evacuation in the near future.
By incorporating such real-time mobile app consumption data into our analysis, we aim to
enhance our ability to detect and respond to emerging mobility patterns, ultimately enabling
more effective management of transportation systems in response to unforeseen events.

Finally, a more ambitious research direction will focus on integrating traditional traffic
flow models with neural-network-based forecasting algorithms towards increased accuracy,
interpretability and adaptability to unseen traffic scenarios. This approach follows the prin-
ciples of physics-guided machine learning, a methodology recently investigated in traffic fore-
casting with promising results [266]. Traffic flow models typically encapsulate conservation
laws that are similar to those observed in fluid dynamics. The incorporation of these phys-
ical laws into machine learning frameworks serves to enhance not only predictive accuracy
but also the model's level of interpretability. This effectively mitigates a significant limitation often associated with purely data-driven methods. When calibrated with empirically derived travel demand data and substantiated by realistic mode choice models, these physics-informed machine learning models should be expected to faithfully reproduce complex traffic scenarios. Such models can be, therefore, leveraged to swiftly detect discrepancies between real-time data and projected traffic patterns. This capability should enable more informed and timely anomaly detection, contributing to more resilient traffic management strategies. Concurrently, it may be worthwhile to investigate emergent learning frameworks such as Liquid Neural Networks (LNN), initially introduced in [267]. LNN constitute a category of brain-inspired, continuous-time neural-network models that demonstrate capability to capture causal relationships and adaptability under changing conditions, as evidenced in [268], where they were successfully deployed for robust flight navigation using autonomous drones in unpredictable and out-of-distribution scenarios. The integration of biologically inspired priors into LNN results in flexible, robust, and interpretable forecasting solutions that can be adapted to new environments and conditions. These attributes make LNN particularly well-suited for developing versatile learning frameworks capable of performing forecasting across a diverse array of rapidly evolving and unobserved traffic scenarios.

## 6.1.2 Synthetic Population Generation via Multi-Source Data

Agent-based models in transportation traditionally rely on synthetic travel demand, which involves creating a synthetic population of households and individuals with socio-demographic attributes, along with daily activity patterns in both time and space. This approach allows for studying traffic evolution and evaluating the impact of specific transport supply and travel demand variations. However, for the generated simulations to be realistic and meaningful, it is crucial that the synthetic population accurately represents the actual population and its multi-modal travel practices.

Creating a realistic synthetic population of agents is a challenging task, traditionally involving multiple data sources that are processed according to the following steps [269]: (i) ensuring that the marginal distributions of the socio-economic features of the actual population, such as gender, age, and occupation, known from sources like census data, match those of the synthetic population; (ii) associating daily activity patterns to the synthesised population by statistically matching the trip chains derived from a second smaller-scale input population sample, such as regional household travel surveys.

A critical limitation of current synthetic population generation approaches is that they traditionally match activity chains without synthesising realistic location information for activities due to the difficulty in identifying reliable and high-resolution spatial data sources for agent locations. While the home location is generally known via census data, it is rarer that work/study and secondary activity locations are actually inferred from real data. The main methods used to generate location activities include gravity models for commuting patterns [270], discrete destination choice models for workplace assignments [271], and applying space-time constraints for the selection of discretionary activity locations [272]. However, while these approaches can produce reasonable trips in some circumstances, the reliability of the final synthetic travel demand remains limited compared to real-world-data matching possibilities.

We are currently working on extending the state-of-the-art pipeline proposed in [269] with a novel data-driven approach for the location of activities, leveraging anonymised mobile phone data from the telecommunications provider Orange. As complete trajectories can be computationally heavy to process, and feature a significant privacy risk, our idea consists in using the time-dependent OD matrix from mobile-phone NSD reconstructed trajectories [273], obtained with the approaches described in Sec. 3.2 and Sec. 3.4. OD matrices are lighter and more manageable than whole trajectory datasets. This makes them more readily available, easier to process and safer regarding the privacy of transportation users. More-
over, they are built based on high-penetration and high-spatial-resolution data on human mobility 3.2, which makes them an ideal candidate for realistically reconstructing activity locations.

The main originality of the proposed approach lies in interpreting the OD matrix as a transition probability between two locations, which is used to sample the location chain. Additionally, the total number of trips in the OD matrix for each time of the day is utilised as a target to match the population agendas, achieved through hierarchical re-scaling for individuals and trips.

Currently, multiple alternative approaches are being investigated to draw location in-
formation associated with a whole activity chain of an agent from the OD-based transition
probability matrix. While census data are used to determine a pre-defined home location
for each agent, unlike state-of-the-art methods, we do not rely upon sources other than the
OD matrix to pre-define work, study, or secondary activity places. This decision is justified
by the fact that 45% of trips in the activity chains of the population are commute trips,
while the mobile data contain the distribution of the sum of commute and non-commute
trips. Using them only for the latter purpose would be invalid.

Among the most promising solutions to associate location information with each agent’s activity, an approach based on Markov Chain Monte Carlo (MCMC) sampling is being considered, integrated as an additional step to the base pipeline defined in [269]. Our MCMC-based sampling process begins by fixing the agent’s home location from census data and ensuring that all work-related activities of any single agent are drawn only once from the OD-based transition probability matrix, remaining the same for all activity chains containing them, assuming each agent has only one work or study place. The MCMC sampling then proceeds by assigning known values to all pre-defined states (e.g., the agent’s known home location at the beginning and end of their activity chain) and random values to all other states. Subsequently, a new value is successively drawn for each state with dependency on the other states considered known. Through iterative sampling of the unknown locations a large number of times, MCMC converges to the stationary distribution within an acceptable error.

The proposed methods have been evaluated for an area of 25,600 meters wide centred
around Lyon, divided into 515 distinct zones, representing municipalities or spatial units
(IRIS) common to many statistical analyses in France. We used the census performed by
the French statistics institute INSEE to generate the synthetic population of travellers.
Each of the 487,628 rows for our study zone features socioeconomic variables and a non-
integer scaling coefficient for a total of 1,366,072 persons. The activity chains were taken,
without location information, from the household travel survey (named Enquête Ménage
Déplacements (EMD)) performed by the French agency for urban planning [274], which
details the complete agendas of 25,203 persons in the Lyon region. Finally, we used OD-
matrices from TRANSIT spatio-temporally describing the movement of residents of the
region of Lyon. These data cover three months, from 19th March 2019 to 18th June 2019
and contain information about 49,681,883 distinct flows. As the scope of the other data
sources is restricted to working days, we removed weekends and main holidays in France.

The current validation of the generated synthetic population primarily relies on verifying
the agreement between the distributions from the input data (i.e., census data, mobile phone
OD matrix, and household travel surveys) to those derived from the synthetic population.
Naturally, the algorithms aim to match the distributions derived from these sources, which
should ideally correspond closely to those related to the synthetic population. However,
discrepancies in the data sources pose challenges in effectively and jointly matching them.
Among other limitations of the data sources, in particular, the socio-economic population
sample is small and thus not fully representative of the entire population. In Fig. 6.4(a), we
illustrate how the socioeconomic composition of our population matches the official census
better than the EMD survey. Each red point represents the share of the population with
a given socio-economic feature value, in the case of the actual one (on the x-axis) versus
<figure><img src="image_82.png" /><figcaption>**Figure 6.4:** Average deviance of each mode per cluster.</figcaption></figure>

the synthetic one (on the y-axis). The same comparison is reported with green points for the EMD population. As the EMD represents only a small number of people, it is normal that the totals have a higher variance. Our synthetic population can then be seen as a version of the EMD that fully agrees with the census on the socioeconomic distribution of the population. Fig. 6.4(c) graphically presents the distribution of trips during a typical day as measured from mobile data (in orange), our synthetic data (in red), and the survey EMD (in green). In this case, we considered the observed volumes from the mobile data to be the ground truth. With no surprise, it can be noted that our population perfectly fits the ground truth<fnref n="2" />. On the contrary, the EMD seems to overestimate the volumes of the morning, noon, and evening peaks while underestimating the volumes during the rest of the day.

Finally, we compared the flows obtained with our approaches with the ones observed from the mobile data. A first natural way of comparing the flows is to compare the distributions of destinations $P(d|o, t)$. In Fig. 6.4(b), we illustrate how the probability of destinations given origin and hour fits the probability table derived from the OD matrix. Each point corresponds to a combination of destination, origin, and time of the day. Its x-coordinate is the probability $P(d|o, t)$ observed in the OD matrix, while its y-coordinate is the same probability as observed in the trips of our synthetic population. We observe a bi-modality in the ground truth, and a generally good fit of the synthetic population flows with the mobile phone OD matrix.

The reported results represent a preliminary step in developing the complete methodology, but further improvements and assessments must be considered. While it can be

<footnote n="2">²Here, we consider the mobile phone data trip volumes as ground truth.</footnote>
expected that our synthetic demand data matches well the distribution $P(D|O, T)$, we have
observed that the distribution of trips $P(O, D|T)$ is not the same as in the input OD ma-
trix. This is because although we sample the destinations using $P(D|O, T)$ as a transition
probability, we have no mechanism to ensure that the correct number of agents leave each
origin $O$ at each time step $T$. This problem can be addressed by decomposing each time
slice of the transition matrix into a sum of transition matrices depending on $O, T$, and on
the time step of the next trip of the agent. This new decomposition of the OD matrix in
our future work amounts to adding explanatory variables to the mobile data. By carefully
choosing such a decomposition, we can make sure the agents taking a trip before time step
$T$ are assigned to destinations such that the map of agents leaving for a trip at time step
$T$ corresponds to the map of origins of trips in the OD matrix. Additional future research
directions include refining the structure of the OD transition matrix by distinguishing recur-
ring commuting trips from non-recurring ones, as well as introducing additional constraints
on trip distance and transport modes for more accurately matching location activities, by
possibly including additional data sources in the pipeline. We also aim to perform validation
using, as ground truth, external data sources that were not considered within the population
generation pipeline. For instance, we plan to incorporate aggregate statistics on mobility in
the considered region from different data providers to validate, e.g., the distributions of trip
distances, transport mode shares, and the joint distribution of activity locations based on
socio-economic features. This validation will help us gain further confidence in the reliability
of our methodology and outcomes.

Additionally, we aim to conduct multi-modal traffic simulations with in-house macro-
scopic simulators, properly calibrated with external data sources, including multi-modal
travel demand data, loop detector flow data, FCD speed data, and ticket validation data.
The resulting simulated traffic states will be compared to those obtained using a multi-
agent multi-modal traffic simulator, such as MATSim [275], using our generated synthetic
population as the main input to run simulations.

## 6.2 Longer-term Perspectives

In the long term, my primary research objective is to investigate and enhance the resilience
characteristics of urban networks within a *multi-domain context* and according to a *holis-
tic* vision. The ultimate goal is to contribute to developing an integrated, data-driven
digital framework enabling *by-design* to assess and promote the joint resilience of diverse
inter-related urban technical networks, such as those related to mobility, energy, telecom-
munications, and more.

These long-term perspectives find their reason in the urgency of addressing climate
change, whose effects are especially critical in large metropolitan areas due to the pres-
ence of heterogeneous yet highly interdependent infrastructures [276, 277]. For example,
electricity pylons collapsing from the consequences of extreme weather-related events can
disrupt traffic by impeding vehicle evacuations, power outages can compromise mobility by
affecting traffic lights and electrified transport means, and a lack of communication networks
may hinder future connected and autonomous vehicle mobility.

Considering the fragility of these interconnected infrastructures, their susceptibility to
climate change and the cascading effects of weather-related disruptions, it appears unavoid-
able to rethink the traditional *in-silo* approach to resilience modelling according to a new
holistic approach aimed to investigate, model, monitor, and optimise as a *whole* the various
networks constituting the urban infrastructure ecosystem, while prioritising energy efficiency,
sustainability and user-privacy requirements.

The vision of resilience as a crucial urban commodity – involving the synergistic in-
teraction of multiple stakeholders, resources and infrastructures – is becoming increasingly
significant in each of the mentioned fields [278, 279, 280] and ties into the notion of Resilience-
as-a-Service (RaaS). The term was initially introduced in the context of high-performance computing [281] to provide cloud solutions with a highly available, distributed, and scalable fault-tolerant service. Very recently, the term has been introduced in the field of transportation [282], defining a conceptual and operational framework where the available resources of various service providers could be integrated to manage and improve over time the resilience of a system faced by hard-to-predict disruptions and extreme events.

Despite some isolated attempts to apply this concept in a multi-domain urban context [283, 284], significant research is still needed to implement and generalise this vision on a large scale, encompassing diverse infrastructural networks within a city. In particular, by considering again the parallel between the urban critical infrastructures and the cloud computing paradigm, I believe this vision should further emphasise the need for *targeted and energy-efficient data collection and processing, network-level sensing, continuous and explainable learning*, as well as *real-time network monitoring and control*. These functionalities should all be conceived as pivotal elements of an integrated urban network digital architecture for the realistic and large-scale implementation of RaaS.

To achieve this, three key macro-objectives are described in the following as the pillars of my future research.

## The City Digital Twin

A key research objective involves the joint modelling of urban, multi-domain, and multi-modal networks as an integrated and autonomous digital system that mirrors physical reality, and which could be named the *city digital twin*. This objective aims to leverage real-time, multi-source data – spanning malfunctions, service usages, and contextual information from various service networks – as primary resources. These data would enable automated solutions to identify standard operational patterns, detect anomalies, and autonomously implement mitigation, adaptation, and emergency-response strategies. This vision is gaining more concrete ground for real-world implementation from emerging network infrastructures, including 6G telecommunication networks [285], IoT-enabled smart grids and micro-grids [286], fleets of autonomous connected vehicles for public transport [287], which are expected to integrate sensing and actuating mechanisms as standard features. The objective aligns closely with the expanding scope of *urban digital twins* [38, 288], a concept referring to the widespread adoption of digital twins technology across the urban ecosystem. This includes buildings, energy systems, vehicles, and telecommunication facilities and devices. A *city digital twin* shall push this vision further, thus functioning as an accurate digital replica of the physical urban environment as a whole. The city digital twin shall be endowed with advanced sensing capabilities, and autonomous interactions with the physical urban infrastructures, and provide self-diagnostic and self-repair functions that are aware of the inter-connections and inter-dependencies among the different technical components of the urban environment.

## Coupled Vulnerability Assessment

A pivotal focus will be the identification and systematic analysis of critical interdependencies across a range of modelled networks. This effort aims to cultivate a modular and extensible platform specifically engineered for conducting coupled vulnerability assessments. The platform shall be capable of modelling, simulating, and quantifying the impacts of a diverse set of disruptive scenarios on various urban networks. A particular emphasis will be placed on understanding cascading effects, which result from failures in one network influencing others [289]. Central to this aim will be the design of the platform's core as middleware, crafted to facilitate seamless integration among a variety of simulation and modelling tools. This core will serve as the backbone for a unified framework wherein dynamic models for, e.g., flood, communication, electrical, and traffic simulation, can coexist and dynamically
interact. Such interaction will enable these models to respond collaboratively to simulated
stressors like heavy rainfall, heatwaves and other kinds of shocks. This integrated approach
will provide holistic insights into the vulnerabilities inherent in the urban environment and
allow for the assessment of the effectiveness of various resilience strategies.

Multi-domain Resilience Strategies

The last fundamental objective I would like to develop concerns the exploration of decision-
making tools and strategies that promote the synergistic integration of various urban net-
works and services. The aim is to enhance the sustainability, robustness, adaptability, and
mitigation capabilities of the urban system in the face of extreme events and network failures.
For instance, optimal scheduling of mitigation and repair actions – such as reinforcing dis-
tribution poles or towers, pre-installing battery backup systems at signalised intersections,
deploying emergency vehicles, and rerouting traffic to facilitate access to disrupted facilities
– serves as a straightforward example aimed at bolstering the resilience of the integrated
transport, electricity, and telecommunications systems. More ambitious initiatives could
focus on the joint optimisation and emergency operation of the integrated urban networks.
These can be enabled by the monitoring capabilities of the city digital replica and by the in-
sights derived from the joint vulnerability assessments, in line with approaches very recently
proposed in the literature [289]. For example, Electric and electrified Vehicles (EVs) can
serve as emergency mobile resilience assets during prolonged power outages, supplying power
to critical buildings like hospitals. Similarly, a fleet of EVs could power nearby facilities or
even assist in restarting substation transformers following a natural disaster [290].

Energy efficiency, reduced carbon footprint and cybernetic resilience shall also be con-
sidered central requirements in this vision for a resilient and sustainable urban landscape.
For instance, I plan to explore data-frugal processing frameworks, such as continual learn-
ing [291], which aim to optimise the capabilities of machine learning models. The objective
here is to store key information from previous data observations, minimise the impact of
catastrophic forgetting [292], and adapt swiftly to changing data distributions. Achieving
this allows for a reduction in data footprint. Specifically, it reduces the computational and
storage resources needed for data processing and analysis.

These perspectives and objectives will serve as the backbone for a European Research Coun-
cil (ERC) consolidator project that I plan to submit in the coming years. Additionally, I
am actively contributing to the development of a large-scale project on this topic with the
Gustave Eiffel Foundation. My involvement includes participating in a dedicated working
group on urban resilience that focuses on the study of the vulnerabilities of city technical
networks, and the definition of proper adaptation strategies.
Bibliography

[1] L. Bevilacqua, A. Furno, V. S. Di Carlo, and E. Zimeo, "A tool for automatic generation of ws-bpel compositions from owl-s described services," in *2011 5th International Conference on Software, Knowledge Information, Industrial Management and Applications (SKIMA) Proceedings*. IEEE, 2011, pp. 1-8. 2

[2] A. Furno and E. Zimeo, "Context-aware composition of semantic web services," *Mobile Networks and Applications*, vol. 19, no. 2, pp. 235-248, 2014. 2

[3] —, "Self-scaling cooperative discovery of service compositions in unstructured p2p networks," *Journal of Parallel and Distributed Computing*, vol. 74, no. 10, pp. 2994-3025, 2014. 2

[4] A. Furno, "Scalable service composition in autonomic computing," PhD dissertation, University of Sannio, Department of Engineering, 2014. [Online]. Available: https://people.licit-lyon.eu/furno/documents/phd_thesis_angelo_furno.pdf 2

[5] P. Suppa and E. Zimeo, "A context-aware mashup recommender based on social networks data mining and user activities," in *2016 IEEE International Conference on Smart Computing (SMARTCOMP)*. IEEE, 2016, pp. 1-6. 2

[6] A. Furno, D. Naboulsi, R. Stanica, and M. Fiore, "Mobile demand profiling for cellular cognitive networking," *IEEE Transactions on Mobile Computing*, vol. 16, no. 3, pp. 772-786, 2016. 3, 26, 27

[7] D. Naboulsi, R. Stanica, and M. Fiore, "Classifying call profiles in large-scale mobile traffic datasets," in *IEEE INFOCOM 2014-IEEE conference on computer communications*. IEEE, 2014, pp. 1806-1814. 3, 26, 29

[8] D. Bega, M. Gramaglia, M. Fiore, A. Banchs, and X. Costa-Perez, "Deepcog: Cognitive network management in sliced 5g networks with deep learning," in *IEEE INFOCOM 2019-IEEE conference on computer communications*. IEEE, 2019, pp. 280-288. 3

[9] A. Furno, M. Fiore, R. Stanica, C. Ziemlicki, and Z. Smoreda, "A tale of ten cities: Characterizing signatures of mobile traffic in urban areas," *IEEE Transactions on Mobile Computing*, vol. 16, no. 10, pp. 2682-2696, 2016. 3, 10, 29, 31

[10] A. Furno, M. Fiore, and R. Stanica, "Joint spatial and temporal classification of mobile traffic demands," in *IEEE INFOCOM 2017-IEEE Conference on Computer Communications*. IEEE, 2017, pp. 1-9. 3, 10, 34

[11] C. Kennedy, S. Pincetl, and P. Bunje, "The study of urban metabolism and its applications to urban planning and design," *Environmental pollution*, vol. 159, no. 8-9, pp. 1965-1973, 2011. 3

[12] D. N. Bristow and C. A. Kennedy, "Urban metabolism and the energy stored in cities: Implications for resilience," *Journal of Industrial Ecology*, vol. 17, no. 5, pp. 656-667, 2013. 3
[13] E. Hollnagel, D. D. Woods, and N. Leveson, *Resilience engineering: Concepts and precepts*. Ashgate Publishing, Ltd., 2006. 5, 14, 16, 19

[14] E. Bellini and P. Nesi, “Exploiting smart technologies to build smart resilient cities,” in *Routledge handbook of sustainable and resilient infrastructure*. Routledge, 2018, pp. 685–705. 5, 17

[15] V. Mahajan, N. Kuehnel, A. Intzevidou, G. Cantelmo, R. Moeckel, and C. Antoniou, “Data to the people: a review of public and proprietary data for transport models,” *Transport reviews*, vol. 42, no. 4, pp. 415–440, 2022. 7

[16] C. Wang and D. B. Hess, “Role of urban big data in travel behavior research,” *Transportation research record*, vol. 2675, no. 4, pp. 222–233, 2021. 7

[17] X. Kong, M. Li, K. Ma, K. Tian, M. Wang, Z. Ning, and F. Xia, “Big trajectory data: A survey of applications and services,” *IEEE Access*, vol. 6, pp. 58295–58306, 2018. 7

[18] P. R. Stopher and S. P. Greaves, “Household travel surveys: Where are we going?” *Transportation Research Part A: Policy and Practice*, vol. 41, no. 5, pp. 367–381, 2007. 8, 35

[19] M. Lin and W.-J. Hsu, “Mining gps data for mobility patterns: A survey,” *Pervasive and mobile computing*, vol. 12, pp. 1–16, 2014. 9

[20] O. Cats, “Identifying human mobility patterns using smart card data,” *arXiv preprint arXiv:2208.05352*, 2022. 9

[21] S. Bührmann, “New seamless mobility services: Public bicycles,” *NICHES Policy Notes*, 2007. 9

[22] S. A. Shaheen, S. Guzman, and H. Zhang, “Bikesharing in europe, the americas, and asia: past, present, and future,” *Transportation research record*, vol. 2143, no. 1, pp. 159–167, 2010. 9

[23] J. Song, L. Zhang, Z. Qin, and M. A. Ramli, “Spatiotemporal evolving patterns of bike-share mobility networks and their associations with land-use conditions before and after the covid-19 outbreak,” *Physica A: Statistical Mechanics and its Applications*, vol. 592, p. 126819, 2022. 9

[24] J. Zhou, Y. Guo, J. Sun, E. Yu, and R. Wang, “Review of bike-sharing system studies using bibliometrics method,” *Journal of traffic and transportation engineering (English edition)*, vol. 9, no. 4, pp. 608–630, 2022. 9

[25] T. H. Silva, A. C. Viana, F. Benevenuto, L. Villas, J. Salles, A. Loureiro, and D. Quer CIA, “Urban computing leveraging location-based social network data: a survey,” *ACM Computing Surveys (CSUR)*, vol. 52, no. 1, pp. 1–39, 2019. 10

[26] V. D. Blondel, A. Decuyper, and G. Krings, “A survey of results on mobile phone datasets analysis,” *EPJ data science*, vol. 4, pp. 1–55, 2015. 10

[27] M. C. González, C. A. Hidalgo, and A.-L. Barabási, “Understanding individual human mobility patterns,” *Nature*, vol. 453, no. 7196, pp. 779–782, 6 2008. [Online]. Available: http://www.nature.com/doifinder/10.1038/nature06958 10

[28] M. S. Iqbal, C. F. Choudhury, P. Wang, and M. C. González, “Development of origin–destination matrices using mobile phone call data,” *Transportation Research Part C: Emerging Technologies*, vol. 40, pp. 63–74, 3 2014. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S0968090X14000059/10
[29] R. W. Douglass, D. A. Meyer, M. Ram, D. Rideout, and D. Song, "High resolution population estimates from telecommunications data," *EPJ Data Science*, vol. 4, no. 1, p. 4, 12 2015. [Online]. Available: http://www.epjdatascience.com/content/4/1/4

[30] A. Wesolowski, C. O. Buckee, K. Engø-Monsen, and C. J. E. Metcalf, "Connecting mobility to infectious diseases: the promise and limits of mobile phone data," *The Journal of infectious diseases*, vol. 214, no. suppl_4, pp. S414-S420, 2016. 10

[31] D. Naboulsi, M. Fiore, S. Ribot, and R. Stanica, "Large-scale mobile traffic analysis: a survey," *IEEE Communications Surveys & Tutorials*, vol. 18, no. 1, pp. 124-161, 2015. 10

[32] Q. Xu, A. Gerber, Z. M. Mao, and J. Pang, "AccuLoc: Practical localization of performance measurements in 3G networks," in *MobiSys'11 - Compilation Proceedings of the 9th International Conference on Mobile Systems, Applications and Services and Co-located Workshops*. New York, New York, USA: ACM Press, 2011, pp. 183-195. [Online]. Available: http://portal.acm.org/citation.cfm?doid=1999995.2000013

[33] G. Chen, A. C. Viana, M. Fiore, and C. Sarraute, "Complete trajectory reconstruction from sparse mobile phone data," *EPJ Data Science*, vol. 8, no. 1, pp. 1-24, 2019. 10

[34] G. Ranjan, H. Zang, Z.-L. Zhang, and J. Bolot, "Are call detail records biased for sampling human mobility?" *ACM SIGMOBILE Mobile Computing and Communications Review*, vol. 16, no. 3, pp. 33-44, 12 2012. [Online]. Available: https://dl.acm.org/doi/10.1145/2412096.2412101

[35] D. Bachir, V. Gauthier, M. E. Yacoubi, and G. Khodabandelou, "Using mobile phone data analysis for the estimation of daily urban dynamics," in *2017 IEEE 20th International Conference on Intelligent Transportation Systems (ITSC)*. IEEE, 10 2017, pp. 626-632. [Online]. Available: http://ieeexplore.ieee.org/document/8317956/ 13

[36] H. R. Meredith, V. Andreani, H. R. Ma, A. J. Lopatkin, A. J. Lee, D. J. Anderson, G. Batt, and L. You, "Applying ecological resistance and resilience to dissect bacterial antibiotic responses," *Science Advances*, vol. 4, no. 12, p. eaau1873, 2018. 14

[37] T. Parsons, E. A. Shils, and N. J. Smelser, "The social system," in *Toward a General Theory of Action*. Routledge, 2017, pp. 190-233. 14

[38] X. Liu, D. Li, M. Ma, B. K. Szymanski, H. E. Stanley, and J. Gao, "Network resilience," *Physics Reports*, vol. 971, pp. 1-108, 2022. 14, 15, 19, 21, 22, 81, 157

[39] R. Billinton and R. N. Allan, *Reliability evaluation of engineering systems*. Springer, 1992, vol. 792. 14

[40] D. D. Woods, "Four concepts for resilience and the implications for the future of resilience engineering," *Reliability Engineering & System Safety*, vol. 141, pp. 5-9, 2015. 14, 15

[41] L. Fisher, "More than 70 ways to show resilience," *Nature*, vol. 518, no. 7537, pp. 35-35, 2015. 15

[42] M. Bruneau, S. E. Chang, R. T. Eguchi, G. C. Lee, T. D. O'Rourke, A. M. Reinhorn, M. Shinozuka, K. Tierney, W. A. Wallace, and D. Von Winterfeldt, "A framework to quantitatively assess and enhance the seismic resilience of communities," *Earthquake spectra*, vol. 19, no. 4, pp. 733-752, 2003. 15, 20

[43] D. Henry and J. E. Ramirez-Marquez, "Generic metrics and quantitative approaches for system resilience as a function of time," *Reliability Engineering & System Safety*, vol. 99, pp. 114-122, 2012. 15, 19, 20
[44] R. Steen and T. Aven, "A risk perspective suitable for resilience engineering," *Safety science*, vol. 49, no. 2, pp. 292-297, 2011. 15

[45] T. Aven, "Risk assessment and risk management: Review of recent advances on their foundation," *European Journal of Operational Research*, vol. 253, no. 1, pp. 1-13, 2016. 15

[46] L.-G. Mattsson and E. Jenelius, "Vulnerability and resilience of transport systems- a discussion of recent research," *Transportation research part A: policy and practice*, vol. 81, pp. 16-34, 2015. 15, 19, 20, 83

[47] J. S. Cañavera-Herrera, J. Tang, T. Nochta, and J. M. Schooling, "On the relation between 'resilience' and 'smartness': A critical review," *International Journal of Disaster Risk Reduction*, p. 102970, 2022. 17, 19

[48] K. Barker, J. H. Lambert, C. W. Zobel, A. H. Tapia, J. E. Ramirez-Marquez, L. Albert, C. D. Nicholson, and C. Caragea, "Defining resilience analytics for interdependent cyber-physical-social networks," *Sustainable and Resilient Infrastructure*, vol. 2, no. 2, pp. 59-67, 2017. 17

[49] M. Z. Serdar, M. Koç, and S. G. Al-Ghamdi, "Urban transportation networks resilience: indicators, disturbances, and assessment methods," *Sustainable Cities and Society*, vol. 76, p. 103452, 2022. 17, 19

[50] E. Henry, "Resilience modeling of urban multimodal transport networks," Ph.D. dissertation, Université de Lyon, 2021. 19

[51] Y. Cheng, E. A. Elsayed, and Z. Huang, "Systems resilience assessments: a review, framework and metrics," *International Journal of Production Research*, vol. 60, no. 2, pp. 595-622, 2022. 19

[52] B. Rathnayaka, C. Siriwardana, D. Robert, D. Amaratunga, and S. Setunge, "Improving the resilience of critical infrastructure: Evidence-based insights from a systematic literature review," *International Journal of Disaster Risk Reduction*, p. 103123, 2022. 19

[53] J. Douglas, "Physical vulnerability modelling in natural hazard risk assessment," *Natural Hazards and Earth System Sciences*, vol. 7, no. 2, pp. 283-288, 2007. 19

[54] V. Latora and M. Marchiori, "Efficient behavior of small-world networks," *Physical review letters*, vol. 87, no. 19, p. 198701, 2001. 20

[55] P. Gauthier, A. Furno, and N.-E. El Faouzi, "Road network resilience: how to identify critical links subject to day-to-day disruptions," *Transportation research record*, vol. 2672, no. 1, pp. 54-65, 2018. 20, 83

[56] M. Kivelä, A. Arenas, M. Barthelemy, J. P. Gleeson, Y. Moreno, and M. A. Porter, "Multilayer networks," *Journal of complex networks*, vol. 2, no. 3, pp. 203-271, 2014. 21, 69

[57] A. Aleta, S. Meloni, and Y. Moreno, "A multilayer perspective for the analysis of urban transportation systems," *Scientific reports*, vol. 7, no. 1, p. 44359, 2017. 21, 23

[58] J. Lin and Y. Ban, "Complex network topology of transportation systems," *Transport reviews*, vol. 33, no. 6, pp. 658-685, 2013. 21

[59] B. Jiang and C. Claramunt, "Topological analysis of urban street networks," *Environment and Planning B: Planning and design*, vol. 31, no. 1, pp. 151-162, 2004. 21
[60] S. V. Buldyrev, R. Parshani, G. Paul, H. E. Stanley, and S. Havlin, "Catastrophic cascade of failures in interdependent networks," *Nature*, vol. 464, no. 7291, pp. 1025-1028, 2010. 21

[61] Y. Duan and F. Lu, "Robustness of city road networks at different granularities," *Physica A: Statistical Mechanics and its Applications*, vol. 411, pp. 21–34, 2014. 21, 83

[62] E. Estrada, *The structure of complex networks: theory and applications*. Oxford University Press, 2012. 21

[63] J. Wang, "Resilience of self-organised and top-down planned cities—a case study on London and Beijing street networks," *PloS one*, vol. 10, no. 12, p. e0141736, 2015. 21

[64] M. Akbarzadeh, S. Memarmontazerin, S. Derrible, and S. F. Salehi Reihani, "The role of travel demand and network centrality on the connectivity and resilience of an urban street system," *Transportation*, vol. 46, pp. 1127–1141, 2019. 21

[65] A. Nair and J. M. Vidal, "Supply network topology and robustness against disruptions—an investigation using multi-agent model," *International Journal of Production Research*, vol. 49, no. 5, pp. 1391–1404, 2011. 21

[66] E. Jenelius and L.-G. Mattsson, "Resilience of transport systems," in *Encyclopedia of Transportation*. Elsevier Amsterdam, Netherlands, 2020. 21, 81

[67] S. Dong, A. Mostafizi, H. Wang, J. Gao, and X. Li, "Measuring the topological robustness of transportation networks to disaster-induced failures: A percolation approach," *Journal of Infrastructure Systems*, vol. 26, no. 2, p. 04020009, 2020. 21

[68] J. Wu, Z. Gao, H. Sun, and H. Huang, "Urban transit system as a scale-free network," *Modern Physics Letters B*, vol. 18, no. 19n20, pp. 1043–1049, 2004. 21

[69] C. von Ferber, B. Berche, T. Holovatch, and Y. Holovatch, "A tale of two cities: Vulnerabilities of the London and Paris transit networks," *Journal of Transportation Security*, vol. 5, pp. 199–216, 2012. 21

[70] Y. Meng, X. Tian, Z. Li, W. Zhou, Z. Zhou, and M. Zhong, "Comparison analysis on complex topological network models of urban rail transit: A case study of shenzhen metro in China," *Physica A: Statistical Mechanics and Its Applications*, vol. 559, p. 125031, 2020. 21

[71] C. Von Ferber, T. Holovatch, Y. Holovatch, and V. Palchykov, "Public transport networks: empirical analysis and modeling," *The European Physical Journal B*, vol. 68, pp. 261–275, 2009. 21

[72] C. Han and L. Liu, "Topological vulnerability of subway networks in China," in *2009 International Conference on Management and Service Science*. IEEE, 2009, pp. 1–4. 21

[73] M. E. O’Kelly, "Network hub structure and resilience," *Networks and Spatial Economics*, vol. 15, pp. 235–251, 2015. 22

[74] K. Adjetey-Bahun, J.-L. Planchet, B. Birregah, and E. Châtelet, "Railway transportation system's resilience: Integration of operating conditions into topological indicators," in *NOMS 2016-2016 IEEE/IFIP Network Operations and Management Symposium*. IEEE, 2016, pp. 1163–1168. 22

[75] Q.-C. Lu, "Modeling network resilience of rail transit under operational incidents," *Transportation Research Part A: Policy and Practice*, vol. 117, pp. 227–237, 2018. 22
[76] M. A. Taylor, “Travel through time: the story of research on travel time reliability,” *Transportmetrica B: transport dynamics*, vol. 1, no. 3, pp. 174–194, 2013. 22

[77] R. D. Connors and D. P. Watling, “Assessing the demand vulnerability of equilibrium traffic networks via network aggregation,” *Networks and Spatial Economics*, vol. 15, pp. 367–395, 2015. 22

[78] E. Jenelius, T. Petersen, and L.-G. Mattsson, “Importance and exposure in road network vulnerability analysis,” *Transportation Research Part A: Policy and Practice*, vol. 40, no. 7, pp. 537–560, 2006. 22, 23, 84

[79] D. M. Scott, D. C. Novak, L. Aultman-Hall, and F. Guo, “Network robustness index: A new method for identifying critical links and evaluating the performance of transportation networks,” *Journal of Transport Geography*, vol. 14, no. 3, pp. 215–227, 2006. 22, 23

[80] A. Nagurney and Q. Qiang, “Fragile networks: identifying vulnerabilities and synergies in an uncertain age,” *International Transactions in Operational Research*, vol. 19, no. 1-2, pp. 123–160, 2012. 22

[81] E. Jenelius, “User inequity implications of road network vulnerability,” *Journal of Transport and Land Use*, vol. 2, no. 3/4, pp. 57–73, 2010. 22

[82] J. L. Sullivan, D. C. Novak, L. Aultman-Hall, and D. M. Scott, “Identifying critical road segments and measuring system-wide robustness in transportation networks with isolating links: A link-based capacity-reduction approach,” *Transportation Research Part A: Policy and Practice*, vol. 44, no. 5, pp. 323–336, 2010. 22, 83, 85

[83] E. Jenelius and L.-G. Mattsson, “Road network vulnerability analysis of area-covering disruptions: A grid-based approach with case study,” *Transportation research part A: policy and practice*, vol. 46, no. 5, pp. 746–760, 2012. 22, 91

[84] E. Mitsakis, A. Papanikolaou, G. Ayfadopoulou, J. Salanova, C. Doll, G. Giannopoulos, and C. Zerefos, “An integrated framework for linking climate change impacts to emergency adaptation strategies for transport networks,” *European Transport Research Review*, vol. 6, no. 2, pp. 103–111, 2014. 22

[85] E. Jenelius, “Network structure and travel patterns: explaining the geographical disparities of road network vulnerability,” *Journal of Transport Geography*, vol. 17, no. 3, pp. 234–244, 2009. 22, 23

[86] M. A. Taylor, S. V. Sekhar, and G. M. D'Este, “Application of accessibility based methods for vulnerability analysis of strategic road networks,” *Networks and Spatial Economics*, vol. 6, pp. 267–291, 2006. 23

[87] M. A. Taylor et al., “Remoteness and accessibility in the vulnerability analysis of regional road networks,” *Transportation research part A: policy and practice*, vol. 46, no. 5, pp. 761–771, 2012. 23

[88] R. Nyberg and M. Johansson, “Indicators of road network vulnerability to storm-felled trees,” *Natural hazards*, vol. 69, pp. 185–199, 2013. 23

[89] V. L. Knoop, M. Snelder, H. J. van Zuylen, and S. P. Hoogendoorn, “Link-level vulnerability indicators for real-world networks,” *Transportation Research Part A: Policy and Practice*, vol. 46, no. 5, pp. 843–854, 2012. 23

[90] X. He and H. X. Liu, “Modeling the day-to-day traffic evolution process after an unexpected network disruption,” *Transportation Research Part B: Methodological*, vol. 46, no. 1, pp. 50–71, 2012. 23
[91] J. Yates and S. Sanjeevi, "A length-based, multiple-resource formulation for shortest path network interdiction problems in the transportation sector," *International Journal of Critical Infrastructure Protection*, vol. 6, no. 2, pp. 107-119, 2013. 23

[92] H. Ho, A. Sumalee, W. H. Lam, and W. Szeto, "A continuum modeling approach for network vulnerability analysis at regional scale," *Procedia-Social and Behavioral Sciences*, vol. 80, pp. 846-859, 2013. 23

[93] E. Rodríguez-Núñez and J. C. García-Palomares, "Measuring the vulnerability of public transport networks," *Journal of transport geography*, vol. 35, pp. 50-63, 2014. 23

[94] O. Cats and E. Jenelius, "Dynamic vulnerability analysis of public transport networks: mitigation effects of real-time information," *Networks and Spatial Economics*, vol. 14, pp. 435-463, 2014. 23

[95] ——, "Planning for the unexpected: The value of reserve capacity for public transport network robustness," *Transportation Research Part A: Policy and Practice*, vol. 81, pp. 47-61, 2015. 23

[96] S. K. Peterson and R. L. Church, "A framework for modeling rail transport vulnerability," *Growth and Change*, vol. 39, no. 4, pp. 617-641, 2008. 23

[97] L. Hong, M. Ouyang, S. Peeta, X. He, and Y. Yan, "Vulnerability assessment and mitigation for the chinese railway system under floods," *Reliability Engineering & System Safety*, vol. 137, pp. 58-68, 2015. 23

[98] M. Janić, "Modelling the resilience, friability and costs of an air transport network affected by a large-scale disruptive event," *Transportation Research Part A: Policy and Practice*, vol. 81, pp. 77-92, 2015. 23

[99] D. Yin, W. Huang, B. Shuai, H. Liu, and Y. Zhang, "Structural characteristics analysis and cascading failure impact analysis of urban rail transit network: From the perspective of multi-layer network," *Reliability Engineering & System Safety*, vol. 218, p. 108161, 2022. 23

[100] L. Bellocchi, V. Latora, and N. Geroliminis, "Dynamical efficiency for multimodal time-varying transportation networks," *Scientific reports*, vol. 11, no. 1, pp. 1-14, 2021. 24

[101] "Ericsson mobility report," Nov 2022. [Online]. Available: https://www.ericsson.com/4ae28d/assets/local/reports-papers/mobility-report/documents/2022/ericsson-mobility-report-november-2022.pdf 25

[102] Key technological challenges of the EC H2020 5G Infrastructure PPP. [Online]. Available: http://5g-ppp.eu/ 25

[103] V. D. Blondel, M. Esch, C. Chan, F. Clerot, P. Deville, E. Huens, F. Morlot, Z. Smoreda, and C. Ziemlicki, "Data for development: the d4d challenge on mobile phone data," 2012. [Online]. Available: https://arxiv.org/abs/1210.0137 26, 77

[104] Y.-A. de Montjoye, Z. Smoreda, R. Trinquart, C. Ziemlicki, and V. D. Blondel, "D4d-senegal: The second mobile phone data for development challenge," *ArXiv*, vol. abs/1407.4885, 2014. 26, 77

[105] T. Italia, "Telecom Italia big data challenge," URL https://dandelion.eu/datamine/open-big-data, 2015. 26

[106] A. Ceselli, M. Fiore, A. Furno, M. Premoli, S. Secci, and R. Stanica, "Prescriptive analytics for mec orchestration," in *2018 IFIP Networking Conference (IFIP Networking) and Workshops*. IEEE, 2018, pp. 1-9. 28
[107] U. Paul, A. P. Subramanian, M. M. Buddhikot, and S. R. Das, "Understanding traffic dynamics in cellular data networks," in *2011 Proceedings IEEE INFOCOM*. IEEE, 2011, pp. 882-890. 29

[108] M. Z. Shafiq, L. Ji, A. X. Liu, J. Pang, S. Venkataraman, and J. Wang, "A first look at cellular network performance during crowded events," *ACM SIGMETRICS performance evaluation review*, vol. 41, no. 1, pp. 17-28, 2013. 29

[109] R. Ahas, A. Aasa, Y. Yuan, M. Raubal, Z. Smoreda, Y. Liu, C. Ziemlicki, M. Tiru, and M. Zook, "Everyday space-time geographies: using mobile phone-based sensor data to monitor urban activity in harbin, paris, and tallinn," *International Journal of Geographical Information Science*, vol. 29, no. 11, pp. 2017-2039, 2015. 29

[110] M. Z. Shafiq, L. Ji, A. X. Liu, J. Pang, and J. Wang, "Large-scale measurement and characterization of cellular machine-to-machine traffic," *IEEE/ACM transactions on Networking*, vol. 21, no. 6, pp. 1960-1973, 2013. 29

[111] S. Almeida, J. Queijo, and L. M. Correia, "Spatial and temporal traffic distribution models for gsm," in *Gateway to 21st Century Communications Village. VTC 1999-Fall. IEEE VTS 50th Vehicular Technology Conference (Cat. No. 99CH36324)*, vol. 1. IEEE, 1999, pp. 131-135. 29

[112] M. R. Vieira, V. Frias-Martinez, N. Oliver, and E. Frias-Martinez, "Characterizing dense urban areas from mobile phone-call data: Discovery and social dynamics," in *2010 IEEE Second International Conference on Social Computing*. IEEE, 2010, pp. 241-248. 29

[113] M. De Nadai, J. Staiano, R. Larcher, N. Sebe, D. Quercia, and B. Lepri, "The death and life of great italian cities: a mobile phone data perspective," in *Proceedings of the 25th international conference on world wide web*, 2016, pp. 413-423. 29

[114] A. Furno, R. Stanica, and M. Fiore, "A comparative evaluation of urban fabric detection techniques based on mobile traffic data," in *Proceedings of the 2015 IEEE/ACM international conference on advances in social networks analysis and mining 2015*, 2015, pp. 689-696. 30

[115] V. Soto and E. Frías-Martínez, "Automated land use identification using cell-phone records," in *Proceedings of the 3rd ACM international workshop on MobiArch*, 2011, pp. 17-22. 30

[116] S. Grauwin, S. Sobolevsky, S. Moritz, I. Gódor, and C. Ratti, "Towards a comparative science of cities: Using mobile traffic records in new york, london, and hong kong," *Computational approaches for urban environments*, pp. 363-387, 2015. 30

[117] B. Cici, M. Gjoka, A. Markopoulou, and C. T. Butts, "On the decomposition of cell phone activity patterns and their connection with urban ecology," in *Proceedings of the 16th ACM International Symposium on Mobile Ad Hoc Networking and Computing*, 2015, pp. 317-326. 30

[118] A. Furno, N.-E. El Faouzi, M. Fiore, and R. Stanica, "Fusing GPS probe and mobile phone data for enhanced land-use detection," in *2017 5th IEEE International Conference on Models and Technologies for Intelligent Transportation Systems (MT-ITS)*. IEEE, 2017, pp. 693-698. 33

[119] I. T. Joliffe and B. Morgan, "Principal component analysis and exploratory factor analysis," *Statistical methods in medical research*, vol. 1, no. 1, pp. 69-95, 1992. 34

[120] T. Arentze, H. Timmermans, F. Hofman, and N. Kalfs, "Data needs, data collection, and data quality requirements of activity-based transport demand models," *Transportation research circular*, no. E-C008, pp. 30-p, 2000. 35
[121] F. Giannotti and D. Pedreschi, Mobility, data mining and privacy: Geographic knowledge discovery. Springer Science & Business Media, 2008. 35

[122] P. Bonnel, “Postal, telephone, and face-to-face surveys: How comparable are they?” in Transport survey quality and innovation. Emerald Group Publishing Limited, 2003, pp. 215–237. 35

[123] J. Wolf, M. Oliveira, and M. Thompson, “Impact of underreporting on mileage and travel time estimates: Results from global positioning system-enhanced household travel survey,” Transportation research record, vol. 1854, no. 1, pp. 189–198, 2003. 35

[124] F. Wang and C. Chen, “On data processing required to derive mobility patterns from passively-generated mobile phone data,” Transportation Research Part C: Emerging Technologies, vol. 87, pp. 58–74, 2018. 35

[125] M. Fekih, T. Bellemans, Z. Smoreda, P. Bonnel, A. Furno, and S. Galland, “A data-driven approach for origin-destination matrix construction from cellular network signaling data: a case study of Lyon region (France),” Transportation, vol. 48, pp. 1671–1702, 2021. 36

[126] l. e. l. c. p. F. Centre d'études sur les réseaux, les transports, L'enquête ménages, déplacements" standard Certu". Ministère de l'écologie, de l'énergie, du développement durable et de l..., 2008. 37

[127] M. Shafiei, M. Nazemi, and S. Seyedabrishami, “Estimating time-dependent origin-destination demand from traffic counts: extended gradient method,” Transportation Letters, vol. 7, no. 4, pp. 210–218, 2015. 39

[128] A. Kuppam, R. Copperman, J. Lemp, T. Rossi, V. Livshits, L. Vallabhaneni, K. Jeon, and E. Brown, “Special events travel surveys and model development,” Transportation Letters, vol. 5, no. 2, pp. 67–82, 2013. 39

[129] M.-H. Wang, S. D. Schrock, N. Vander Broek, and T. Mulinazzi, “Estimating dynamic origin-destination data and travel demand using cell phone network data,” International Journal of Intelligent Transportation Systems Research, vol. 11, pp. 76–86, 2013. 39

[130] Z. Zhao, S.-L. Shaw, Y. Xu, F. Lu, J. Chen, and L. Yin, “Understanding the bias of cell detail records in human mobility research,” International Journal of Geographical Information Science, vol. 30, no. 9, pp. 1738–1762, 2016. 39

[131] M. Fekih, L. Bonnetain, A. Furno, P. Bonnel, Z. Smoreda, S. Galland, and T. Bellemans, “Potential of cellular signaling data for time-of-day estimation and spatial classification of travel demand: a large-scale comparative study with travel survey and land use data,” Transportation Letters, vol. 14, no. 7, pp. 787–805, 2022. 39, 58

[132] M. Seppecher, L. Leclercq, A. Furno, D. Lejri, and T. V. da Rocha, “Estimation of urban zonal speed dynamics from user-activity-dependent positioning data and regional paths,” Transportation Research Part C: Emerging Technologies, vol. 129, p. 103183, 2021. 39, 49, 51

[133] M. Janzen, M. Vanhoof, Z. Smoreda, and K. W. Axhausen, “Closer to the total? Long-distance travel of French mobile phone users,” Travel Behaviour and Society, vol. 11, pp. 31–42, 4 2018. 42

[134] M. Saberi, H. S. Mahmassani, T. Hou, and A. Zockaie, “Estimating network fundamental diagram using three-dimensional vehicle trajectories: extending die’s definitions of traffic flow variables to networks,” Transportation Research Record, vol. 2422, no. 1, pp. 12–20, 2014. 45, 47
[135] S. Uppoor, C. Ziemlicki, S. Secci, and Z. Smoreda, "On mobile traffic distribution over cellular backhauling network nodes," in *2016 13th IEEE Annual Consumer Communications & Networking Conference (CCNC)*. IEEE, 2016, pp. 726-731. 45

[136] N.-E. El Faouzi, A. Furno, L. Bonnetain, and M. Fiore, "Détermination de trajectoires à partir de données de téléphonie mobile," January 2023, FR3125197, EP4117318. 46

[137] M. Gramaglia, M. Fiore, A. Furno, and R. Stanica, "Glove: towards privacy-preserving publishing of record-level-truthful mobile phone trajectories," *ACM/IMS Transactions on Data Science (TDS)*, vol. 2, no. 3, pp. 1-36, 2021. 46, 74, 76, 77

[138] J. Zhang, F.-Y. Wang, K. Wang, W.-H. Lin, X. Xu, and C. Chen, "Data-driven intelligent transportation systems: A survey," *IEEE Transactions on Intelligent Transportation Systems*, vol. 12, no. 4, pp. 1624-1639, 2011. 47

[139] Y. Zheng, F. Liu, and H.-P. Hsieh, "U-air: When urban air quality inference meets big data," in *Proceedings of the 19th ACM SIGKDD international conference on Knowledge discovery and data mining*, 2013, pp. 1436-1444. 47

[140] J. Shang, Y. Zheng, W. Tong, E. Chang, and Y. Yu, "Inferring gas consumption and pollution emission of vehicles throughout a city," in *Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining*, 2014, pp. 1027-1036. 47

[141] L. Ntziachristos, D. Gkatzoflias, C. Kouridis, and Z. Samaras, "Copert: a european road transport emission inventory model," in *Information Technologies in Environmental Engineering: Proceedings of the 4th International ICSC Symposium Thessaloniki, Greece, May 28-29, 2009*. Springer, 2009, pp. 491-504. 47

[142] M. Seppecher, L. Leclercq, A. Furno, T. Vieira da Rocha, J.-M. André, and J. Boutang, "Identification of aggregate urban mobility patterns of nonregular travellers from mobile phone data," *Future Transportation*, vol. 3, no. 1, pp. 254-273, 2023. 47

[143] M. Seppecher, "Mining call detail records to reconstruct global urban mobility patterns for large scale emissions calculation," Ph.D. dissertation, Université de Lyon, 2022. 47

[144] C. F. Daganzo, "Urban gridlock: Macroscopic modeling and mitigation approaches," *Transportation Research Part B: Methodological*, vol. 41, no. 1, pp. 49-62, 2007. 48, 52

[145] M. Yildirimoglu and N. Geroliminis, "Approximating dynamic equilibrium conditions with macroscopic fundamental diagrams," *Transportation Research Part B: Methodological*, vol. 70, pp. 186-200, 2014. 48

[146] S. Batista, L. Leclercq, and N. Geroliminis, "Estimation of regional trip length distributions for the calibration of the aggregated network traffic models," *Transportation Research Part B: Methodological*, vol. 122, pp. 192-217, 2019. 48, 49

[147] G. Mariotte, L. Leclercq, S. Batista, J. Krug, and M. Paipuri, "Calibration and validation of multi-reservoir mfd models: A case study in Lyon," *Transportation Research Part B: Methodological*, vol. 136, pp. 62-86, 2020. 49, 54

[148] S. Batista, L. Leclercq, and M. Menendez, "Dynamic traffic assignment for regional networks with traffic-dependent trip lengths and regional paths," *Transportation Research Part C: Emerging Technologies*, vol. 127, p. 103076, 2021. 49

[149] G. Chen, S. Hoteit, A. C. Viana, M. Fiore, and C. Sarraute, "Enriching sparse mobility information in call detail records," *Computer Communications*, vol. 122, pp. 44-58, 2018. 50
[150] N. Geroliminis and C. F. Daganzo, "Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings," *Transportation Research Part B: Methodological*, vol. 42, no. 9, pp. 759-770, 2008. 52

[151] L. Bonnetain, A. Furno, J. Krug, and N.-E. E. Faouzi, "Can we map-match individual cellular network signaling trajectories in urban environments? data-driven study," *Transportation Research Record*, vol. 2673, no. 7, pp. 74-88, 2019. 55, 60, 62, 69, 70

[152] W. Wu, Y. Wang, J. B. Gomes, D. T. Anh, S. Antonatos, M. Xue, P. Yang, G. E. Yap, X. Li, S. Krishnaswamy et al., "Oscillation resolution for mobile phone cellular tower data to enable mobility modelling," in *2014 IEEE 15th International Conference on Mobile Data Management*, vol. 1. IEEE, 2014, pp. 321-328. 56

[153] B. C. Csáji, A. Browet, V. A. Traag, J.-C. Delvenne, E. Huens, P. Van Dooren, Z. Smoreda, and V. D. Blondel, "Exploring the mobility of mobile phone users," *Physica A: statistical mechanics and its applications*, vol. 392, no. 6, pp. 1459-1473, 2013. 56

[154] D. Bachir, G. Khodabandelou, V. Gauthier, M. El Yacoubi, and J. Puchinger, "Inferring dynamic origin-destination flows by transport mode using mobile phone data," *Transportation Research Part C: Emerging Technologies*, vol. 101, pp. 254-275, 2019. 56

[155] S. Jiang, J. Ferreira, and M. C. Gonzalez, "Activity-based human mobility patterns inferred from mobile phone data: A case study of singapore," *IEEE Transactions on Big Data*, vol. 3, no. 2, pp. 208-219, 2017. 58

[156] P. Katsikouli, M. Fiore, A. Furno, and R. Stanica, "Characterizing and removing oscillations in mobile phone location data," in *2019 IEEE 20th International Symposium on* "A World of Wireless, Mobile and Multimedia Networks" (WoWMoM). IEEE, 2019, pp. 1-9. 58

[157] A. A. Taha and A. Hanbury, "An efficient algorithm for calculating the exact hausdorff distance," *IEEE transactions on pattern analysis and machine intelligence*, vol. 37, no. 11, pp. 2153-2163, 2015. 60

[158] E. Akopyan, A. Furno, N.-E. El Faouzi, and E. Gaume, "Unsupervised real-time anomaly detection for multivariate mobile phone traffic series." in *ESANN*, 2021. 67

[159] L. Bonnetain, "Unlocking the potential of mobile phone data for large scale urban mobility estimation," Ph.D. dissertation, Université de Lyon, 2022. 69

[160] L. Bonnetain, A. Furno, and N.-E. E. Faouzi, "Multi-modal fine-grained map-matching of mobile phone network signaling data in urban area," in *101st Annual Meeting of the Transportation Research Board (TRB)*, 2022. [Online]. Available: https://people.licit-lyon.eu/furno/documents/trb_2022_bonnetain.pdf 69, 70, 71, 72

[161] Z. Shen, W. Du, X. Zhao, and J. Zou, "Dmm: Fast map matching for cellular data," in *Proceedings of the 26th annual international conference on mobile computing and networking*, 2020, pp. 1-14. 69, 70

[162] R. Mohamed, H. Aly, and M. Youssef, "Accurate real-time map matching for challenging environments," *IEEE Transactions on Intelligent Transportation Systems*, vol. 18, no. 4, pp. 847-857, 2016. 70

[163] E. Algizawy, T. Ogawa, and A. El-Mahdy, "Real-time large-scale map matching using mobile phone data," *ACM Transactions on Knowledge Discovery from Data (TKDD)*, vol. 11, no. 4, pp. 1-38, 2017. 70, 71