WSN Wireless Sensor Network ADMM Alternating Direction of Multipliers Method MSE Mean Square Error FFT Fast Fourier Transform SRG Strongly Regular Graph

SREL Success rate of existing links SRNL Success rate of non-existing links vii Acknownledgement First and foremost, I would like to express my appreciation to my supervisor, Dr. Alain Y. Kibangou, for his guidance, advice, kindness and support during my study in Necsteam, Gipsa-Lab of the

by means of the method of Lagrange multipliers. The second formulation (indirect approach) is obtained after an adequate re-parametrization. The problem is then convex and solved by using the distributed sub-gradient algorithm and the Alternating Direction of Multipliers Method (ADMM). Furthermore, three cases are considered: the final average value is perfectly known, approximated, or completely unknown. We also provide a way for computing the multiplicities of the estimated eigenvalues by means of integer programming.

In this spectral approach, given the Laplacian spectrum, the network topology can be reconstructed through estimation of Laplacian eigenvector. In particular, we study the reconstruction of the network topology in the presence of anonymous nodes.

The efficiency of the proposed solutions is evaluated by means of simulations. However, in several cases, convergence of the proposed algorithms is slow and needs to be improved in future works. In addition, the indirect approach is not scalable to very large graphs since it involves the computation of roots of a polynomial with degree equal to the size of the network. However, instead of estimating all the spectrum, it can be possible to recover only a few number of eigenvalues and then deduce some significant bounds on robustness indices.

Résumé

Durant la dernière décennie, une attention considérable a été consacrée au consensus des systèmes multi-agents. Le consensus est un processus coopératif dans lequel les agents interagissent afin de parvenir à un accord. La plupart des études ont été orientées vers l'analyse en régime permanent de ce processus d'agrément. Or, durant le régime transitoire, les algorithmes de consensus à convergence asymptotique génèrent un grand nombre de données. Dans cette thèse, notre objectif est d'exploiter ces données afin de concevoir des protocoles de consensus moyenne en temps fini, évaluer la robustesse du graphique, et éventuellement reconstituer la topologie du réseau de manière distribuée.

Le consensus de moyenne en temps fini garantit un temps d'exécution minimal qui peut assurer l'efficacité et la précision des algorithmes distribués complexes dans lesquels il est utilisé comme une sous-composante. Nous nous concentrons d'abord sur l'étape de configuration consacrée à la conception de protocoles de consensus qui garantissent la convergence vers la moyenne de la condition initiale dans un nombre donné d'étapes. En considérant des réseaux d'agents modélisés avec des graphes non orientés connectés, nous formulons ce problème comme étant un problème de factorisation de la matrice de moyenne J N = 1 N 11 T et étudions des solutions distribuées à ce problème. Puisque, les objets communicants doivent apprendre leur environnement avant d'établir des liens de communication, nous suggérons l'utilisation de séquences d'apprentissage afin de résoudre le problème de la factorisation. Ensuite, un algorithme semblable à l'algorithme de rétro-propagation du gradient est proposé pour résoudre un problème d'optimisation non convexe sous contrainte. Nous montrons que tout minimum local de la fonction de coût donne une factorisation exacte de la matrice de moyenne.

En contraignant les matrices facteur à être comme les matrices de consensus basées sur la matrice laplacienne, il est maintenant bien connu que la factorisation de la matrice de moyenne est entièrement caractérisée par les valeurs propres non nulles de la matrice Laplacienne du graphe. Par conséquent, la résolution de la factorisation de la matrice de moyenne de manière distribuée, avec une telle contrainte, permet d'estimer le spectre de la matrice Laplacienne. Puisque, spectre peut être utilisé pour calculer des indices de robustesse (Nombre d'arbres couvrant et résistance effective du graphe, aussi connu comme l'indice de Kirchhoff), la deuxième partie de cette thèse est consacrée à l'évaluation de la robustesse du réseau à travers l'estimation distribuée du spectre du Laplacien. Le problème est posé comme un problème de consensus sous contrainte formulé de deux façons différentes. La première formulation (approche directe) conduit à un problème d'optimisation non-convexe résolu de manière distribuée au moyen de la méthode des multiplicateurs de Lagrange. La seconde formulation (approche indirecte) est obtenue après une reparamétrisation adaptée. Le problème devient alors convexe et est résolu en utilisant l'algorithme du sous-gradient distribué et la méthode de direction alternée de multiplicateurs (ADMM). En outre, trois cas sont considérés: la valeur moyenne est parfaitement connue, elle est approximativement connue, ou complètement inconnue. Nous fournissons également une technique de calcul des multiplicités des valeurs propres estimées au moyen d'une programmation linéaire en nombres entiers.

Dans cette approche spectrale, compte tenu du spectre du Laplacien, la topologie du réseau peut être reconstruite à travers l'estimation des vecteurs propres du Laplacien. Nous proposons une nouvelle approche permettant à chaque noeud du réseau de pouvoir reconstruire la topologie du réseau même en cas de présence de noeuds anonymes.

L'efficacité des solutions proposées est évaluée au moyen de simulations. Cependant, dans plusieurs cas, la convergence des algorithmes proposés est lente et doit être améliorée dans les travaux futurs. La passage à l'échelle des méthodes développées est une question crucial qui mériterait d'être approfondie. Cependant, au lieu d'estimer tout le spectre, il peut être plus judicieux de ne restreindre l'estimation qu'à quelques valeurs propres et de pouvoir en déduire des bornes significatives sur les mesures de robustesse.

List of Notations

Context of the thesis 1.Multi-agent systems

Multi-agent systems (MASs) have received a growing interest in the last decades. They are developed for the demand of flexibility, robustness, and re-configuration features that appear in various applications domains including manufacturing, logistics, smart power grids, building automation, disaster relief, intelligent transportation systems, surveillance, environmental monitoring and exploration, infrastructure security and protection, etc.. A MAS is a system composed of multiple interacting intelligent agents (sensors, plants, vehicles, robots, etc.) and their environment as shown in Figure 1.1. An intelligent agent possesses the following key characteristics [START_REF] Ferber | Multi-Agent Systems: An Introduction to Distributed Artificial Intelligence[END_REF]:

• Autonomy: the ability to partially operate without the intervention of a human.

• Reactivity: the ability to react to changes in its environment and to behave correctly in order to satisfy its goal.

• Pro-activeness: the ability to take a responsibility for setting its own goals.

• Social ability: the ability to interact with humans and other agents.

• Local views: No agent has a full global view of the system, or the system is too complex for an agent to make practical use of such a knowledge.

A MAS can deal with tasks that are difficult or even impossible to be accomplished by an individual agent. During the recent decades, MASs have gained a widespread interest in many disciplines such as mathematics, physics, biology, computer science, social science. An increasing range of research topics in MASs includes cooperation and coordination, distributed computation, automatic control, wireless communication networks, etc..

In a typical centralized structure, a Fusion Center (FC) collects all measurements from the agents and then makes the final computations. However, due to the high information flow to FC, congestion can arise. Such a structure is vulnerable to FC failure. Also, the hardware requirements to build wireless communications can be one of reasons for an increase in the cost of the devices and thus, a higher overall cost of the network. For these reasons, a centralized structure can be inefficient. Hence, the research trend of MASs have shifted to decentralized structure where the interaction between agents is implemented locally without global knowledge. A good example is Wireless Sensor Network (WSN)s, which find broad application domains such as military applications (battlefield surveillance, monitoring friendly forces, equipment and ammunition,etc.), environment applications (forest fire detection, food detection,etc.), health applications (telemonitoring of human physiological data,etc.), home automation, formation control, etc.. Figure 1.2 depicts a WSN that collects data for the air quality, light intensity, sound volume, heat, precipitation and wind. (Source:http://mostepicstuff.com/chicago-city-installing-smart-sensors-city-wide-to-monitoreverything)

Motivation

There are three points that attract a great attention of the research community:

(A) Based on local information and interactions between agents, how can all agents reach an agreement?

This problem is called consensus problem, which is to design a network protocol based on the local information obtained by each agent such that all agents finally reach an agreement on certain quantities of interest.

Consensus problems of MASs have received tremendous attention from various research communities due to their broad applications in many areas including multisensors data fusion [START_REF] Olfati-Saber | Consensus filters for sensor networks and distributed sensor fusion[END_REF], flocking behavior of swarms [START_REF] Blondel | Convergence in multiagent coordination, consensus, and flocking[END_REF][START_REF] Jadbabaie | Coordination of groups of mobile autonomous agents using nearest neighbor rules[END_REF][START_REF] Olfati-Saber | Flocking for multi-agent dynamic systems: algorithms and theory[END_REF][START_REF] Olfati-Saber | Swarms on sphere: A programmable swarm with synchronous behaviors like oscillator networks[END_REF], multi-vehicle formation control [START_REF] Dimarogonas | Formation control and collision avoidance for multi-agent systems and a connection between formation infeasibility and flocking behavior[END_REF][START_REF] Moore | Formation control via distributed optimization of alignment error[END_REF][START_REF] Olfati-Saber | Distributed cooperative control of multiple vehicle formations using structural potential functions[END_REF], distributed computation [START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF][START_REF] Borkar | Asymptotic agreement in distributed estimation[END_REF][START_REF] Lynch | Distributed Algorithms[END_REF], rendez-vous problems [START_REF] Cortes | Robust rendezvous for mobile autonomous agents via proximity graphs in arbitrary dimensions[END_REF] and so on. More specifically, average consensus algorithms (i.e. the consensus value corresponds to the average of the initial states) are commonly used as building block for distributed control, estimation or inference algorithms.

In the recent literature, one can find average consensus algorithms embedded in the distributed Kalman filter, [START_REF] Olfati-Saber | Distributed Kalman filtering for sensor networks[END_REF]; Distributed Least Squares algorithms, [START_REF] Bolognani | Distributed sensor calibration and least-squares parameter identification in WSNs using consensus al-gorithms[END_REF]; Distributed Alternating Least Squares for tensors factorization, [START_REF] Kibangou | Distributed PARAFAC based DS-CDMA blind receiver for wireless sensor networks[END_REF]; Distributed Principal Component Analysis, [START_REF] Macua | Consensus-based distributed principal component analysis in wireless sensor networks[END_REF]; or distributed joint input and state estimation, [START_REF] Esna-Ashari | Distributed input and state estimation for linear discrete-time systems[END_REF] to cite few. However, the asymptotic convergence of the consensus algorithms is not suitable for these kinds of sophisticated distributed algorithms. A slow asymptotic convergence can not ensure the efficiency and the accuracy of the algorithms, which can lead to other unexpected effects. For example, regarding to the WSNs, a reduction in the total number of iterations until convergence can lead to a reduction in the total amount of energy consumption of the network, which is essential to guarantee a longer life time for the entire network. On the other hand, the protocols that guarantee a minimal execution time are much more appealing than those ensuring asymptotic convergence. For this purpose, several contributions dedicated to finite-time consensus have been recently published in the literature, meaning that, consensus is obtained in a finite number of iterations.

(B) Can we assess the robustness of a network?

Generally, we are surrounded by networks of different kinds (social networks, sensor networks, power networks, etc.). Therefore, the most important thing is that these networks are robust. Meaning that, they can undergo through damage or failure. Hence, the terminology robustness rapidly invades the research community. Definition 1

Robustness is the ability of a network to continue performing well when it is subject to failures or attacks [START_REF] Ellens | Effective graph resistance[END_REF].

In order to assess whether the network is robust, measuring the network robustness is needed [START_REF] Ellens | Effective graph resistance[END_REF]. There are various proposed approaches to measure network robustness metrics. However, the most popular approach in the literature is based on the analysis of the graph. From data produced during consensus protocol can we infer network robustness?

(C) Can we estimate the network topology?

Network topology is usually a schematic description of the arrangement of a network, including its nodes and connecting edges. Network topology may describe how data is transferred between these nodes.

Network topology identification refers to detecting and identifying the interested network elements and the relationship between elements in this network and represents the topology construction in an appropriate form. Therefore, identification of networks of systems becomes a growing attractive task for solving many problems in different science and engineering domains.

Theoretically, there are some matrices that are associated with the graph representing a network. Therefore, one can find a method to identify the network topology via the state matrix in a state space representation of linear dynamic systems, which is related to the network topology. Eigenvalue decomposition of a matrix shows that the structure of this matrix is linked to both eigenvalues and eigenvectors of that matrix. In the other words, estimating both the eigenvalues and the eigenvectors of network matrix can obviously infer the network structure.

Goals of the thesis

It is well-known that consensus protocols being iterative, a huge amount of data is produced. However, most of studies focus only on convergence properties of such algorithms. With data collected from each iteration of the asymptotic consensus protocol in a distributed way, we aim to:

• Design finite-time average consensus protocols.

• Assess the network robustness.

• Reconstruct the network topology.

The consensus problem

Consensus issue in networks of autonomous agents has been widely investigated in various fields, including computer science and engineering. In such networks, according to an a priori specified rule, also called protocol, each agent updates its state based on the information received from its neighbors with the aim of reaching an agreement to a common value. When the common value corresponds to the average of the initial states, average consensus is to be achieved.

Preliminaries

Let us consider a network modelled as a graph. In what follows, we first give some basic notations and definitions from graph theory [START_REF] Godsil | Algebraic graph theory[END_REF][START_REF] Mesbahi | Graph theoretic methods in multiagent networks[END_REF].

Given a connected undirected graph G(V, E), where V = {1, 2, . . . , N} is the set of vertices (agents or nodes), and E ⊆ V × V is the set of edges (links) between agents.

The set of neighbors of node i and its degree are denoted by N i = {j|(i, j) ∈ E} and

d i = |N i | respectively.
It is common to resort to some matrices for characterizing graphs. That is the case of adjacency matrix A, degree matrix D of the graph, which has vertex degrees d i , i ∈ V on its diagonal and zeroes elsewhere, Laplacian Matrix L = D -A, clearly defined in Chapter 2.

Through out this thesis, an undirected graph G(V, E) is mainly considered.

Example:

Consider an arbitrary network of 5 agents communicating with each other as described in Figure 1. 3. Each agent has an initial value. A consensus protocol is an interaction rule that specifies the information exchange between an agent and all of its neighbors on the network to reach an agreement regarding a certain quantity of interest that depends on the state of all agents. Informally, despite the initial values of all agents, they converge to the common value (in this case, the average of initial values).

Consensus algorithms

The literature on consensus protocols can be organized as in Figure 1.4.

In this thesis, we concentrate on the consensus problem for discrete-time systems with fixed topology. Given a graph G(V, E), each node has an associated value x i defined as the state of node i. Let x(0) = [x 1 (0), x 2 (0), . . . , x N (0)] T be the vector of initial states of the given network. In general, given the initial values at each node x i (0), i ∈ V , the main task is to compute the final consensus value using distributed linear iterations. Each iteration involves local communication between nodes. In particular, each node repeatedly updates its value as a linear combination of its own value and those of its neighbors. The main benefit of using a linear-iteration scheme is that, at each time-step, each node only has to transmit a single value to each of its neighbors. The linear iteration-based consensus update equation is expressed as:

x i (k + 1) = w ii (k)x i (k) + j∈N i w ij (k)x j (k), i = 1, 2, . . . , N
or equivalently in matrix form:

x(k + 1) = W(k)x(k),
where: W(k) is the matrix with entries w ij (k) = 0 if (i, j) / ∈ E and

j∈N i ∪{i} w ij (k) = 1.
The asymptotic consensus is reached if lim k→∞ x(k) = µ1, meaning that all the nodes agreed on the value µ. When µ is equal to the average of the initial values, i.e.

µ = 1 N N i=1
x i (0), it is called average consensus problem.

The distributed average consensus algorithm is described as follows:

Algorithm 1 (Distributed average consensus algorithm)

1. Input:

• Network topology modelled with a graph G(V, E),

• Initial values x(0),

• Maximum number of iterations k max .

2. Initialization: Define weight consensus matrix W(k), k = 0, 1, . . ..

3.

Set k = 0.

While k < k max do

• Update:

x i (k + 1) = w ii (k)x i (k) + j∈N i w ij (k)x j (k), i = 1, . . . , N. • k = k + 1;
5. end 6. Output: x(k max ) = µ1.

Convergence conditions

According to the classification of consensus algorithms in Figure 1.4, the convergence analysis is also studied for specified categories such as static topology and dynamic topology, directed and undirected topologies. In this thesis, we focus on distributed discrete-time average consensus protocols for a network modelled by an undirected graph with fixed topologies. Meaning that, the evolution of the x(k) can be written in matrix form as follows:

x(k + 1) = Wx(k), k = 0, 1, . . . (1.2.1) or equivalently, x(k) = W k x(0). (1.2.2)
Generally speaking, the system is said to achieve distributed consensus asymptotically if

lim k→∞ x(k) = lim k→∞ W k x(0) = 1c c c T x(0), (1.2.3)
where 1 1 1 ∈ R N ×1 is the all-ones column vector and c c c T is some constant row vector. Therefore, it is good to view the conditions on the weight matrix W ensuring the convergence of the consensus protocol as well as the average consensus protocol:

lim k→∞ W k = 1c 1c 1c T . (1.2.4)
The convergence conditions are described as follows:

Theorem 1.2.1 ([74])
Consider the linear iteration protocol (1.2.1), distributed consensus is achieved if and only if the weight consensus matrix W satisfies the following conditions:

1. W1 = 1 1 1 2. ρ(W -1c 1c 1c T ) < 1,
where c c c is chosen so that 1 1 1 T c c c = 1 and ρ(W -1c 1c 1c T ) is the spectral radius of W -1c 1c 1c T .

To sum up, the weight matrix has row-sum equal to 1, 1 is a simple eigenvalue of W and all other eigenvalues are strictly less than one in magnitude. It means that the weight matrix W is a row-stochastic matrix.

When c c c = 1 N 1 1 1, the system is said to achieve average consensus (i.e the consensus value is the average of all initial values):

lim k→∞ W k = 1 N 11 11 11 T . (1.2.5)
In order to achieve average consensus, beside the conditions in Theorem 1.2.1, the left and right eigenvector of W corresponding to eigenvalue 1 are c T and 1, respectively.

Theorem 1.2.2

The equation (1.2.5) holds if and only if [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF]:

1.

1 T W = 1 1 1 2. W1 = 1 1 1 3. ρ(W -1 N 11 11 11 T ) < 1
A necessary and sufficient condition to ensure the convergence of the average consensus protocol for undirected graphs is that the weight consensus matrix W is a doubly stochastic matrix.

1.2.4 Design of the weight matrix.

In the literature, there are some works devoted to the design of the weight matrix W that satisfies the convergence conditions of the average consensus algorithms. For instance, in [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF][START_REF] Xiao | Distributed average consensus with least-mean-square deviation[END_REF], it has been shown that the following weight matrices satisfy the convergence conditions pointed out above.

1. Maximum-degree weights:

An approach to design the weight matrix W in a graph with fixed topology (timeinvariant topology) consists in assigning a weight on each edge equal to the inverse of the maximum degree of the network, i.e.,

w ij =            1 dmax+1 if j ∈ N i 1 -d i dmax+1 if i = j 0 otherwise (1.2.6)
where

d max = max i∈V d i ≤ N.
We can see that before running the consensus algorithm, the maximum degree d max should be known in a priori. The required information on the graph topology is global, but can be obtained in a distributed way by running a max consensus protocol [START_REF] Olfati-Saber | Consensus problems in networks of agents with switching topology and time-delays[END_REF].

Metropolis hasting weights:

The metropolis weights for a graph with a time-invariant topology are proposed in [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF]:

w ij =              1 1+max{d i ,d j } , if j ∈ N i 1 - j∈N i 1 1+max{d i ,d j } , if i = j 0, otherwise (1.2.7)
Here, the nodes do not need any global information of the communication graph, not even the number of the vertices N. They only need to know their own degree and that of their neighbors to compute the weights. The metropolis weights are very simple to compute and suited to a distributed implementation.

Constant edge weights:

By a result in [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF], the simplest approach is to set all the edge weights for neighboring nodes equal to a constant α. This is the most widely applied model for the weight matrix in both time-variant and time-invariant topologies. The weight matrix is defined as follows:

w ij =            α, if j ∈ N i 1 -α|N i |, if i = j 0, otherwise (1.2.8)
where | • | denotes the cardinality of the set. The weight matrix can be expressed in matrix form as follows:

W = I N -αL, (1.2.9)
where L, I N are the Laplacian matrix of the graph G and the N × N identity matrix respectively. In general, α > 0 for the convergence condition 1.2.2 to hold. Additionally, we can express the eigenvalues of W in terms of those of L:

λ i (W) = 1 -αλ i (L), i = 1, . . . , N,
where λ i (.) denotes the i-th largest eigenvalue of L. Therefore,

ρ(W - 1 N 11 T ) = max{λ 2 (W), -λ N (W)} = max{1 -αλ 2 (L), αλ N (L) -1}. (1.2.10)
The value of α that can ensure the convergence condition of the consensus protocols is 0 < α < 2 λ N (L) . There are some simple choices that do not require the knowledge of the Laplacian spectrum, for example, the upper bound of the largest Laplacian eigenvalue λ N ≤ 2d max . Hence, we can have a simple choice of α as 0 < α < 1 dmax . Moreover, following the result in [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF], the best possible constant edge weight in terms of convergence speed is given by the following α:

α = 2 λ N (L) + λ 2 (L) . (1.2.11)
In the analysis of consensus problems, convergence rate is an important index to evaluate the performance of a consensus protocol. Therefore, there are several works dealing with accelerating the rate of convergence of the consensus protocol by solving some optimization problems in a centralized way or using Chebyshev polynomials.

Since the conditions of the weight consensus matrix W obviously influence on the convergence of consensus algorithms, the main direction of research is now shifted to the computation of W in order to improve the convergence rate, [START_REF] Kokiopoulou | Polynomial filtering for fast convergence in distributed consensus[END_REF][START_REF] Montijano | Fast distributed consensus with Chebyshev polynomials[END_REF][START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF].

A. [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF] proposed an optimization method to obtain the optimum weight matrix W achieving average consensus in linear time-invariant topologies as the solution of a semi-definite convex programming.

min W ρ(W -1 N 11 T ) subject to W ∈ S G , 1 T W = 1 T , W1 = 1.
The condition W ∈ S G expresses the constraint on the sparsity pattern of the matrix W with the set S G defined as follows:

S G = {W ∈ R N ×N |w ij = 0 if (i, j) / ∈ E and i = j}.
B. In [START_REF] Kokiopoulou | Polynomial filtering for fast convergence in distributed consensus[END_REF], the author proposed an interesting method to accelerate the speed of convergence through polynomial filtering. The main idea is to apply a polynomial filter on the consensus matrix W that will shape its spectrum in order to increase the convergence rate.

Denote by p k (λ(W)) the polynomial filter of degree k that is applied on the spectrum of W,

p k (λ(W)) = α 0 + α 1 λ(W) + α 2 λ 2 (W) + . . . + α k λ k (W). (1.2.12)
The matrix polynomial is defined as:

p k (W) = α 0 I N + α 1 W + α 2 W 2 + . . . + α k W k . (1.2.13)
where α 0 , . . . , α k are the coefficients of the polynomial.

For a given matrix W, and a certain degree k, instead of finding weight matrix W to minimize the spectral radius ρ(W -1 N 11 T ), the problem is now to find the polynomial coefficients α i , i = 0, . . . , k such that min

α α α∈R k+1 ρ( k i=0 α i W i -1 N 11 T ) subject to ( k i=0 α i W i )1 = 1.
After obtaining the coefficients α i , one can run the polynomial filtered distributed consensus proposed in [START_REF] Kokiopoulou | Polynomial filtering for fast convergence in distributed consensus[END_REF]:

Algorithm 2 (Polynomial filtered distributed consensus) (a) Inputs:

• polynomial coefficients α i , i = 0, . . . , k + 1, k degree of the polynomial filter.

• tolerance ε.

(b) Initialization: Consensus Matrix W, x(0);

(c) Set t=1; (d) Repeat • if mod(t, k + 1) == 0 then x i [t] = α 0 x i [t -k -1] + α 1 x i [t -k] + α 2 x i [t -k + 1] + . . . + α k x i [t -1]. x i [t] = w ii x i [t] + j∈N i w ij x j [t], i = 1, 2, . . . , N else x i [t] = w ii x i [t -1] + j∈N i w ij x j [t -1]
• end if

• t = t + 1; • x[t] = x i [t]; (e) until x[t] -x[t -1] < ε. (f) Outputs x[t]
C. Chebyschev polynomials for distributed consensus, [START_REF] Montijano | Fast distributed consensus with Chebyshev polynomials[END_REF]:

Denote the Chebyshev polynomial of degree k by T k (x). These polynomials satisfy:

T k (x) = cos(k arccos x), ∀x ∈ [-1, 1],
and |T k (x)| > 1 when |x| > 1, ∀k ∈ N. These polynomials are defined by recurrence:

• T 0 (x) = 1, • T 1 (x) = x, • T k (x) = 2xT k-1 (x) -T k-2 (x), k ≥ 2.
Chebyshev polynomials have the following properties:

• T k (x) = cos(k arccos x) ≤ 1, for x ∈ [-1, 1]. • max x∈[-1,1] |T k (x)| = 1, T k (-1) = (-1) k and |T k (x)| > 1, ∀|x| > 1.
On the other hand, the direct expression of T k (x) is characterized by:

T k (x) = 1 + τ 2k 2τ k , τ = τ (x) = x - √ x 2 -1.
In order to speed up the consensus, the main idea consists in designing a distributed linear iteration such that the execution of a fixed number of n steps is equivalent to the evaluation of some polynomial, P k (x), in the fixed matrix W. The polynomial must satisfy that P k (1) = 1 and

|P k (x)| < 1 if |x| < 1.
With two real coefficients κ, υ (-1 < κ < υ < 1), the polynomial is defined as:

P k (x) = T k (cx -d) T k (c -d) , with c = 2 κ -υ , d = κ + υ κ -υ , ∀k
which has the following properties:

• if x ∈ [κ, υ], then cx -d ∈ [-1, 1]. • P k (1) = 1 and P k (κ + υ -1) = (-1) k . • |P k (x)| < 1, ∀x ∈ (κ + υ -1, 1) and |P k (x)| > 1 otherwise.
The polynomial satisfies the recurrence:

P k (x) = 2 T k-1 (c -d) T k (c -d) (cx -s)P k-1 (x) - T k-2 (c -d) T k (c -d) P k-2 (x).
The consensus rule is then:

x(k) = P k (W)x(0), (1.2.14)
where,

x(1) = P 1 (W)x(0) = 1 T 1 (c -d) (cW -dI N )x(0), and, x(k) = P k (W)x(0) = 2 T k-1 (c-d) T k (c-d) (cW -dI N )P k-1 (W) -T k-2 (c-d) T k (c-d) (cW -dI N )P k-2 (W) x(0) = 2 T k-1 (c-d) T k (c-d) (cW -dI N )x(k -1) -T k-2 (c-d) T k (c-d) (cW -dI N )x(k -2), k ≥ 2.
The choice of the polynomial determine the convergence speed of the algorithm, given by max

λ i (W) |P k (λ i (W))|, where λ i (W) is the i -th eigenvalue of W.
The convergence rate is given by: max

λ i (W) =1 |T k (cλ i -d)| T k (c -d) = max{ |T k (cλ N (W) -d)| T k (c -d) , |T k (cλ 2 (W) -d)| T k (c -d) }
Using the Chebyshev polynomials means that we are minimizing max λ∈[-1,1] P k (λ), therefore, getting high chances to obtain a good convergence rate. By an optimal selection of κ, υ, we can maximize the convergence speed. Following the Theorem 4.2 in [START_REF] Montijano | Fast distributed consensus with Chebyshev polynomials[END_REF], the optimal parameters are κ = λ N (W) and υ = λ 2 (W).

Let compare the convergence speed between the standard consensus protocol (1.2.1) and the consensus protocol using the Chebyshev polynomials (1.2.14).

Denoting by v i , λ i (W), i = 1, . . . , N the eigenvectors and eigenvalues of W, (1.2.1) can be expressed as:

x(k) = W k x(0) = v 1 + λ k 2 (W)v 2 + . . . + λ k N (W)v N . Theorem 1.2.3 ([53]) Let λ = max(|λ 2 (W)|, |λ N (W)|)
be the convergence rate of the standard consensus protocol (1.2.1). For any 0 < υ < 2λ λ 2 +1 and κ = -υ, P k (λ) goes to zero faster than λ k when k goes to infinity. Therefore, the consensus protocol using Chebyshev polynomial converges faster than the standard one.

Since the average consensus protocols can be embedded in more sophisticated distributed algorithms, protocols that guarantee a minimal execution time are more appealing than those ensuring asymptotic convergence. For this purpose, several contributions dedicated to finite-time consensus have been recently published in the literature, meaning that consensus is obtained in a finite number of steps. The main contribution of this thesis is to address the design of finite-time average consensus problem by making use of distributed optimization methods.

Network robustness

The aim of network robustness research is to find a robustness measure to evaluate the performance of a network, [START_REF] Dekker | Network robustness and graph topology[END_REF][START_REF] Ellens | Effective graph resistance[END_REF]. Furthermore, understanding whether a network is robust can protect and improve the performance of the network efficiently. By the way, it is also used to design new networks that are able to perform well when facing with failures or attacks.

Several robustness measures have been proposed in the literature. However, we are only interested in topological measures, since it is well known that certain topologies exhibit high robustness against failures. Therefore, they can be used to change the network topology and decrease failures.

Vertex (Edge) connectivity

The vertex (edge) connectivity of an incomplete graph K v (K e ) respectively represents the minimal number of vertices (edges) that has to be removed to disconnect the graph. Hence, K v (K e ) depends on the least connected part of the graph. It has been shown in [START_REF] Godsil | Algebraic graph theory[END_REF] that K v ≤ K e ≤ d min , where d min is the minimum degree of the network. For a complete graph of N nodes, K v = K e = N -1.

Remark 1

• The measure value K v (K e ) is an integer,

• The larger K v (K e ), the harder it is to disconnect the graph, thus the graph is more robust.

• In order to compute K v (K e ), the network topology should be known at the first sight.

Algebraic connectivity

The algebraic connectivity is the second smallest eigenvalue of the Laplacian matrix. If λ 2 (L) = 0, the graph is disconnected. The algebraic connectivity of an incomplete graph is not greater than the vertex connectivity

K v : 0 ≤ λ 2 (L) ≤ K v ≤ K e ≤ d min
Remark 2

• A higher value of λ 2 (L) means that the graph is more robust.

• Sometimes, there is a problem that the algebraic connectivity is not strictly increasing when an edge is added, [START_REF] Ellens | Effective graph resistance[END_REF].

Number of spanning trees ξ

A spanning tree is a sub-graph containing N -1 edges, all N vertices, and no cycles. The number of spanning trees is defined as a function of the Laplacian eigenvalues [START_REF] Ellens | Effective graph resistance[END_REF]

ξ = 1 N N i=2 λ i (L) (1.3.1)
The number of spanning trees can be used to assess the robustness of the network, [4].

Remark 3

The greater ξ, more robust is the graph.

Effective graph resistance R

Assume the graph is seen as an electrical circuit, where an edge (i, j) corresponds to a resistor of unit resistance R ij = 1 Ohm. Informally, the effective resistance between two vertices of a network (when a voltage is applied across them) can be calculated by series and parallel manipulations. The effective graph resistance is the sum of the effective resistances over all pairs of vertices, [START_REF] Ellens | Effective graph resistance[END_REF].

R = 1≤i<j≤N R ij = N N i=2 1 λ i (L) (1.3.2)
The effective graph resistance is the sum of the inverse of nonzero Laplacian eigenvalues.

Remark 4

The smaller R, more robust is the graph.

Outline of this Thesis

For the general picture, the contributions of this thesis can be viewed as three-fold one as shown in the following scheme. The first task is dedicated to the design of the finite-time average consensus protocols, while the second task focuses on network robustness assessment. The remaining goal is dedicated to network topology reconstruction. Hereinafter, we present the outline of this thesis as follows:

Distributed Estimation of Laplacian Eigenvalues

• Chapter 2: Graph Theory.

In this chapter, we give general definitions of graph theory.

• Chapter 3: Distributed Design of Finite-time Average Consensus Protocols. This is one of the main contributions of this thesis, where the design of the finitetime average consensus is solved in a distributed way. The problem is formulated as a factorization of the averaging matrix J N = 11 T N . The solution is obtained by solving a constrained optimization problem in a distributed way. For the simple understanding, each linear iteration k, k = 0, . . . , D with D being the number of steps is illustrated as a layer of a network (including D + 1 layers). The main idea is to compute the partial derivatives of the error between the actual output and desired output, then propagate them back to each iteration k in order to update the weights. Then, the proposed algorithm is based on the back-propagation method and an associated gradient descent method. The evaluation of the proposed method is assessed by numerical results.

• Chapter 4: Distributed network robustness assessment.

In this chapter, we make use of the metrics which are the number of spanning tree ξ and the effective graph resistance R to assess the robustness of a given network. By the nature of these metrics, they are all functions of the Laplacian eigenvalues λ i (L), hence, the Laplacian eigenvalues estimation is now the main purpose of this task. In this chapter, instead of estimating the Laplacian eigenvalues directly, we solve this problem as a constrained optimization problem with respect to the inverse of the nonzero distinct Laplacian eigenvalues. As the result, by applying the gradient descent method, we finally get the set of nonzero distinct Laplacian eigenvalues. Furthermore, by using an integer programming algorithm, we can estimate the multiplicities corresponding to these eigenvalues. Then, we deduce the whole spectrum of the Laplacian matrix associated with the graph G(V, E). Thanks to these results, the robustness measures can be estimated to evaluate the robustness of the network.

For the task of Laplacian spectrum estimation, all methods proposed in this chapter are divided into two categories: non-convex and convex formulations. According to non-convex optimization problem (it is also called direct method), the problem is viewed as a constrained optimization problem with respect to the stepsizes α k of the Laplacian-based consensus protocol, whose inverses are equal to the nonzero distinct Laplacian eigenvalues.

Solving the non-convex optimization problem, the main issue is the slowness of the convergence speed. Therefore, we try to make the speed of the estimation of the Laplacian eigenvalues as fast as possible. The main point is to transform the nonconvex problem into a convex one. We solve the problem by two algorithms using the distributed projection method, [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF] and the Alternating Direction of Multipliers Method (ADMM), [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF]. We consider three cases: final average value perfectly known, noisy and completely unknown.

• Chapter 5: Network topology reconstruction.

Since, the components used to infer the network structure are the eigenvalues and eigenvectors of the weight consensus matrix, we first propose an algorithm to obtain these eigenvectors. Besides that, the corresponding eigenvalues can be obtained by using the algorithms proposed in the previous chapters. Finally, we can make use of received eigenvectors and eigenvalues to reconstruct the network topology.

• Chapter 6: Conclusion and Future Works. This chapter summarises the contributions of this thesis and also outlines possible future directions.

Publications list

The main contributions of this thesis are the principal subject of the following publications:

1.5.1 International conference papers with proceedings According to the communication policy, the graph G(V, E) can be undirected or directed.

• Undirected Graph:

The communication topology of the network of N agents is represented by a graph G(V, E), where V = {1, 2, . . . , N} is the set of vertices (agents or nodes), and E ⊆ V × V is the set of edges (links between agents).

If there is no direction assigned to the edges, then both edges (i, j) and (j, i) are included in the set of edges E. The graph is called undirected graph.

If agent i and j are connected, then the link between i and j is included in E, (i, j) ∈ E and i and j are called neighbors. The set of neighbors of agent i is denoted by N i and its degree is denoted by d i = |N i |, where |.| stands for the cardinality.

• Directed Graph:

If a direction is assigned to the edges, the relations are asymmetric and the graph is called a directed graph (or a digraph). For a directed graph as shown in Figure 2.1, for a directed edge (i, j), i is called the head and j is called the tail.

A vertex i is connected to j by a directed edge, or that j is a neighbor of i if (i, j) ∈ E. The edge (i, j) is then an outgoing edge for i and an ingoing edge for j.

The out-degree d out i of a vertex i is its number of outgoing edges, or number of vertices to which it is connected. The in-degree d in i is its number of ingoing edges or number of vertices that are connected to it. • a 6-node undirected graph with set of vertices V = {1, 2, . . . , 6}, and set of edges E = {(1, 2), [START_REF] Acharya | On the hardness of topology inference[END_REF]3), [START_REF] Acharya | On the hardness of topology inference[END_REF]4), [START_REF] Acharya | On the hardness of topology inference[END_REF][START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF], [START_REF] Aderhold | Hierarchical Bayesian models in ecology: Reconstructing species interaction networks from non-homogeneous species abundance data[END_REF]3), [START_REF] Aderhold | Hierarchical Bayesian models in ecology: Reconstructing species interaction networks from non-homogeneous species abundance data[END_REF][START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF], (3,4), (4,[START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF], (4,[START_REF] Blondel | Convergence in multiagent coordination, consensus, and flocking[END_REF], [START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF][START_REF] Blondel | Convergence in multiagent coordination, consensus, and flocking[END_REF])}

i∈V d out i = i∈V d in i = |E|.
• a 5-node directed graph with set of vertices V = {1, 2, . . . , 5}, and set of edges E = {(1, 5), (1, 3), (2, 1), (2, 4), (3, 2), (3,[START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF], (4, 3), (4, 5)}.

Connectivity of a graph

A path from a vertex i to a vertex j is a sequence of distinct vertices starting with vertex i and ending with vertex j such that consecutive vertices are adjacent. A simple path is a path with no repeated vertices.

• In an undirected graph G, two vertices i and j are connected if there is a path from i to j. An undirected graph G is connected if for any two vertices in G there is a path between them, see Figure 2.2. Conversely, two vertices i and j in G are disconnected if there is no path from i to j. An undirected graph G is disconnected if we can partition its vertices into two nonempty sets χ and Γ such that no vertex in χ is adjacent to a vertex in Γ, see Figure 2.3.

• A directed graph is strongly connected if between every pair of distinct vertices (i, j) in G, there is a directed path that begins at i and ends at j. It is called weakly connected if replacing all of its directed edges with undirected edges produces a connected undirected graph.

• A graph is said to be complete (fully-connected) if every pair of vertices has an edge connecting them, meaning that the number of neighbors of each vertex is equal to N -1, see Figure 2.4.

Algebraic graph properties

• Two vertices joined by an edge are called the endpoints of the edge. If vertex i and vertex j are endpoints of the same edge, then i and j are said to be adjacent to each other. In an undirected graph, vertices that are adjacent to a vertex i are called the neighbors of i. The set of all neighbors of a vertex i is defined as

N i = {j ∈ V : (i, j) ∈ E}.
• Given two vertices i and j, the distance dist(i, j) is the length of the shortest simple path between i and j.

• The eccentricity ǫ i of a vertex i is the greatest distance between i and any other vertex j ∈ V .

• The radius r(G) of a graph G is the minimum eccentricity of any vertex.

• The diameter d(G) of a graph is the maximum eccentricity of any vertex in the graph, d(G) = max i,j dist(i, j).

• The structure of a graph with N nodes is described by means of an N × N matrix. The adjacency matrix A is the matrix with entries a i,j given by

a ij =    1, if (i, j) ∈ E 0, otherwise (2.2.1)
meaning that, the (i, j) entry of A is 1 only if vertex j is a neighbor of vertex i.

The diagonal entries of A are all equal to 0. If G is an undirected graph, a ij = a ji , i.e., A is symmetric. In case of a directed graph, A is asymmetric.

For instance, take the undirected graph in Figure 2.1 as an example, we have the adjacency matrix defined as follows:

A =             0 1 1 1 1 0 1 0 1 0 1 0 1 1 0 1 0 0 1 0 1 0 1 1 1 1 0 1 0 1 0 0 0 1 1 0            
An adjacency matrix for a directed graph is defined as follows:

A =             0 0 1 1 1 0 1 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 1 0 0 0 1 0 0 0 1 0 0            
• The in-degree and out-degree of a vertex i are defined by the sums of the weights of the outgoing and the incoming edges respectively, i.e.,

d in i = N j=1 a ij and d out i = N j=1 a ji .
A vertex i is said to be balanced if its in-degree and out-degree are equal, i.e., d in i = d out i . Therefore, all undirected graphs are balanced graphs.

• A digraph G is called balanced if N j =i a ij = N j =i a ji .
The degree matrix D of G is the N × N diagonal matrix with (i, j) entry given by:

D ij =    d out i , if i = j 0, otherwise (2.2.2)
or equivalently, D = diag(A1).

(2.2.3)

Spectral graph properties

• The Laplacian matrix is used for mathematical convenience to describe the connectivity in a more compact form. The graph Laplacian L is defined as the matrix with entries l ij given by:

l ij =      N k=1,k =i a ik if i = j -a ij if i = j (2.3.1)
The Laplacian matrix L can be expressed in matrix form as follows:

L = D -A. (2.3.2)
Again, we take the 6-node graph in Figure 2.1 as an example:

1 2 3 4 5 6 Undirected graph L =             4 -1 -1 -1 -1 0 -1 3 -1 0 -1 0 -1 -1 3 -1 0 0 -1 0 -1 4 -1 -1 -1 -1 0 -1 4 -1 0 0 0 -1 -1 2            
• Some important properties of the Laplacian matrix are:

1. As can be seen in the definition of L, the row-sum of L is zero. Since L1 = 0, where 0 is a zero vector of length N, L has at least one eigenvalue zero with associated right eigenvector 1.

If the graph G(V, E) is undirected, then:

2. The spectrum of L, that is the set of all eigenvalues of the Laplacian matrix, is denoted as

sp(L) = {λ m 1 1 , λ m 2 2 , . . . , λ m D+1
D+1 , }, where m i stands for the multiplicity of the ith eigenvalue and D is the number of nonzero distinct Laplacian eigenvalues. The multiplicity of the zero eigenvalue represents the number of connected components of the graph. For undirected graphs, L is symmetric and positive semi-definite, hence 0 = λ 1 < λ 2 < . . . < λ D+1 < 2d max , and 3. λ 2 is the second smallest eigenvalue, which is known as the algebraic connectivity and reflects the degree of connectivity of the graph G, [START_REF] Friedler | Algebraic connectivity of graphs[END_REF]. It is well-known that the speed of convergence of a consensus protocol depends on this second smallest Laplacian eigenvalue λ 2 , [START_REF] Godsil | Algebraic graph theory[END_REF].

min x =0,1 T x=0 x T Lx x 2 = λ 2 (L). (2.3.3)
From the literature, the second smallest Laplacian eigenvalue is upper bounded by max

j∈N i {d i + d j }.
4. In [START_REF] Merris | Laplacian matrices of a graph: a survey[END_REF], another interesting property has been pointed out:

N i=1 λ i (L) = N i=1 d i .

Standard classes of graphs

1. A graph is called weighted if a weight is associated with every edge according to a proper map W : E → R, such that if (i, j) ∈ E, then w ij = 0, otherwise w ij = 0. In other words, a weighted graph is a graph whose edges have been labelled with real numbers. The length of a path in a weighted graph is the sum of the weights of the edges in the path. V,E) with diameter d(G) is said to be a distance-regular graph if for any vertices i and j of V and any integers v, w = 0, 1, . . . , d(G) (where d(G) is diameter of graph G), the number of vertices at distance v from i and distance w from j depends only on v, w, and the graph distance between i and j, independently of the choice of i and j. A graph G is said to be strongly regular, SRG(N,κ,a,c), if it is neither complete nor empty and there are integers κ, a, and c such that:

• G is regular with valency κ;

• any two adjacent vertices have exactly a common neighbors;

• any two distinct non-adjacent vertices have exactly c common neighbors.

These parameters are linked as follows:

(Nκ -1)c = κ(κa -1).

(2.4.1)

Examples of strongly regular graphs are:

• Petersen graph (SRG(10, 3, 0, 1)),

• Clebsch graph (SRG(16, 5, 0, 2)),

• Payley graphs (SRG(q, q-1 2 , q-5 4 , q-1 4 ), with q congruent 1 (mod 4)), • n-dimensional Hamming graphs (SRG(n 2 , 2n -2, n -2, 2), n ≥ 2),

• Cycle graphs with N < 6.

The most obvious property of the strongly regular graphs is that the diameter is equal to two (d(G) = 2).

Introduction

Consider an undirected connected graph G(V, E) with x i (k) being the state of node i at time-step k and define the initial state of the given network as x(0) = [x 1 (0), x 2 (0), . . . , x N (0)] T where N = |V |. A novel method is proposed in this chapter to design finite-time average consensus protocols for MASs or WSNs in a distributed way. In the discrete-time framework, a linear iteration scheme expresses its benefits for a distributed algorithm since each node repeatedly updates its value as a weighted linear combination of its own value and those of its neighbors.

x i (k + 1) = w ii (k)x i (k) + j∈N i w ij (k)x j (k), (3.1.1)
where

w ij =    0, if (i, j) / ∈ E w ij , if (i, j) ∈ E
or equivalently in matrix form:

x(k + 1) = W(k)x(k), (3.1.2)
W(k) ∈ S G where S G is the set of matrices that can be factorized as

W(k) = Q(k) • (I N + A)
where Q(k) stands for an arbitrary square matrix, I N as the N × N identity matrix, and • denotes the Hadamard matrix product that corresponds to an entry-wise matrix product.

As presented in the Chapter 1.2, average consensus is reached asymptotically if

lim k→∞ x(k) = 1 N 11 T x(0), (3.1.3)
meaning that

lim k→∞ W k = 1 N 11 T = J N . (3.1.4)
However, the foremost thing to note is that in order to run an average consensus algorithm, two main steps are required: the configuration step (also called design step) and the execution step. During the configuration step, a task can be achieved through a selfconfiguration algorithm instead of resorting to a network manager. Self-configuration can include graph discovering and distributed decision on some parameters. For instance, if the protocol is the maximum degree weights one ( 1.2.6), each agent first computes the number of its neighbors before running a max-consensus algorithm for computing the maximum degree of the underlying graph. In the case of the Metropolis-Hasting based protocol ( 1.2.7), each agent compares its degree with that of its neighbors in order to compute the weights of the average consensus protocol. One commonly used protocol is the constant edge weights, or graph Laplacian based average consensus protocol ( 1.2.8), where a common stepsize is used by all agents. Asymptotic convergence is hence guaranteed if the stepsize is strictly positive and lower than 2 λ N , where λ N stands for the largest Laplacian eigenvalue. Even through there are some simple bounds that give choices for the stepsize without requiring exact knowledge of the Laplacian spectrum, agents have to agree on an adequate stepsize.

To the best of our knowledge, there is no published work dealing with self-configuration protocols for the constant edge weights based average consensus protocol.

In real systems, the execution time is getting more and more impact. Therefore, the purpose is now to design a finite-time average consensus algorithm allowing all nodes to reach the average consensus value in a finite number of steps D for self-configuration protocols, i.e.

x(D) = 1 N 11 T x(0) = J N x(0). (3.1.5)

Literature review

In the literature, the proposed finite-time consensus protocols can be classified by two main approaches, namely: the minimal polynomial concept based approach and matrix factorization based approach.

Minimal polynomial concept based approach

In [START_REF] Sundaram | Finite-time distributed consensus in graphs with time-invariant topologies[END_REF], a finite-time consensus algorithm for linear time-invariant topologies based on the minimal polynomial of the weight consensus matrix W was proposed. Specifically, the minimal polynomial of the matrix W is the unique monic polynomial q(W) of the smallest degree such that q(W) = 0.

Suppose that the minimal polynomial of W has degree D + 1 ≤ N. For some constants α 0 , α 1 , . . . , α D , which are coefficients of the minimal polynomial q(W), since q(W) = 0, we get:

W D+1 + α D W D + . . . + α 1 W + α 0 I N = 0.
By multiplying the two sides of this equation with x(k), we obtain:

W D+1 x(k) + α D W D x(k) + . . . + α 1 Wx(k) + α 0 x(k) = 0.
or equivalently:

x(k + D + 1) + α D x(k + D) + . . . + α 1 x(k + 1) + α 0 x(k) = 0.
Therefore, ∀k ≥ 0 and 1 ≤ i ≤ N, x i (k) satisfies a linear difference equation of the form:

x i (k + D + 1) + α D x i (k + D) + . . . + α 1 x i (k + 1) + α 0 x i (k) = 0. (3.2.1)
It was shown that:

xi = [x i (D) x i (D -1) . . . x i (0)]s [1 1 . . . 1]s , where s =            1 1 + α D 1 + α D-1 . . . 1 + D j=1 α j           
The nodes need a sufficient observation time which is the degree of the minimal polynomial to reconstruct the initial states. However, the design of this protocol requires a strong knowledge on the underlying network topology. In fact, each node needs to store N(N -1) samples with N being the network size. Therefore, a decentralized calculation of the minimal polynomial was proposed although it is not entirely local, in the sense that a local calculation is repeated over N independent iterations (where N is the total number of nodes of the network) and at each iteration, it requires each node to store its own values for N + 1 time steps in order to obtain the minimal polynomial.

An improved method which requires much less data storage has been proposed in [START_REF] Yuan | Decentralized minimimtime consensus[END_REF].

The authors have proposed algorithms to compute the consensus value using the minimal number of observations of an arbitrary chosen node in a network using 2(D i +1) successive discrete-time steps, where D i + 1 is the degree of the monic minimal polynomial.

From equation (3.2.1), the measurable value at node i is determined by the minimal of polynomial of the corresponding matrix observability pair [W, e T i ], where e T i = [0, . . . , 0, 1 i , 0, . . . , 0] ∈ IR 1×N .

Definition 1 ([78])

The minimal polynomial associated with the matrix pair [W, e T i ] denoted by q i (t) t D i +1 + D i k=0 α i,k t k is the unique monic polynomial of minimal degree D i + 1 that satisfies e T i q i (W) = 0.

Hence, the minimal integer value D i +1 necessary for (3.2.1) to hold for almost any initial condition x(0) is given by the degree of the minimal polynomial of the observability pair [W, e T i ] [START_REF] Yuan | Decentralised final value theorem for discrete-time LTI systems with application to minimal-time distributed consensus[END_REF]. Consider the vector of difference between successive discrete-time values at node i, Xi (0, 1, 2 . . . , 2k) = {x i (1)x i (0), x i (2)x i (1), . . . , x i (2k + 1)x i (2k)}, k ∈ Z and its associated Hankel matrix:

Γ{ Xi (0, 1, 2 . . . , 2k)}      x i (1) -x i (0) x i (2) -x i (1) • • • x i (k + 1) -x i (k) x i (2) -x i (1) x i (3) -x i (2) • • • x i (k + 2) -x i (k + 1) . . . . . . . . . . . . x i (k + 1) -x i (k) x i (k + 2) -x i (k + 1) • • • x i (2k + 1) -x i (2k)      . (3.2.2)
The minimal number of steps is computed by checking a rank condition of a Hankel matrix constructed from local output observations in a distributed way.

Increase the dimension k of the square Hankel matrix Γ{ Xi (0, 1, 2 . . . , 2k)} until it loses rank and store the first defective Hankel matrix. In [START_REF] Yuan | Decentralized minimimtime consensus[END_REF], the authors have shown that by using a Vandermonde factorization of the Hankel matrix, one can see that Γ{ Xi (0, 1, 2 . . . ,

2k)} is defective if k ≥ D i + 1.
Then, the normalized kernel β β β = [β 0 . . . β D i +1 1] T of the first defective Hankel matrix is estimated.

Based on linear iteration (3.2.1), the following regression equation denoted by the output observations of node i in minimal time is defined as:

x i (k + D i + 1) + α i,D i x i (k + D i ) + . . . + α i,1 x i (k + 1) + α i,0 x i (k) = 0. (3.2.3)
Consider the Z-transform of x i (k):

X(z) = D i +1 j=1 α j j-1 l=0 x i (l)z j-l q i (z) H(z) q i (z) .
In order to estimate the final consensus value at node i, we define the following polynomial:

p i (z) q i (z) z -1 D i j=0 β j z j .
Then, the final consensus value is defined:

φ = lim z→1 (z -1)X(z) = H(1) p i (1) = y T i (D i )β β β 1 T β β β , where y T i (D i ) = [x i (0) x i (1) . . . x i (D i )] and β β β = [β 0 . . . β D i -1 1]
T is the vector of coefficients of the polynomial p i (z).

In [START_REF] Yuan | Decentralized minimimtime consensus[END_REF], the authors show that for an arbitrary initial condition, except for a set of initial conditions with Lebesgue measure zero [START_REF] Blondel | Continuous-time averagepreserving opinion dynamics with opinion-dependent communications[END_REF], the consensus value can be obtained from local observations in a minimal number of steps that does not depend on the total size of the graph. However, the computation of the rank of a given matrix and that of its kernel are needed. Therefore, the computational cost is still the weakness of these methods.

Matrix factorization based approach

This approach is based on the factorization of the averaging matrix J N . More precisely, we desire to find a finite sequence of matrices {W(1),

W(2), . . . , W(D)}, W(k) ∈ S G such that x(D) = 1 k=D W(k)x(0) = J N x(0) for all x(0) ∈ IR N , (3.2.4) with 1 k=D W(k) = J N . (3.2.5)
Finding this set of matrices is equivalent to solve a multi-variable polynomial system of equations, which can be solved using Grobner basis theory. However, the computational complexity is also burdens. In general, studying the existence of solutions to such a system of equations can be intractable.

Existence

The existence issue has been deeply considered in [START_REF] Georgopoulos | Definitive consensus for distributed data inference[END_REF] and [START_REF] Hendrickx | Graph diameter, eigenvalues, and minimum-time consensus[END_REF]. It has been pointed out that no solution exists if the factor matrices W(k) are all equal except if the graph is complete. Since the diameter of the graph characterizes the time necessary for a given information to reach all the nodes in the network, the number of factor matrices cannot be lower than the diameter d(G).

In [START_REF] Hendrickx | Graph diameter, eigenvalues, and minimum-time consensus[END_REF], an upper-bound of the number of factor matrices is given by 2r(G), with r(G) being the radius of the given graph (trees and graphs with minimum diameter spanning tree). Furthermore, for more general graphs, according to the results in [START_REF] Kibangou | Graph Laplacian based matrix design for finite-time distributed average consensus[END_REF], the number D is upper bounded by N -1. Hence, the number of steps D is bounded by:

d(G) ≤ D ≤ N -1. (3.2.6)
The solution of the factorization problem is not unique and no distributed solution exists.

In [START_REF] Hendrickx | Finite-time consensus using stochastic matrices with positive diagonals[END_REF], the authors have investigated the possibility of reaching (average) consensus in finite time by making use of a linear iterations scheme with a finite sequence of T consensus matrices, which are restricted to be stochastic with positive diagonals and consistent with a connected undirected graph structure G(V, E) where N = |V |. In this paper, the upper bound on the number of matrices was also given by N (N -1)

2

. On the other hand, for directed graphs, the authors have proven that finite-time consensus is reached by a sequence of stochastic matrices with positive diagonals only if the graph is strongly connected and contains a simple directed cycle with even length. However, if the graph G is a simple directed cycle, then there is no sequence of stochastic matrices with positive diagonals achieving consensus.

Solutions for particular graphs

For distance regular graphs, a closed form expression of the matrix factorization exists [START_REF] Kibangou | Step-size sequence design for finite-time average consensus in secure wireless sensor networks[END_REF]. The factorization of the averaging matrix has been proposed based on parameters of the intersection array characterizing the underlying graph constrained as a distance regular graph.

Definition 2 ([12]

)

A graph G(V, E) with diameter d(G) is said to be distance regular if there exist integers a k , b k , c k , k = 0, 1, . . . , d(G) such that for any two vertices i, j ∈ V and distance k = dist(i, j), there are exactly a k neighbors of j in V k (i), b k neighbors of j in V k+1 (i), and c k neighbors of j in V k-1 (i), where V k (i) is the set of vertices j of G(V, E) with dist(i, j) = k.
The parameters of a distance regular graph with valency K are linked as

a i + b i + c i = K, with a 0 = c 0 = b d(G) = 0 and c 1 = 1. Therefore, it is usual to characterize a distance regular graph by its intersection array {b 0 , b 1 , . . . , b d(G) ; c 1 , c 2 , . . . , c d(G) } where b 0 ≥ b 1 ≥ . . . ≥ b d(G) and 0 ≤ c 1 ≤ c 2 ≤ . . . ≤ c d(g) . Let us define G k (V, E k ) with (i, j) ∈ E k if and only if dist(i, j) = k and A k be its adjacency matrix. d(G) k=0 A k = J N .
One property can be noted as follows:

AA k = b k-1 A k-1 + a k A k + c k+1 A k+1 ,
where A = A 1 stands for the graph adjacency matrix.

Example : The Hamming graph H(d(H), n) with d(H) being the diameter of the graph, and n d(H) being the number of vertices. Therefore, its intersection array is given by {4, 3, 2, 1; 1, 2, 3, 4}.

Additionally, in [START_REF] Kibangou | Step-size sequence design for finite-time average consensus in secure wireless sensor networks[END_REF], the author has presented the factorization of the averaging matrix for distance regular graphs. Herein, we review two main theorems that help to determine the sequence of the consensus matrices. Theorem 3.2.1 ( [START_REF] Hoffman | The polynomial of a directed graph[END_REF])

Consider a graph G with N vertices and adjacency matrix A

1. There exist a polynomial H(ν) such that H(A) = J N if and only if G is a connected regular graphs.

2. For a connected regular graph with valency K, the unique polynomial of least degree satisfying

H(A) = J N is H(ν) = N g(K) g(ν) where (ν -K)g(ν) is the minimal polynomial of A.
The polynomial H(A) is called the Hoffman polynomial of G. We can note that the degree of this polynomial is unknown, except for distance regular graph. Lemma 1 ( [START_REF] Kibangou | Step-size sequence design for finite-time average consensus in secure wireless sensor networks[END_REF])

Consider a distance regular graph G(V, E) with |V | = N, valency K, diameter d(G) and intersection array {b 0 , b 1 , . . . , b d(G) ; c 1 , c 2 , . . . , c d(G) }. The averaging matrix 1 N J N is a d(G)
-th order polynomial of the graph adjacency matrix A:

1 N J N = P G (A) = d(G) q=0 γ q A q , where γ q = d(G) k=0 β k,q
N , and the coefficients β k,q depend on the parameters of the intersection array as follows:

β k+1,0 = - b k-1 β k-1,0 + a k β k,0 c k+1 , β k+1,q = β k,q-1 -b k-1 β k-1,q + a k β k,q c k+1 , 1 ≤ q -1, β k+1,k = β k,k-1 -a k β k,k c k+1 , β k+1,k+1 = β k,k c k+1 , (3.2.7
)

with a k = K -b k -c k , β 0,0 = 1 and β k,q = 0 if k < q.
Here, P G (A) stands for the polynomial of A.

Since P G (A) = 1 N H(A), the polynomial P G (ν) = d(G) q=0
γ q ν q will be called the averaging Hoffman polynomial of

G (P G (K) = 1 N H(K) = 1)
. The roots of P G (ν) are also those of the minimal polynomial of A. Theorem 3.2.

([39])

Let P G be the averaging Hoffman polynomial of the distance regular graph G(V, E) of valency K, diameter d(G), adjacency matrix A, and |V | = N. The averaging matrix 1 N J N can be factorized as:

1 N J N = d(G) t=1 ((1 -α t K)I N + α t A), α t = 1 K -ǫ t , (3.2.8)
where ǫ t ∈ IR \ {K} stands for the roots of P G .

Let take the Hamming graph H(4, 2) depicted in Figure 3.1 with valency K = 4, diameter d(H) = 4 and N = 2 4 = 16. From ( 3.2.7), we obtain the coefficients as follows:

β 0,0 = 1; β 1,0 = 0; β 1,1 = 1; β 2,0 = -2; β 2,1 = 0; β 2,2 = 1 2 ; β 3,0 = 0 β 3,1 = - 5 3 ; β 3,2 = 0; β 3,3 = 1 6 ; β 4,0 = 1; β 4,1 = 0; β 4,2 = - 2 3 ; β 4,3 = 0; β 4,4 = 1 24 .
As the result, we compute

γ q = d(H) k=0 β k,q
N as follows:

γ 0 = 0; γ 1 = - 2 3 × N ; γ 2 = - 1 6 × N ; γ 3 = 1 6 × N ; γ 4 = 1 24 × N .
Hence, the averaging polynomial is given by P

H (ν) = 1 N (-2 3 ν -1 6 ν 2 + 1 6 ν 3 + 1 24 ν 4
). Solving this polynomial, we get the roots ǫ = {-4, -2, 0, 2}. Furthermore, from (3.2.8), the corresponding stepsizes are obtained as α = { 1 8 , 1 6 , 1 4 , 1 2 }, which are exactly the nonzero inverses of the Laplacian eigenvalues of the graph. Therefore, the sequence of the consensus matrices are:

W(1) = 1 2 I N - 1 8 A; W(2) = 1 3 I N - 1 6 A; W(3) = A; W(4) = -I N + 1 2 A.

General solutions

In [START_REF] Ko | On matrix factorization and scheduling for finite-time average consensus[END_REF], the resulting solution yields a link scheduling on the complete graph to achieve finite-time consensus. Such scheduling is to be controlled by a central node. They provide necessary and sufficient conditions for finite-time consensus and compute the minimum consensus time on the Boolean hypercube. In order to prove the existence of a finite factorization of J N , the mechanism is that starting with a spanning tree T of graph G and arbitrarily designate a node as the root. The basic idea is that all nodes pass all their goods to the chosen root. Then, the root propagates the appropriate amount of goods back into the network so that every node has an equal amount at the end. The proposed method includes two actions. The first action is described that from vertices farthest from the root of a spanning tree T , the proposed algorithm traverses upwards to give all goods to their parents (they can be children of other parents), then remove themself. This action terminates when only a single vertex (the root) remains. At this point, the remain node (the root) contains the sum of all initial values. Then, the second action is to propagate back down the tree while re-distributing the values to achieve the average consensus value at termination.

For more general graphs, solutions based on graph Laplacian have been recently introduced in closed-form, [START_REF] Kibangou | Finite-time average consensus based protocol for distributed estimation over awgn channels[END_REF][START_REF] Kibangou | Graph Laplacian based matrix design for finite-time distributed average consensus[END_REF].

Lemma 2 ([41])

The matrices

W(k) = α k I N + β Â, k = 1, .
. . , D are jointly diagonalizable, consistent with the network topology, and admit 1 as eigenvector.

Where,  is defined as  = d max I N -L whose properties are:

(a) Â is symmetric;

(b) Its eigenvalues λ k ( Â) = d max -λ k (L)
, where λ k (L), d max are the eigenvalues of the Laplacian matrix L and the maximum degree respectively. In particular, d max is a simple eigenvalue.

(c) The eigenvectors of the Laplacian matrix are also the eigenvectors of Â. Meaning that, 1 is the eigenvector associated with the eigenvalue d max . Theorem 3.

([41])

Given a connected undirected graph associated with the Laplacian matrix L, the set of matrices W(k) = (α k + d max β)I N -βL, k = 1, . . . , D, β = 0, allows reaching the average consensus in D steps if:

(a) D + 1 is the number of distinct eigenvalues of the Laplacian matrix;

(b) the parameter β and α k are defined as:

β = 1 D+1 i=2 D λ i (L) α k = λ k+1 (L) -d max D+1 i=2 D λ i (L) , k = 1, . . . , D.
with λ i (L), i = 2, . . . , D + 1 being the nonzero distinct Laplacian eigenvalues.

More precisely, in [START_REF] Kibangou | Graph Laplacian based matrix design for finite-time distributed average consensus[END_REF], the solution was given by α k = 1 and

β k = -1 λ k+1 (L) , λ k (L) being a nonzero Laplacian eigenvalue. Remark 5
• The number of factor matrices is the number of distinct nonzero Laplacian eigenvalues.

• The proposed solutions make use of the Laplacian spectrum. For this purpose, we propose distributed methods for Laplacian spectrum estimation in Chapter 4.

Comparison between the minimal polynomial approach and the matrix factorization approach

To summarise the difference between the two approaches, we focus on two terms that are design time and execution time.

• Design time:

(a) For the minimum polynomial approach, each node chooses weights to generate the weight consensus matrix W that satisfies the convergence conditions 1.2.2.

(b) The matrix factorization approach requires the storage and computation resources for global information (network topology or Laplacian spectrum) to design the weight consensus matrix W.

• Execution time:

(a) The minimum polynomial approach mainly needs data storage and computation resources: 2(D r + 1) successive steps, [START_REF] Yuan | Decentralized minimimtime consensus[END_REF].

(b) At this stage, the factorization of averaging matrix does not require storage of past data. Hence, the execution time is as simple as a standard consensus protocol.

The purpose of this thesis is to study a new algorithm for self-configuration for finitetime average consensus where the weight matrices are not necessarily based on the graph Laplacian. The aim is to design protocols that allow achieving average consensus in the fastest possible time, possibly as fast as diameter of the underlying graph. More precisely, we solve a matrix factorization problem in a distributed way.

Distributed solution to the matrix factorization problem.

For a connected graph G(V, E), let x(0) = [x 1 (0), x 2 (0), . . . , x N (0)] T be the state vector of the N-agent network associated with the graph G. For the sake of simplicity, N agents (nodes) of the network (graph G) are arranged in the first layer that corresponds to the initial state x(0) as described in the Figure 3.2 illustrating the linear iteration scheme (3.1.2) in space and time.

As presented in Chapter 1.2, the linear iteration scheme is an useful tool that allows the system to reach the finite-time average consensus in D steps. Therefore, according to each iteration, the updated state vectors x(k), k = 1, . . . , D -1 are assigned to be consecutive layers respectively. Finally, the last iteration is dedicated to the output values x(D).

The weight matrices W(k), k = 1, . . . , D have the entries w ij (k) are correspondent to the links between the nodes of two consecutive layers at time-step k. For instance, W(1) is the weight matrix assigned to the link between the first layer and the second one. x(0)

x 1 (0)

x 2 (0)

x N (0) 1 2 N x(1) 1 2 N x(D -1) 1 2 N x(D)
x 1 (D)

x 2 (D)

x N (D) W 1 W D w 1 11 w 1 N N w 1 2N w 1 1N w D 1N w D 2N w D N N w D

Remark 6

The weight matrix W(k) ∈ S G is associated with the graph G. Therefore, the entries w ij (k) are designed following the rule defined as:

w ij (k) = w ij (k), if j ∈ N i

0, otherwise

To simplify the notation, we denote W(k) as W k .

The linear iteration scheme in space and time can be viewed as a multilayer neural network. The selection of the weights can be analyzed through the scope of training. By using a set of P learning sequences, the mechanism is divided into two main steps, namely, forward step and backward step. The main idea is to compute the partial derivatives of the error between the actual output and desired output, then propagate them back to each layer of the network in order to update the weights. The aim is to reduce the error criteria E efficiently.

Remark 7

Indeed, in most systems where communication is involved, learning sequences are used for communication channel identification or for synchronization. These sequences are used during the mandatory configuration step before transmitting informative data, i.e. running average consensus in our case. We assume that all the agents know the average of the learning sequence. For instance, after running a standard consensus protocol (Algorithm 1) with the weight consensus matrix as mentioned in Subsection 1.2.4 or the minimal time consensus algorithm proposed in [START_REF] Yuan | Decentralized minimimtime consensus[END_REF]. x i,p (0). Our aim is to estimate the factor matrices W k , k = 1, • • • , D, by minimizing the quadratic error

Let {x i,p (0), y i,p }, i = 1, • • • , N, p = 1, • • • , P ,
E(W 1 , • • • , W D ) = 1 2 N i=1 P p=1 (x i,p (D) -y p ) 2 , (3.3.1) with x i,p (k) = j∈N i ∪{i} w ij (k)x j,p (k -1), w ij (k) being the entries of matrices W k ∈ S G .
We can rewrite the cost function (3.3.1) as

E(W 1 , • • • , W D ) = 1 2 1 k=D W k X(0) -Y 2 F , s.t W k ∈ S G , (3.3.2)
where . F stands for the Frobenius norm, Y = J N X(0), Y and X(0) being N × P matrices with y i,p and x i,p as entries, respectively.

We assume that X(0)X(0) T = I N which means that the input vector is orthogonal. For instance, vectors of the canonical basis of IR N can be used as inputs. Hence, we can

note that E(W 1 , • • • , W D ) = 1 2 1 k=D W k -J N 2 F
, meaning that minimizing (3.3.2) is equivalent to solving the factorization problem (3.2.5):

{W * k }| k=1,...,D = arg min {W k } t=1,...,D 1 2
P p=1 tr ǫ p (W)ǫ T p (W) , (3.3.3)
with

ǫ p (W) = 1 k=D W k x p (0) -y p , (3.3.4)
where tr(.) denotes the trace operator and y p = J N x p (0).

Denoting E p (W) = 1 2 tr ǫ p (W)ǫ T p (W)
, the solution of this optimization problem can then be obtained iteratively by means of a gradient descent method:

W k := W k -α P p=1 ∂E p (W) ∂W k
or a stochastic gradient one:

W k := W k -α ∂E p (W) ∂W k ,
where the gradient of the cost function is approximated by the gradient at a single inputoutput sequence. In what follows, we resort to a stochastic gradient method. For this purpose, we first state the following technical lemma: Lemma 3

The derivatives of the cost function E p (W) = 1 2 tr ǫ p (W)ǫ T p (W) with ǫ p (W) defined in (3.3.4) can be computed as follows: x i,p (0); and

∂E p (W) ∂W D = δ δ δ D,p x T p (D -1) (3.3.5) ∂E p (W) ∂W k = δ δ δ k,p x T p (k -1), k = 1, . . . , D -1, (3.3
δ δ δ k-1,p = W T k δ δ δ k,p , k = 1, . . . , D. (3.3.

7)

Proof : The consensus network being a linear system we know that:

x p (k) = W k x p (k -1),
Therefore, we can explicitly write the output according to the weighting matrix of interest, i.e.

x p (D) = W D x p (D -1) = k+1 j=D W j W k x p (k -1), k = 1, • • • , D -1.
Equivalently, by defining

Z k+1 = k+1 j=D W j , we get x p (D) = Z k+1 W k x p (k -1)
. The cost function can be written as

E p (W) = 1 2 tr((W D x p (D -1) -xp )(W D x p (D -1) -xp ) T ).
We can easily deduce that ∂Ep(W)

∂W D = δ δ δ D,p x T p (D -1) with δ δ δ D,p = x p (D) -xp .

Now, we can express the cost function according to any matrix W

t , t = 1, • • • , D -1 as E p (W) = 1 2 tr((Z k+1 W k x p (k -1) -xp )(Z k+1 W k x p (k -1) -xp ) T ).
Expanding the above expression and taking into account the linearity of the trace operator yield

E p (W) = 1 2 [tr(Z k+1 W k x p (k -1)x p (k -1) T W T k Z T k+1 ) -tr(Z k+1 W k x p (k -1)x T p ) -tr(x p x p (k -1) T W T k Z T k+1 ) -tr(x p xT p )].
Now, computing the derivative, we get:

∂E p (W) ∂W k = 1 2 [2 × Z T k+1 Z k+1 W k x p (k -1)x p (k -1) T -2 × Z T k+1 xp x p (k -1) T ] = Z T k+1 (Z k+1 W t x p (k -1) -x)x p (k -1) T = Z T k+1 (x p (D) -xp ) δ δ δ D,p x p (k -1) T = Z T k+1 δ δ δ D,p x p (k -1) T = W T k+1 W T k+2 . . . W T D δ δ δ D δ δ δ D-1,p x p (k -1) T = W T k+1 δ δ δ k+1,p x p (k -1) T = δ δ δ k,p x(k -1) T .
Applying the results of Lemma 3, the updating scheme of the optimization algorithm is as follows:

W k [t + 1] = W k [t] -α ∂E p(t) (W) ∂W k = W k [t] -αδ δ δ k,p(t) x T p(t) (k -1), (3.3.8)
where p(t) ∈ {1, 2, • • • , P }, and t stands for the t th iteration of the optimization process. Then, we project the updated weights on S G that gives:

W k [t + 1] = W k [t] -α ∂E p(t) (W) ∂W k = W k [t] -α δ δ δ k,p(t) x T p(t) (k -1) • (I N + A) , (3.3.9)
The gradient descent algorithm (3.3.9) acts by alternating the following steps: the learning sequence is first propagated forward. Then the error between the targeted output and x(D) is computed and then propagated . Proposition 1

Let (W * 1 , . . . , W * D ) be a local minimum, of E(W 1 , . . . , W D ), then 1 k=D W * k = J N 60 . Proof : Let W * 1 , . . . , W * D be the local minimum of the problem (3.3.2), then the derivatives ∂E p (W * ) ∂W * k | k=1,...,D = 0.
Following Lemma 3,

δ δ δ k,p x T p (k -1) = 0, k = 1, . . . , D.
Thus,

δ δ δ D,p x T p (D -1) = 0 ⇔ δ δ δ D,p = 0 ⇔ x p (D) = xp , p = 1, . . . , P.
It can be rewritten in matrix form:

1 k=D W * k X(0) = J N X(0) ⇔ 1 k=D W * k -J N X(0) = 0.
Since X(0)X T (0) = I N , then

1 k=D W * k -J N = 0 ⇔ 1 k=D W * k = J N .
This mechanism is similar to the gradient back-propagation algorithm. Its convergence has been well studied in the literature, see [START_REF] Chong | An Introduction to Optimization[END_REF] and [START_REF] Mangasarian | Serial and parallel backpropagation convergence via nonmonotone perturbed minimization[END_REF][START_REF] Wu | Convergence of gradient method with momentum for back-propagation neural networks[END_REF] for instance. Convergence towards a local minimum is guaranteed if the step-size α is appropriately chosen inside the interval (0, 1). If the step-size α is set too large, the error E tends to oscillate, never reaching a minimum. In contrast, too small α prevents the optimization from making reasonable progress. Therefore, choosing an appropriate α is also an encountered difficulty in practice. For a particular problem, choosing a suitable stepsize α involves experimenting with different values to see how the convergence reacts. Such a step-size also influences on speed of convergence of the consensus algorithm. Several rules have been proposed for accelerating the convergence speed such as using the adaptive learning rate α [START_REF] Vrahatis | Improving the convergence of the backpropagation algorithm using learning rate adaptation methods[END_REF] or adaptive momentum [START_REF] Shao | Convergence analysis of a back-propagation algorithm with adaptive momentum[END_REF] . However, these techniques are not consistent with the distributed framework of the proposed algorithm. One trick to speed-up the convergence is to add a regularization term in the cost function to be minimized.

{W * k }| k=1,...,D = arg min {W k } k=1,...,D P p=1 E p(t) (W) + 1 2 D k=1 β W k [t] -W k [t -1] 2 .
By minimizing such a cost function, the update equation is given by:

W k [t + 1] = W k [t] -α δ δ δ k,p(t) x T p(t) (k -1) • (I N + A) + β(W k [t] -W k [t -1]).
Entry-wise, each agent updates its entries as follows:

w k ij [t + 1] := w k ij [t] -α[t](δ i,k x j (k -1)â ij ) + β(w k ij [t] -w (k-1) ij [t -1]),
where:

• δ i,k is the i-th entry of δ δ δ k and x j (k -1) the j-th entry of x(k -1). Here, âij -the entry of matrix (I N + A) is defined as:

âij =    a ij if i = j. 1 if i = j.
• α and β are the learning rate and the momentum rate respectively. They can be chosen as constants or as time varying parameters. In our study, the learning rate α[t] varies at each optimization iteration step t while the momentum rate is fixed.

Remark 8

In order to avoid the confusion between the notation of the time-step k of the finitetime average consensus protocol and the optimization iteration t, we denote w ij (k) as w k ij .

For a clear understanding, Figure 3.3 shows the mechanism of the proposed algorithm.

In detail, for each learning sequence, the algorithm first executes the forward step, where each node updates repeatedly its value. At the end of this step, each layer has stored the updated values, then the partial derivatives of the error δ D between x(D) and the desired average consensus values x are estimated. The backward step starts by propagating back the estimated partial derivatives δ D to each layer until the input layer receives δ 1 . Now, all weights can be updated. The regime continues with the next learning sequence. The mechanism is repeated until the error criterion E satisfies a given threshold θ, which is a significantly small number. The distributed algorithm is then described as follows: Algorithm 3

1 2 N x p (0) 1 2 N x p (1) 1 2 N x p (D -1) 1 2 N x p (D) W 1 W D w 1 11 w 1 N N w 1 2N w 1 1N w D 1N w D 2N w D N N w D
1. Initialization:

• Number of steps D, number of patterns P

• Learning sequence {x i,p (0), y p }, where i = 1, • • • , N, p = 1, • • • , P , with y p = 1 N N i=1
x i,p (0).

• Random initial weight matrices W k [0], k = 1, . . . , D, and W k [-1] = 0 • Learning rate: 0 < α[0] < 1;
• Momentum terms: 0 < β < 1;

• Select a threshold θ

• Set t = 0 2. Set p = 0;
(a) Set p := p + 1, (b) Select the corresponding input-output sequence x i (0) = x i,p (0), x = y p .

(c) Learning sequence propagation:

x i (k) = j∈N i ∪{i} w k ij [t]x j (k -1), k = 1, • • • , D.
(d) Error computation:

δ i,D = x i (D) -x; e i,p = δ 2 i,D .
(e) Error propagation:

δ i,k-1 = j∈N i ∪{i} w k ji [t]δ j,k , k = D, • • • , 2. (f) Matrices updating: for k = 1, • • • , D, i = 1, • • • , N, and j ∈ N i ∪ {i}, w k ij [t + 1] = w k ij [t] -α[t]δ i,k x j (k -1)â ij + β(w k ij [t] -w k ij [t -1]).
(g) Increment t.

(h) If p = P , compute the mean square error

E i = 1 N P p=1
e i,p , else return to 2a.

(i) If E i < θ stop the learning process, else return to 2.

Numerical results

Example 1

Let us consider the graph depicted in Figure 3.4. Its diameter is d(G) = 3, hence the bounds of the number of iterations D are d(G) < D < N -1 = 5. The Laplacian matrix is defined as follows:

L =             2 -1 0 -1 0 0 -1 2 -1 0 0 0 0 -1 3 -1 0 -1 -1 0 -1 3 -1 0 0 0 0 -1 2 -1 0 0 -1 0 -1 2            
The eigenvalues of the Laplacian matrix L are: 0, 1, 2, 3 and 5. According to [START_REF] Kibangou | Graph Laplacian based matrix design for finite-time distributed average consensus[END_REF], the factorization of the average consensus matrix by means of graph Laplacian-based consensus matrices needs D factor matrices, D being the number of nonzero distinct eigenvalues of the Laplacian matrix, this number being, in general, greater than or equal to the diameter of the graph. Therefore, we get the following factorization:

(I N -L)(I N - 1 2 L)(I N - 1 3 L)(I N - 1 5 L) = J 6 ,
meaning that consensus is achieved in 4 steps.

J 6 =             0.
           
Since the optimal number of steps is not a priori known, matrices proposed in [START_REF] Kibangou | Finite-time average consensus based protocol for distributed estimation over awgn channels[END_REF][START_REF] Kibangou | Graph Laplacian based matrix design for finite-time distributed average consensus[END_REF]. The weight matrices are:

W 1 =          -0.
         W 2 =          0.
         W 3 =          0.
        
We can easily check that:

W 3 W 2 W 1 = J 6 .
Now, let us consider the finite-time consensus protocol in its execution step, i.e. after the configuration step. For arbitrary initial values, the trajectories of the state of the network are depicted in Figure 3.6(a). Exact average consensus is achieved in 3 steps. When using the optimal constant edge weights protocol W = I N -αL where α = are actually needed to reach consensus. The Laplacian matrix is defined as:

W α=0.3333 =             0.
            0 2 4 6 8
L =                      3 -1 0 -1 0 0 0 0 -1 0 -1 4 -1 -1 -1 0 0 0 0 0 0 -1 3 0 -1 0 0 0 0 -1 -1 -1 0 4 0 -1 -1 0 0 0 0 -1 -1 0 4 0 -1 -1 0 0 0 0 0 -1 0 3 -1 0 0 -1 0 0 0 0 -1 0 -1 3 -1 0 -1 0 0 0 0 0 0 -1 3 -1 0 0 -1 0 0 -1 0 0 -1 3                     
The nonzero corresponding eigenvalues are {2, 5, 6}. Hence, if the Laplacian spectrum is known, then we can find the consensus matrices based on the Laplacian spectrum. That are:

W k = I N -1 λ k L, k = 1, 2, 3. Therefore, (I N - 1 2 L)(I N - 1 5 L)(I N - 1 6 L) = J 10 .
Running the proposed algorithm, the performance in Figure 3.8 tells us that this network cannot achieve average consensus in 2 steps, which is the diameter of a graph as suggested in [START_REF] Georgopoulos | Definitive consensus for distributed data inference[END_REF]. Such a result has been demonstrated analytically by [START_REF] Hendrickx | Graph diameter, eigenvalues, and minimum-time consensus[END_REF]. However, by running the algorithm as well as varying the value of D in the interval [d(G), 2r(G)], we can pick the optimal D in this case. The final MSE was obtained by considering a sufficient large number of iterations. For this topology, D = 3 gives us the best solution for finite-time average consensus, see Figure 3.8. Also, the trajectory of an arbitrary state vector is depicted in Figure 3.9 after getting the sequence of weight matrices. The sequence of consensus matrices are:

W 1 =                   0.
                  W 2 =                   -0.
                  W 3 =                   0.
                 
As the result, the product of consensus matrices gives out the averaging matrix.

W 3 W 2 W 1 = J 10 = 1 10 11 T .

Conclusion and discussion

In this chapter, we have proposed a way for distributively designing the finite-time average consensus protocol. By using a learning sequence, we have shown how to solve the matrix factorization problem ( 3.2.5) in a fully distributed way. The method is based on the gradient back-propagation method. The factorization gives rise to factor matrices that are not necessarily symmetric or stochastic. Given the diameter of the graph d(G), we can find out the optimal value of the number of steps necessary for reaching average consensus by varying the value of D in the interval (d(G) ≤ D ≤ N -1). However, it is obviously inconvenient for large-scale networks.

On the other hand, we can see that, the obtained consensus matrices are not doubly stochastic matrices.

Finally, the speed of convergence of the proposed algorithm is still an open issue. For this purpose, future works encompass optimal step-size for the gradient descent based method and other optimization methods with a particular focus on large size graphs.

In addition, the dependency of the convergence speed on the learning sequence is to be studied. Hence, another interesting issue will be the design of optimal learning sequences.

Beside that, the robustness of the finite-time average consensus should be definitely paid a great attention.

Introduction

The efficiency of a network is evaluated through its functionality and its robustness.

According to the latter property, some questions may arise. For instance, if there is an accidental event, how does the network react? Will the network survive? In order to answer these questions, the study of the robustness of networks has increasingly attracted the attention of the scientific community.

In the literature, there is quite a lot of robustness measures (see Section 1.3), [START_REF] Ellens | Graph measures and network robustness. e-version[END_REF][START_REF] Ellens | Effective graph resistance[END_REF]. For instance, the earliest measure is based on the connectivity of the graph representing the network [START_REF] Dekker | Network robustness and graph topology[END_REF]. In the other words, a network is robust if there exists the presence of alternate paths, which ensure the possibility of communication in spite of damage [START_REF] Dekker | Network robustness and graph topology[END_REF]. There are two classical concepts of connectivity for a graph which can be used to model network robustness that are the vertex and edge connectivities. The vertex (edge) connectivity K v (K e ) is the minimal number of vertices (edges) to be removed in order to disconnect the given graph (see Subsection 1.3.1). By this analytical approach, in order to deduce the robustness index, we have to remove each edge to check for the connectivity of the graph. Even if this approach is simple, it is also a burden if we analyse reasonable larger complex networks.

In [START_REF] Wu | Spectral measure of structural robustness in complex networks[END_REF], the authors have proposed the concept of natural connectivity as a spectral measure of robustness in complex networks. The natural connectivity is expressed in mathematical form as the average eigenvalue of the adjacency matrix of the graph representing the network topology. The advantage of this method is that the proposed measure works in both connected and disconnected networks. On the other hand, the Laplacian spectrum sp(L) can be used to compute the robustness indices that evaluate the performance of the given network. For instance, the second smallest eigenvalue of the Laplacian spectrum λ 2 (L), which is also known as algebraic connectivity [START_REF] Friedler | Algebraic connectivity of graphs[END_REF], is well-known as a critical parameter that influences on the performance and robustness of dynamical systems due to its main role in the connectivity of the graph [START_REF] Ellens | Graph measures and network robustness. e-version[END_REF]. Hence, there is a lot of researches on this algebraic connectivity in the literature [3,[START_REF] Merris | Laplacian matrices of a graph: a survey[END_REF][START_REF] Yang | Decentralized estimation and control of graph connectivity in mobile sensor networks[END_REF]. Furthermore, there are other remarkable measures for getting the desired robustness indices: the effective graph resistance R and the number of spanning tree ξ (see Section 1.3), which are based on the whole non-zero Laplacian spectrum.

We propose distributed methods to estimate both the effective graph resistance R and the number of spanning trees ξ. Figure 4.1 shows that the distributed estimation of the network robustness is divided into two main stages including Laplacian spectrum sp(L) estimation, effective graph resistance and number of spanning trees (R, ξ) computation. Firstly, we propose consensus-based methods to estimate the whole spectrum of the Laplacian matrix in a distributed way. Primarily, we compute the average consensus value x by using the standard consensus protocols (Algorithm 1) or the decentralized minimum-time consensus value computation algorithm [START_REF] Yuan | Decentralized minimimtime consensus[END_REF]. Next, the estimation of the Laplacian spectrum can be divided in two steps that are the estimation of the distinct Laplacian eigenvalues, and then evaluation of the corresponding multiplicities. Finally, from the estimated sp(L) the robustness indices can be computed. Our focus is therefore on the estimation of the Laplacian spectrum sp(L).

Standard consensus protocols

x i (k), k = 0, . . . , h = N -1.

Laplacian spectrum estimation:

Distinct Laplacian eigenvalues estimation

Laplacian eigenvalues multiplicities estimation

Computation of network robustness indices

x q i (k), k = 0, . . . , h = N -1 λ 2 , λ 3 , . . . , λ D+1 sp(L) R, ξ

Literature review

The estimation of the Laplacian spectrum sp(L) can be executed by means of centralized algorithms where a global knowledge about the network topology is available. However, in order to circumvent the inherent vulnerability of centralized schemes, various studies have been carried out on decentralized algorithms for estimating the Laplacian eigenvalues. Most of them are restricted to the estimation of the second smallest Laplacian eigenvalue λ 2 (L). For instance, in [START_REF] Yang | Decentralized estimation and control of graph connectivity for mobile sensor networks[END_REF], the second smallest Laplacian eigenvalue was estimated by resorting to a continuous-time decentralized power iteration method that enables for each node i to estimate x i , where x = [x 1 . . . x N ] T being the estimate of the eigenvector v 2 , corresponding to λ 2 (L). In [3], the power iteration method consists of transforming the Laplacian matrix L associated with the undirected and connected graph G into a new matrix C, which is the deflated version of the Perron matrix of L, where the algebraic connectivity λ 2 (L) becomes the leading eigenvalue, i.e. the spectral radius ρ(C) is an expression that depends on λ 2 (L). First, nodes compute the powers of a deflated Laplacian matrix C, C k in a distributed way. Then, at each time step k for each node i, the sum c i (k) of the absolute values of the entries of C k is defined to execute the max-consensus on c i (k), whose maximum value converges to the spectral radius ρ(C) when k tends to infinity. Once ρ(C) is computed, the estimated algebraic connectivity is finally obtained by reversing the transformation.

For computing the whole Laplacian spectrum, several approaches have been considered in the recent literature: Fast Fourier Transform (FFT) based methods [START_REF] Franceschelli | Decentralized estimation of Laplacian eigenvalues in multi-agent systems[END_REF][START_REF] Sahai | Hearing the cluster of a graph: A distributed algorithm[END_REF], local eigenvalue decomposition of given observability based matrices [START_REF] Kibangou | Decentralized Laplacian eigenvalues estimation and collaborative network topology identification[END_REF].

For FFT-based methods, the main idea is to make the state of each agent oscillate only at frequencies corresponding to the eigenvalues of the network topology. The problem is then mapped into a signal processing one that can be efficiently and independently solved by each agent by means of FFT. In [START_REF] Franceschelli | Decentralized estimation of Laplacian eigenvalues in multi-agent systems[END_REF], the amplitude and phase of the oscillations have been characterized as function of the eigenvectors of the Laplacian and the initial conditions. In this paper, a connected undirected graph G(V, E) with N = |V | has been considered. From an initial condition x i (0), z i (0) generated from a uniform distribution in [-1, 1], each node i simulates the following local interaction rule with its neighbors

N i (t):    ẋi (t) = z i (t) + j∈N i (z i (t) -z j (t)), żi (t) = -x i (t) -j∈N i (x i (t) -x j (t)),
which can be expressed in vector form as:

ẋ(t) ż(t) = M x(t) z(t) , where M = 0 I N -L I N + L 0 is skew symmetric and 0 is the null N × N matrix.
To any eigenvalue λ(L) it corresponds a couple of complex and conjugate eigenvalues λ(M), λ(M), that is:

λ(M) = j(1 + λ(L)), λ(M) = -j(1 + λ(L)),
while the corresponding eigenvectors v M of M are function of the eigenvectors v L of L:

v M = [v T L jv T L ] T , vM = [v T L -jv T L ] T .
It follows that each state of each node x i (t) follows an oscillating trajectory which is a linear combination of sinusoids, whose amplitudes and phases shifted are function of the initial conditions and of the graph topology, oscillating at and only at frequencies corresponding to the eigenvalues of the network at time t. Also, the state trajectory of x i (t) and z i (t) can be viewed as a function of Laplacian eigenvectors and the initial conditions:

x i (t) = N j=1 [v Lj (i)(cos((1 + λ j (L))t)v T Lj x(0) + sin((1 + λ j (L))t)v T Lj z(0))], z i (t) = N j=1 [v Lj (i)(-sin((1 + λ j (L))t)v T Lj x(0) + cos((1 + λ j (L))t)v T Lj z(0))].
In detail, the module of the Fourier transform of the i-th state components x i (t) and z i (t), i = 1, . . . , N can be written as:

|F [x i (t)]| = |X i (f )| = m j=1 a j,i 2 δ f ± 1 + λ j (L) 2π , |F [z i (t)]| = |Z i (f )| = m j=1 b j,i 2 δ f ± 1 + λ j (L) 2π ,
where f is the frequency domain variable; m is the number of distinct Laplacian eigenvalues; and a j,i , b j,i are coefficients given by:

• For λ 1 (L) = 0 :

a 1,i = v 1 (i)v T 1 x(0) = 1 T N x(0) N , b 1,i = v 1 (i)v T 1 z(0) = 1 T N z(0) N ,
where v 1 is the unitary norm eigenvector corresponding to λ 1 (L) = 0.

• For λ j (L) > 0 :

a j,i = b j,i = ν j k=1 v (k) j (i)v (k) j T x(0) 2 + ν j k=1 v (k) j (i)v (k) j T z(0) 2 ,
where v (k) j , k = 1, . . . , ν j is the unitary norm eigenvector associated to λ j (L) > 0, with ν j be the algebraic multiplicity.

In a time window of length T , node i estimates the frequencies of the sinusoids of which signal x i (t) is composed. The frequencies of the spectrum of x i (t) and z i (t) are given by

f j = ± 1+λ j (L) 2π
. Hence, λ j (L) = 2πf j ± 1. However, this approach inherits the limitations of the FFT algorithm. For instance, since △ f = fs M where M is the number of samples to be recorded, the resolution of the estimated eigenvalues is strongly dependent on that of the FFT method and the accuracy depends on the amount of stored data (M). On the other hand, the node can only estimate the eigenvalues associated with the mode that is observable. We know that the system is not observable from a single node if there exists an eigenvalue with multiplicity greater than 2. In addition, some nodes may observe only a subset of the eigenvalues. Then, they need to execute additional coordination rules for propagating their data.

In contrast, in [START_REF] Kibangou | Decentralized Laplacian eigenvalues estimation and collaborative network topology identification[END_REF], an algebraic method using the observability properties of the network was proposed. With this method, the eigenvalues of the network matrix can be recovered by solving a local eigenvalue decomposition problem on an appropriately constructed matrix of observed data. In this paper, a connected undirected graph G(V, E) with N = |V | has been also considered. Suppose that each node i has some initial values organized in an M-length row vector x i (0), M ≥ N. The dynamics of the network is represented as:

X(k + 1) = WX(k), Y i (k) = E i X(k),
where:

• X(k) ∈ IR N ×M is a matrix with x i (k), i = 1, . . . , N as rows;

• Y i (k) are the outputs or node values that are seen by node i at the k th time step;

• E i ∈ IR di ×N is the row selection matrix with di bounded by [0, d i + 1]; • W = I N -ǫL is the consensus matrix. Its eigenvalue decomposition is W = UDU T ,
where the eigenvectors and eigenvalues are respectively organized in the orthogonal matrix U and the diagonal matrix D.

Since L = U∆ ∆ ∆U T with ∆ ∆ ∆ being diagonal matrix of λ i (L), then D = I N -ǫ∆ ∆ ∆.
The purpose of [START_REF] Kibangou | Decentralized Laplacian eigenvalues estimation and collaborative network topology identification[END_REF] is to estimate the eigenvalues λ i (L), i = 1, . . . , N from the observations Y i (k), k = 0, 1, . . . , K i + 1 where K i stands for the observability index of the pair (W, E i ). From the available observations Y i (k), k = 0, 1, . . . , K i + 1, two matrices are generated:

Ȳi =        Y i (0) Y i (1) . . . Y i (K i )        ; Ȳi =        Y i (1) Y i (2) . . . Y i (K i + 1)       
.

Define the block matrix

Ȳi Ȳi = ŨΣ Σ Σ ṼT in a singular value decomposition form. From here, the matrix Ũ of left singular vector of the matrix ( ȲT

i ȲT i ) T is computed.
Then, we partition Ũ as Ũ = ( ŨT

1 ŨT 2 )
T , the two blocks having the same number of rows. [START_REF] Kibangou | Decentralized Laplacian eigenvalues estimation and collaborative network topology identification[END_REF] points out that the consensus matrix and the matrix R built from the observations have exactly the same spectrum. Therefore, D results from the eigenvalues decomposition of the matrix ( ŨT 1 Ũ2 )( ŨT 1 Ũ1 ) -1 . Finally, the Laplacian eigenvalues can be obtained by the equation ∆ ∆ ∆ = 1 ǫ (I N -D).

Computing the matrices

R 1 = ŨT 1 Ũ1 , R 2 = ŨT 1 Ũ2 and R = R 2 R -1 1 ,
However, this method is only applicable to networks having nodes with sufficient storage and computation capabilities. Actually, data to be stored equals to di M(K i + 2).

In this thesis, the approach is based on the fact that the averaging matrix can be viewed as a product of Laplacian-based consensus matrices with step-sizes given by the inverses of the Laplacian eigenvalues [START_REF] Kibangou | Graph Laplacian based matrix design for finite-time distributed average consensus[END_REF]. Therefore, by computing a factorization of the averaging matrix in a distributed way, a distributed estimation of the Laplacian eigenvalues is achieved.

Remark 9

Hereafter, we denote the Laplacian eigenvalue λ i (L) as λ i .

G(V, E). Its state x(k) = [x 1 (k), x 2 (k), . . . , x N (k)]
T evolves as:

x(k + 1) = (I N -αL)x(k) = (I N -αL) k x(0), (4.3.1)
where α is selected such as 0 < α < 1 dmax to ensure asymptotic convergence [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF]. From the state x(k), one can define q(k) = L k x(0) as:

q(k) = (q 1,k , . . . , q N,k ) T = 1 (-α) k x(k) - k-1 i=0 k i (-α) i q(i) , (4.3.2)
with q(0) = x(0).

Assumption 1:

The number N of nodes is known.

2. Assumption 2: Each node i stores its initial value x i (0) and the consecutive values q i,k , k = 0, . . . , h.

The problem under study can therefore be formulated as follows:

Given an arbitrary input-output pair {x(0), x = x1}, where x = 11 T N x(0), and the intermediate observations q(k), k = 0, . . . , h, estimate in a distributed way the Laplacian eigenvalues of the graph representing the network.

However, instead of computing directly the Laplacian eigenvalues, we will compute the inverse of the eigenvalues. In addition, since the number D of distinct nonzero Laplacian eigenvalues is unknown, we will re-parametrize the factorization. Let us first state the following proposition: Proposition 2

Consider a connected graph with Laplacian matrix L and distinct nonzero eigenvalues λ 2 , • • • , λ D+1 . Assume that the local values vector x(0) is not orthogonal to any eigenvector of L, the cost function

E(α α α) = x(h) -x 2 = 1 k=h (I N -α k L)x(0) -x 2 , (4.3.3
)

with α α α = (α 1 • • • α h ) T , h ≥ D, is minimal if and only if {1/λ 2 , • • • , 1/λ D+1 } ⊆ {α 1 , α 2 , • • • , α h }.
Proof : Consider the eigenvalues decomposition of the Laplacian matrix: L = U∆ ∆ ∆U T , where ∆ ∆ ∆ = diag(λ 1 , λ 2 , . . . , λ N ), U being the orthogonal matrix of the eigenvectors. Without loss of generality, we assume that the (D + 1) first entries of the diagonal matrix ∆ ∆ ∆ are the (D + 1) distinct graph Laplacian eigenvalues. The graph being connected, we know that λ 1 = 0. Therefore, the corresponding eigenvector (first column of U) is equal to 1 √ N 1. Let e i be the ith vector of the canonical basis of IR N . Then the averaging matrix J N is equal to Udiag(e 1 )U T . We can therefore rewrite the cost function E(α α α) as:

E(α α α) = h k=1 I N -α k U∆ ∆ ∆U T x(0) -J N x(0) 2 = U h k=1 (I N -α k ∆ ∆ ∆) -diag(e 1 ) U T x(0) 2 = Udiag (ǫ ǫ ǫ) U T x(0) 2 , with ǫ ǫ ǫ = (0, h k=1 (1 -α k λ 2 ), • • • , h k=1 (1 -α k λ N )) T . Now, let denote by w i , i = 1, • • • , N , the
entries of U T x(0). Then:

E(α α α) = N i=2 h k=1 (1 -α k λ i ) 2 w 2 i .
We can note that this function is a sum of square terms. Therefore, it vanishes if and only if each From Proposition 2, we can conclude that the set S 1 = {α k , k = 1, . . . , h}, of solutions of (4.3.3) contains the inverse of Laplacian eigenvalues S 2 = { 1 λ i , i = 2, . . . , D+1}: S 2 ⊂ S 1 . Therefore, the estimation procedure is to be divided in two parts. We first get S 1 by solving (4.3.3) in a distributed way, then we retrieve the set S 2 from S 1 and deduce the Laplacian eigenvalues locally.

term h k=1 (1-α k λ i ) 2 w 2 i , i = 2, • • • , N is equal to zero. Since, x (
We know that if {α k } h k=1 contains D inverses of the nonzero Laplacian eigenvalues, then

x = h k=1 (I -α k L)x(0), or equivalently, x = h k=0 c k L k x(0) = h k=0 c k q(k) = Qc, (4.3.4)
where

Q = [q(0) q(1) • • • q(h)] ∈ IR N ×(h+1) and c = [c 0 c 1 • • • c h ] T ∈ IR (h+1)×1
, the coefficient c i being defined as:

c k =            1, if k = 0. (-1) k i<j<•••<k α i α j • • • α k , if k = 1, . . . , h -1. (-1) h h k=1 α k , if k = h. (4.3.5)
It is obvious that α k are the roots of the polynomial with coefficients c k , and c has an alternating sign pattern.

We can state the following proposition which is a corollary of Proposition 2: Proposition 3

Let P(c, :) be the h th degree polynomial with coefficients c = [c 0 c 1 . . . c h ] T ∈ IR (h+1)×1 built from the set of stepsizes S 1 = {α 1 , α 2 , . . . , α h } by using (4.3.5). For any node i and intermediate measurements q i = [q i,0 , q i,1 , . . . , q i,h ] T ∈ IR (h+1)×1 of the consensus protocol (4.3.1), we get

x = P(c, q i ) = q T i c (4.3.6)
if and only if the D inverses of the distinct Laplacian eigenvalues are contained in S 1 .

Given an upper-bound h for the number D of distinct Laplacian eigenvalues, the idea is to minimize (4.3.3) and then reduce step by step the obtained set of coefficients c k built from S 1 until it is minimal. After all, the set of the inverses of distinct Laplacian eigenvalues S 2 is achieved. In order to obtain the whole Laplacian spectrum sp(L), the corresponding multiplicities m i , i = 1, . . . , D + 1 of the distinct Laplacian eigenvalues Λ = {λ 1 , . . . , λ D+1 } should be estimated (it will be considered in Section 4.3.2).

For a simple understanding of the structure of this Chapter, Figure 4.2 illustrates its organization, showing the way to estimate the network robustness.

Standard consensus protocols

x i (k), k = 1, . . . , N -1
Distinct Laplacian eigenvalues estimation:

Direct approach Indirect approach x q i (k), k = 0, . . . , N -1

α 1 , α 2 , . . . , α N -1 c 0 , c 1 , . . . , c N -1 λ 1 , λ 2 , . . . , λ N -1
Laplacian eigenvalues corresponding multiplicities estimation Network robustness indices sp(L) In fact, for distributed Laplacian eigenvalues estimation, three cases have been considered:

1. The case of perfect knowledge of the average consensus value x, 2. The case of imperfect knowledge of the average consensus value x, 3. The case of unknown average consensus value x. In this way, the proposed direct method with respect to the stepsizes α α α consists in solving a non-convex optimization (4.3.3) through the Lagrange method of multipliers.

Method 1: Direct distributed factorization of the averaging matrix.

Equation (4.3.3) shows that α k are global parameters. For distributively carrying out the factorization of the averaging matrix as factors of Laplacian-based consensus matrices, we can first note that the cost function (4.3.3) can be rewritten as:

E(α α α) = x(h) -x 2 = 1 k=h (I N -diag(α α α k )L)x(0) -x 2 , (4.3.7) with α α α k = [α 1,k , α 2,k , . . . , α N,k ] T = α k 1, k = 1, 2, . . . , h = N -1.
The idea of the proposed method is to minimize the disagreement between neighbors on the value of α k while ensuring that the factorization of the averaging matrix is achieved. Such a factorization is assessed by constraining the values of the nodes after h iterations of the consensus algorithm to be equal to the average of the initial values: Any solution of (4.3.9) contains the inverse of Laplacian eigenvalues according to Proposition 2.

min α α α k ∈IR N×1 ,k=1,2,...,h 1 2 h k=1 i∈V j∈N i (α j,k -α i,k ) 2 , (4.3.8)
The constrained optimization problem (4.3.9) can then be performed as an unconstrained optimization problem by means of a Lagrange method with a Lagrange function defined as follows:

H(α α α 1 , α α α 2 , . . . , α α α h , y) = 1 2 h k=1 α α α T k Lα α α k + y T (x(h) -x) (4.3.10)
where y ∈ IR N ×1 is the vector of Lagrange multipliers.

The minimization of the Lagrange function (4.3.10) can be viewed as h consensus problems to be solved simultaneously with a constraint that provides a kind of external reference to be tracked.

To avoid misunderstanding between time-step k and optimization iteration t, we denote x(k), W(k) as x k , W k respectively.

In order to derive the corresponding distributed algorithm, we now state the following technical lemma: Lemma 4

The derivatives of the cost function H(α α α k , y) defined in (4.3.9) can be computed as follows:

∂H(α α α k , y) ∂α α α k = Lα α α k -diag -1 (α α α k )diag(x k-1 -x k )δ δ δ k , (4.3.11)
where

δ δ δ h = y and δ δ δ k-1 = W k δ δ δ k , k = 1, . . . , h. (4.3.

12)

Proof : The consensus network being a linear system we know that x k = W k x k-1 . Therefore, we can explicitly write the output according to the weighting matrix of interest, i.e.

x h = W h x h-1 and x h = k+1 j=h W j W k x k-1 , k = 1, • • • , h -1.
The cost function can be written as:

H(α α α 1 , α α α 2 , . . . , α α α h , y) = 1 2 D k=1 α α α T k Lα α α k + y T ( k+1 i=h W i W k x k-1 -x) = 1 2 h k=1 α α α T k Lα α α k + y T ( k+1 i=h W i (I N -diag(α α α k )Lx k-1 -x) = 1 2 h k=1 α α α T k Lα α α k + y T ( k+1 i=h W i x k-1 -x) -y T k+1 i=h W i diag(α α α k )Lx k-1 ).
Now, using the property (A.0.1) of the Khatri-Rao product, we get:

H(α α α 1 , α α α 2 , . . . , α α α h , y) = 1 2 h k=1 α α α T k Lα α α k + y T ( k+1 i=h W i x k-1 -(x T k-1 L T ⊙ k+1 i=h W i )α α α k -x).
Therefore, we can easily deduce that

∂H(α α α k , y) ∂α α α k = Lα α α k -(x T k-1 L T ⊙ k+1 i=h W i ) T y.
Taking the symmetry of the consensus matrices into account, we get:

∂H(α α α k , y) ∂α α α k = Lα α α k -diag(x T k-1 L T ) T h i=k+1 W i y = Lα α α k -diag(Lx k-1 ) h-2 i=k+1 W i W h-1 W h y δ δ δ h-1 = Lα α α k -diag(Lx k-1 ) D-2 i=k+1 W i W h-1 δ δ δ h-1 δ δ δ h-2 = Lα α α k -diag(Lx k-1 )δ δ δ k . Since x k = (I N -diag(α α α k )L)x k-1 , we can deduce that Lx k-1 = diag -1 (α α α k )(x k-1 -x k ). Hence, diag(Lx k-1 ) = diag -1 (α α α k )diag(x k-1 -x k ). As a consequence, ∂H(α α α k , y) ∂α α α k = Lα α α k -diag -1 (α α α k )diag(x k-1 -x k )δ δ δ k .
Applying the results of Lemma 4, the updating scheme of the optimization algorithm is as follows:

α α α k [t + 1] = α α α k [t] -β ∂H(α α α k , y) ∂α α α k [t] = (I N -βL)α α α k [t] + βdiag -1 (α α α k )diag(x k-1 [t] -x k [t])δ δ δ k ) (4.3.13) y[t + 1] = y[t] + µ(x h [t] -x). (4.3.14)
The proposed distributed algorithm can then be described as follows:

Algorithm 4 (Lagrange multipliers method for distributed averaging matrix factorization)

Inputs:

• An upper-bound h of the number of distinct Laplacian eigenvalues h = N -1,

• learning rates 0 < µ < 1 and 0 < β < 1,

• initial input-output values {x i (0), xi }, i = 1, . . . , N with x = 1 N N i=1
x i (0) obtained from standard consensus protocols.

1. Initialization: Each node i, i = 1, . . . , N generates random steps-sizes α i,k [0], k = 1, . . . , h and Lagrange multipliers y i [0] and then set t = 0.

Propagation of Lagrange multipliers:

(a) Set:

δ i,h [t] = y i [t] (b) Propagate the Lagrange multipliers for k = h, • • • , 2: δ i,k-1 [t] = δ i,k [t] + j∈N i (δ j,k [t] -δ i,k [t]).
3. Update of the factor matrices step-sizes: (c) After receiving the messages m j,k [t] from its neighbors j ∈ N i , each node i carries out the following updates:

(a) x i,0 [t] = x i (0). (b) For k = 1, 2, • • • , h,
x i,k = x i,k-1 + α i,k j∈N i (x j,k-1 -x i,k-1 ) α i,k [t + 1] = α i,k [t] -β( j∈N i (α j,k [t] -α k,i [t])) -β j∈N i (x j,k-1 [t] -x i,k-1 [t])δ i,k [t]
4. Update the Lagrange multiplier using (4.3.13)

5. Exit if a stopping criterion is reached else set t := t + 1 and return to 2.

Outputs: Set S 1,i of step-size α i,k [t].
Each node deduces the Laplacian eigenvalues as

λ i,k = 1 α i,k [t+1] .
This algorithm is similar to a gradient back-propagation algorithm. According to [START_REF] Chong | An Introduction to Optimization[END_REF], selecting β, µ in the interval of (0, 1) ensures convergence to a local minimum as long as the constraint x(h) = x and the constraint on α k are achieved. The obtained solution contains the inverse of the Laplacian eigenvalues. The speed of convergence depends on the choice of these parameters.

In order to accelerate the speed of convergence, we can re-parameterize the problem (4.3.1) by means of the equation (4.3.4).

x =

h k=0 c k L k x(0) = h k=0 c k q(k) = Qc.
Instead of the minimization with respect to the step-sizes α k , the new optimization with respect to the polynomial coefficients c k , k = 0, . . . , h can be applied. That is called indirect distributed factorization of the averaging matrix.

Method 2: Indirect distributed factorization of the averaging matrix.

The aim of this method is to find the polynomial coefficients c k in a distributed way by solving the problem: min

c x -Qc 2 , with Q = [q(0) q(1) • • • q(h)] ∈ IR N ×(h+1) .
All the initial conditions x(0) should not be eigenvectors of Laplacian matrix L and not be identical, meaning that x(0) = δ1, δ ∈ IR, otherwise all the columns of Q are zero except the first one. Since the initial condition is randomly generated, almost surely it cannot be an eigenvector of Laplacian matrix L. Lemma 5

Consider a connected undirected graph G(V, E) with Laplacian spectrum sp(L) = {λ 1 , λ 2 , . . . , λ N }. Herein, sp(L) has D + 1 distinct Laplacian eigenvalues. If x(0) is not an eigenvector of L , then the matrix

Q = Udiag(U T x(0))      1 λ 1 λ 2 1 • • • λ h 1 1 λ 2 λ 2 2 • • • λ h 2 . . . . . . . 1 λ N λ 2 N • • • λ h N     
has rank min(D + 1, N), where U is the matrix of Laplacian eigenvectors.

Proof : Q = [q(0) q(1) . . . q(h)] ∈ IR N ×(h+1) , where q(k) = L k x(0).

Equivalently,

Q = [I N x(0) Lx(0) L 2 x(0) . . . L h x(0)] = [I N L L 2 . . . L h ]diag(x(0), . . . , x (0) h+1 times

]

Knowing that L = U∆ ∆ ∆U T , where ∆ ∆ ∆ = diag(λ 1 , λ 2 , . . . , λ N ), U are the diagonal matrix of Laplacian eigenvalues and the orthogonal matrix of the eigenvectors respectively. Hence,

[I N L L 2 . . . L h ] = U[I N ∆ ∆ ∆ ∆ ∆ ∆ 2 . . . ∆ ∆ ∆ h ]diag(U T , . . . , U T h+1 times ) Therefore, Q = U[I N ∆ ∆ ∆ ∆ ∆ ∆ 2 . . . ∆ ∆ ∆ h ]diag(U T x(0), . . . , U T x(0) h+1 times
)

where diag(U T x(0), . . . , U T x(0)) ∈ IR ((h+1)×N )×(h+1) .

Denote by Λ Λ Λ = [λ 1 , . . . , λ N ] T the vector of Laplacian eigenvalues. Therefore,

Q = U        1 N • U T x(0) Λ Λ Λ • U T x(0) Λ Λ Λ2 • U T x(0) . . . . . . Λ Λ Λh • U T x(0)        T = U        diag(U T x(0))1 N ) diag(U T x(0)) Λ Λ Λ diag(U T x(0)) Λ Λ Λ2 . . . . . . . . diag(U T x(0)) Λ Λ Λh        T = Udiag(U T x(0))      1 λ 1 λ 2 1 • • • λ h 1 1 λ 2 λ 2 2 • • • λ h 2 . . . . . . . . 1 λ N λ 2 N • • • λ h N      = UFB
where • denotes the Hamadard product.

Besides that, h = N -1, hence the Vandermonde matrix B is a square matrix. Therefore, Rank(Q) = min(rank(F), rank(B)).

Here, if U T x(0) = 0, or equivalently x(0) is not an eigenvector of the Laplacian matrix, then rank(F) = rank(diag(U T x(0))) = N . On the other hand, since B is a Vandermonde matrix, one property can be reviewed is that rank of the Vandermonde is equal to the number of distinct generators (Laplacian eigenvalues), meaning that rank(B) = D + 1. Thus, rank(Q) = D + 1.

From the other side, if there exist elements of the vector U T x(0) equal to 0, then rank(Q) = min(D + 1, l) with l being the number of nonzero elements of U T x(0).

Remark 11

rank(Q) < N if there exists λ i ∈ sp(L) with multiplicity higher than 1, ∀x(0).

In a centralized way, in order to find c from the equation Qc = x, Q must be a full rank matrix. Then,

c = Q -1 x
In constrast, if Q is not a full rank matrix, then there are several solutions for c. To deal with this situation, we make use of a least-square solution to find the most common solution S 1 inside the set of solutions (So 1 , So 2 , . . . , So m ) of c as described in The problem now becomes:

min c c 2 subject to Qc = x
where Q is a full-row-rank matrix.

In order to solve this problem in a distributed way, we solve the problem locally according to c i , i = 1, . . . , N. Denoting by c i the local version of c and Q i ∈ IR (d i +1)×(h+1) the submatrix of Q constituted with the rows associated with node i and its d i neighbors, the solution of the previous problem can be obtained in a distributed way by solving min

c i N i=1 j∈N i c i -c j 2 s.t. Q i c i = x1, i = 1, • • • , N, (4.3.15)
provided that Q i is full row rank.

Let us call Se i a row selection matrix with respect to node i. Se is a matrix, where Se j,j = 1 if j ∈ (i ∪ N i ). Q i can be expressed as:

Q i = Se i Q = Se i Udiag(U T x(0))        1 λ 1 λ 2 1 • • • λ h 1 1 λ 2 λ 2 2 • • • λ h 2 . . . . . . . 1 λ N λ 2 N • • • λ h N        Obviously, rank(Se i U) = d i + 1.
As shown in (4.3.15), Q i should be a full row rank matrix to ensure the constraint Q i c i = x1 to be satisfied. Therefore, we can construct

Q i = q i q j
, j ∈ N i by selecting the q j , j ∈ N i until the rank(Q i ) cannot increase more.

There is a lot of rules that can be applied for the choice of q j . However, in our study, we select q j based on the condition number of the established Q i . Meaning that, q j is selected so that the condition number of Q i is the smallest.

From (4.3.5), knowing that the coefficients c k have an alternating sign pattern, the problem can be reformulated as follows:

min c i ∈IR (h+1) ,i=1,...,N 1 2 N i=1 j∈N i c i -c j 2 . (4.3.16) subject to c i ∈ C i = {c i ∈ IR (h+1) : Q i c i = x1}, i = 1, 2, • • • , N. (4.3.17) c i,k = sign(k)c i,k , (4.3.18)
where sign(k) =

   -1 if k is odd, 1 if k is even.
The problem here is to optimize the sum of convex objective functions corresponding to the nodes in the network. In contrast to the previous direct method, the optimization problem (4.3.16) is convex. To solve it we will resort to a projected subgradient method inspired from [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF] and an alternating direction method of multipliers [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF][START_REF] Erseghe | A distributed and scalable processing method based upon ADMM[END_REF].

(A) Distributed projected subgradient method

The main idea of the distributed projected subgradient method is to allow each node i to update its estimate following two steps.

Firstly, a subgradient step is taken by minimizing the local objective function

1 2 j∈N i c i -c j 2 .
Then, the intermediate result is projected on the constraint set C i . Therefore, the first step yields:

ĉi [t + 1] = c i [t] -γ[t] j∈N i (c i [t] -c j [t]) , (4.3.19)
where γ[t] > 0 is a stepsize, while the second step requires the computation of the projection

c i [t+1] = Ω C i [ĉ i [t+1]]
. This projection step is equivalent to solve a constrained optimization problem formulated as:

min 1 2 c i [t + 1] -ĉi [t + 1] 2 s.t. Q i c i [t + 1] = x1.
To solve this constrained optimization problem, generally, the Lagrangian augmented function is defined as follows:

L(c i [t + 1], y) = 1 2 c i [t + 1] -ĉi [t + 1] 2 + y T (Q i c i [t + 1] -x1).
Following the KKT (Karush-Kuhn-Tucker) conditions, we have to solve:

∇ c i [t+1] L(c i [t + 1], y) = 0. (4.3.20) ∇ y L(c i [t + 1], y) = 0. (4.3.21)
Let solve (4.3.20) first:

c i [t + 1] -ĉi [t + 1] + Q T i y = 0 ⇔ c i [t + 1] = ĉi [t + 1] -Q T i y. (4.3.22)
Then, we substitue c i into (4.3.21):

Q i c i [t + 1] -x1 = 0 ⇔ Q i (ĉ i [t + 1] -Q T i y) = x1 ⇔ Q i ĉi [t + 1] -Q i Q T i y = x1 ⇔ y = (Q i Q T i ) -1 (Q i ĉi [t + 1] -x1). ( 4
c i [t + 1] = ĉi [t + 1] -Q T i (Q i Q T i ) -1 (Q i ĉi [t + 1] -x1) ⇔ c i [t + 1] = ĉi [t + 1] -Q T i (Q i Q T i ) -1 Q i ĉi [t + 1] + Q T i (Q i Q T i ) -1 x1 ⇔ c i [t + 1] = (I h+1 -Q T i (Q i Q T i ) -1 Q i )ĉ i [t + 1] + Q T i (Q i Q T i ) -1 x1. (4.3.24)
or equivalently,

c i [t + 1] = Qi xi + (I h+1 -Qi Q i )ĉ i [t + 1], (4.3.25) with Qi = Q T i (Q i Q T i ) -1
, provided h ≥ d i and Q i full row rank. Then, we project c i [t+1] on the sign constraint as:

c i,k [t + 1] =    max(0, c i,k [t + 1]) if mod(k, 2) = 0. min(0, c i,k [t + 1]) if mod(k, 2) = 0. (4.3.26)
The corresponding distributed projected subgradient algorithm is then described in Algorithm 5:

Algorithm 5 (Projected subgradient algorithm)

Given h = N -1, where N stands for the number of nodes, and a learning rate 0 < γ 0 < 1, the following steps are to be carried out:

1. Initialization:

• Initial input-output values {x i (0), xi }, i = 1, • • • , N, with x = 1 N N i=1 x i (0)
obtained from a standard average consensus algorithm.

• Each node i forms the vector q i , and the matrices Q i and Qi .

• Each node i generates random initial values for c i [0].

Set t = 0;

(a) Set t := t + 1;

(b) Update the learning rate γ[t] by using the diminishing rule: (d) Return to 2a or stops the iterations if a stopping criterion is reached.

γ[t] = γ 0 √ t (c)
3. Each node builds a polynomial p i with the entries of c i [t] and compute the set of stepsizes {α i,k } as the roots of the polynomial p i .

In [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF], the convergence of this method has been well studied. With an appropriate choice of the stepsize γ[0], Algorithm 5 will converge to a global minimum. The speed of convergence depends on the choice of the gradient stepsize.

(B)Alternating direction method of multipliers (ADMM).

ADMM is relatively easy to implement and its convergence properties have been well studied in recent works such that in [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF][START_REF] Erseghe | Fast consensus by the alternating direction multipliers method[END_REF]. The pros of this method is that it is guaranteed to converge for any choice of the tunning parameters under mild conditions [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF] while other techniques such as subgradient projection or dual decomposition depend on the choice of the stepsize for the gradient updates. Unlike other methods, ADMM has only one single penalty parameter ρ that can influence on the speed of convergence. ADMM assures a very good convergence speed when its parameters are appropriately chosen. There exist some works dealing with penalty parameter selection for accelerating the convergence rate of ADMM [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Shames | Optimal parameter selection for the alternating direction method of multipliers ADMM: quadratic problems[END_REF][START_REF] Ghadimi | On the optimal step-size selection for the alternating direction method of multipliers[END_REF][START_REF] Raghunathan | Optimal step-size selection in alternating direction method of multipliers for convex quadratic programs and model predictive control[END_REF][START_REF] Teixeira | Optimal scaling of the ADMM algorithm for distributed quadratic programming[END_REF]. The convergence analysis is also studied in [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF][START_REF] Shames | Optimal parameter selection for the alternating direction method of multipliers ADMM: quadratic problems[END_REF][START_REF] Erseghe | Fast consensus by the alternating direction multipliers method[END_REF][START_REF] Ph | Explicit convergence rate of a distributed alternating direction method of multipliers[END_REF][START_REF] Raghunathan | Optimal step-size selection in alternating direction method of multipliers for convex quadratic programs and model predictive control[END_REF]. However, in this study, our penalty parameter is a fixed constant.

The ADMM algorithm acts by introducing some auxiliary variables z ij , so that a strictly equivalent optimization problem can be formulated as follows:

min c i ∈IR (h+1) ,i=1,...,N 1 2 N i=1 j∈N i c i -c j 2 . (4.3.27) subject to c i = z ij , j ∈ N i , z ji = z ij , z k ij = sign(k)z k ij , c i ∈ C i = {c i ∈ IR h+1 : Q i c i = x1}, i = 1, 2, • • • , N. where sign(k) =    1 if k is even, -1 if k is odd .
The corresponding augmented Lagrangian function is given by:

L ρ (c, z, y) = N i=1 L (i) ρ (c i , z, y), with
L (i) ρ (c i , z, y) = j∈N i 1 2 c i -c j 2 + y T ij (c i -z ij ) + ρ 2 c i -z ij 2 , (4.3.28)
where ρ is a penalty coefficient. In contrast to the method of Lagrange multipliers, in ADMM the minimization of the cost function according to the primal variables c i and z ij is carried out sequentially. So, given the primal variables and the dual variables y ij , the cost function is first minimized according to c i and the obtained solution projected on the constrained set (c-minimization). Then z ij are obtained from the given dual variables and the updated values of c i (z-minimization). The third step corresponds to the update of the dual variables from the updated primal ones (y-minimization):

c i [t + 1] = Ω C i [argmin L ρ (ĉ i [t + 1], z ij [t], y ij [t])]; (4.3.29) z ij [t + 1] = argmin L ρ (c i [t + 1], z ij , y ij [t])]; (4.3.30) y ij [t + 1] = argmin L ρ (c i [t + 1], z ij [t + 1], y ij [t])]. (4.3.31)
These iterations indicate that the method is particularly useful when the c-and zminimizations can be carried out efficiently. First, we solve the c-minimization:

∇ c i L ρ (ĉ i [t + 1], z ij [t], y ij [t]) = 0. ⇔ j∈N i (ĉ i [t + 1] -c j [t]) + j∈N i y ij [t] + ρ(ĉ i [t + 1] -z ij [t]) = 0 ⇔ (d i (1 + ρ))ĉ i [t + 1] = j∈N i c j [t] + ρ j∈N i z ij [t] - j∈N i y ij [t] ⇔ ĉi [t + 1] = (d i (1 + ρ)) -1 ( j∈N i c j [t] + ρ j∈N i z ij [t] - j∈N i y ij [t]),
where d i is the degree of node i.

For the cost function considered herein, we can note that the c-minimization problem (4.3.29) is solved by: ĉi

[t + 1] = (d i (1 + ρ)) -1 ( j∈N i c j [t] + ρ j∈N i z ij [t] - j∈N i y ij [t]), (4.3.32)
d i being the degree of node i. Then, as for the projected subgradient case we get the same expression for c i [t + 1] as in (4.3.25).

Next, we solve the sub-optimization (4.3.30) analogically, meaning that

∇ z ij L(c i [t + 1], z ij [t + 1], y i [t]
) = 0 with respect to the term that z ij = z ji . Rewrite the augmented Lagrange function:

L ρ (c i , z, y) = j∈N i 1 2 c i -c j 2 + y T ij (c i -z ij ) + ρ 2 c i -z ij 2 + y T ji (c j -z ji ) + ρ 2 c j -z ji 2 , Since ∇ z ij L(c i [t + 1], z ij [t + 1], y i [t]) = 0 with z ij = z ji , it is equivalent to: -(y ij [t] + y ji [t]) -ρ(c i [t + 1] + c j [t + 1]) + 2ρz ij [t + 1] = 0 ⇔ z ij [t + 1] = 1 2 (c i [t + 1] + c j [t + 1]) + 1 2ρ (y ij [t] + y ji [t]). (4.3.33)
Then, we project z ij [t + 1] on the sign constraint as:

z k ij [t + 1] =    max(0, z k ij [t + 1]) if mod(k, 2) = 0. min(0, z k ij [t + 1]) if mod(k, 2) = 0. (4.3.34)
The Lagrange multipliers are updated as:

y ij [t + 1] = y ij [t] + ρ(c i [t + 1] -z ij [t + 1]). (4.3.35)
This fully decentralized algorithm is summarized as follows:

Algorithm 6 (Distributed ADMM algorithm)

Given h = N -1, where N stands for the number of nodes N, and a penalty parameter ρ, the following steps are to be carried out:

1. Initialization:

• Initial input-output values {x i (0), xi }, i = 1, • • • , N, with x = 1 N N i=1 x i (0)
obtained from a standard average consensus algorithm,

• Each node i forms the vector q i , and the matrices Q i and Qi .

• Each node i generates random values for

c i [0], z ij [0] and y ij [0], j ∈ N i .
2. Set t = 0;

(a) Set t := t + 1;

(b) Update equations: each node i

• updates its primary variable c i [t] using (4.3.32) and (4.3.25).

• updates the primary variables z ij [t] associated with its neighbors using (4.3.33) and (4.3.34).

• update the dual variables y ij using (4.3.35).

• Return to 2a until a stopping criterion is reached.

(c) Each node builds a polynomial p i with the entries of c i [t] and compute the set of the stepsizes {α i,k } as the roots of the polynomial p i .

The convergence analysis for distributed quadratic programming has been studied deeply in [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Ph | Explicit convergence rate of a distributed alternating direction method of multipliers[END_REF].

Laplacian-based factorization of the averaging matrix with imperfect knowledge of x

For all methods proposed in previous case 4.3.1.2, the consensus value x is the foremost required parameter. We now relax this requirement by means of a novel method that deals with Laplacian eigenvalues estimation when an exact average consensus value is unknown.

In what follows, we assume that the nodes run a standard average consensus protocol ( 1.2.1), which is stopped after M iterations, and store their local values x i,M . In this study, we want to show the performance of the proposed algorithm with respect to some values of M.

Since x is unknown exactly, for a sufficiently large value of M, x(M) can be viewed as its rational approximation. Therefore, our aim is to find the coefficients c k that minimize the quadratic error on the final consensus value:

E(c) = h k=0 c k q(k) -x(M) 2 . (4.3.36)
Then, we can deduce the set S 1 = {α α α k } k=1,••• ,h and that of the inverse of the Laplacian eigenvalues S 2 = { 1 λ i , i = 2, . . . , D}, which belongs to S 1 : S 2 ⊂ S 1 . In order to solve the problem in a distributed way, it can be reformulated as follows:

min c c c k ∈IR N×1 ,k=0,...,h 1 2 N i=1 ( h k=0 c i,k q i,k -xi,M ) 2 , (4.3.37)
subject to c i,k = c j,k , i = 1, . . . , N; j ∈ N i ; k = 0, . . . , h.

c i = (c i,0 , c i,1 , • • • , c i,h ) T ∈ C i ,
where c i,k is the k th local polynomial coefficient and C i stands for the constrained set with respect to node i. It is defined as

C i = {c i ∈ IR h+1 , i = 1, . . . , N : q T l c i = xi,M with q l ∈ IR h+1 , l ∈ N i ∪ i} or equivalently as C i = {c i ∈ IR h+1 , i = 1, . . . , N : Q i c i = xi,M 1, Q i ∈ IR d i +1×h+1
} where xi,M being the local version of x at iteration M and Q i = [q(0) q(1) . . . q(h)], being full row rank by construction. This constraint set allows to enforce not only the local coefficients vector c i to be equal in a given neighborhood but also to enforce the scalar product q T i c i to be globally equal. In the literature, there exist several methods dealing with distributed consensus optimization (4.3.37), including distributed descent algorithms [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF], dual averaging methods [START_REF] Duchi | Dual averaging for distributed optimization: Convergence analysis and network scaling[END_REF], and the ADMM [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF][START_REF] Erseghe | Fast consensus by the alternating direction multipliers method[END_REF]. However, we still make use of ADMM, which has become a popular approach for solving convex minimization problems by means of parallelization.

One can note that by introducing the auxiliary variables {z ij }, the following constrained convex optimization problem is absolutely equivalent to problem (4.3.37):

min c c c k ∈IR N×1 ,k=0,...,h 1 2 N i=1 ( h r=0 c i,k q i,k -xi,M ) 2 . (4.3.38) subject to c i,k = z k ij , i = 1, . . . , N; j ∈ N i z k ji = z k ij , k = 0, . . . , h. (4.3.39) z k ij = sign(k)z k ij , c i ∈ C i , where sign(k) =    1 if k is even, -1 if k is odd
. Since the graph is assumed to be connected, the constraints (4.3.39) force c i = c j , j ∈ N i . We can then write the associated augmented Lagrangian:

L ρ (c, z, y) = N i=1 1 2 h k=0 (c i,k q i,k -xi,M ) 2 + h k=0 j∈N i y k ij (c i,k -z ij,k ) + j∈N i ρ 2 c i -z ij 2 .
or in vector form:

L ρ (c, z, y) = N i=1 1 2 (q T i c i -xi,M ) 2 + j∈N i y T ij (c i -z ij ) + j∈N i ρ 2 c i -z ij 2 .
The ADMM solution to this problem acts in three steps:

• First, minimization with respect to c i :

c i [t + 1] = Ω C i [argmin L(c i , z ij [t], y i [t])]. (4.3.40)
where Ω C i [•] stands for the projection in the constraints set of the vector in the argument.

• Second, minimization with respect to z ij with the constraint z ji = z ij :

z ij [t + 1] = argmin L(c i [t + 1], z ij , y i [t]). (4.3.41)
• Third, Lagrange multiplier update:

y ij [t + 1] = y ij [t] + ρ(c i [t + 1] -z ij [t + 1]). ( 4

.3.42)

Solving the sub-optimization problem (4.3.40) through ∂Lρ(c,z,y) ∂c i = 0, we get:

(q i q T i -xi,M )q i + j∈N i y T ij (c i -z ij ) + ρ 2 j∈N i (c i -z ij ) 2 = 0. ⇔ (q i q T i + ρd i I h+1 )c i = xi,M q i + ρ j∈N i z ij - j∈N i y ij . Then, ĉi [t + 1] = (q i q T i + ρd i I h+1 ) -1 (x i,M q i + ρ j∈N i z ij [t] - j∈N i y ij [t]).
Then, projecting the obtained solution in the constraints set

c i [t + 1] = Ω C i [ĉ i [t + 1]],
we get:

c i [t + 1] = Qi xi,M + (I h+1 -Qi Q)ĉ i [t + 1], (4.3.43) with Qi = Q T i (Q i Q T i ) -1
and Q i be a full-row rank matrix.

Next, solving the sub-optimization problem (4.3.41) yields

z ij [t + 1] = 1 2 (c i [t + 1] + c j [t + 1]) + 1 2ρ (y ij [t] + y ji [t]). (4.3.44)
Then, we project z ij [t + 1] on the sign constraint as:

z k ij [t + 1] =    max(0, z k ij [t + 1]) if mod(k, 2) = 0. min(0, z k ij [t + 1]) if mod(k, 2) = 0. (4.3.45)
The corresponding ADMM algorithm is then described as follows:

Algorithm 7 (ADMM-based distributed estimation of stepsize α k in the case of inexact x)

1. Inputs:

• Given the number of agents N,

• Define h = N -1, the initial input {x i (0)}, i = 1, • • • , N, and xi,M = x i (M)
for some M,

• Each node forms the vector q i , and the matrices Q i and Qi .

Initialization:

• Penalty parameter ρ,

• Random initial values of the coefficients {c i,k },k = 0, . . . , h, at each node i, i = 1, . . . , N.

• Random initial values of z ij [0] and y ij [0].

• Set t = 0;

3. Update Process:

(a) Set t := t + 1; At each iteration, node i sends a message including its local coefficient c i,k .

After receiving the information from its neighbors, each node carries out the update process.

(b) Update equations:

• Compute c i [t + 1] using (4.3.43).

• Compute z ij [t + 1] using (4.3.44), (4.3.45).

• Update the Lagrange multipliers y ij [t + 1] using (4.3.42).

(c) Return to 3a or stop the iterations if a stopping criterion is reached.

4. Each agent builds a polynomial p i with the set of coefficients c i,k [t + 1].

5. Each agent computes the set of the stepsizes {α i,k } as the roots of the polynomial p i .

When the coefficients c k are approximately obtained, we can make use of Algorithms 9 and 10 that will be presented later in Subsection 4.3.2 to estimate the whole spectrum of Laplacian matrix L, that can be used to estimate the network robustness indices.

Laplacian-based factorization of the averaging matrix with unknown consensus value x

According to the previous section, one way for estimating the Laplacian eigenvalues is:

• solve the average consensus problem.

• solve the averaging matrix factorization problem (4.3.15), Algorithms 5 and 6.

• retrieve the Laplacian eigenvalues using Algorithm 9 and 10.

Instead of following these steps, in this section, we propose a method to merge the first two steps. For this purpose, we propose to solve a convex combination of two convex functions, i.e.

min xi ,c i ,i=1,••• ,N θ 2 N i=1 j∈N i c i -c j 2 + 1 -θ 2 N i=1 (x i -x i (0)) 2 , s.t. Q i c i = xi 1, i = 1, . . . , N, (4.3.46)
with θ-an adjustable scalar parameter (0 ≤ θ ≤ 1), xi -a local estimate of the average value, and x i (0)-a fixed initial value of node i.

It can be seen that, the objective function (4.3.46) is a combination of two convex functions. The first one takes into account the estimation of the polynomial coefficients whose roots are the inverse of the nonzero distinct Laplacian eigenvalues while the second one concerns the average consensus problem. And this convex optimization problem can be solved efficiently using the ADMM approach.

In order to derive an ADMM algorithm for the problem (4.3.46), we introduce the auxiliary variables {z ij , µ ij }. Therefore, the following constrained convex optimization problem is absolutely equivalent to the problem (4.3.46):

min xi ,c i ,i=1,••• ,N θ 2 N i=1 j∈N i c i -c j 2 + 1 -θ 2 N i=1 (x i -x i (0)) 2 , s.t xi = µ ij , i = 1, . . . , N, j ∈ N i . c i = z ij . z k ij = sign(k)z k ij , c i ∈ C i (x i ).
where

• sign(k) =    1 if k is even, -1 if k is odd . • C i (x i
) stands for the constrained set with respect to node i. It is defined as

C i = {c ∈ IR h+1 | Q i c = xi 1}.
Now, introducing the dual variables υ υ υ and y, the augmented Lagrangian is defined as :

L ρ 1 ρ 2 (x, c, µ µ µ, z, υ υ υ, y) = L ρ 1 (c, z, y) + L ρ 2 (x, µ µ µ, υ υ υ).
where:

L ρ 1 (c, z, y) = θ 2 N i=1 j∈N i c i -c j 2 + j∈N i ρ 1 2 c i -z ij 2 + j∈N i y T ij (c i -z ij ), L ρ 2 (x, µ µ µ, υ υ υ) = 1 -θ 2 N i=1 (x i -x i (0)) 2 + j∈N i ρ 2 2 (x i -µ ij ) 2 + j∈N i υ ij (x i -µ ij ),
where ρ 1 , ρ 2 are penalty terms assumed to be constant.

ADMM provides an iterative solution that acts following the five steps below:

• Minimization with respect to the average consensus value xi :

xi [t + 1] = arg min L ρ 2 (x i [t], µ ij [t], υ ij [t]). (4.3.47)
• Minimization with respect to the polynomial coefficients c i :

ĉi = argmin L ρ 1 (c i [t], z ij [t], y ij [t]), (4.3.48) c i [t + 1] = Ω C i (x i [t+1]) [ĉ i ] (4.3.49)
where

Ω C i (x i [t+1]) [•]
stands for the projection into the constraints set C i (x i [t + 1]) of the vector in the argument.

• Minimization with respect to the auxiliary variables µ ij with the constraint µ ji = µ ij :

µ ij [t + 1] = argmin L ρ 2 (x i [t + 1], µ ij [t], υ ij [t]). ( 4

.3.50)

• Minimization with respect to the auxiliary variables z ij :

z ij [t + 1] = argmin L ρ 1 (c i [t + 1], z ij [t], y ij [t]). (4.3.51)
• Lagrange multipliers update:

υ ij [t + 1] = υ ij [t] + ρ 2 (x i [t + 1] -µ ij [t + 1]
), (4.3.52)

y ij [t + 1] = y ij [t] + ρ 1 (c i [t + 1] -z ij [t + 1]). (4.3.53)
Solving the sub-optimization problems (4.3.47) and (4.3.48) is equivalent to solve the following equations:

1. ∂ Lρ 1 (c,z,y) ∂c i = 0 ∂L ρ 1 (c, z, y) ∂c i = θ j∈N i (c i -c j ) + j∈N i y ij + j∈N i ρ 1 (c i -z ij ) = 0.
We get:

⇔ ĉi [t + 1] = θ j∈N i c j [t] + ρ 1 j∈N i z ij [t] -j∈N i y ij [t] d i (θ + ρ 1 ) . 2. ∂ Lρ 2 (x,µ µ µ,υ υ υ) ∂ xi = 0 ∂L ρ 2 (x, µ µ µ, υ υ υ) ∂ xi = (1 -θ)( xi -x i (0)) + j∈N i υ ij + j∈N i ρ 2 ( xi -µ ij ) = 0 ⇔ (1 -θ + ρ 2 d i ) xi = (1 -θ)x i (0) + ρ j∈N i µ ij - j∈N i υ ij ⇔ xi [t + 1] = (1 -θ)x i (0) + ρ 2 j∈N i µ ij [t] -j∈N i υ ij [t] 1 -θ + ρ 2 d i . (4.3.54)
Then, projecting the obtained solution in the constraints set

c i [t + 1] = Ω C i (x i [t+1]) [ĉ i [t + 1]],
we get:

c i [t + 1] = Qi xi [t + 1]1 + (I h+1 -Qi Q)ĉ i [t + 1], (4.3.55) with Qi = Q T i (Q i Q T i ) -1
, Q i being full row rank by construction. Next, solving the sub-optimization problems (4.3.50) and (4.3.51), in the same way, yields 1. The derivative of L ρ 2 (x, µ µ µ, υ υ υ) with respect to µ ij with respect to the constraint

µ ij = µ ji : ∂L ρ 2 (x, µ µ µ, υ υ υ) ∂µ ij = -υ ij -υ ji -ρ 2 (x i + x j ) + 2ρ 2 µ ij = 0 ⇔ µ ij [t + 1] = xi [t + 1] + xj [t + 1] 2 + υ ij [t] + υ ji [t] 2ρ 2 . (4.3.56)
2. The derivative of L ρ 1 (c, z, y) with respect to

z ij ∂ Lρ 1 (c i ,z ij ,y ij ) z ij = 0 ⇔ -y ij [t] -ρ(c i [t + 1] -z ij [t + 1]) = 0 ⇔ z ij [t + 1] = c i [t + 1] + 1 ρ 1 y ij [t]. (4.3.57)
Then, we project z ij [t + 1] on the sign constraint as:

z k ij [t + 1] =    max(0, z k ij [t + 1]) if mod(k, 2) = 0. min(0, z k ij [t + 1]) if mod(k, 2) = 0. (4.3.58)
The parameter θ allows weighting the two combined cost functions. Since the constraint set C i (x i ) depends on the consensus value xi . The optimization procedure can initially give more weight to the average consensus problem and then to the factorization one. For this purpose, instead of keeping θ constant, we adopt a time varying coefficient:

θ(t) = 1 -e -βt
1 + e -βt , 0 < β < 1, (4.3.59) to keep it in the range [0, 1] and to increase the importance of the factorization progressively while decreasing that of the average consensus.

Remark 12

Practically, we can see that the average consensus term converges much faster than the term relative to the averaging matrix factorization problem. Hence, θ can be chosen close to 1.

The corresponding ADMM algorithm is then described as follows:

Algorithm 8 (ADMM-based joint distributed average consensus and estimation of Laplacian eigenvalues)

1. Inputs:

• Given the number of nodes N, h = N -1,

• The initial measurements values {x i (0)}, i = 1, • • • , N,

Initialization:

• Penalty parameters ρ 1 , ρ 2 , and β.

• Random initial values of the coefficients {c i,k },k = 0, . . . , h, at each node i, i = 1, . . . , N.

• Random initial values of z ij , y ij and µ ij , υ ij .

• Set t = 0;

3. Update Process:

(a) Set t := t + 1; At each iteration, node i sends a message including its local value xi , and local coefficient c i,k .

After receiving the information from its neighbors, each node i carries out the update process.

(b) Update equations:

• Compute xi [t + 1] using (4.3.54)
• Compute c i [t + 1] using (4.3.55).

• Compute µ ij [t + 1] using (4.3.56).

• Compute z ij [t + 1] using (4.3.57), (4.3.58).

• Update the Lagrange multipliers υ ij , y ij [t + 1] using (4.3.52), (4.3.53).

• Update θ(t) using (4.3.59).

(c) Return to 3a or stop the iterations if a stopping criterion is reached.

4. Each agent build a polynomial P i with the set of coefficients c i [t + 1].

5. Each agent compute the set of the stepsizes {α i,k } as the roots of the polynomial P i , the set of Laplacian eigenvalues using Algorithm 1.

The average consensus value is xi

[t + 1].

Laplacian spectrum retrieving

At each node, the set S 1 of step-sizes obtained by means of the proposed Algorithms (Algorithms 4, 5, 6, 7, 8) contains the set S 2 of the inverse of nonzero distinct Laplacian eigenvalues.

Let xi be the final consensus value reconstructed by xi = h k=0 c k q i,k = P i (c, q i ) with the coefficients obtained by means of Algorithms 4, 5, 6, 7 and 8.

Following the Proposition 3, the idea is to reduce, step-by-step, the degree of the polynomial P i (c, q i ) by removing one element of S 1 at each step. We know that if the removed element is also an element of S 2 then xi = Pi ( c, q i ) where Pi is the polynomial with reduced degree with respect to node i. Algorithm 9 describes the proposed procedure. Algorithm 9 (Laplacian eigenvalues retrieving )

Given the set of step size S 1 , the final consensus value xi , the measurements vector q i and the threshold ǫ 1. Initialization: S 2 = ∅ and S 3 = S 1 .

2. while S 3 = ∅, select an element α j of S 3 .

3. From the set S 3 \{α j } compute the coefficients c i using (4.3.5) and then form ci .

4. If | xi -cT q i | > ǫ, include α j in S 2 so that S 2 := S 2 {α j }, set S 3 := S 3 /{α j }, and return to 2.

5. If | xi -cT q i | < ǫ, set S 3 := S 3 /{α j }, and return to 2.

6. If S 3 = ∅, deduce the Laplacian eigenvalues as the inverses of the elements in S 2 .

We can note that this algorithm is to be run strictly locally.

As a result, we obtain the set of nonzero distinct Laplacian eigenvalues Λ = {λ 2 , . . . , λ D+1 }. In order to obtain the whole spectrum sp(L) = {λ m 1 1 , λ m 2 2 , . . . , λ m D+1 D+1 }, the multiplicities m i should be defined. Since all of them are integers, the problem now resorts to an integer programming problem.

Given the set of nonzero distinct Laplacian eigenvalues Λ, d i the degree of each node i, inspired from an interesting property N i=1 λ i (L) = N i=1 d i (see in Section 2.3), the optimization is now described as follows: where D = |Λ| and m is the vector of eigenvalues multiplicities.

min m∈N D Λ T m - N i=1 d i , ( 4
In the literature, the branch and bound method is the basic workhorse technique for solving integer programming problems. The linear programming (LP) relaxation (with no integer restrictions) is adopted to estimate the optimal solution of an integer programming [START_REF] Chinneck | Practical optimization: A gentle introduction[END_REF]. Therefore, problem (4.3.60) is now rewritten strictly equivalently in linear programming form as follows:

min m∈N D Λ T m, (4.3.61) subject to D j=1 m j = N -1; Λ T m = N i=1 d i , (M) m j ≥ 1,
where M is the constraint set without integer restriction. Hence, M can be expressed as (c)

M = {m ∈ IR D | D j=1 m j = N -1, Λ T m = N i=1 d i
m 1 = argmin Λ T m, s.t m ∈ M 1 = M ∪ {m j ≤ ⌊ mj ⌋} . (d) m 2 = argmin Λ T m, s.t m ∈ M 2 = M ∪ {m j ≥ ⌈ mj ⌉} • If m 1 (or m 2 ) is integer, S = m 1 (or m 2 ) then stop.
• If there is an infeasible solution, then prun the corresponding node.

• Compute UB 1 (or UB 2 ).

(e) If UB 1 (UB 2 ) < UB 2 (UB 1 ), then pick m 1 (or m 2 ) is the next node to create 2 new nodes and V = V ∪ m 2 (m 1 ), and M = M 1 (M 2 ).

Return to step 2a.

The policy of selecting a variable to be branched on in this algorithm is that the variable with the greatest fractional part is the first order then the variable with the smallest fractional part is the next. After that, the policy is repeated for the rest of unselected variables.

Since, the Laplacian spectrum sp(L) is obtained, we can compute the effective graph resistance R and the number of spanning trees ξ to assess the robustness of a given network using ( 1.3.2) and ( 1.3.1).

Simulation results

In this section, several examples will be carried out to test the efficiency of the proposed methods to assess the network robustness (via R, ξ). As presented obviously in previous Section, the procedure of all proposed methods are devided into three categories according to the consensus value x that are: perfect consensus value (Algorithms 4, 5 and 6); imperfect consensus value (Algorithm 7) and unknown consensus value (Algorithm 8).

In order to evaluate the network robustness, let us take two networks which have the same nodes but different number of edges as examples. The corresponding Laplacian matrices are given by: Assume that the final average consensus value x = 0.6884 is obtained perfectly by Algorithm 1.

L =             2 -1 0 -1 0 0 -1 2 -1 0 0 0 0 -1 3 -1 0 -1 -1 0 -1 3 -1 0 0 0 0 -1 2 -1 0 0 -1 0 -1 2             with sp(L) {0, 1, 2, 3, 3, 5} R 14.2 ξ 15 L =             2 -1 0 0 0 -1 -1 2 -1 0 0 0 0 -1 2 -1 0 0 0 0 -1 2 -1 0 0 0 0 -1 2 -1 -1 0 0 0 -1 2             with sp(L) {0, 1, 1, 3, 3, 4} R
Applying Algorithms proposed in Subsection 4.3.1.2, we can see the difference between them via their performances.

Lagrange multiplier method for distributed averaging matrix factorization

We make use of Algorithms 4, 9 to find the stepsizes α k = {1, 0.5, 0.3333, 0.2} and deduce the nonzero distinct eigenvalues {1, 2, 3, 5}. Furthermore, by using Algorithm 10, the multiplicities of these nonzero distinct Laplacian eigenvalues are achieved to get the whole Laplacian spectrum sp(L).

In practice, the parameters β, µ are chosen to ensure the convergence of the Algorithm 4: β = 0.0222, µ = 0.2. After running Algorithm 4, we obtain the set of stepsizes {0.9993, 0.2, 0.5, 0.0895, 0.3333}. The trajectory of these stepsizes is illustrated in Figure 4.5.

Here, running Algorithm 9, the real set of the nonzero distinct Laplacian eigenvalues can be inferred Λ = {1, 5, 2, 3}. Finally, by Algorithm 10, the corresponding multiplicities are obtained m = {1, 1, 1, 2}. Therefore, we get the whole Laplacian spectrum sp(L).

Due to the fact that the proposed method is based on the gradient concept of a nonconvex optimization problem, the convergence speed is obviously slow. Besides that, the choice of the parameters β, µ plays an important role. By using the equations ( 1.3.1) and ( 1.3.2), the estimated robustness indices can be obtained:

• Number of spanning trees ξ:

ξ = 1 N N i=2 λ i (L) = 15.0135. • Effective graph resistance R R = N N i=2 1 λ i (L) = 14.1954.

Distributed projected subgradient method

This convex optimization is implemented with respect to the coefficients c k , k = 1, . . . , 6. Therefore, the performance is evaluated via the mean square error (E(c)) defined as:

E(c) = 1 N N i=1 h k=0 j∈N i (c i,k -c j,k ) 2 .
We now try to find the coefficients c k , k = 0, . . . , 6 by using the Algorithm 5 and 6, which deduces the stepsizes that are the inverses of the nonzero distinct Laplacian eigenvalues. Then, via Algorithms 9 and 10, the whole spectrum of the Laplacian matrix is obtained. In comparison with the Lagrange-based method (Algorithm 4), the speed of convergence is improved.

In this implementation, we have constrained the sign of the elements of the coefficient vector c. Therefore, by c i,5 , i = 1, . . . , N is set to be 0 due to the fact that c i,5 > 0, i = 1, . . . , N. Furthermore, by solving the roots of the polynomial of coefficients c k , we obtain the stepsizes α k = {0.9998, 0.4996, 0.3342, 0.1993} that are definitely close to {1, 0.5, 0.3333, 0.2}. By making use of Algorithm 10, the whole Laplacian spectrum sp(L) = {1, 2, 3, 3, 5} is obtained. Using Algorithm 9, we get S 2 = {1, 0.5, 0.3333, 0.2} from which we deduce Λ = {1, 2, 3, 5}.

ADMM-based method

Through Algorithm 10, the whole spectrum sp(L) is achieved.

Another interesting note is the performance of the proposed algorithm with respect to the penalty parameter ρ. Hereafter, we describe a comparison with different choices of ρ. In fact, the main point is that for the distributed projected subgradient method, which is still based on gradient method, the convergence depends on the choice of the stepsize γ. On the other hand, the ADMM-based method seems to be the most attracting method when the only parameter needs to be noticed is the penalty parameter ρ to regulate the speed of convergence.

Figure 4.10 shows that when ρ increases, the speed of convergence is faster. But, the error oscilates as well. However, we still obtain the set of the inverses of the Laplacian eigenvalues.

Since both methods (distributed projected subgradient method and ADMM-based method) obtain the exact Laplacian spectrum sp(L), we can deduce the estimated robustness indices by using equations ( 1.3.1) and ( 1.3.2):

• Number of spanning trees ξ:

ξ = 1 N N i=2 λ i (L) = 15. (4.4.1)
• Effective graph resistance R:

R = 1≤i,j≤N R ij = N N i=2 1 λ i (L) = 14.2 (4.4.2)
The last thing can be noted here is that the solution of the coefficients c k is not unique since the number of distinct Laplacian eigenvalues (D = 4) is less than h = N -1 = 5. The considered network is modelled by a 6-node cycle graph in Figure 4.4(b) with the initial values for the consensus protocol:

x(0) = [0.3074, 0.4561, 0.1017, 0.9954, 0.3321, 0.2973] T .

Moreover, we know that the number of distinct nonzero Laplacian eigenvalue D is equal to the diameter of this graph (d(G) = 3). Therefore, we set h = D = 3 and start to run the algorithm to assess the network robustness.

Lagrange multiplier method for distributed averaging matrix factorization

Picking µ = 0.15, β = 0.003, and running Algorithm 4, we obtain 3 distinct inverses of the Laplacian eigenvalues Λ = {0.3333, 1, 0.25}. Then, by applying the Algorithm 10, the corresponding multiplicities are estimated m = {2, 2, 1}. Therefore, the set of distinct Laplacian eigenvalues is S 2 = {1, 0.3333, 0.25} which deduces Λ = {1, 1, 3, 3, 4} through Algorithm 10.

ADMM-based method

R ij = N N i=2 1 λ i (L) = 17.5 (4.4.4)
To conclude, we can see that, network (a), which has a larger number of edges, is more robust than the network (b).

Imperfect consensus value x

The performance is evaluated by means of the mean square error (MSE) between the estimated Laplacian eigenvalues λj,i and the actual ones.

MSE = 1 N N i=1 D j=1 (λ j -λj,i ) 2 . (4.4.5)

Example 1

Firstly, we run an average consensus protocol using the constant edge weights consensus protocol 1.2.9 with α = 0.2 as a stepsize, and x(0) = [0.4401, 0.5271, 0.4574, 0.8754, 0.5181, 0.9436] T as an initial condition. We can note that a reasonable agreement is obtained after iteration 22. Then, Algorithm 7 for computing the factorization of the averaging matrix was performed for different values for the number of iterations M of the average consensus algorithm.

For M ≤ 22, x M was too far from x. As a consequence, the obtained polynomial gave rise to negative or complex-valued roots meaning that the coefficients c created from these roots do not satisfy xi = cT i q i in Algorithm 9. For M > 22, Algorithm 9 was successfully applied for retrieving the Laplacian eigenvalues approximately. As expected, we can note that the closer x M is to x = 0.6269, the higher is the precision on the computation of the Laplacian eigenvalues.

In addition, the tight bounds around the MSE inform us about the agreement on the estimated Laplacian eigenvalues.

We can note that when the consensus value is close to the actual average, the Laplacian eigenvalues are very well estimated and the agreement on the estimated eigenvalues is perfect. Figure 4.19 shows that when the consensus value is far from the actual average the nodes fail to compute the Laplacian eigenvalues. However, some of them are able to carry out such a computation while the other fail. For instance, at iteration M = 22, there are only 3 nodes estimating the inverse of the Laplacian eigenvalues approximately.

RE R = 1R -R R R 2 R 2 , RE ξ = 1ξ -ξ ξ ξ 2 ξ 2 .
Where R R R, ξ ξ ξ are the estimated vectors of local estimated Effective Graph Resistance R i , i = 1, . . . , N and locally estimated the number of spanning trees ξ i , i = 1, . . . , N respectively. Figure 4.22 shows that the effective graph resistance R is more sensitive to the MSE than the number of spanning trees ξ. Let see the estimation of the polynomial coefficients c k at M = 24 via Figure 4.23. We can see that, for some imperfect value of x close enough to the actual consensus value x, we still assess approximately the network robustness.

Example 2

For this network, the diameter d(G) = 3 can be known in a priori. Therefore, we can

set h = d(G) = 3.
As depicted in Figure 4.24, we first run an average consensus protocol using the constant edge weights consensus protocol 1.2.9 with α = 0.2 as a stepsize and set x(0) = [0.2785, 0.5469, 0.9575, 0.9649, 0.1576, 0.9706] T as an initial condition. Figure 4.24 depicts the trajectory of the state of each node. We can note that a reasonable agreement is obtained from iteration 7.

For instance, Figure 4.25 depicts the trajectories of the estimated polynomial coefficients. We can note the good agreement between nodes on the estimation of these polynomial coefficients, which are then used to build a 3rd-order polynomial whose roots contain the inverse of the Laplacian eigenvalues.

For M > 7, Algorithm 9 was successfully performed for retrieving the Laplacian eigenvalues. Figure 4.26 depicts the MSE on the estimation of the Laplacian eigenvalues and the corresponding standard deviation. As expected, we can note that the closer x M is to x = 0.646, the higher is the precision on the computation of the Laplacian eigenvalues. In addition, the tight bounds around the MSE inform us about the agreement on the estimated Laplacian eigenvalues. We can note that when the consensus value is close to the actual average, the Laplacian eigenvalues are very well estimated and the agreement on the estimated eigenvalues is perfect. Figure 4.27 shows that when the consensus value is far from the actual average (M < 7) the nodes fail to compute the Laplacian eigenvalues. However, some of them are able to carry out such a computation while the other fail (M = 7, . . . , 17). be far from the actual value x, the convergence of Algorithm 7 is not enough to deduce the exact polynomial of c k that can be used to infer the Laplacian spectrum. To assess the independence of the estimates on the initial condition x(0), a Monte-Carlo simulation has been carried out. 100 initial vectors x(0) drawn independently from a uniform distribution have been generated.

MSE c = 1 N N i=1 c i -c 2 c 2 , and MSE x = 1 N N i=1 (x i -x) 2 x2 ,
where x, c stand for the actual average consensus value and actual polynomial coefficients respectively. As can be seen from the figure, when θ(t) is near to 0, the convergence speed of the average consensus problem is faster. However, when θ(t) tends to 1, the convergence speed of MSE c becomes faster. On the other hand, the choice of θ(t) also influences on the convergence speed of the Algorithm 8. Moreover, with the same choice of penalty parameter ρ = 0.1 and threshold (stopping criterion), regarding to the accuracy of the estimation of the coefficients c t , the choice of value β = 0.01 is the best one. Besides that, all of the solutions of the polynomial coefficients c k give the same assessment of the network robustness. Regarding to the speed of the convergence, β = 0.06 seems to be the good choice in this example.

Example 2

We first generate the initial values x(0) = [0.6531, 0.0403, 0.5047, 0.8945, 0.3857, 0.2921] T .

According to the proposed Algorithm 8, the nodes have to achieve agreement on the values of the polynomial coefficients c k and that of the average consensus value x. As can be seen from this Table, despite of the different choices of β, the proposed method gives exactly the same polynomial coefficients c k , which can be used to infer the same assessment of the network robustness. However, with β = 0.1, the speed of convergence is the fastest.

Conclusion

In this Chapter, we have proposed some methods to assess the robustness of a given network in a distributed way. Formally, the robustness metrics called the effective graph resistance R and the number of spanning trees ξ are functions of the Laplacian spectrum sp(L). Therefore, these proposed methods are distinguished by the way to estimate the Laplacian spectrum. In detail, three cases have been considered with respect to the consensus value x, that are: perfect consensus value x, imperfect consensus value xM and unknown consensus value x. Based on recent results, the average consensus matrix can be factored in D Laplacian based consensus matrices, where D stands for the number of nonzero distinct Laplacian eigenvalues. For this purpose, all proposed methods are characterized into two approaches, which are the direct approach and indirect approach in the case of perfect consensus value x. A direct approach was first presented. It gives rise to a non-convex optimization problem with respect to the stepsizes α k , k = 1, . . . , D = N -1 that are the inverse of the nonzero distinct Laplacian eigenvalues. For this approach, since it is based on the gradient descent method, the main issue is that it suffers from slow convergence towards a local minimum. To improve that, a convex optimization problem has been considered to estimate the optimal coefficients c k , k = 0, . . . , N -1 of a given polynomial whose roots are the stepsizes α k , k = 1, . . . , N -1 that allow to factorize the averaging matrix. In this way, we have applied two well-known methods that are the projected subgradient method and the ADMM to solve the convex problem.

After that, the Laplacian eigenvalues are thus estimated by taking the inverse of these obtained stepsizes α k , k = 1, . . . , N -1. Then the real set of nonzero distinct Laplacian eigenvalues Λ = {λ i , i = 2, . . . , D} is captured by a retrieving algorithm. Finally, the whole Laplacian spectrum is obtained via defining the corresponding multiplicities of these eigenvalues. In what follows, the network robustness measure can be easily deduced.

The drawback of the proposed methods is their scalability. In fact, for large graphs, the whole Laplacian spectrum estimation becomes a significant burden. Therefore, the estimation can be restricted to that of the most significant eigenvalues used to estimate bound of the network robustness index. For example, λ 2 (L) is a lower-bound of the vertex connectivity (K). Furthermore, we may encounter ill-conditioned problem in implementation of the indirect proposed methods. Hence, we can investigate new method for preconditioning the problem.

In the case of imperfect consensus value, we assume that a standard consensus protocol has been run and stores some intermediate values of xM according to the iteration M. Then, we make use of the ADMM-based Algorithm 7 to see how the network robustness reacts. Finally, it claims that when the intermediate value of x is as close as possible to the real one, the network robustness index can be taken approximately.

In the case of unknown consensus value x, we consider the convex optimization problem which resorts to a convex combination of two convex functions. The first one takes into account the average consensus problem while the second one concerns the estimation of the Laplacian matrix. It is an indirect method where parameters of the factorization of the averaging matrix are computed first and then the Laplacian eigenvalues are estimated as the roots of a given polynomial. For large graphs, such a solution is clearly intractable. Therefore, the proposed Algorithms 5, 6, 7, 8 are efficient and well adapted for networks with medium size.

The ADMM seems to be the best method compared with all the methods proposed. Therefore, the study of the convergence rate must be more interesting. In the near future, we will investigate the convergence rate of the ADMM for the quadratic programming based on the work in [START_REF] Ph | Explicit convergence rate of a distributed alternating direction method of multipliers[END_REF].

Chapter 5

Network topology reconstruction

Introduction

Network topology is usually a schematic description of the arrangement of a network, including its nodes and connecting edges. Briefly speaking, network topology may describe how the data is transferred between these nodes. 1Network topology identification refers to detecting and identifying the interested network elements and the relationship between elements in target network and represents the topology construction in an appropriate form. Therefore, identification of networks of systems becomes a growing attractive task for solving many problems in different science and engineering domains, for instance, in biology (biochemical, neural and ecological networks [START_REF] Aderhold | Hierarchical Bayesian models in ecology: Reconstructing species interaction networks from non-homogeneous species abundance data[END_REF]), finance, computer science (Internet and World Wide Web), transportation (delivery and distribution network), and electrical engineering.

For example, the architecture of an overlay network -how it allocates addresses, etc.may be significantly optimized by knowledge of the distribution and connectivity of the nodes on the underlay network that actually carries the traffic. Several important systems, such as P4P (stands for provider portal for Peer-to-Peer (P2P) applications) and RMTP (reliable multicast transport protocol), utilize information about the topology of the underlay network for optimization as well as management [START_REF] Acharya | On the hardness of topology inference[END_REF].

In the Internet, for instance, the usual mechanism for generating the topology of a network is by the use of Traceroute. Traceroute is executed on a node, called the source, by specifying the address of a destination node. This execution produces a sequence of identifiers, called a trace, corresponding to the route taken by packets travelling from the source to the destination. A trace set T is generated by repeatedly executing Traceroute over a network, varying the terminal nodes, i.e. the source and destination. If T contains traces that identify every instance when an edge is incident on a node, it is possible to reconstruct the network exactly [START_REF] Acharya | On the hardness of topology inference[END_REF].

Following our study, from a general consensus network and consensus measurements, we deal with the problem of inferring the correct network topology in the presence of anonymous nodes (or unlabelled nodes). In other words, is it possible for an arbitrary node j to reconstruct the correct network from average consensus measurements?

Literature review

In the literature, there are various methods proposed for reconstructing the network topology.

For centralized framework, in [START_REF] Nabi-Abdolyousefi | Network identification via node knockout[END_REF], the authors introduced a network identification scheme which involves the excitation and observation of nodes running consensus-type protocols.

The problem is to identify the interaction geometry among a known number of nodes, adopting a (weighted) consensus-type algorithm for their coordination. In the considered procedure, the node broadcasts a zero state to its neighbors without being removed from the network. This node is called node knockout, which is essentially a grounding procedure. Starting with the number of vertices in the network as a known parameter as well as the controllability and observability of the resulting steered-and-observed network, define the characteristic polynomial of the grounded consensus at assumed known node i, which leads to deduce the matrix Ψ Ψ Ψ as the adjugate of the matrix (I N + L) in order to determine the entries of L. The advantage of this paper is the ramifications for exact identification from boundary nodes of networks that have an embedded consensus-based distributed algorithm used for flocking and distributed estimation.

In [START_REF] Shahrampour | Topology identification of directed dynamical networks via power spectral analysis[END_REF], the problem of identifying the topology of an unknown directed network of linear time-invariant (LTI) systems stimulated by a wide-sense stationary noise of unknown power spectral density has been considered. The authors have proposed several reconstruction algorithms based on the power spectral properties of the network response to the noise. In [START_REF] Hayden | Network reconstruction from intrinsic noise: Minimum-phase systems[END_REF], the identifiability of the structure and dynamics of an unknown network driven by unknown noise has been assessed based on factorizations of the output spectral density. In this paper, the authors have claimed that for networks with closedloop transfer matrices that are minimum phase, the network reconstruction problem can have a unique solution from its output spectral density.

In [START_REF] Ren | Noise bridges dynamical correlation and topology in coupled oscillator networks[END_REF], the reconstruction concepts are based on paradigmatic dynamical models such as phase oscillators. The authors have discovered a general connection between dynamical correlation among oscillators and the underlying topology in the presence of noise. In particular, the correspondence between the dynamical correlation matrix of nodal time series and the connection matrix of structures, due to the presence of noise.

In order to bridge dynamical correlation and topology, consider N non-identical coupled oscillators x i ∈ IR m , i = 1, . . . , N, which denotes the m-dimensional state variables of i-th oscillator. In the presence of noise, the dynamics of the whole coupled-oscillator system can be expressed as:

ẋi = F i (x i ) -γ N j=1 l ij P(x j ) + η η η i , (5.2.1)
where:

• F i : IR m → IR m is intrinsic dynamics of the i-th oscillator,

• P : IR m → IR m is the coupling function of oscillators,

• l ij is an element at i-th row and j-th column of the Laplacian matrix L,

• γ > 0 is the coupling strength,

• η η η = [η 1 , . . . , η N ] T is the noise vector with covariance matrix D.

Let xi be the counterpart of x i in the absence of noise and assume a small perturbation ξ i , then x i = xi + ξ i . Substituting this into (5.

0 = d dt   1 T T 0 ξ ξ ξξ ξ ξ T dt   = -MC -CM T + 1 T T 0 η η ηξ ξ ξ T dt + 1 T T 0 ξ ξ ξη η η T dt, (5.2.2)
where M = -(DF(x) + γL ⊗ DP(x)). After a few basic manipulations, (5.2.2) can be rewritten as:

MC + CM T = D. (5.2.3)
This equation illustrates a general relationship between the dynamical correlation C and the Laplacian L in the presence of noise as characterized by D. If we consider onedimensional state variable and linear coupling such that DP = 1, with Gaussian white noise D = σ 2 I mN , and further regard the intrinsic dynamics DF as small perturbations.

Then, (5.2.3) can be simplified to LC + CL T = σ 2 I mN γ . For an undirected network with symmetric coupling matrix, the solution of C can be expressed as

C = σ 2 2γ L + , (5.2.4)
where L + denotes the pseudo inverse of the Laplacian matrix.

Hence, the dynamic correlation matrix C is closely related to the network connection matrix L, which can be used to infer network structures when no knowledge about nodal dynamics is available.

From (5.2.4), the Laplacian L can be extracted from the knowledge of noise-corrupted time series:

L = σ 2 2γ C + , (5.2.5)
where C + is the pseudo inverse. The correlation matrix is inversely proportional to the Laplacian matrix that contains full information about the network structure. The reconstruction of the network topology based on time series then becomes possible for undirected networks. In case of directed networks, there is no unique solution for L from C, because the asymmetric L has twofold degree of freedom as that of symmetric of C. Thus, the global structure of directed networks cannot be inferred solely depending on the correlation, but we can still infer the local structure.

Based from this method, in [START_REF] Morbidi | A distributed solution to the network reconstruction problem[END_REF] the authors have presented a fully distributed algorithm for reconstructing the topology of an undirected consensus network by using a distributed least-squares algorithm.

Assume the number of nodes N is known, the noisy consensus dynamics is given by: ẋ = -γLx + η η η, (5.2.6)

(5.2.7) where x = [x 1 , . . . , x N ] T ∈ IR N and η η η ∈ IR N .

From (5.2.5), we obtain:

L = σ 2 2γ C + = σ 2 2γ C + J N -1 -J N , (5.2.8)
where J N is the averaging matrix. And C is obtained as:

C = 1 T T 0 (x(t) -x1)(x -x1) T dt where x = 1 N N i=1 x i .
Therefore,

L = σ 2 2γ (I N -J N )C R (I N -J N ) + J N -1 -J N , (5.2.9)
where

C R = 1 T T 0 x(t)x T (t)dt.
Assuming s discrete observations in the time interval [0, T ], C R can be expressed approximately as:

C R ≃ 1 s s-1 k=0 x(k)x T (k).
In order to estimate the C R , the measurement equation for node i is introduced:

y i (k) = E i x(k)
where E i ∈ {0, 1} (|N i |+1)×N is a selection matrix. By collecting the measurements of all nodes at time k,

y(k) = Ex(k), k ∈ {0, 1, . . . , S -1} where y = [y T 1 . . . y T N ] T and E = [E T 1 . . . E T N ] T . If E T E is invertible, then the least-squares estimate of x(k) is x(k) = E + y(k) = (E T E) -1 E T y(k) = 1 N N i=1 E T i E i -1 1 N N i=1 E T i y i (k)
Let assume node i maintains a matrix B i ∈ IR N ×N and a vector z i ∈ IR N and executes the following iterations for h ∈ {0, 1, . . . , D -1}, D > 1 (5.2.11) where

B i (h + 1) = B i (h) + ζ j∈N i (B j (h) -B i (h)), (5.2.10)
z i (h + 1) = z i (h) + ζ j∈N i (z j (h) -z i (h)),
B i (0) = E T i E i and z i (0) = E T i y i (k), i = 1, . . . , N. If ζ ∈ (0, 1) satisfies the condition ζ < 2 λ N (L) , then lim h→∞ B i (h) = 1 N N i=1 E T i E i (5.2.12) lim h→∞ z i (h) = 1 N N i=1 E T i y i (k), (5.2.13)
and node i asymptotically obtains an estimate xi of vector x at time k:

xi (k) = lim h→∞ B -1 i (h)z i (h).
Node i can then compute its local estimate C R,i of C R as:

ĈR,i = 1 s s-1 k=0 xi (k)x T i (k).
and eventually obtains the Laplacian matrix through (5.2.9).

In [START_REF] Kibangou | Decentralized Laplacian eigenvalues estimation and collaborative network topology identification[END_REF], the authors have proposed a decentralized reconstruction algorithm using an algebraic method using observability properties of a network. Consider a connected undirected graph G(V, E) with the dynamics presented as:

X(k + 1) = WX(k). Y i (k) = E i X(k).
Where:

• X(k) ∈ IR N ×M is a matrix with x i (k), i = 1, . . . , N as rows;

• Y i (k) is the outputs or node values that are seen by node i at the k th time step;

• E i ∈ IR di ×N is the row selection matrix with di in the interval [0, d i + 1];

• W = I N -ǫL is the consensus matrix. Its eigenvalue decomposition is W = UDU T , where the eigenvector and eigenvalues are respectively organized in the orthogonal matrix U and the diagonal matrix D.

The identification of the network topology is implemented by estimating both eigenvalues (in Chapter 4.2) and eigenvectors of the network matrix. From the observability theory, the pair (W, E i ) is observable if and only if the observability matrix O i,K i ∈ IR (K i +1) di ×N defined as:

O i,K i =           E i E i W E i W 2 .
. .

E i W K i           =        V i V i D . . . V i D K i        U T = Ōi,K i U T with V i = E i U is full column rank, i.e (O i ) = N.
From the available observations Y i (k), k = 0, 1, . . . , K i + 1, and the two generated matrices:

Ȳi =        Y i (0) Y i (1)
. . .

Y i (K i )        ; Ȳi =        Y i (1) Y i (2)
. . .

Y i (K i + 1)        , we can show that Ȳi = Ōi,K i C T Ȳi = Ōi,K i DC T , with C T = U T X(0).
As presented in Chapter 4.2, Ũ1 , Ũ2 are respectively the top and the bottom parts of Ũ of left singular vectors of matrix ( ȲT

i ȲT i ) T corresponding to Ũ1 = Ōi,K i T T and Ũ2 = Ōi,K i DT T . By defining G = T ŌT i,K i Ōi,K i , we get R 1 = GT T and R 2 = GDT T . Then, with R = R 2 R -1
1 , we get RG = GD, which is a standard eigenvalue problem. Then, D, G can be deduced.

Then, we can estimate

T T = G -1 R 1 and Ōi,K i = R -1 1 G Ũ1 . As a result, C T = G -1 R 1 Ȳ1 .
In order to infer the structure of the network, assume that each node knows E i , C, D, the authors try to compute U.

At first, build the matrices Φ Φ Φ i and the vectors Ψ Ψ Ψ i defined as:

Φ Φ Φ i =        C ⊗ E i CD ⊗ E i . . . CD K i ⊗ E i        , Ψ Ψ Ψ i =        vec(Y i (0)) vec(Y i (1)) . . . vec(Y i (K i ))       
Then, compute the local quantities Φ Φ Φ T i Φ Φ Φ i and Φ Φ Φ T i Ψ Ψ Ψ i with an average consensus algorithm to get:

Φ Φ Φ = 1 N N i=1 Φ Φ Φ T i Φ Φ Φ i Ψ Ψ Ψ = 1 N N i=1 Φ Φ Φ T i Ψ Ψ Ψ Ψ Ψ Ψ Ψ Ψ Ψ i .
Finally, by solving the least-squares problem Φ Φ Φvec(U) = Ψ Ψ Ψ, the matrix U can be computed.

At the end, each node can estimate the consensus matrix as W = UDU T , and this is the network topology.

In all these methods, the presence of anonymous nodes is not considered.

Problem statement

From the previous chapter we now know that node i can obtain the overall Laplacian spectrum from average consensus measurements. We state the following assumption:

Assumption 1: All the nonzero Laplacian eigenvalues are distinct.

By considering a network represented by a graph G(V, E), we recall that at time instant k, the state of this network is given by the vector

x(k) = Wx(k -1)
or equivalently

x(k) = W k x(0)
Let us consider the following eigenvalue decomposition of the consensus matrix W:

W = UDU T , (5.3.1)
with D = I N -ǫ∆ ∆ ∆ and ∆ ∆ ∆ is the diagonal matrix of Laplacian eigenvalues. ǫ can be chosen as

2 λ 2 (L)+λ N (L [74], then λ i (W) = 1 -ǫλ i (L), i = 1, . . . , N.
Since the Laplacian eigenvectors form a basis of IR N , the initial condition x(0) can be expanded such that x(0) = N i=1

β i U .i = Ub, where b = β 1 β 2 • • • β N T contains
the expansion coefficients, and U .i stands for the i-th eigenvector, i.e. the i-th column of U. Therefore:

x(k) = UD k U T Ub = UD k b = Udiag(b)vecd(D k ),
where vecd(•) is the column vector built with the diagonal entries of the matrix in the argument.

We can equivalently write the above equation as:

x(k) = Ũvecd(D k ), with Ũ = Udiag(b)
Therefore, Ũ ŨT = Udiag(b. 2 )U T and ŨT Ũ = diag(b. 2 ). We can conclude that: At node j, the state at instant k is then given by:

W = ŨDdiag(b. 2 ) -1 ŨT . ( 5
x j (k) = N i=1
λ k i (W)β i u j,i , u j,i being the jth entry of the eigenvector U .i . By stacking N consecutive measurements node j obtains:

      
x j (0)

x j (1) . . .

x j (N -1)        =        1 1 • • • 1 1 λ 2 (W) • • • λ N (W) . . . . . . . . . . . . 1 λ N -1 2 (W) • • • λ N -1 N (W)               ũj,1 ũj,2 . . . ũj,N        , ũj,i = u j,i β i or equivalently x j = Υ Υ Υ ŨT j.
(5.4.1)

From Assumption 1, we can conclude that Υ Υ Υ is a full rank matrix. Therefore node j can efficiently solve the above problem and get:

ŨT j. = Υ Υ Υ -1 x j . (5.4.2)
Now, the node only knows the jth row of Ũ.

Let us consider a message passing scheme where the nodes send to their neighbors the rows they have. Then after a number of message exchanges, at least equal to the diameter of the graph, node j has all the rows of Ũ. However, since some of the nodes are assumed to be anonymous, they do not send the label of the row they have. Two kinds of messages are transmitted in the network. The rows associated to non-anonymous nodes are correctly labelled while those associated to anonymous nodes are not labelled. As a consequence, node j only receives Ũ up to rows permutation, i.e. it receives Û = Π Π Π Ũ, Π Π Π being a permutation matrix.

Network topology reconstruction

We can note that ÛT Û = diag(b 2 ). Therefore, from (5.3.2), we get:

W = Π Π Π ŴΠ Π Π T , with Ŵ = ÛD( ÛT Û) -1 ÛT . (5.4.3)
We can conclude that from the estimated matrix Ŵ, we obtain a graph that is isomorphic to the original graph. In what follows we state two sufficient conditions to ensure correct reconstruction of the graph. We first recall the following notions of observability:

A network is said to be [START_REF] Kibangou | Observability in connected strongly regular graphs and distance regular graphs[END_REF]:

• node-observable from a given node if that node is able to reconstruct the entire network state from its own measurements. This issue has been studied for instance in [START_REF] O'clery | Observability and coarse graining of consensus dynamics through the external equitable partition[END_REF] and [START_REF] Parlangeli | On the reachability and observability of path and cycle graphs[END_REF] where it has been stated that a network with a consensus matrix having at least one non-simple eigenvalue is not node-observable.

• neighborhood observable from a given node if that node can reconstruct the entire network state from its own measurements and those of its neighbors. This issue was studied in [START_REF] Kibangou | Algebraic characterization of observability in distance-regular consensus networks[END_REF].

• globally observable if it is neighborhood observable from any node.

Proposition 4

Assume that the graph is node-observable or neighborhood observable from node j. If all the entries of the initial condition are distinct then the network topology can be exactly reconstructed from node j.

Proof: From Ũ and diag(b 2 ) the initial condition x(0) can be reconstructed as follows:

x(0) = Ũdiag(b 2 )β β β, with β β β = 1/β 2 1 1/β 2 2 • • • 1/β 2 N T
. Since Ũ is known only up to rows permutation through Û, the reconstructed initial condition x(0) = Ûdiag(b 2 )β β β is related to the actual initial condition as follows: x(0) = Π Π Πx(0). The graph being node-observable or neighborhood-observable from node j, node j can reconstruct the initial condition from its consecutive measurements. It can therefore compare x(0) with x(0). If all the entries of x(0) are distinct then the permutation matrix can be retrieved and the matrix W exactly recovered from Ŵ.

Proposition 5

Assume that the network contains a set A of anonymous nodes and the graph is nodeobservable or neighborhood observable from node j. If all the entries of the initial condition associated with the anonymous nodes are all distinct then the network topology can be exactly reconstructed from node j.

The proof of this proposition is similar to the previous one. Here, the permutation ambiguity is restricted to the entries associated with the anonymous nodes.

We can also note that if the initial condition is driven from a continuous joint probability distribution, then the reconstruction of the network topology is almost sure.

The proposed network topology reconstruction algorithm is described as followed:

Algorithm 11 (Network topology reconstruction)

1. Inputs:

• Given the Laplacian spectrum sp(L) and the number of nodes N.

• Given a stepsize ǫ of consensus protocol with matrix W = I N -ǫL.

• Given the first N -1 consensus measurements x i = [x i (0), x i (1), . . . , x i (N -1)] T , i = 1, . . . , N.

Initialization:

• Set of eigenvalues λ i (W) = 1ǫλ i (L) for i = 1, . . . , N.

• Matrix Υ Υ Υ in equation (5.4.1).

3. Each node j estimates the j-th row of Ũ using (5.4.2).

4. Each node sends to its neighbors the rows it has. After a number of message exchanges, each node j obtains Û = Π Π Π Ũ.

5. Since diag(b 2 ) = ÛT Û is a diagonal matrix, then each node can then estimate the matrix topology as W = Π Π Π ÛDdiag( ÛT Û) -1 ÛT Π Π Π T 6. Define the estimated Laplacian matrix L = 1 ǫ (I N -W). Then take the rounded estimated Laplacian matrix to see the effectiveness of the proposed method.

Numerical result

In this section, we first consider the case of a network as depicted in Figure 5.1. The Laplacian matrix is defined as:

L =             3 -1 0 0 -1 -1 -1 4 0 -1 -1 -1 0 0 3 -1 -1 -1 0 -1 -1 2 0 0 -1 -1 -1 0 4 -1 -1 -1 -1 0 -1 4            
The corresponding Laplacian spectrum of this graph is sp(L) = {0, 1.7857, 3, 4.5392, 5, 5.6751}.

We consider two cases: perfect Laplacian eigenvalues sp(L) and two imperfect Laplacian spectrum sp(L) to see how Algorithm 11 reacts.

At first, we run the Algorithms 7, 9 and 10 to obtain the perfect Laplacian spectrum sp(L) and the two arbitrary imperfect Laplacian spectrums sp 1 (L), sp 2 (L). Hence, we can deduce the eigenvalues of consensus matrix W corresponding to both imperfect cases:

sp(W) sp 1 (W) sp 2 (W) 1
Υ Υ Υ =             1
            Υ Υ Υ 1 =             1 1 1 1 1
            Υ Υ Υ 2 =             1 1 1 1 1
           
Next, we take the first N -1 consecutive measurements x i that can be expressed in matrix form:

X =             0.
            Û2 =             0.
            Ŵ1 =             0
            Ŵ2 =             0.
           
Define the estimated rounded Laplacian matrices L, L1 , L2 as:

L = L1 = L2 =             3 -1 0 0 -1 -1 -1 4 0 -1 -1 -1 0 0 3 -1 -1 -1 0 -1 -1 2 0 0 -1 -1 -1 0 4 -1 -1 -1 -1 0 -1 4            
We can see that, even with the inexact Laplacian spectrum (with a relevant tolerance of the error between the actual Laplacian spectrum ( sp(L))and the real Laplacian spectrum sp(L)), the proposed method can infer the network topology. Now, to achieve a more general evaluation of the proposed algorithm, we generate 100 random graphs with different number of nodes. We assess the quality of network reconstruction for different values of Laplacian eigenvalues estimation defined as

Error = sp(L) -sp(L) .
The performance is evaluated using the Success rate of existing links (SREL) and the Success rate of non-existing links (SRNL), [START_REF] Ren | Noise bridges dynamical correlation and topology in coupled oscillator networks[END_REF]. SREL(SRNL)is defined as the ratio of the number successfully predicted existent (non-existent) links to the total number of existent (nonexistent) links. As can be seen, whatever the size of the graph, SREL and SRNL decrease with the quality of estimation of Laplacian eigenvalues. In addition, for a given estimation error, the performance also decreases when the size of the graph grows. This is due to the quality of estimation of the eigenvectors. Indeed, by increasing the size of the graph, the rows at the bottom of matrix Υ Υ Υ becomes smaller and smaller and leads to loose the numerical rank. Therefore the quality of estimation can be significantly lost. It is then necessary to resort to more efficient procedure to solve the eigenvector estimation problem for large graphs.

Conclusion

In this chapter, we have proposed a way for reconstructing the network topology in the presence of anonymous nodes from average consensus measurements. The consensus matrix is determined from its eigenstructure estimated in distributed way. Since the Laplacian spectrum can be obtained by means of the algorithms proposed in the previous chapter, each node can deduce the matrix Û up to rows permutation because of anonymous nodes. Therefore, the obtained graph is isomorphic to the original graph. In the case that the graph is node-observable or neighborhood observable from node j, if all the entries of x(0) are distinct then node j can exactly reconstruct the network topology.

Algorithm 11 works with the assumption that all eigenvalues are distinct. Therefore, future works encompass the design of network reconstruction protocols that deal with spectrum in which the multiplicities of the corresponding eigenvalues can be higher than 1. In addition, numerical issues for large graphs are to be considered for making the proposed method scalable.

Chapter 6

General conclusions and future works To conclude this dissertation, we give a brief summary and discussion on future works.

Review of contributions and conclusions

Finite-time consensus problems have received a vast attention from the research community due to their benefits in many application domains, especially in distributed computation. Briefly speaking, the finite-time average consensus algorithms (for which consensus can be obtained in a finite number of steps) are commonly used as building blocks for more sophisticated distributed control, estimation or inference algorithms. Due to the fact that they can guarantee a minimal execution time, they are much more appealing than those ensuring asymptotic convergence. This thesis impacts with the design of finite-time average consensus protocols regardless of the Laplacian matrix of a given graph in a distributed way and its application to network robustness assessment.

First, we proposed a new distributed algorithm for self-configuration of finite-time average consensus where matrices are not necessarily based on the Laplacian matrix of a given graph. More precisely, we solved a matrix factorization problem in a distributed way by using a learning sequence. The method is based on the back-propagation method and an associated gradient descent method. Conceptually, each linear iteration k, k = 0, . . . , D is illustrated as a layer of a network. The main idea is to compute the partial derivatives of the error between the actual output and desired output, and then propagate them back to each iteration k in order to update the weights.

Obviously, this method is impractical for large-scale networks. However, one interesting point can be noted in this work is that the obtained consensus matrices are not doubly stochastic matrices and allow to reach consensus in a minimum number of steps.

Second, the study on the robustness of networks has increasingly attracted the attention of the research community. Therefore, in this thesis, the assessment of the network robustness via the effective graph resistance (R) and the number of spanning trees (ξ) in a distributed way, obviously becomes a promising application domain. Since both robustness metrics are functions of the Laplacian spectrum sp(L), thus the main task is now to estimate the Laplacian spectrum. More precisely, three cases have been considered: perfectly known x, approximated, and unknown x.

The average consensus matrix can be factored in D Laplacian based consensus matrices, where D stands for the number of nonzero distinct Laplacian eigenvalues. For the first case (perfectly known x), all methods are characterized into two approaches, which are direct approach and indirect approach. The direct method gives rise to a non-convex optimization problem with respect to the stepsizes α k , k = 1, . . . , N -1 that are the inverse of the nonzero distinct Laplacian eigenvalues. Since, this is a gradient based method, it suffers the slow convergence. To improve that, a convex optimization problem is introduced to estimate the coefficients c k , k = 0, . . . , N -1, of a given polynomial whose roots are the stepsizes α k , k = 1, . . . , N -1, that allows to factorize the averaging matrix J N = 11 T N . In this way, two famous methods (the subgradient projected method and ADMM) are applied to solve this convex optimization problem. By means of these proposed methods, the distinct Laplacian eigenvalues can be deduced. Then, by means of an integer programming, the whole spectrum sp(L) is obtained to compute the robustness indices.

The drawback of these methods is scalability. For large-scale graphs, the whole spectrum sp(L) estimation is significantly difficult to obtain. Therefore, the estimation can be restricted to that of the most significant eigenvalues to estimate some bounds of the network robustness indices.

Projected subgradient method is very applicable due to its simplicity. However, the drawback is that it can be slow. Therefore, the choice of the stepsize γ(t) is typical critical. In contrast, under mild conditions, the ADMM is guaranteed to converge for all tunning parameters. It has a single penalty parameter ρ that influences on the speed of convergence. Therefore, in the left cases, ADMM is to be the best choice.

In the case of imperfect consensus value x, we assumed that a standard consensus protocol has been run and stored some intermediate values of xM according to the iteration M.

Then, we made use of the ADMM-based Algorithm 7 to evaluate the network robustness. Finally, it claims that when the intermediate value of x is as close as possible to the real one, the network robustness index can be estimated with a very good accuracy.

In the case of unknown consensus value x, we considered a convex optimization problem which resorts to a convex combination of two convex functions. The first one takes into account the average consensus problem while the second one concerns the estimation of the Laplacian eigenvalues. For large graphs, numerical problems can appear. Therefore, the current formulations of the proposed Algorithms 5, 6, 7, 8 are efficient and well adapted for networks of medium size.

Finally, network topology identification refers to detecting and identifying the interested network elements and the relationship between elements in target network and represents the topology construction in an appropriate form. Therefore, identification of networks of systems becomes a growing attractive task for solving many problems in different science and engineering domains.

Theoretically, the consensus matrix W admits the eigenvalues decomposition W = UDU T , where D = I N -ǫ∆ ∆ ∆, with ∆ ∆ ∆ be diagonal matrix of λ i (L), stands for the diagonal matrix of the eigenvalues and U the matrix of eigenvectors. Since Laplacian spectrum sp(L) can be achieved by means of proposed methods in Chapter 4, we proposed a method to reconstruct the network topology via the estimation of eigenvectors of the matrix W.

Ongoing and future works

In what follows, we state some issues that are still open and worth of investigations:

1. Design of finite-time average consensus protocols:

(a) The dependence of the convergence speed on the learning sequence is to be studied. Thus, the design of optimal learning sequences can be dig in.

(b) Researching on the robustness of the finite-time average consensus.

Network robustness assessment:

(a) Novel methods that can be applied to large-scale networks.

(b) For indirect method, in practice, it may encounter with the ill-conditioned problem when establishing submatrices Q i , which leads to loose the numerical rank of these submatrices. Therefore, one can investigate in the preconditioning of these matrices to keep the accuracy of the proposed methods.

• Vues locales: Aucun agent a une vue globale du système, ou le système est trop complexe pour un agent de faire un usage pratique de ces connaissances.

Un MAS peut traiter des tâches difficiles, voire impossibles à accomplir par un agent individuel. Au cours des dernières décennies, les MASs gagnent un intérêt généralisé dans de nombreuses disciplines telles que les mathématiques, la physique, la biologie, l'informatique, les sciences sociales. Un nombre croissant de sujets de recherche des MASs comprend la coopération et la coordination, calcul distribué, contrôle automatique, les réseaux de communication sans fil, etc.

Dans la structure centralisée typique , un centre de fusion (FC) collecte toutes les mesures des agents, puis effectue les calculs finaux. Néanmoins, en raison de l'information à haut débit au FC, la congestion peut être créé. Une telle structure est vulnérable à la défaillance du FC. En outre, les besoins en matériel pour construire les communications sans fil peut engendrer d'une augmentation du coût des dispositifs et, par conséquent, un coût global du réseau plus élevé . Pour ces raisons, une structure centralisée peut devenir inefficace. Ainsi, la tendance de la recherche de la masse déplacée vers MASs décentralisées où l'interaction entre agents est implémenté localement sans connaissance globale. Un bon exemple est les réseaux de capteurs sans fil (WSN), qui trouvent de larges domaines d'application telles que les applications militaires (la surveillance de bataille, le monitering forces amies, de l'équipement et des munitions, etc.), les applications de l'environnement (la détection des feux de forêt, la détection de la nourriture, etc.), les applications de santé (de télésurveillance des données physiologiques humaines, etc. ), domotiques, le contrôle de la formation, et cetera. Figure 2 dépeint un WSN qui collecte des données concernant la qualité de l'air, l'intensité lumineuse, le volume du son, la chaleur, les précipitations et le vent. La première tâche est dédiée à la conception des protocoles de consensus de moyenne en temps fini, tandis que la seconde tâche se concentre sur l'évaluation de la robustesse du réseau. Nous présentons les grandes lignes de cette thèse comme suivant:

• Chapitre 2: Théorie des graphes.

Dans ce chapitre, nous donnons des définitions générales sur la théorie des graphes.

• Chapitre 3: Conception distribuée des Protocoles de Consensus de Moyenne en temps fini.

Ce chapitre est une contribution principale de cette thèse, où la conception du consensus de moyenne en temps fini est résolu de manière distribuée. Le problème est formulé comme une factorisation de la matrice moyenne J N = 11 T N . La solution est obtenue en résolvant un problème d'optimisation sous contrainte d'une manière distribuée. Pour la simple compréhension, chaque itération linéaire k, k = 0, . . . , D avec D soit le nombre d'étapes est illustrée comme une couche d'un réseau (y compris D + 1 couches). L'idée principale est de calculer les dérivées partielles de l'erreur entre la sortie réelle et la sortie désirée, puis propager les ramener à chaque itération k pour mettre à jour les poids. Ensuite, l'algorithme proposé est basé sur la méthode de 'back-propagation' et une méthode de descente de gradient associée. L'évaluation de la méthode proposée est évaluée par les résultats numériques.

• Chapitre 4: Évaluation distribuée de la robustesse du réseau.

Dans ce chapitre, nous utilisons des mesures qui sont le nombre d'arbres couvrants ξ et la résistance effective du graphe R pour évaluer la robustesse d'un réseau. prise distribué sur certains paramètres. Par exemple, si le protocole a les poids de degré maximum ( 1.2.6), chaque agent calcule d'abord le nombre de ses voisins avant d'exécuter un algorithme de Max-consensus pour le calcul du degré maximum de le graphe. Dans le cas du protocole basé sur Metropolis-Hasting ( 1.2.7), chaque agent compare son degré avec celle de ses voisins afin de calculer les poids du protocole de consensus de moyenne. Un protocole couramment est le poids des bords constantes, ou un protocole de consensus moyen basé sur le graphe de Laplace ( 1.2.8), où un stepsize commun est utilisé par tous les agents. Convergence asymptotique est donc garantie si le stepsize est strictement positif et inférieur à 2 λ N , où λ N représente la plus grande valeur propre de Laplace. Même à travers, il y a quelques bornes simples qui donnent des choix pour le stepsize sans nécessiter une connaissance exacte du spectre de Laplace, les agents doivent s'accorder dans un stepsize adéquate.

Au meilleur de notre connaissance, il n'y a pas de travail publié traitant avec protocoles d'auto-configuration pour le protocole de consensus de moyenne sur la base des poids des bords constantes. Dans les systèmes complexes réels, le temps d'exécution devient de plus en plus impactant. Conséquent, le but est maintenant de concevoir un algorithme de consensus de moyenne en temps fini permettant à tous les noeuds d'atteindre la valeur moyenne du consensus en un nombre fini d'étapes D pour les protocoles d'autoconfiguration, c'est-à-dire, Les matrices de poids W(k), k = 1, . . . , D qui ont les entrées w k ij sont correspondant aux liens entre les noeuds de deux couches consécutives à étape temporelle k.

Remark 13

Les matrices de poids W(k) ∈ S G est associées à un graphe G. Par conséquent, les Dans la littérature, il y a beaucoup de mesures de robustesse proposé (voir la Section 1.3), [START_REF] Ellens | Graph measures and network robustness. e-version[END_REF][START_REF] Ellens | Effective graph resistance[END_REF]. Par exemple, la première mesure est basée sur la connexité du graphe représentant le réseau [START_REF] Dekker | Network robustness and graph topology[END_REF]. Par les autres mots, le réseau est robuste s'il existe la présence des chemins alternatif, qui assurent la possibilité de communication en dépit des dommages [START_REF] Dekker | Network robustness and graph topology[END_REF]. Il y a deux concepts classiques de connectivité pour un graphe qui peut être utilisé pour modéliser la robustesse du réseau qui sont les sommet et bord connectivités. Le sommet (bord) connectivité K v (K e ) est le nombre minimal de sommets (de bords) d'être retiré afin de déconnecter le graphe donné (voir Subsection 1.3.1). Par cette approche analytique, afin d'en déduire l'indice de la robustesse, nous devons enlever chaque bord pour vérifier la connectivité du graphe. Même si cette approche est simple, c'est aussi un fardeau si nous analysons des réseaux plus grands et plus complexes raisonnables.

Dans [START_REF] Wu | Spectral measure of structural robustness in complex networks[END_REF], les auteurs ont proposé le concept de la connectivité naturelle comme une mesure spectrale de robustesse dans les réseaux complexes. La connectivité naturel est exprimée sous forme mathématique comme la moyenne des valeurs propres de la matrice d'adjacence du graphe représentant la topologie du réseau. L'avantage de cette méthode est que la mesure proposée fonctionne dans les deux réseaux connectés et déconnectés. D'autre part, le spectre de Laplace sp(L) peut être utilisé pour calculer les indices de robustesse qui évalue la performance du réseau donné. Par exemple, la deuxième plus petite valeur propre du spectre de Laplace λ 2 (L), défini également comme la connectivité algébrique [START_REF] Friedler | Algebraic connectivity of graphs[END_REF], est bien connu comme un paramètre critique que influence sur la performance et robustesse des systèmes dynamiques en raison de son rôle principal dans la connectivité du graphe [START_REF] Ellens | Graph measures and network robustness. e-version[END_REF]. Par conséquent, il y a beaucoup de recherches sur cette connectivité algébrique dans la littérature [3,[START_REF] Merris | Laplacian matrices of a graph: a survey[END_REF][START_REF] Yang | Decentralized estimation and control of graph connectivity in mobile sensor networks[END_REF]. En outre, il existe des mesures remarquables pour obtenir les indices de robustesse souhaités: la résistance effective du graphe R et Nombre d'arbres couvrants ξ (voir Section 1.3), qui sont basées sur la tout le spectre de Laplace différents de zéro.

Nous proposons des méthodes pour estimer la résistance effective du graphe R et Nombre d'arbres couvrants ξ. Figure 5 montre que l'estimation distribué de la robustesse du cette méthode est qu'elle garantie de converger pour tous les paramètres [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF] tandis que d'autres techniques telles que la projection de sous-gradient ou double décomposition dépend du choix du stepsize pour les mises à jour de gradient. Contrairement à d'autres méthodes, ADMM a seulement d'un seul paramètre de pénalité ρ qui peut être influencée sur la vitesse de convergence. ADMM assure une très bonne vitesse de convergence lorsque ses paramètres sont choisis de manière appropriée. Il existe quelques travaux traitant de la sélection des paramètres de pénalité pour accélérer la vitesse de ADMM [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Ghadimi | On the optimal step-size selection for the alternating direction method of multipliers[END_REF][START_REF] Teixeira | Optimal scaling of the ADMM algorithm for distributed quadratic programming[END_REF]. L'analyse de convergence est également étudié dans [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF][START_REF] Erseghe | Fast consensus by the alternating direction multipliers method[END_REF]. Cependant, dans cette étude, notre paramètre de pénalité est une constante fixe.

L'algorithme de ADMM actes en introduisant certaines variables auxiliaires z ij , de sorte qu'un problème d'optimisation strictement équivalente peut être formulé comme suit: subject to c i = z ij , j ∈ N i ,

z ji = z ij , z ij = sign(k)z k ij , c i ∈ C i = {c ∈ IR (h+1) : Q i c i = x1}, i = 1, 2, • • • , N.

B.5.1 Estimation des vecteurs propres de la matrice de consensus basée matrice de Laplace

Au noeud j, l'état à l'instant k est alors donnée par:

x j (k) = N i=1 λ k i (W)β i u j,i ,
En empilant N mesures consécutives noeud j obtient:

      
x j (0)

x j (1) . . . Maintenant, le noeud ne connaît que le j-ere rangée de Ũ.

x j (N -1)        =        1 1 • • • 1 1 λ 2 (W) • • • λ N (W) . . . . . . . . . . . . 1 λ N -1 2 (W) • • • λ N -1 N (W)               ũj,
Considérons un message passant schéma où les noeuds envoient à leurs voisins les rangées qu'ils ont.

Puis, un certain nombre d'échanges des messages, au moins égale au diamètre du graphe, noeud j a toutes les rangées de Ũ. Cependant, étant donné que certains des noeuds sont supposés être anonymes, ils ne envoient pas le label de la rangée qu'ils ont. Deux types de messages sont transmis dans le réseau. les rangées associées à des noeuds non-anonymes sont correctement étiquetés tandis que ceux qui sont associés à des noeuds anonymes ne sont pas étiquetés. En conséquence, noeud j ne reçoit Ũ allant jusqu'à rangées permutation, i.e. il reçoit les Û = Π Π Π Ũ, Π Π Π étant une matrice de permutation. L'inconvénient de ces méthodes est l'evolutivité. Pour les graphes de grande échelle, l'estimation de tout le spectre sp(L) est significativement difficile à exécuter. Par conséquent, l'estimation peut être limitée à celle des valeurs propres les plus significatifs pour estimer la limite de l'indice de la robustesse du réseau (λ 2 (L)).

Méthode de sous-gradient projeté est très applicable en raison de sa simplicité. Cependant, l'inconvénient est qu'il peut être lente. Conséquent, il faut faire attention au choix du stepsize γ(t). En revanche, l'ADMM est garanti à converger pour tous les paramètres. Il a seulement un seul paramètre de pénalité ρ qui peut être influencée sur la vitesse de convergence. Ainsi, dans les cas restants, ADMM est le meilleur choix.

Dans le cas de la valeur de consensus imparfaite x,nous avons supposé un protocole standard de consensus a été exécuté et conservé certaines valeurs intermédiaires de xM selon l'itération M.Ensuite, nous avons utilisé la basée algorithme basée sur ADMM 7 de voir comment la robustesse du réseau réagit. Enfin, il prétend que quand la valeur intermédiaire de x est aussi proche que possible de la valeur réelle, l'indice de la robustesse du réseau peut être pris environ.

Dans le cas de la valeur de consensus inconnue x, nous avons considéré le problème d'optimisation convexe qui recourt à une combinaison convexe de deux fonctions convexes. Le premier prend en compte le problème du consensus de moyenne, tandis que le second concerne l'estimation des valeurs propres du laplacien. Pour les grands graphes, une telle solution est évidemment insolubles. Conséquent, les Algorithmes 5, 6, 7, 8 sont efficaces et bien adaptées pour les réseaux de taille moyenn.

L'identification de la topologie du réseau se réfère à détecter et identifier les éléments du réseau intéressés et la relation entre les éléments du réseau destinataire et représente la construction de la topologie sous une forme appropriée. Conséquent, l'identification des réseaux de systèmes devient une tâche de plus en plus attrayante pour résoudre de nombreux problèmes dans différents domaines scientifiques et techniques.

En théorie, la matrice de consensus W admet la décomposition des valeurs propres W = UΛ Λ ΛU T , où Λ Λ Λ = I N -ǫ∆ ∆ ∆, avec ∆ ∆ ∆ est matrice diagonale de λ i (L), représente la matrice diagonale des valeurs propres et U la matrice de vecteurs propres. Nous

4 . 1 1 . 1 1 . 4

 411114 The i th element of the vector x The achievement of the proposed algorithms with respect to the variation of β . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 4.2 The achievement of the proposed algorithms with respect to the variation of β . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129List of Figures 1.1 General architecture of a MAS. . . . . . . . . . . . . . . . . . . . . . . . 1.2 A network of wireless sensors on the light poles all over the city of Chicago. (Source:http://mostepicstuff.com/chicago-city-installing-smartsensors-city-wide-to-monitor-everything) . . . . . . . . . . . . . . . . . . 1.3 Average consensus in a network: initial condition (left) and steady state (right). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Classification of Consensus protocols. . . . . . . . . . . . . . . . . . . . . 1.5 Scheme of the thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1 6-node undirected and 5-node directed networks. . . . . . . . . . . . . . . 2.2 A connected undirected graph. . . . . . . . . . . . . . . . . . . . . . . . . 2.3 A disconnected graph. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4 A complete graph. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.5 A weighted graph. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.6 A 3-regular graph . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.7 A distance-regular graph . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.8 A 5-regular Clebsch graph (SRG(16, 5, 0, 2)) . . . . . . . . . . . . . . . . 3.1 A H(4, 2) Hamming graph. . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Linear iteration scheme in space and time. . . . . . . . . . . . . . . . . . 3.3 Mechanism of the proposed back-propagation based method. . . . . . . . 3.4 6-node graph . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.5 MSE comparison for the two possible numbers of factors . . . . . . . . . Context of the thesis . . . . . . . . . . . . . . . . . . . . . . . 1.1.1 Multi-agent systems . . . . . . . . . . . . . . . . . . . . . . . . 1.1.2 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.3 Goals of the thesis . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 The consensus problem . . . . . . . . . . . . . . . . . . . . . . 1.2.1 Preliminaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.2 Consensus algorithms . . . . . . . . . . . . . . . . . . . . . . . 1.2.3 Convergence conditions . . . . . . . . . . . . . . . . . . . . . . 1.2.4 Design of the weight matrix. . . . . . . . . . . . . . . . . . . . 1.3 Network robustness . . . . . . . . . . . . . . . . . . . . . . . . 1.3.1 Vertex (Edge) connectivity . . . . . . . . . . . . . . . . . . . . 1.3.2 Algebraic connectivity . . . . . . . . . . . . . . . . . . . . . . . 1.3.3 Number of spanning trees ξ . . . . . . . . . . . . . . . . . . . . 1.3.4 Effective graph resistance R . . . . . . . . . . . . . . . . . . . Outline of this Thesis . . . . . . . . . . . . . . . . . . . . . . . 1.5 Publications list . . . . . . . . . . . . . . . . . . . . . . . . . . 1.5.1 International conference papers with proceedings . . . . . . . . 1.5.2 Journals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Figure 1 . 1 :

 11 Figure 1.1: General architecture of a MAS.

Figure 1 . 2 :

 12 Figure 1.2: A network of wireless sensors on the light poles all over the city of Chicago. (Source:http://mostepicstuff.com/chicago-city-installing-smart-sensors-city-wide-to-monitoreverything)

FinishFigure 1 . 3 :

 13 Figure 1.3: Average consensus in a network: initial condition (left) and steady state (right).

Figure 1 . 4 :

 14 Figure 1.4: Classification of Consensus protocols.

Figure 1 . 5 :

 15 Figure 1.5: Scheme of the thesis

Figure 2 . 1 : 6 -

 216 Figure 2.1: 6-node undirected and 5-node directed networks.

Figure 2 .

 2 Figure 2.1 shows:

Figure 2 . 2 :

 22 Figure 2.2: A connected undirected graph.

Figure 2 . 3 :

 23 Figure 2.3: A disconnected graph.

Figure 2 . 4 :

 24 Figure 2.4: A complete graph.

m

  i = N with m 1 = 1, where d max is the maximum degree of the graph, i.e., d max = max i d i . We denote by Λ = {λ 2 , . . . , λ D+1 } the set of nonzero distinct Laplacian eigenvalues.

Figure 2 . 5 :

 25 Figure 2.5: A weighted graph.

Figure 2 . 6 :

 26 Figure 2.6: A 3-regular graph

Figure 2 . 7 :

 27 Figure 2.7: A distance-regular graph

Figure 2 . 8 :

 28 Figure 2.8: A 5-regular Clebsch graph (SRG(16, 5, 0, 2))

Figure 3 . 1 :

 31 Figure 3.1: A H(4, 2) Hamming graph.

Figure 3 . 2 :

 32 Figure 3.2: Linear iteration scheme in space and time.

  be the input-output pairs defining the learning sequences, with y i,p = y p = 1 N N i=1

. 6 )

 6 where δ δ δ D,p = x p (D)xp is the difference between the actual output and the desired output with xp = y p = 1 N N i=1

Figure 3 . 3 :

 33 Figure 3.3: Mechanism of the proposed back-propagation based method.

Figure 3 . 4 :

 34 Figure 3.4: 6-node graph

Figure 3 . 5 :

 35 Figure 3.5: MSE comparison for the two possible numbers of factors

Figure 3 . 6 :

 36 Figure 3.6: (a)Trajectory of the proposed finite-time consensus protocol (b)Trajectory of an asymptotic consensus protocol by using optimal constant edge weights

Figure 3 . 8 :

 38 Figure 3.8: Final MSE comparison for different values of D

Figure 3 . 9 :

 39 Figure 3.9: Trajectory of an arbitrary state vector for a 10-node network

Figure 4 . 1 :

 41 Figure 4.1: Diagram of distributed estimation of network robustness.

4. 3

 3 Distributed solutions to the network robustness 4.3.1 Distributed Laplacian eigenvalues estimation 4.3.1.1 Problem statement Consider a network modelled by a connected undirected graph

( 1 -( 1 -

 11 0) is not orthogonal to a Laplacian eigenvector, then w i = 0, ∀i. Therefore, E(α α α) vanishes if and only ifh k=1 α k λ i ) 2 = 0, i = 2, • • • , N .Due to the multiplicities of the Laplacian eigenvalues, we have in fact D distinct equations h k=1 α k λ i ) 2 = 0, i = 2, • • • , D + 1. These equations are fulfilled if and only if there are D coefficients α k respectively equal to the inverse of the D distinct Laplacian eigenvalues.

Figure 4 . 2 :

 42 Figure 4.2: Mechanism of network robustness estimation

4. 3 . 1 . 2

 312 Laplacian-based factorization of the averaging matrix with perfect knowledge of x Given an arbitrary input-output pair {x(0), x} and the intermediate observations q(k), k = 0, 1, . . . , h = N -1, we can minimize the cost function (4.3.3) to get the set of the stepsizes α α α = (α 1 , . . . , α N -1 ), whose inverses are the Laplacian eigenvalues {λ 1 , . . . , λ N -1 }.

Figure 4 1 Figure 4 . 3 :

 4143 Figure 4.3: Several solutions of c in the case of a not-a-full-rank matrix Q.

  N -1, j = 1, . . . , D. m j ≥ 1, integer.

  and m ≥ 1}. The problem (4.3.61) without integer restriction is called linear programming (LP) relaxation. The branch-and-bound-based algorithm is described as: Algorithm 10 (Brand-and-bound-based multiplicities estimation) 1. Initialization: (a) Number of nodes N, d i , Λ. (b) t = 0; (c) S = ∅ set of integer solutions. (d) F best = ∅ function value of feasible solution. (e) V = ∅ set of nodes, that the program will visit. (f) Compute m[0] = argmin Λ T m, s.t m ∈ M • If the stopping criteria is satisfied, then stop. • Compute the upper bound (UB), which is the function value of the rounded-up feasible solution, F best = UB. • S 4 = m[0].(g) Order of variables in S to be branched on.2. While S 4 = ∅, select an ordered variable mj to be branched on and S 4 := S 4 \{ mj }.

(a) t = t + 1 ;

 1 (b) Create 2 new nodes for 2 new sub-problems.

Figure 4 . 4 :

 44 Figure 4.4: Two graphs with the same number of nodes but different number of edges.

17 .5 ξ 6 4. 4 . 1

 17641 Cases of perfect consensus value x 4.4.1.1 Example 1 (Figure 4.4 (a)) Consider the graph in Figure 4.4(a) with initial vector of state: x(0) = [0.5832, 0.7400, 0.2348, 0.7350, 0.9706, 0.8669] T .

5 Figure 4 . 5 :

 545 Figure 4.5: Trajectory of the stepsizes α k .

Figure 4 .Figure 4 . 6 :

 446 Figure 4.6: Convergence of the Algorithm decribed via the cost function H(α).

Figure 4 .

 4 Figure 4.7 shows the trajectory of the coefficients c k obtained by means of Algorithm 5.

5 Figure 4 . 7 :Figure 4 . 8 :

 54748 Figure 4.7: Trajectory of the estimated coefficients c k at each node obtained by means of distributed projected subgradient method

Figures 4 .

 4 Figures 4.9 shows the trajectory of the coefficients c k obtained corresponding to the value of penalty parameter ρ = 0.003 by means of Algorithms 6 . The obtained coefficients are c = [1, -7.0065, 11.4787, -7.1633, 1.8568, -0.1658 T ] and then S 1 = {4.9732, 1, 0.5, 0.3333, 0.2}.

Figure 4 .

 4 10 illustrates the difference in convergence rate corresponding to different values of penalty parameter ρ

5 Figure 4 . 9 :

 549 Figure 4.9: Trajectory of the estimated coefficients c k at each node obtained by means of ADMM

Figure 4 . 10 :

 410 Figure 4.10: Convergence rate of the ADMM-based method with different values of ρ.

4. 4 . 1 . 2

 412 Example 2 (Figure 4.4 (b))

Figure 4 .

 4 Figure 4.11 shows the trajectory of each node on the stepsize α k , k = 1, . . . , 3.

3 Figure 4 . 11 :

 3411 Figure 4.11: Zoom in the trajectory of each node on the stepsize α k

Figure 4 .

 4 Figure 4.12 shows the convergence of Algorithm 4 in this example.

Figure 4 . 12 : 2 .Figure 4 .

 41224 Figure 4.12: Convergence of Algorithm 4 performed through the cost function H(α)

3 Figure 4 . 13 :

 3413 Figure 4.13: Trajectory of the estimated coefficients c k at each node obtained by means of distributed projected subgradient method

Figure 4 . 14 :

 414 Figure 4.14: Convergence speed of the distributed subgradient projected method via E(c)

Figure 4 . 3 Figure 4 . 15 :

 43415 Figure 4.15 shows the trajectory of the coefficients c k obtained by means of Algorithms 6.

Figure 4 . 16 :

 416 Figure 4.16: Convergence of the distributed ADMM-based method via E(c)

Figure 4 .

 4 Figure 4.17 depicts the trajectory of the state of each node.

Figure 4 . 17 :

 417 Figure 4.17: Trajectory of the network state during average consensus protocol.

Figure 4 . 18 :

 418 Figure 4.18: Performance of Laplacian eigenvalues estimation for different values of the number of iterations of the standard average consensus algorithm.

Figure 4 . 19 :

 419 Figure 4.19: Rate of nodes succeeded to compute the Laplacian eigenvalues.

Figure 4 .

 4 Figure 4.18 depicts the MSE on the estimation of the Laplacian eigenvalues and the corresponding standard deviation. As expected, we can note that the closer x M is to x = 0.6269, the higher is the precision on the computation of the Laplacian eigenvalues. In addition, the tight bounds around the MSE inform us about the agreement on the estimated Laplacian eigenvalues.

Figures 4 .

 4 Figures 4.20and 4.21 show the trajectories of each node on the network robustness indices (R, ξ) corresponding to the iterations of the standard consensus protocol (M), started from M = 23.

Figure 4 . 20 :

 420 Figure 4.20: Trajectory of each node on the effective graph resistance R

Figure 4 . 21 :

 421 Figure 4.21: Trajectory of each node on the number of spanning trees ξ From these two Figures, the nodes estimate different values of R and ξ until M = 35. Then, they estimate almost the same value of the Robustness Indices.

Figure 4 .

 4 Figure 4.22 describes the sensitiveness of the effective graph resistance (R) and the number spanning trees (ξ) to the MSE via the relative errors which are defined respectively as:

Figure 4 . 22 :

 422 Figure 4.22: Relative errors RE R and RE ξ versus MSE.

5 Figure 4 . 23 :

 5423 Figure 4.23: Estimation of the polynomial coefficients c k at M = 24.

Figure 4 . 24 :

 424 Figure 4.24: Trajectory of the network state during average consensus protocol.

3 Figure 4 . 25 :

 3425 Figure 4.25: Estimation of the polynomial coefficients for (M = 17).

Figure 4 . 26 :

 426 Figure 4.26: Performance of Laplacian eigenvalues estimation for different values of the number of iterations of the standard average consensus algorithm.

Figure 4 . 27 :

 427 Figure 4.27: Rate of the nodes succeed to compute the Laplacian eigenvalues.

Figure 4 . 28 :

 428 Figure 4.28: Trajectory of each node on the effective graph resistance R

Figure 4 . 29 :

 429 Figure 4.29: Trajectory of each node on the number of spanning trees ξ Furthermore, Figure 4.30 describes the comparison between two robustness metrics on the sensitiveness to the MSE.Figure 4.30 shows that the effective graph resistance R is more sensitive to the MSE than the number of spanning trees ξ.

Figure 4 .Figure 4 . 30 :

 4430 Figure 4.30: Relative Errors RE R and RE ξ versus MSE.

Figure 4 . 31 :

 431 Figure 4.31: Convergence of Algorithm 7 for distinct Laplacian estimation for M = {7, 17, 40}.

  5356, 0.1649, 0.8834, 0.6665, 0.8477, 0.7627] T .According to the proposed Algorithm 8, the nodes have to achieve agreement on the values of the polynomial coefficients c k and that of the average consensus value x simultaneously.

Figure 4 . 5 Figure 4 . 32 :

 45432 Figure 4.32 and Figure 4.33 depict the trajectories of the nodes for the estimation of the coefficients c i,k , k = 0, . . . , N -1, and that of the consensus value x. Note that, the average consensus value is obtained much more faster than the common values of the coefficients c k .

Figure 4 . 33 :

 433 Figure 4.33: Nodes trajectories converging to the average of the initial condition.

Figure 4 .

 4 Figure 4.34 describes the performance in terms of the MSE between the estimated values and the actual values for the polynomial coefficients c i and the average consensus values xi as:

Figure 4 . 34 :

 434 Figure 4.34: Mean square error between the estimated values and actual values with respect to c i (M SE c ) and xi (M SE x ).

Figure 4 .

 4 Figure 4.34 indicates the performances of MSE c , MSE x with respect to the variation of θ(t). As can be seen from the figure, when θ(t) is near to 0, the convergence speed of

Figure 4 .

 4 Figure 4.35 and Figure4.36 depict the trajectories of the nodes for the estimation of the coefficients c i,k , k = 0, . . . , N -1, and that of the consensus value x. All the nodes agree on the actual values of polynomial coefficients c i,k and average consensus value x, and convergence occurs in less than 100 iterations. Moreover, the average consensus value is obtained much more faster than the common values of the coefficients c k .

3 Figure 4 . 35 :Figure 4 . 36 :

 3435436 Figure 4.35: Nodes trajectories for the estimation of the coefficients c k .

Contents 5 . 1

 51 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 5.2 Literature review . . . . . . . . . . . . . . . . . . . . . . . . . 132 5.3 Problem statement . . . . . . . . . . . . . . . . . . . . . . . . 138 5.4 Distributed solution to the network topology reconstruction problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 5.4.1 Estimation of the eigenvectors of the Laplacian-based consensus matrix . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 5.4.2 Network topology reconstruction . . . . . . . . . . . . . . . . . 140 5.5 Numerical result . . . . . . . . . . . . . . . . . . . . . . . . . . 142 5.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 147

2 . 1 )

 21 , the variational equation can be obtained:ξ ξ ξ = (DF(x) -γL ⊗ DP(x)) ξ ξ ξ + η η η,where ξ ξ ξ = [ξ 1 , . . . , ξ N ] T denotes the deviation vector, DF(x) = diag (DF 1 (x 1 ), . . . , DF N (x N )) and DF i , DP are m × m Jacobian matrices of F i and P respectively. Denoting C ξ ξ T dt the dynamical correlation of oscillators, with [0, T ] is the time interval in which the time average is performed, we have:

.3. 2 ) 5 . 4

 254 Distributed solution to the network topology reconstruction problem 5.4.1 Estimation of the eigenvectors of the Laplacian-based consensus matrix

Figure 5 . 1 :

 51 Figure 5.1: An generic network with 6 nodes.

Figures 5 .Figure 5 . 2 :

 552 Figures 5.2 and 5.3 describe the average SREL and SRNL for graphs of different sizes with respect to the norm of the estimation error of Laplacian eigenvalues.

Figure 5 . 3 :

 53 Figure 5.3: Success rate of non-existing links for 4 graphs with different sizes on the different values of norm of estimation error of Laplacian eigenvalues.

Contents

6. 1

 1 Review of contributions and conclusions . . . . . . . . . . . 149 6.2 Ongoing and future works . . . . . . . . . . . . . . . . . . . . 151

Figure 2 :

 2 Figure 2: Un réseau de capteurs sans fil sur les lampadaires dans toute la ville. (Source:http://mostepicstuff.com/chicago-city-installing-smart-sensors-city-wide-to-monitoreverything)

Figure 3 :

 3 Figure 3: Schéma de la thèse

x(D) = 1 N 11 T

 111 x(0) = J N x(0). (B.2.5) B.2.2 La solution distribuée du problème de la factorisation de la matrice Pour un graphe connexe G(V, E), soit x(0) = [x 1 (0), x 2 (0), . . . , x N (0)] T le vecteur d'état du reseau N-agent associé à un graphe G. Comme présenté dans le chapitre 1.2, schéma d'itération linéaire est un outil utile qui permet au système d'atteindre le consensus de moyenne en temps fini dans D étapes. Par conséquent, conformément à chaque itération, les vecteurs d'état mises à jour x(k), k = 1, . . . , D -1 sont assignés à être des couches consécutives respectivement. Enfin, la dernière itération est dédiée à des valeurs output x(D).

Figure 4 :

 4 Figure 4: Mécanisme de la méthode proposée basée sur le concept de rétropropagation.

B. 3 L

 3 'évaluation distribuée de la robustesse du réseauB.3.1 IntroductionL'efficacité d'un réseau est évaluée par sa fonctionnalité et sa robustesse. Selon la propriété plus tard, certaines questions peuvent se poser. Par exemple, s'il y a un événement accidentel, comment réagit-il le réseau? Sera le réseau survécu? Afin de répondre à ces questions, l'étude sur la robustesse des réseaux de plus en plus attiré une attention de la communauté scientifique.

2 N

 2 min c i ∈IR (h+1) ,i=1,...,N 1 i=1 j∈N i c ic j 2 . (B.3.22)

- 1 m

 1 if k is odd . La fonction Lagrangien augmenté correspondante est donnée par:L ρ (c, z, y) = N i=1 L (i) ρ (c i , z, y),avecL (i) ρ (c i , z, y) = j∈N i 1 2 c ic j 2 + y T ij (c iz ij ) + ρ 2 c iz ij 2 , (B.3.23)avec ρ-un coefficient de pénalité. Contrairement à la méthode des multiplicateurs de Lagrange, dans ADMM la minimisation de la fonction de coût selon les variables primaires c i et z ij est effectuée de façon séquentielle. Donc, étant donné les variables primales et les variables duales y ij , la fonction de coût est d'abord minimisé selon à c i et la solution obtenue projetée sur l'ensemble des contraintes (c-minimisation). Ensuite z ij sont obtenus à partir des variables duales données et les valeurs mises à jour de c i (zminimisation). La troisième étape correspond à la mise à jour des variables duales de celles primaires mises à jour (y-minimisation):c i [t + 1] = Ω C i [argmin L ρ (ĉ i [t + 1], z ij [t], y ij [t])]; (B.3.24)enlevé est également un élément de S 2 puis xi = Pi ( c, q i ) où Pi est le polynôme de degré réduit en ce qui concern noeud i. Algorithme 9 décrit la procédure proposée.Ensuite, nous obtenons l'ensemble des valeurs propres du laplacien distincts différents de zéro Λ = {λ 2 , . . . , λ D+1 }. Afin d'obtenir tout le spectreSp(L) = {λ m 1 1 , λ m 2 2 , . . . , λ m D+1D+1 }, les multiplicités m i devraient être définis. Puisque, elles sont toutes entier, le problème a recours actuellement à un problème de la programmation en nombres entiers.Étant donné l'ensemble des valeurs propres du laplacien distinctes différents de zéro Λ,d i le degré de chaque noeud i, inspiré d'une propriété intéressante N i=1 λ i (L) = N i=1 d i (voirdans Section 2.3), l'optimisation est maintenant décrit comme suit: j = N -1, j = 1, . . . , D. m j ≥ 1, integer. où D = |Λ| et m est le vecteur des multiplicités des valeurs propres. Dans la littérature, la méthode de séparation et d'évaluation progressive est la technique basique pour résoudre le problème de programmation entier. La programmation linéaire (LP) de la relaxation (sans restriction de nombre entier) est adoptée pour estimer la solution optimale d'une programmation en nombres entiers [13]. Conséquent, le problème (B.3.53) est maintenant réécrit strictement équivalente sous la forme de programmation linéaire comme suit: m j ≥ 1, integer. Où M est l'ensemble des contraintes sans restriction en nombres entiers. De ce fait, M peut être exprimée comme M = {m ∈ IR D | D j=1 m j = N -1, Λ T m = N i=1 d i and m ≥ 1}. Le problème (B.3.54) sans restriction en nombres entiers est appelé la programmation linéaire (LP) de la relaxation. La Algorithme de séparation et d'évaluation progressive est décrit comme dans l'algorithme 10.

,

  ũj,i = u j,i β i ou de façon équivalentex j = Υ Υ Υ ŨT j. (B.5.3) De Hypothèse 1, nous pouvons conclure que Υ Υ Υ est une matrice de rang plein. Par conséquent, noeud j peut résoudre le problème ci-dessus efficacement et obtenir: ŨT j. = Υ Υ Υ -1 x j . (B.5.4)

B. 5 . 2

 52 Network topology reconstructionNous pouvons noter que ÛT Û = diag(b 2 ). Par conséquent, de (B.5.2), nous obtenons:W = Π Π Π ŴΠ Π Π T , with Ŵ = ÛD( ÛT Û) -1 ÛT .(B.5.5) différents de zéro. Puisque, c'est une méthode basée sur gradient, il souffre de convergence lent. Pour améliorer la vitesse de convergence, un problème d'optimisation convexe est introduit pour estimer les coefficients c k , k = 0, . . . , N -1 d'un polynôme donné dont les racines sont les stepsizes α k , k = 1, . . . , N -1 qui permet de factoriser la matrice de moyennage J N = 11 T N . De cette façon, deux méthodes fameux (les méthode de sousgradient projeté et ADMM) sont appliquées pour résoudre ce problème d'optimisation convexe. Par ces méthodes proposées, les valeurs propres du laplacien distinctes peuvent être déduites. Ensuite, au moyen d'une programmation en nombres entiers, tout le spectre sp(L) est obtenu pour calculer les indices de la robustesse.

  2. A graph is called a regular graph if all the vertices have the same number of neighbors. The graph is written as n-regular if the number of neighboring vertices is n for all vertices. n is called valency of the n-regular graph.

  Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 3.2 Literature review . . . . . . . . . . . . . . . . . . . . . . . . . 47

	Chapter 3
	Distributed design of finite-time
	average consensus protocols
	Contents
	3.1 3.2.1 Minimal polynomial concept based approach . . . . . . . . . . 47
	3.2.2 Matrix factorization based approach . . . . . . . . . . . . . . . 49
	3.2.3 Comparison between the minimal polynomial approach and the
	matrix factorization approach . . . . . . . . . . . . . . . . . . . 55
	3.3 Distributed solution to the matrix factorization problem. . 56
	3.4 Numerical results . . . . . . . . . . . . . . . . . . . . . . . . . 64
	3.4.1 Example 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
	3.4.2 Example 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
	3.5 Conclusion and discussion . . . . . . . . . . . . . . . . . . . . 70

  each node sends a message m i,k [t] containing its current local value x i,k [t], and the local step-size α k [t].

Table 4 .

 4 1 describes the the achievement of the ( 4.3.46) with respect to different values of β. As can be seen from this Table, all of them is shown to obtain the nonzero distinct Laplacian eigenvalues.

Table 4 . 1 :

 41 The achievement of the proposed algorithms with respect to the variation of β

	β	0.1	0.01	0.06
	Estimated Laplacian eigenvalues Estimated coefficient c k	0.9995	{1, 2, 3, 5} 1	0.9998
		-2.8378 -3.0864 -2.7880
		3.0039 3.5080 2.9015
		-1.4674 -1.8059 -1.3984
		0.3287 0.4195 0.3102
		-0.0269 -0.0351 -0.0252
	Number of iterations	301	1170	259
	R ξ		14.2 15

Table 4 .

 4 [START_REF] Aderhold | Hierarchical Bayesian models in ecology: Reconstructing species interaction networks from non-homogeneous species abundance data[END_REF] shows the performance of the proposed method (Algorithm 8) with respect to different choice of the β.

Table 4 . 2 :

 42 The achievement of the proposed algorithms with respect to the variation of β

	β	0.1 0.01 0.06
	Estimated Laplacian eigenvalues Estimated coefficient c k	{1, 3, 4} 1
		-1.5833
		0.6667
		-0.0833
	Number of iterations	151 1188 232
	R ξ	17.5 6

  • α et β sont le taux d'apprentissage et le taux de l'élan respectivement. Ils peuvent être choisis comme constantes ou des paramètres variables. Dans notre étude, le taux d'apprentissage α[t] varie à chaque étape d'itération d'optimisation t tandis que le taux de l'élan est fixé. Pour éclairer le problème, Figure 4 montre le mécanisme de l'algorithme proposé. Dans le détail, pour chaque séquence d'apprentissage, l'algorithme exécute d'abord l'étape en avant, où chaque noeud met à jour à plusieurs reprises sa valeur. À la fin de cette étape, chaque couche stocke les valeurs mises à jour, puis les dérivées partielles de l'erreur δ D entre x(D) et les valeurs moyennes souhaitées de consensus x sont estimés. L'étape en arrière commence à propager retour les dérivés partiels estimés δ D à chaque couche jusqu'à ce que la couche d'entrée reçoit le δ 1 . Maintenant, tous les poids peuvent être mis à jour. Le régime se poursuit avec la séquence d'apprentissage suivante. Le mécanisme est répété jusqu'à ce que les critères d'erreur E satisfait le seuil θ, qui est un petit nombre.

λ 2 (L)+λ N (L)proposed in[START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF], we get the trajectories in Figure

3.6(b). Showing that more iterations

Source: http://whatis.techtarget.com/definition/network-topology

(c) Investigating in the convergence rate analysis based on the study in [START_REF] Ph | Explicit convergence rate of a distributed alternating direction method of multipliers[END_REF].

Reconstruction of the network topology:

(a) Design a new method that deals with a Laplacian spectrum sp(L), whose eigenvalues are not all distinct.

(b) Assess the robustness of these methods when the data exchange are prone to imperfect communications.

A. Khatri-Rao Product

Given two matrices A ∈ IR I×F , B ∈ IR J×F with the same number of columns, the Khatri-Rao product A ⊙ B ∈ IR IJ×F is defined as:

where a f is the f th column of A, similarly for b f , and ⊗ denotes the Kronecker product of two column vectors, whereas diag i (A) stands for the diagonal matrix built with the ith row of A. The Khatri-Rao product can be viewed as a column-wise Kronecker product. Its main property, which will be used herein, is related to the vectorization operation of a given matrix: vec(Adiag(d)B T ) = (B ⊙ A)d.

where the vec(.) operator stacks the columns of its matrix argument and d ∈ R F ×1 . Therefore we can deduce the following identity: Les systèmes multi-agents (MAS) ont reçu un intérêt croissant au cours des dernières décennies. Ils sont développés pour la demande des caractéristiques de flexibilité, robustesse, et reconfiguration qui apparaissent dans divers domaines d'applications, notamment la fabrication, la logistique, les réseaux électriques intelligents, l'automatisation des bâtiments, les opérations de secours en cas de catastrophe, les systèmes de transport intelligents, la surveillance, la surveillance de l'environnement et l'exploration des infrastructures de sécurité et protection, et cetera. Un MAS est un système composé de plusieurs agents interagissent intelligents (capteurs, les plantes, véhicules, robots, etc.) et de leur environnement, comme indiqué dans la Figure 1.

Un agent intelligent possède les caractéristiques clés suivantes:

• L'autonome: la capacité de travailler sans intervention directe, jusqu'à un point défini par l'utilisateur.

• La Réactivité: la capacité de percevoir son environnement et réagir à ses changements, que ce soit la modification des objectifs de l'utilisateur ou des ressources disponibles.

• Proactive: la capacité de prendre la responsabilité de fixer ses propres objectifs.

• La capacité sociale: la capacité de communiquer avec les humains et d'autres agents

B.1.1.2 Motivation

Il y a deux points qui attirent une grande attention de la communauté de recherche:

(A) Sur la base de l'information locale et les interactions entre les agents, comment peuvent tous les agents parvenir à un accord?

Ce problème est appelé problème du consensus, ce qui est pour concevoir un protocole de réseau basé sur les informations locales obtenus par chaque agent de telle sorte que tous les agents atteignent finalement un accord sur certaines quantités d'intérêt Problèmes de consensus ont reçu une attention considérable de diverses communitities de recherche en raison de ses larges applications dans de nombreux domaines, y compris multi-capteurs fusion de données [START_REF] Olfati-Saber | Consensus filters for sensor networks and distributed sensor fusion[END_REF], le flocage comportement des essaims [START_REF] Blondel | Convergence in multiagent coordination, consensus, and flocking[END_REF][START_REF] Jadbabaie | Coordination of groups of mobile autonomous agents using nearest neighbor rules[END_REF][START_REF] Olfati-Saber | Flocking for multi-agent dynamic systems: algorithms and theory[END_REF][START_REF] Olfati-Saber | Swarms on sphere: A programmable swarm with synchronous behaviors like oscillator networks[END_REF], multi-véhicule de contrôle de la formation [START_REF] Dimarogonas | Formation control and collision avoidance for multi-agent systems and a connection between formation infeasibility and flocking behavior[END_REF][START_REF] Moore | Formation control via distributed optimization of alignment error[END_REF][START_REF] Olfati-Saber | Distributed cooperative control of multiple vehicle formations using structural potential functions[END_REF], distribué calcul [START_REF] Bertsekas | Parallel and Distributed Computation: Numerical Methods[END_REF][START_REF] Borkar | Asymptotic agreement in distributed estimation[END_REF][START_REF] Lynch | Distributed Algorithms[END_REF], rendez-vous des problèmes [START_REF] Cortes | Robust rendezvous for mobile autonomous agents via proximity graphs in arbitrary dimensions[END_REF], etc.. Plus particulièrement, les algorithmes de consensus de moyenne (c'est à dire l'accord correspond à la moyenne des valeurs initiales) sont couramment utilisés comme la composante pour plusieurs algorithmes de contrôle, estimation ou inférence distribués.

Dans la littérature récente, on peut trouver des algorithmes de consensus moyens embarqués dans le filtre de Kalman distribué, [START_REF] Olfati-Saber | Distributed Kalman filtering for sensor networks[END_REF]; Distribué algorithme des moindres carrés, [START_REF] Bolognani | Distributed sensor calibration and least-squares parameter identification in WSNs using consensus al-gorithms[END_REF]; Distribué alternatif méthode des moindres carrés pour la factorisation des tenseurs, [START_REF] Kibangou | Distributed PARAFAC based DS-CDMA blind receiver for wireless sensor networks[END_REF]; Principal Distributed Component Analysis, [START_REF] Macua | Consensus-based distributed principal component analysis in wireless sensor networks[END_REF]; ou distribué entrée commune et estimation d'état, [START_REF] Esna-Ashari | Distributed input and state estimation for linear discrete-time systems[END_REF] pour citer quelques-uns. Néanmoins, la convergence asymptotique des algorithmes de consensus ne sont pas adaptés à ces types d'algorithmes distribués sophistiqués. La convergence faible asymptotique ne peut pas assurer l'efficacité et l'exactitude des algorithmes, ce qui peut conduire à d'autres effets inattendus. Par exemple, en ce qui concerne les réseaux de capteurs sans fil (WSNs), une réduction du nombre total d'itérations jusqu'à ce que la convergence peut conduire à une réduction de la quantité totale de consommation d'énergie du réseau, ce qui est indispensable pour garantir une durée de vie plus longue pour l'ensemble du réseau. D'autre part, les protocoles qui garantissent un temps d'exécution minimum sont plus attrayants que ceux assurant la convergence asymptotique. Dans ce but, plusieurs contributions dédiées à consensus en temps fini ont récemment été publiés dans la littérature. Ce qui signifie que le consensus est obtenu dans un nombre fini d'étapes.

(B) Pouvons-nous évaluer la robustesse de ce réseau?

En général, nous sommes entourés par des réseaux de différents types (réseaux sociaux, réseaux de capteurs, réseaux d'énergie, et cetera). Conséquent, la chose la plus importante est que ces réseaux sont robustes. Ce qui signifie que, ils peuvent subir des dommages ou des défaillances travers. Par conséquent, la terminologie robustesse envahit rapidement le communauté de la recherche.

Definition 2

La robustesse est la capacité d'un réseau de continuer à bien fonctionner quand il est sujette à des défaillances ou des attaques, [START_REF] Ellens | Effective graph resistance[END_REF].

Pour évaluer si le réseau est robuste, la mesure de la robustesse du réseau est nécessaire [START_REF] Ellens | Effective graph resistance[END_REF]. Il existe différentes approches proposées pour mesurer les paramètres de robustesse du réseau. Cependant, l'approche la plus populaire dans la littérature est basée sur l'analyse du graphe. Avec les données produites durant protocole de consensus pouvons-nous déduire des mesures de robustesse du réseau?

(C) Pouvons-nous reconstruit la topologie du réseau? L'identification de la topologie du réseau se réfère à détecter et identifier les éléments du réseau intéressés et la relation entre les éléments du réseau destinataire et représente la construction de la topologie sous une forme appropriée. Conséquent, l'identification des réseaux de systèmes devient une tâche de plus en plus attrayante pour résoudre de nombreux problèmes dans différents domaines scientifiques et techniques.

B.1.1.3 Objectifs de la thèse

Il est bien connu que les protocoles de consensus étant itérative, de nombreuses données sont produites. Cependant, la plupart des études se concentrent uniquement sur les propriétés de la convergence de ces algorithmes. Avec les données collectées auprès chaque itération du protocole de consensus asymptotique de manière distribuée, nous visons à:

• Concevoir des protocoles de consensus de moyenne en temps fini pour les accords subséquents.

• Évaluer la robustesse du réseau.

• Reconstruit la topologie du réseau.

B.1.2 Structure du document

Afin d'avoir une vue globale, les contributions de cette thèse peuvent être considérés triplé comme indiqué dans le schéma suivant.

Par la nature de ces mesures, elles sont toutes les fonctions des valeurs propres du laplacien λ i (L), par conséquent, l'estimation des valeurs propres du laplacien est le but principal de cette tâche. Dans ce chapitre, au lieu d'estimer directement des valeurs propres du laplacien, nous résolvons le problème comme un problème d'optimisation sous contrainte en ce qui concerne à l'inverse des valeurs propres du laplacien distinctes différents de zéro. Par conséquent, en appliquant la méthode de descente de gradient, nous obtenons finalement l'ensemble des valeurs propres du laplacien distinctes différents de zéro. En outre, en utilisant un algorithme de programmation en nombres entiers, on peut estimer les multiplicités correspondant à ces valeurs propres. Ensuite, on en déduit l'ensemble du spectre de la matrice de Laplace associée au graphe G(V, E). A partir de ce moment, les mesures de robustesse peuvent être estimés à évaluer la robustesse du réseau.

Pour la tâche d'estimation du spectre de Laplace, toutes les méthodes proposées dans ce chapitre sont divisées en deux catégories qui sont des formulations nonconvexes et convexes. Selon problème d'optimisation non convexe (appelle aussi méthode directe), le problème est vue comme un problème d'optimisation sous contrainte en ce qui concerne la stepsizes α k de protocole basé sur la Laplace matrice de consensus, dont les inverses sont égales aux valeurs propres du laplacien distinctes différents de zéro.

Résoudre le problème d'optimisation non convexe, le principal problème est la lenteur de la vitesse de convergence. Nous essayons donc d'accélérer la vitesse de l'estimation des valeurs propres du laplacien le plus vite possible. Le point principal est de tranformer le problème non convexe dans l'une convexe. Nous résolvons le problème en deux algorithmes utilisant la méthode de projection distribué, [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF] et la Direction alternatif de la Méthode de multiplicateurs (ADMM), [START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF]. Nous considérons trois cas: la valeur moyenne finale parfaitement connue, bruyant et totalement inconnu.

ou de façon équivalente sous forme de matrice:

W(k) ∈ S G où S G est l'ensemble des matrices qui peuvent être factorisées comme

Cependant, il faut noter que pour exécuter l'algorithme de consensus de moyenne, deux étapes principales sont nécessaires: l'étape de configuration (également appelé étape de la conception) et l'étape de l'exécution. Au cours de l'étape de configuration, une tâche peut être réalisée par un algorithme d'auto-configuration au lieu de recourir à un gestionnaire de réseau. L'auto-configuration peut inclure des graphes découvertes et la entrées w ij (k) sont conçus suivant la règle définie comme:

w k ij x j,p (k -1), w k ij étant les entrées de matrices W k ∈ S G . Nous pouvons réécrire la fonction de coût (B.2.6) comme

où . F supports pour la norme de Frobenius, Y = J N X(0), Y et X(0) étant N × P matrices avec y i,p and x i,p comme entrées, respectivement.

Nous supposons que X(0)X(0) T = I N ce qui signifie que le vecteur d'entrée est orthogonal. Par exemple, les vecteurs de la base canonique de IR N peuvent être utilisés comme 'inputs'. Par conséquent, nous pouvons noter que où tr(.) désigne l'opérateur de trace et y p = J N x p (0).

En dénotant E p (W) = 1 2 tr ǫ p (W)ǫ T p (W) , la solution de ce problème d'optimisation peut alors être obtenue par itération au moyen d'une méthode de descente de gradient:

ou une méthode de gradient stochastique:

où le gradient de la fonction coût est approximé par le gradient à une séquence d'entréesortie unique. Dans ce qui suit, nous avons recours à une méthode de gradient stochastique. Dans ce but, nous affirmons d'abord le lemme technique suivant: Lemma 6

Les dérivés de la fonction de coût E p (W) = 1 2 tr ǫ p (W)ǫ T p (W) avec ǫ p (W) défini dans la (B.2.9) peuvent être calculées comme suit:

où δ δ δ D,p = x p (D)xp est la différence entre la sortie réelle et la sortie désirée avec

En appliquant les résultats de Lemme 6, le schéma de mise à jour de l'algorithme d'optimisation est la suivante:

Ensuite, nous projetons les poids mis à jour sur S G qui donne:

Ce mécanisme est similaire à l'algorithme de rétro-propagation de gradient. Sa convergence a été bien étudié dans la littérature, voir [START_REF] Chong | An Introduction to Optimization[END_REF] et [START_REF] Mangasarian | Serial and parallel backpropagation convergence via nonmonotone perturbed minimization[END_REF] par exemple. Convergence vers un minimum local est garanti si le step-size α est choisi de manière appropriée dans l'intervalle (0 < α < 1). Si le step-size α est trop grande, l'erreur E tend à osciller, n'atteint jamais un minimum. En revanche, trop petite α empêche l'optimisation de faire des progrès raisonnables. Conséquent, le choix d'un appropriée α est aussi une difficulté. Une telle step-size influe également sur la vitesse de convergence de l'algorithme. Une astuce pour accélérer la convergence est d'ajouter un terme de régularisation dans la fonction de coût à minimiser.

{W * k }| k=1,...,D = arg min

En minimisant une fonction de coût, la mise à jour équation est donnée par:

Élément par élément, chaque agent met à jour les entrées de la façon suivante:

), où:

Protocoles standards de Consensus

L'estimation du spectre de Laplace:

L'estimation des valeurs propres distinctes de Laplace L'estimation des multiplicits des valeurs propres de Laplace

Calcul du index de la robustesse du rseau réseau est divisé en deux phases principales y compris le spectre de Laplace sp(L), la computation de la résistance effective du graphe et Nombre d'arbres couvrants (R, ξ). Tout d'abord, nous proposons des méthodes basées sur le consensus pour estimer tout le spectre de la matrice de Laplace de manière distribuée. Principalement, nous calculons la valeur moyenne de consensus x en utilisant les protocoles standards de consensus (Algorithme 1) ou algorithme décentralisée de calcul de valeur de consensus en temps minimum [START_REF] Yuan | Decentralized minimimtime consensus[END_REF]. Ensuite, l'estimation du spectre de Laplace peut être divisée en deux étapes qui sont l'estimation des valeurs propres du laplacien distinctes, puis l'évaluation des multiplicités correspondantes. Enfin, à partir de sp(L) estimé les indices de robustesse peuvent être calculées. Notre objectif est donc sur l'estimation du spectre de Laplace sp(L).

B.3.2 Les solutions distribuées de la robustesse du réseau B.3.2.1 L'estimation Distribuée des valeurs propres du laplacien B.3.2.2 Formulation du problème

Considérons un réseau modélisé par un graphe non orienté connecté G(V, E). Son état

T évolue comme suit:

où α est choisi tel que 0 < α < 1 dmax pour assurer la convergence asymptotique [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF]. De l'état x(k), on peut définir q(k) = L k x(0) comme:

avec q(0) = x(0).

1. Hypothèse 1: Le nombre N de noeuds est connue.

2. Hypothèse 2: Chaque noeud i sauve sa valeur initiale x i (0) et les valeurs consécutives q i,k , k = 0, . . . , h.

Le problème étude peut donc être formulé comme suit:

Étant donné une paire d'entrée-sortie arbitraire {x(0), x = x1}, où x = 11 T N x(0), et les observations intermédiaires q(k), k = 0, . . . , h, estimer de manière distribuée les valeurs propres du laplacien du graphe représentant le réseau.

Cependant, au lieu de calculer directement les valeurs propres du laplacien, nous allons calculer l'inverse des valeurs propres du laplacien. En outre, puisque le nombre D de valeurs propres du laplacien distinctes différents de zéro est inconnu, la factorisation sera ré-paramétré. Énonçons d'abord la proposition suivante: Proposition 7

Considérons un graphe connexe avec la matrice de Laplace L et ses valeurs propres distinctes différents de zéro λ 2 , • • • , λ D+1 . Supposons que le vecteur des valeurs locales x(0) n'est pas orthogonal à n'importe quel vecteur propre de L, la fonction de coût

De la Proposition 7, nous pouvons conclure que l'ensemble

Il est évident que α k sont les racines du polynôme avec des coefficients c k . Et, c a un motif de signe alternatif.

Nous pouvons affirmer la proposition suivante qui est un corollaire de Proposition 7: Proposition 8

Soit P(c, :) le polynôme de

Pour ne importe quel noeud i et des mesures intermédiaires q i = [q i,0 , q i,1 , . . . , q i,h ] T ∈ IR (h+1)×1 du protocole de consensus (B.3.1), nous obtenons Afin de mieux comprendre la structure de ce Chapitre, Figure 6 illustre sa construction en montrant la façon d'estimer la robustesse du réseau.

En fait, pour estimation distribué des valeurs propres du laplacien, trois cas ont été profondément considéré:

1. Le cas de la parfaite connaissance de la valeur moyenne du consensus x, Protocoles Standards de Consensus

Estimation des valeurs propres distinctes de Laplace:

Mthode Directe Mthodes indirectes

Estimation des multiplicits correspondantes des valeurs propres de Laplace Indices de la robustesse du rseau sp(L)

L'idée de la méthode proposée est de minimiser le désaccord entre voisins sur la valeur de α k tout en veillant que la factorisation de la matrice moyenne est obtenue. Une telle factorisation est évaluée en contraignant les valeurs des noeuds après h itérations de l'algorithme de consensus pour être égale à la moyenne des valeurs initiales:

Remark 14

Importe quelle solution de (B.3.9) contient l'inverse des valeurs propres du laplacien selon à Proposition 7.

Le problème d'optimisation sous contrainte (B.3.9) peut ensuite être réalisée comme un problème d'optimisation non contrainte au moyen d'une méthode de Lagrange avec une fonction de Lagrange définis comme suit:

où y ∈ IR N ×1 est le vecteur des multiplicateurs de Lagrange.

La minimisation de la fonction de Lagrange (B.3.10) peut être considéré comme h des problèmes de consensus à résoudre simultanément avec une contrainte qui donne une sorte de référence externe à suivre.

Pour éviter tout malentendu entre l'étape temporelle k et l'itération de l'optimisation t, nous notons x(k), W(k) comme x k , W k respectivement.

Afin de tirer l'algorithme distribué correspondant, nous énonçons maintenant le lemme technique suivant:

Lemma 7

Les dérivés de la fonction de coût H(α α α k , y) défini en (B.3.9) peut être calculée comme suit:

En appliquant les résultats de Lemme 7, le régime de la mise à jour de l'algorithme d'optimisation est la suivante:

L'algorithme correspondante est alors décrit comme dans l'algorithm 4.

Afin d'accélérer la vitesse de convergence, on peut re-paramétrer le problème (B.3.1) au moyen de l'équation (B.3.4).

x =

Au lieu de la minimisation les step-sizes α k , la nouvelle optimisation les coefficients du polynôme c k , k = 0, . . . , h peut être appliquée. Cela s'appelle factorisation distribué indirect de la matrice moyenne.

Méthode 2: La factorisation indirecte distribuée de la matrice moyenne.

Le but de cette méthode est de trouver le polynôme coefficents c k de manière distribuée en résolvant le problème: min

Toutes les conditions initiales x i (0) ne devrait pas être vecteurs propres de la matrice de Laplace L et ne pas être identiques, qui signifie que x(0) = δ1, δ ∈ IR, sinon toutes les colonnes de Q sont nuls sauf la première. Puisque la condition initiale est générée de façon aléatoire, presque sûrement, il ne peut pas être un vecteur propre de la matrice de Laplace L. Lemma 8

Considérons un graphe non orienté connecté G(V, E) avec le spectre de Laplace sp(L) = {λ 1 , λ 2 , . . . , λ N }. Ici, sp(L) a D + 1 valeurs propres du laplacien distinctes. Si x(0) n'est pas un vecteur propre de la matrice de Laplace L , alors la matrice

a rang min(D + 1, N). Où U est la matrice des vecteurs propres de Laplace.

Remark 15

rank(Q) < N si existe λ i ∈ sp(L) avec multiplicité supérieure à 1, ∀x(0).

D'une manière centralisée, afin de trouver c de l'équation Qc = x, Q doit être une matrice de rang plein . Ensuite,

En contraste, si Q n'est pas une matrice de rang plein, alors il y a plusieurs solutions de c. Pour faire face à cette situation, nous utilisons la solution des moindres carrés permettant de trouver la plus courante solution S 1 à l'intérieur de L'ensemble des solutions (So 1 , So 2 , . . . , So m ) de c comme décrit dans Figure 7. Afin de résoudre ce problème de manière distribuée, nous résolvons le problème localement selon à c i , i = 1, . . . , N. En désignant par c i la version locale de c et Q i ∈ IR (d i +1)×(h+1) la sous-matrice de Q constitué avec les lignes associées à noeud i et ses d i voisins, la solution du problème précédent peut être obtenue d'une manière distribuée en résolvant min

pourvu que Q i est rangée rang plein.

Appelons Se i une matrice de sélection rangée qui concerne au noeud i. Se est une matrice, où Se jj = 1 si j ∈ (i ∪ N i ). Q i peut être exprimée comme:

Comme montré dans (B.3.15), Q i devrait être une matrice de rang rangée plein pour assurer la contrainte Q i c i = x1 pour être satisfaits.

Conséquent, Q i = q i q j , j ∈ N i peut être établie en sélectionnant le q j , j ∈ N i jusqu'à ce que le rank(Q i ) ne peut plus augmenter.

Il y a beaucoup de règles qui peuvent être appliquées pour le choix de q j . Cependant, dans notre étude, nous choisissons q j basé sur le conditionnement établi Q i . Ce qui signifie que, q j est choisi de manière que le conditionnement établi Q i est le plus petit.

De (B.3.5), sachant que les coefficients c k ont un motif de signe alternatif, le problème peut être reformulé comme suit:

Pour résoudre ce problème d'optimisation sous contrainte, de manière générale, la fonction Lagrangien augmenté est défini comme suit:

Chaque noeud met à jour son vecteur de paramètres en utilisant:

De [START_REF] Nedic | Constrained consensus and optimization in multi-agent networks[END_REF], l'analyse de la convergence de cette méthode a été bien étudié. Avec un choix approprié du stepsize γ[0], l'algorithme 5 converge vers un minimum global. La vitesse de convergence dépend encore le choix du stepsize gradient.

(B):Méthode de la direction alternée des multiplicateurs (ADMM).

ADMM est relativement facile à mettre en oeuvre et ses propriétés de convergence ont été bien étudiée dans des travaux récents tels que, dans [START_REF] Boley | Local linear convergence of ADMM on quadratic or linear programs[END_REF][START_REF] Boyd | Distributed optimization and statistical learning via the alternating direction method of multipliers[END_REF][START_REF] Erseghe | Fast consensus by the alternating direction multipliers method[END_REF]. Les avantages de

Ces itérations indiquent que la méthode est particulièrement utile lorsque les c-et zminimisations sont réalisée de manière efficace.

Pour la fonction de coût considéré ici, nous pouvons noter que le problème de cminimisation (B.3.24) est résolu par:

) = 0 en ce qui concerne le terme que z ij = z ji . Réécrire la fonction Lagrange augmentée:

Ensuite, nous projetons z ij [t + 1] sur la contrainte de signe comme:

Multiplicateurs de Lagrange sont mis à jour:

L'algorithme correspondante est alors décrit comme dans l'algorithm 6.

B.3.2.4 La factorisation de la matrice moyenne sur la base de Laplace avec une connaissance imparfaite de x

Pour toutes les méthodes proposées dans le cas précédente, la valeur de consensus x est préalablement un paramètre obligatoire.

Nous détendons maintenant cette exigence au moyen d'une nouvelle méthode qui traite l'estimation des valeurs propres du laplacien quand une valeur de consensus de moyenne exacte n'est pas connue.

Dans ce qui suit, nous supposons que les noeuds exécutent un protocole standard de consensus de moyenne ( 1.2.1), qui est arrêté après M itérations, et rangent leurs valeurs locales x i,M . Dans cette étude, nous voulons montrer la performance de l'algorithme proposé en ce qui concerne certaines valeurs de M.

Puisque on ne sait pas exactement x, pour une valeur suffisamment grande de M, x(M) peut être considéré comme son approximation rationnelle. Conséquent, notre objectif est de trouver les coefficients c k qui minimisent l'erreur quadratique sur la valeur finale de consensus:

Ensuite, nous pouvons en déduire l'ensemble S 1 = {α α α k } k=1,••• ,h et celle de l'inverse des valeurs propres du laplacien S 2 = 1 λ i , i = 2, . . . , D, qui appartient à S 1 : S 2 ⊂ S 1 . Afin de résoudre le problème d'une manière répartie, le problème d'optimisation peut être reformulé comme suit:

subject to c i,k = c j,k , i = 1, . . . , N; j ∈ N i ; k = 0, . . . , h.

où c i,k le k th coefficient de polynôme local et C i signifie l'ensemble des contraintes en ce qui concerne noeud i. On le définit comme

q(0) q(1) . . . q(h)], étant rang plein de rangée en construction.

Cet ensemble de contraintes permet de respecter non seulement l'égalité du vecteur des coefficients locaux c i si dans un voisinage donné, mais aussi pour faire appliquer le produit scalaire q T i c i être égal globalement. a noter que, en introduisant les variables auxiliaires {z ij }, le problème d'optimisation convexe sous contrainte suivant est absolument équivalent au problème (B.3.32):

Puisque le graphe est supposé être connecté, les contraintes (B.3.34) forcent c i = c j , j ∈ N i . Nous pouvons alors écrire le Lagrangien augmenté associé:

ou sous forme de vecteur:

La solution à ce problème ADMM agit en trois étapes:

• Tout d'abord, la minimisation de c i :

où Ω C i [•] représente la projection dans l'ensemble des contraintes du vecteur dans l'argumentation.

• Deuxièmement, la minimisation en ce qui concerne à z ji = z ij :

• Troisième, la mise à jour de multiplicateur de Lagrange:

Résoudre le problème sous-optimisation (B.3.35) par ∂Lρ(c,z,y) ∂c i = 0, nous obtenons:

Puis, projeter la solution obtenue à l'ensemble des contraintes

nous obtenons:

L'algorithme correspondante est alors décrit comme dans l'algorithm 7. Quand les coefficients c k sont obtenus approximativement, nous pouvons utiliser des Algorithmes 9 et 10 pour estimer tout le spectre de matrice de Laplace L, qui peut être utilisé pour estimer les indices de la robustesse du réseau.

B.3.2.5 La factorisation de la matrice moyenne sur la base de Laplace avec une valeur de consensus inconnue x

Selon la section précédente, une façon pour l'estimation des valeurs propres du laplacien est:

• résoudre le problème du consensus de moyenne.

• résoudre le problème de la factorisation de la matrice de moyennage (B.3.15), Algorithmes 5 et 6.

• récupérer les valeurs propres du laplacien en utilisant les algorithmes 9 et 10.

Au lieu de suivre ces étapes, dans cette section, nous proposons une méthode pour fusionner les deux premières étapes. Dans ce but, nous proposons de résoudre une combinaison convexe des deux fonctions convexes, c'est-à-dire

avec θ-un paramètre scalaire ajustable (0 ≤ θ ≤ 1), xi -une estimation locale de la valeur moyenne, tandis que x i (0) est une valeur initiale fixe de noeud i.

On voit que, la fonction d'objectif (B.3.40) est une combinaison de deux fonctions convexes. La première prend en compte le problème du consensus de moyenne, tandis que la seconde concerne l'estimation des coefficients du polynôme dont les racines sont l'inverse des valeurs propres du laplacien distinctes différents de zéro. Et ce problème d'optimisation convexe peut être résolu efficacement en utilisant l'approche de l'ADMM.

Maintenant, en introduisant certaines variables duales υ υ υ and y, le Lagrangien augmenté est défini comme:

où:

où ρ 1 , ρ 2 sont des termes de pénalité supposés constants.

ADMM fournit une solution itérative qui agit les cinq étapes suivant ci-dessous:

• Minimisation ec ce qui concerm la valeur moyenne du consensus xi :

• Minimisation en ce qui concern les coefficients polynomiaux c i :

dresse pour la projection dans l'ensemble des contraintes C i (x i [k + 1]) du vecteur dans l'argumentation.

• Minimisation en ce qui concerm les variables auxiliaires µ ij avec la contrainte µ ji = µ ij :

• Minimisation en ce qui concerm les variables auxiliaires z ij :

• La mise à jour Multiplicateurs de Lagrange: Nous obtenons:

Puis, projeter la solution obtenue dans l'ensemble des contraintes:

nous obtenons:

, Q i étant rang plein de rangée en construction. Ensuite, la solution des problèmes sous-optimisation (B.3.44) et (B.3.45), de la même façon, déduit 1. Le dérivé de L ρ 2 (x, µ µ µ, υ υ υ) en ce qui concern µ ij par rapport à la contrainte µ ij = µ ji :

contient les coefficients de expansion, et U .i est le i-ere vecteur propre, i.e. le i-ere colonne de U.

Par conséquent:

où vecd(•) est le vecteur colonne construit avec les éléments diagonaux de la matrice en argument.

Nous pouvons écrire l'équation équivalente ci-dessus:

Nous pouvons conclure que de la matrice estimée Ŵ, on obtient un graphe qui est isomorphe au graphe original. Dans ce qui suit, nous indiquons deux conditions suffisantes pour assurer la reconstruction correcte du graphique. Nous rappelons d'abord les notions suivantes de observabilité:

Un réseau est dit être [START_REF] Kibangou | Observability in connected strongly regular graphs and distance regular graphs[END_REF]:

• noeud observable d'un noeud donné, si ce noeud est en mesure de reconstruire l'état du réseau de ses propres mesures. Cette question a été étudiée par exemple dans [START_REF] O'clery | Observability and coarse graining of consensus dynamics through the external equitable partition[END_REF] and [START_REF] Parlangeli | On the reachability and observability of path and cycle graphs[END_REF] où il a été déclaré que un réseau avec une matrice d'état ayant au moins une valeur propre non simple n'est pas le noeud observable.

• voisinage observable d'un noeud donné si ce noeud peut reconstruire l'état du réseau entier de ses propres mesures et celles de ses voisins. Cette question a été étudiée dans [START_REF] Kibangou | Algebraic characterization of observability in distance-regular consensus networks[END_REF].

• globalement observable si elle est la voisinage observable depuis ne importe quel noeud.

Proposition 9

Supposons que le graphe est noeud-observable ou voisinage observable du noeud j. Si toutes les entrées de la condition initiale sont distinctes alors la topologie du réseau peut être reconstruit exactement du noeud j.

Proposition 10

Supposons que le réseau contient un ensemble A de noeuds anonymes et le graphe est noeud-observable ou voisinage observable du noeud j. Si toutes les entrées de la condition initiale associée aux noeuds anonymes sont toutes distinctes alors la topologie de réseau peut être reconstruit exactement du noeud j.

Ici, l'ambiguïté de permutation est limitée à des entrées associées aux noeuds anonymes.

B.6 Conclusions générales B.6.1 Résumé des contributions et conclusions

Ce chapitre a pour l'objective de conclure la thèse a travers un bref summerization ainsi qu'une perspective pour des travaux futurs.

Les Problèmes de Consensus en temps fini ont été reçus une vaste attention de la communauté de recherche en raison de ses avantages dans de nombreux domaines d'application, en particulier dans le calcul distribué. Brièvement parlant, les algorithmes de consensus de moyenne en temps fini (pour lesquelles un consensus peut être obtenu en un nombre fini d'étapes) sont couramment utilisés comme élément de base pour plusieurs algorithmes de contrôle, estimation ou inférence distribués en raison du fait qu'elles peuvent garantir une exécution en temps minimal, ces algorithmes sont donc beaucoup plus attrayant que ceux qui assurent la convergence asymptotique.

Ceci impact sur la thèse avec la conception de protocole de consensus de moyenne en temps fini indépendamment de la matrice de Laplace d'un graphe donné de manière distribuée et son application dans l'évaluation de la robustesse du réseau. Évidemment, il est incommode pour les réseaux de grande envergure. Cependant, un point intéressant peut être indiquée dans ce travail est que les matrices de consensus obtenus ne sont pas des matrices doublement stochastiques. Ceci est la condition de la convergence pour un protocole de consensus selon [START_REF] Xiao | Fast linear iterations for distributed averaging[END_REF].

Deuxièmement, l'étude sur la robustesse des réseaux de plus en plus attiré un attention de la communauté scientifique. Conséquent, dans cette thèse, l'évaluation de la robustesse du réseau, qui est mis en oeuvre par la Résistance effective du graphe (R) et le nombre d'arbres couvrants (ξ) d'une manière distribuée, devient évidemment une application domaine prometteuse. Puisque deux mesures de la robustesse sont des fonctions du spectre de Laplace sp(L), donc la tâche principale est maintenant d'estimer le spectre de Laplace. Précisément, trois cas ont été considérés: parfaitement connue x, bruyant, et inconnu x.

La matrice de consensus moyen peut être prise en D la matrice basée sur Laplace de consensus, où D représente le nombre de valeurs propres du laplacien distinctes différents de zéro. Pour le premier cas (parfaitement connue x), toutes les méthodes sont caractérisées en deux approches, qui sont approche directe et l'approche indirecte. La méthode directe donne lieu à un problème d'optimisation non convexe en ce qui concern les stepsizes α k , k = 1, . . . , N -1 qui sont l'inverse des valeurs propres du laplacien distinctes pouvons trouver un moyen d'obtenir U. Puisque la spectre de Laplace sp(L) peut être réalisée en utilisant le méthodes proposées dans le Chapitre 4, nous pouvons reconstruire la topologie du réseau.

B.6.2 Travaux en cours et à venir 1. Conception des Protocoles de consensus de moyenne en temps fini:

(a) La dépendance de la vitesse de convergence sur la séquence d'apprentissage reste à étudier. Ainsi, la conception de séquences d'apprentissage optimales peut être étudiée.

(b) La recherche sur la robustesse du consensus de moyenne en temps fini.

2. L'évaluation de la robustesse du réseau:

(a) Nouveaux méthodes qui peuvent être appliquées pour les graphes de grande échelle.

(b) Pour la méthode indirecte, dans la pratique, il peut rencontrer le problème mal conditionné lors de l'établissement sous-matrices Q i , c le problème mal conditionné lors de l'établissement sous-matrices Q i qui conduit à perdre le rang de ces sous-matrices. Conséquent, on peut étudier dans la précondition de ces matrices pour conserver la précision des méthodes proposées.

3. Reconstruire la topologie du réseau:

(a) Concevoir une nouvelle méthode qui traite avec le spectre de Laplace sp(L) dont les valeurs propres ne sont pas toutes distinctes.

(b) D'évaluer la robustesse de ces méthodes lorsque l'échange de données sont sujettes à des communications imparfaites.