Keywords: Graph signal processing, Graph translation, Graph analysis, Stochastic graph signals, Stationarity, Stationary graph signal Kronecker matrix product convolution

When one tries continuously, one ends up succeeding. Thus, the more one fails, the greater the chance that it will work.

-Shadocks' Motto.

To Diane and my parents. . . This dissertation is the result of three years of research within the Dante team in the LIP laboratory

Résumé

Ce manuscrit introduit dans les trois premiers chapitres le domaine du traitement du signal sur graphe en commençant par poser les bases d'algèbre linéaire et de théorie spectrale des graphes. Nous définissons ensuite le traitement du signal sur graphe et donnons des intuitions sur ses forces et faiblesses actuelles comparativement au traitement du signal classique.

Le chapitre 4 cible ensuite l'étude de la structure d'un graphe par l'analyse des signaux temporels via une transformation d'un graphe en une série temporelle. Ce faisant, nous exploitons une approche unifiée d'apprentissage semi-supervisé sur graphe dédiée à la classification pour obtenir une série temporelle lisse. Enfin, nous montrons que cette approche s'apparente à du lissage de signaux sur graphe, comparable à un filtre de Wiener.

Le chapitre 5 introduit un nouvel opérateur de translation sur graphe définit par analogie avec l'opérateur classique de translation en temps et vérifiant la propriété clé d'isométrie. Cet opérateur est comparé aux deux opérateurs de la littérature et son action est décrite empiriquement sur quelques graphes clés.

Le chapitre 6 décrit l'utilisation de l'opérateur ci-dessus pour définir la notion de signal stationnaire sur graphe. Après avoir étudié la caractérisation spectrale de tels signaux, nous donnons plusieurs outils essentiels pour étudier et tester cette propriété sur des signaux réels.

Le dernier chapitre s'attache à décrire la boite à outils GRASP pour MATLAB développée et utilisée tout au long de cette thèse.

Mots clefs :

Traitement du signal sur graphe, Translation sur graphe, Analyse de graphe, Signaux sur graphe aléatoires, Stationnarité, Signaux sur graphe stationnaires.

F

classical Fourier transform for 1D signals.

F

graph Fourier transform matrix.

G n grid graph with n × n vertices.

G n directed grid graph with n × n vertices.

Γ = diag(γ 1 , . . . , γ N -1 ) eigenvalues matrix of A.

g κ heat kernel.

IMF (k) x k th IMF of the signal x.

J k (λ) Jordan block.

K N complete graph with N vertices.

L

Laplacian matrix.

L normalized Laplacian matrix.

Λ = diag(λ 1 , . T n torus graph with n × n vertices.

T n directed torus graph with n × n vertices.

T G graph translation based on L.

T G graph translation based on L .

Θ = diag(θ 1 , . . . , θ N -1 ) eigenvalues matrix of L .

T i generalized translation localized on vertex i .

U = [u 0 • • • u N -1 ] eigenvectors matrix of A.
U N ,U N plane random sampling graph.

W

Molène graph.

Introduction

The advent of new technologies has given us means and opportunities to quantify the world around us by collecting more and more data leading inexorably to the now infamous Big Data trend in research and development. This deluge of data offers plenty of opportunities to develop new data analysis methods and validate, or invalidate, them using these data. In this context, a new and emerging thrilling field aims at the study of the subset of structured data:

The field of Signal Processing on Graphs.

Two centuries of research have given the mother field of Signal Processing the strength to widely spread within all sciences or in our every day life. During this long time scale, the focus has been understanding data structured in time (e. g. sound) or in space (e.g. heat diffusion). These are mainly regular structures with well defined Euclidean metrics. We now are in the days where data is irregularly structured, with unlimited richness.

Can we capitalize on the structure to understand the data?

This grand question can be applied to many types of datasets. Indeed, dealing with large datasets, the data can almost always be split into data on entities. In turn, these entities are structured with relations between them. We now give examples from many fields that show such structures.

With computer networks, and networks of electronic devices in general, questions such as the detection of a faulty device (e.g. a router) or an abnormal behavior (e.g. a Distributed Denial of Service (DDoS) attack) are essential for security reasons. Fair sharing (e.g. load balancing) of network resources for a good Quality of Service (QoS) is another question appearing often in the cloud computing community.

Closely related to computer networks, sensor networks are widely studied in the literature with the ultimate goal of efficiently sensing the real world. Key questions are the placement of the sensors (e.g. ground weather stations), routing strategies (e.g. in content centric networks), transmission of the right amount of data (compressive sensing) [START_REF] Anis | Towards a sampling theorem for signals on arbitrary graphs[END_REF], Shuman et al., 2013a], detection of a faulty sensor [Sandryhaila and Moura, 2014a], or even onboard preprocessing of data [START_REF] Sandryhaila | Discrete Signal Processing on Graphs[END_REF], Shuman et al., 2011, Shi et al., 2015].

Opposite to technological networks lie the social networks. People naturally tend to relate to each other through shared interests, a shared employer, or simply a shared neighborhood. With the increased usage of online social networks, reliable data describing relations between people are easier to obtain. However, one key elements that pertains to these datasets is the irregular nature of the structures, with for example a power law distribution of the number of connections per individual. Such structures are very far from the regular Euclidean spaces studied by classical signal processing. Many works of the past decade have been trying to answer the question of community detection within these networks, such as the recent multiscale approach of [Tremblay and Borgnat, 2014]. These works focus on the structure, but social networks also comprise data directly related to each individual, such as activity in the network, centers of interest, physical or psychological traits to cite only a few. Classical questions are then the inference from data on a few individual [START_REF] Avrachenkov | Generalized Optimization Framework for Graph-based Semi-supervised Learning[END_REF], characterization of individuals from their networks, or non human behavior detection.

Still on the human side, several other kinds of networks lie right within the body: the biological networks. For example, can we characterize the brain activity through its neural network [START_REF] Hammond | Wavelets on graphs via spectral graph theory[END_REF]? or is it possible to characterize the protein-protein interaction network from the dynamics of a few proteins? These questions have a direct impact on medicine and pharmacology. On a related topic, can we model the body as a network of limbs, head and torso to better understand the body language [START_REF] Yang | Gesture dynamics modeling for attitude analysis using graph based transform[END_REF], Kao et al., 2014]?

These are only a few examples of structured data and communications related to the field signal processing on graphs. As of today, the toolbox of this field is growing rapidly, and several successes have appeared. We only cite a few: The graph Fourier transform decomposing data into frequency components [START_REF] Hammond | Wavelets on graphs via spectral graph theory[END_REF][START_REF] Shuman | The Emerging Field of Signal Processing on Graphs: Extending High-Dimensional Data Analysis to Networks and Other Irregular Domains[END_REF], the spectral graph wavelet transform decomposing data into local and global components w.r.t. to each entity [START_REF] Hammond | Wavelets on graphs via spectral graph theory[END_REF], the vertex-frequency decomposition producing similar components but with low and high frequency components instead of local and global [START_REF] Shuman | Vertexfrequency analysis on graphs[END_REF], graph coarsening reducing the size of a graph so as to preserve the meaning of a signal [Shuman et al., 2013a], the Empirical Mode Decomposition (EMD) decomposes a signal into slowly and quickly varying components in a local sense [Tremblay et al., 2014], or sampling theory selecting the minimum number of entities for which data is necessary to reconstruct a signal [START_REF] Anis | Towards a sampling theorem for signals on arbitrary graphs[END_REF].

Introduction

Graph-Based Semi-Supervised Learning (GSSL). In the course, we observe that GSSL is actually a specific smoothing filter on graph signals: a Wiener-like filter.

In the examples of this introduction, in many cases, we are in the presence of many time series relating to each other through a structure. For example, in a sensor network, each sensor collects data over time, in a computer network, routers measure the load over time, and in the brain the activity is a highly dynamic process. The question we have then is whether it is possible to consider these data as a set of graph signals, rather than a set of time series, and study the statistical properties of these graph signals. In the classical setting of time series, the key tool to study statistics of time series is that of stationary signal. Our next contribution to the field of signal processing on graphs extend the framework of stationary signals to the graph setting. In turn, this allows us to properly define Wiener filters for graph signals.

Stationarity being defined as a statistical invariance through a translation operator, we look in chapter 5 into the choice of such an operator. Unfortunately, operators from the literature are missing the key ingredient of energy preservation: If an operator alters the energy of a graph signal, then the statistics of the signal and the altered signal can not be the same. We therefore introduce a new translation operator preserving the energy and call it the graph translation.

We have then all the ingredients to define stationary graph signals. Chapter 6 introduces the key definitions of Wide (or Weak) Sense Stationary (WSS) and Strict Sense Stationary (SSS) graph signals. We show then that WSS graph signals have the same spectral characterization than temporal WSS signals. After a discussion on how to interpret non-stationarity, we apply our results to a dataset of weather readings to study stationarity with a concrete example.

Last but not least, in chapter 7, we present another important contribution to the field: a MATLAB toolbox. This toolbox has been developed during the course of this PhD, and focuses on implementing powerful graph drawing functions, together with graph and graph signals primitives and operators. As we will see in chapter 2 and chapter 3, the field of graph signal processing is heavily relying on matrix representations of graphs. Therefore, we will need results from linear algebra to study matrices and vectors. We review in this section classical definitions and results of the field. We refer the interested reader to [START_REF] Horn | Matrix analysis[END_REF] and [Golub and Van Loan, 1996] for details and in-depth studies of these objects.

Fundamental Definitions and Properties

Let A be a matrix of n rows and m columns. We denote A i j the element at row i and column j of the matrix A. A square matrix is a matrix such that n = m. A column (resp. row) vector is a matrix such that m = 1 (resp. n = 1). An element of C n is therefore a vector. An example of a square matrix is shown on Figure 1. 1. In particular, these figures show two color scales that will be used extensively throughout this dissertation. A random matrix A with entries uniformly drawn in [-1, 1]. These two figures illustrate two color scales that are used throughout this dissertation for both matrices and graphs. The first color scale is essential to illustrate amplitudes with darker tones for higher amplitudes and lighter tones for smaller amplitudes. Yet the sign of the entries can be identified with red (positive) and blue (negative) colors. When all elements of the matrix are positive, a better color scale is the second one with more tones from white to black.

Chapter 1. Matrices, Vectors and Matrix Decompositions

The matrix E i j denotes the matrix having zeros everywhere except for the element of the i th row and j th column that equals 1. The vector δ i denotes the vector having zeros everywhere except for a 1 at its i th component:

E i j = j th         0 0 . . . 1 • • • i th 0 0 δ i =             0 . . . 1 i th . . . 0 .
The set {δ i } i is in particular the canonical basis of the vector space C n . We denote 0 the all zero matrix and 1 the all one matrix. Note that 0 and 1 can also be vectors.

The canonical scalar product between two vectors is denoted 〈x, y〉 = y * x = i x i y * i . We use the notation x * for the complex conjugate x when x is a number or the transpose conjugate if x is a matrix. Without further precision, the scalar product is the canonical scalar product. Its associated norm is the 2 -norm denoted x = 〈x, x〉 = i |x i | 2 . Two vectors are orthogonal if their scalar product equals zero. By extension, a set of vectors is orthogonal if the vectors are pairwise orthogonal. If the vectors are all of unit norm, such a set of vectors is orthonormal.

A matrix A is diagonal if its off-diagonal elements are all zero (∀i = j , A i j = 0). This definition is valid for any matrix but is more relevant to square matrices, in which case we write A = diag(A 11 , . . . , A N N ). By extension, a block diagonal matrix is a matrix that can be written as a diagonal of smaller matrices. A matrix is upper triangular if its non zero elements lie either on its diagonal or on its upper right part. The identity matrix I n is the square diagonal matrix with only ones on its diagonal:

I n = diag(1, . . . , 1) =   1 . . . 1   .
When the value of n is clear from the context, we simply denote the identity matrix I .

We will be using several matrix operations hereafter. The first two are the classical matrix addition and multiplication:

(A + B ) i j = A i j + B i j (AB ) i j = k A i k B k j .
The identity element of the addition is the matrix 0 and that of the multiplication is the identity matrix I . Other matrix operations useful in this dissertation are the transpose operation inverting lines and columns and the transpose conjugate operation taking the complex conjugate of the transpose matrix:

(A T ) i j = A j i A * = A T .
A matrix is then symmetric if A = A T . For complex matrices, the invariance by the complex conjugate is more useful, in which case the matrix is said Hermitian: A = A * .

A matrix A is said invertible if there exists a matrix B = A -1 such that B A = I = AB . A matrix is said unitary if A A * = I = A * A, in which case A -1 = A * . The columns of a unitary

I = A * A =   a * 1 . . . a * n   [ a 1 • • • a n ] =      a * 1 a 1 a * 2 a 1 • • • a * 1 a n a * 1 a 2 a * 2 a 2 . . . . . . . . . . . . a * n a 1 • • • • • • a * n a n      ,
such that ∀i , a * i a i = 1 and ∀i = j , a * i a j = 0. These two properties allow to introduce two equivalence classes between matrices. The matrix A is similar to the matrix B if there exists an invertible matrix V such that A = V BV -1 . The matrix A is unitarily similar to B if there exists a unitary matrix U such that A = U BU * . These two equivalences will be of particular interest when defining a graph Fourier transform.

Building on these equivalences, a matrix is diagonalizable if it is similar to a diagonal matrix. By extension, the diagonalization of a diagonalizable matrix is the process of computing A = V DV -1 with D diagonal and V invertible. A matrix A is unitarily diagonalizable if it is unitarily similar to a diagonal matrix D: A = U DU * . In particular, Hermitian matrices are unitarily diagonalizable with real matrix D [Horn and Johnson, 2013, Theorem 4.1.5]. Note that if A is real, the matrix U may be chosen with real entries.

An eigenvector of a matrix A is a non zero vector x verifying Ax = λx, where λ is its associated eigenvalue. An eigenvector together with an eigenvalue is called an eigenpair. The set of all eigenvalues of the matrix A is called its spectrum and denoted Sp(A) while the set of all its eigenpairs is its spectral decomposition. The characteristic polynomial of the matrix A is the polynomial p A (x) = det(x I -A), i.e. the determinant of the matrix x I -A. The zeros of the polynomial p A are then the eigenvalues of the matrix A, and we have 1 :

p A (x) = q i =1 (x -λ i ) m i ,
with {λ i } 1≤i ≤q the distinct eigenvalues of A. The exponent m i is called the multiplicity of the eigenvalue λ i . Hereafter, we will consider for simplicity the set of eigenvalues {λ i } 1≤i ≤n where the eigenvalues appear as many times as their multiplicity. It is worth noticing that when all eigenvalues are distinct (i.e. of multiplicity 1), the matrix is diagonalizable [Horn and Johnson, 2013, Theorem 1.3.9].

Let A be a diagonalizable (resp. unitarily diagonalizable) matrix with A = V DV -1 (resp. A = U DU * ) and D a diagonal matrix. The columns of V (resp. U ) are then the eigenvectors of A, and the diagonal of D holds its eigenvalues. When A is unitarily diagonalizable, the eigenvectors of A are orthonormal since U is unitary.

An Hermitian matrix is definite positive (resp. semi-definite positive) if and only if ∀x = 0, x * Ax > 0 (resp. ∀x = 0, x * Ax ≥ 0). Equivalently, its eigenvalues are positive (resp. nonnegative). More generally, the quantity x * Ax for an Hermitian matrix A can be used to find the maximum and minimum eigenvalues of an Hermitian matrix using the Rayleigh quotient R A (x) = x * Ax x * x [Horn and Johnson, 2013, Theorem 4.2.2]:

min λ ∈ Sp(A) = min R A (x) : x ∈ C n , x = 0 = min x∈C n x =0
x * Ax x * x Chapter 1. Matrices, Vectors and Matrix Decompositions

max λ ∈ Sp(A) = max R A (x) : x ∈ C n , x = 0 = max x∈C n x =0
x * Ax x * x

We close this section with an alternative definition of matrix product that we will use in chapter 2: the Kronecker matrix product between a matrix A of size n × m and a matrix B of arbitrary size. The Kronecker matrix product of A and B can be written as the following block matrix:

A ⊗ B =   A 11 B • • • A 1m B . . . . . . . . . A n1 B • • • A nm B   .
The Kronecker matrix product is distributive over the matrix addition: A ⊗(B + C ) = A ⊗ B + A ⊗C . Moreover, if Ax = λx and B y = µy, we can show that λµ is an eigenvalue associated to the eigenvector x ⊗ y: A ⊗ B (x ⊗ y) = λµ(x ⊗ y).

( 1.1) The spectral decomposition of A ⊗ B is therefore completely defined by the equation above.

Matrix Decompositions

We introduced the spectral decomposition (λ i , x i ) : Ax i = λ i x i i of a matrix A in the previous section as the most basic decomposition of a matrix. In particular, when the matrix is diagonalizable, this decomposition entirely describe the properties of the matrix A. However, there exists alternative decompositions. We now present several of them that will be useful in this dissertation, in particular to study properties of general matrices.

QR Decomposition

This simple QR decomposition is at the core of several other decomposition. We use it in this dissertation only for the justification of the matrix decompositions in the subsequent sections of this chapter.

Given a matrix A, the QR decomposition computes a unitary matrix Q and an upper triangular matrix R such that A = QR [Horn and Johnson, 2013, Theorem 2.1.14]. The justification is the following. Let a 1 be the first column of A and r 1 = a 1 . Let U 1 be a unitary matrix such that U 1 a 1 = r 1 δ 1 . It follows that:

U 1 A = r 1 0 A 2 ,
with an arbitrary matrix. Applying this scheme recursively gives A 2 = Q 2 R 2 , and:

1 0 0 Q * 2 U 1 A = 1 0 0 Q * 2 r 1 0 A 2 = r 1 0 R 2 = R,
and finally

Q = U * 1 1 0 0 Q 2 .

Jordan Decomposition

We have seen that the spectral decomposition of a diagonalizable matrix leads to the matrix decomposition A = V DV -1 . Unfortunately, when A is not diagonalizable, A is not similar

Matrix Decompositions

to any diagonal matrix. Instead of a diagonal matrix, the Jordan decomposition computes a triangular and block diagonal matrix J similar to A.

These blocks have a very specific form called Jordan blocks. A Jordan block of size k is the following square triangular matrix:

J k (λ) =      λ 1 . . . . . . . . . 1 λ      .
Up to scalar multiplication, δ 1 is then the only eigenvector of J k (λ). It is associated to the eigenvalue λ. The vectors δ 2 to δ k verify the relation J k (λ)δ i = λδ i + δ i -1 . They are called generalized eigenvector of the matrix J k (λ). The canonical basis δ i i is called in this context a Jordan chain.

The spectral decomposition states that a matrix A has a set of independent eigenvectors each associated to an eigenvalue. If A is not diagonalizable, this set is of size less than n, in which case, one or more eigenvectors have an associated Jordan chain of size greater than one. Overall, the set of eigenvectors and generalized eigenvectors has a size of exactly n. Let V be the matrix whose columns are these eigenvectors and generalized eigenvectors, ordered by Jordan chains. We obtain then the Jordan decomposition A = V JV -1 such that:

J =   J k 1 (λ 1 )
. . .

J k q (λ q )   ,
with k i the size of the i th Jordan chain [Horn and Johnson, 2013, Theorem 3.1.11]. In other words, A is similar to the matrix J .

More details on the Jordan decomposition can be found in [Horn and Johnson, 2013, Section 3.1]. In [Golub and Van Loan, 1996, Section 7.6.5], the authors point out that despite the theoretical power of the Jordan decomposition, in practice its computation is not practical for numerical stability reasons. Indeed, a very small perturbation of one or several elements of the matrix, due for example to rounding errors, can change the unitarily similar matrix J at a much bigger scale than said perturbation (i.e. a 1 on the upper diagonal is changed to a 0 or vice versa).

Schur Decomposition

We mentioned that the Jordan decomposition is not numerically stable. We present now a more stable decomposition: the Schur decomposition.

Let λ 1 , . . . , λ n be the eigenvalues of A. The goal of the Schur decomposition is to compute a triangular matrix T which is unitarily similar to the matrix A, while having the diagonal of T equal to the eigenvalues of A. [Horn and Johnson, 2013, Theorem 2.3.1] states that this is always possible, and when the eigenvalues of A are real, then it is possible to have A = U T U * with U having only real elements.

The proof is very simple and intuitive. We begin by choosing an eigenpair λ 1 , x 1 , and a unitary matrix

U 1 = [x 1 u
U * 1 AU 1 =     x * 1 u * 2 . . . u * n     A [ x 1 u 2 • • • u n ] =     λ 1 x * 1 x 1 x * 1 Au 2 • • • x * 1 Au n λ 1 u * 2 x 1 u * 2 Au 2 • • • u * 2 Au n . . . . . . . . . . . . λ 1 u * n x 1 u * n Au 2 • • • u * n Au n     =     λ 1 x * 1 Au 2 • • • x * 1 Au n 0 u * 2 Au 2 • • • u * 2 Au n . . . . . . . . . . . . 0 u * n Au 2 • • • u * n Au n     = λ 1 0 B ,
with representing a matrix of arbitrary elements. We can show that the eigenvalues of B are also eigenvalues of A. Recursively applying this scheme to the matrix B leads to

B = U 2 T B U * 2
with T B triangular. We define then the following matrix U :

U = U 1 1 0 0 U 2 .
Using UU * , we can show that U is unitary. We obtain then:

U * AU = 1 0 0 U 2 * U * 1 AU 1 1 0 0 U 2 = 1 0 0 U 2 * λ 1 0 B 1 0 0 U 2 = λ 1 0 U * 2 BU 2 = λ 1 0 T B = T ,
which is a triangular matrix with the eigenvalues of A on its diagonal, and we have A = U T U * .

A corollary of the Schur decomposition is that every square matrix is block diagonalizable, i.e. for every square matrix A, there exists an invertible matrix V (not necessarily unitary) and a block diagonal matrix D = diag(T 11 , . . . , T kk ) such that A = V DV -1 , every T i i is upper triangular, and the diagonal of T i i is constant and equal to an eigenvalue of A [Horn and Johnson, 2013, Theorem 2.4.6.1].

In practice, the Schur decomposition is efficiently computed using the QR decomposition [Golub and Van Loan, 1996, Chapter 7], in particular, Section 7.5 for the Schur decomposition and Section 7. 6.3 for the block diagonalization.

Cholesky Decomposition

Cholesky decomposition is restricted to a specific type of matrix, namely the semi-definite positive matrices. Let A be such a matrix. There exists a lower triangular matrix L verifying A = LL * . The converse is also true: If A is Hermitian and there exists a lower triangular matrix [Horn and Johnson, 2013, Corollary 7.2.9].

Matrix Decompositions

verifying A = LL * , then A is semi-definite positive. If A is definite positive, L is unique
The decomposition is easy to see using the QR decomposition. If A is semi-definite positive, then it has a square root A 1/2 such that (A 1/2 ) 2 = A. Using the QR decomposition on this square root gives A 1/2 = QR. Since A is Hermitian, A 1/2 is also Hermitian and we have

A = A 1/2 * A 1/2 = R * Q * QR. Since Q is unitary, it follows that A = R * R and L = R * is a lower triangular matrix.

Graphs

CHAPTER 2

Contents

We introduce in this chapter the discrete graph structure. Using tools of linear algebra presented in the previous chapter, we study the basic mathematical objects describing the properties of a graph. In the course of this chapter, we present key results of the field of spectral graph theory that focuses on the spectral decomposition of matrices describing graphs. This chapter closes with example graphs that will be of interest in this dissertation.

Notations

A graph G = (V, E ) is a discrete structure composed of a set of vertices V (G ) linked by edges E (G ) ⊆ V (G ) × V (G ).
When clear from context we simply denote these sets as V and E . Note that a calligraphic symbol such as G will most often denote a graph. We call a loop Chapter 2. Graphs an edge linking a vertex to itself. Throughout this dissertation we only consider graphs without loop. Table 2 .2 (page 24) shows several examples of graphs. The number of vertices is denoted

N = |V |. A subgraph H of a graph G = (V, E ) is a graph with V (H ) ⊆ V (G ), and E (H ) = E (G ) ∩ (V (H ) × V (H )).
A weighted graph is a graph whose edges have weights. The weight function is denoted w G : E → R, and the weighted graph G = (V, E , w). Although complex valued weights or real negative weights could be considered, we only consider real non-negative weights for which strong theoretical results using linear algebra exist. Remark, that a graph without weights is a graph for which all weights can be chosen as 1. By extension, we denote w i j the weight of the edge i j .

A graph is symmetric or undirected if its edges do not have an orientation, i.e. if i j is an edge, then j i is also an edge. On the contrary, if at least one edge has an orientation, the graph is directed or non symmetric. An edge i j such that j i is not an edge is called a directed edge and is represented by an arrow from its source vertex i to its target vertex j . Other edges are represented by a line between the two vertices. We denote a directed edge by i j or by the couple (i , j ). An edge is denoted i j or by the set {i , j }. For simplicity reasons, we will use the shorter notation i j for both edges and directed edges. Note that we will be mostly dealing with symmetric graphs in Part II.

We denote d G i = j w i j the degree of a vertex, i.e. the total weight of the edges connecting i to it neighbors. When the graph is directed this is the out-degree due to the outgoing nature of the directed edges summed. Its counterpart is the in-degree d G i ,in = j w j i . A graph is regular if all its vertices have equal degree. When all edges have equal weight, this is equivalent to all vertices having the same number of neighbors. An example of a regular graph with degree 3 is shown on Figure 2.1(a).

A graph is bipartite if its vertex set can be partitioned into two subsets such that all edges lie between those two sets. Such a graph is advantageously represented by drawing the two subsets on top and on the bottom of a figure. The edges are therefore crossing the figure, and none of them is horizontal. An example of such a graph is shown on Figure 2.1(b). In particular, bipartite graphs appear naturally when the nature of the vertices represent two different objects. For example, people and topics lead to a bipartitioned set of vertices, and edges are drawn whenever a person is interested in the corresponding topic. A path is a succession of vertices connected by edges: p = (i 0 , . . . , i k ) with i j i j +1 a directed edge (or simply an edge). In this context, k is the length of the path. A symmetric graph is connected if there exists a path between any pair of vertices. The equivalent definition for directed graphs is that of strongly connected graph having a path between any (ordered)

Graph Transformations

pair of vertices. Finally, if the symmetrized 1 graph of a directed graph is connected, the graph is said to be weakly connected. By extension, a connected component of a symmetric graph is a connected subgraph. For simplicity, and without loss of generality, we suppose the graphs to be connected.

Graph Transformations

In this section we discuss some transformations of graphs. First, we consider graph transformations as functions modifying the weights of one or several graphs. Let f : R k → R be a multivariate function. By extension, we denote also f the graph transformation:

f (G 1 , . . . , G k ) = k i =1 V (G i ), k i =1 E (G i ), f (w G 1 , . . . , w G k ) (2.1)
These very basic graph transformations are limited to functions of the weights. There exists also transformations at the vertex and edge level. We start with the Cartesian product G 1 G 2 between two graphs G 1 and G 2 [Brouwer and Haemers, 2012, Section 1.4.6]. The vertex set of G 1 G 2 is the Cartesian product V (G 1 ) × V (G 2 ), whereas each edge corresponds to one and only one edge of l ) . Mathematically, it translates to:

G 1 or G 2 . For example, if i is a vertex of G 1 and j a vertex of G 2 , then (i , j ) is a vertex of G 1 G 2 , and if i j is an edge of G 1 and k l is an edge of G 2 , then (i , k) ( j , k) and (i , k) (i , l ) are edges of G 1 G 2 but not (i , k) ( j ,
V (G 1 G 2 ) =V (G 1 ) × V (G 2 ) E (G 1 G 2 ) = (i , j )(i , k) : i ∈ V (G 1 ) , j k ∈ E (G 2 ) (i , k)( j , k) : i j ∈ E (G 1 ) , k ∈ V (G 2 ) .
If the graph is weighted, the weight of an edge of G 1 G 2 can be simply the weight of the corresponding edge in G 1 or G 2 . Intuitively, a path in the resulting graph corresponds to a sequential walk in the sense that it is a succession of paths in G 1 and G 2 . An example of Cartesian product is shown on Figure 2.2(a). The Kronecker product G 1 ⊗ G 2 of G 1 and G 2 is another graph product whose vertex set is also the Cartesian product V (G 1 ) × V (G 2 ). The edge set is however defined as a parallel walk 1 Given a graph G , its symmetrized graph is the graph H on the same set of vertices and such that E (G ) ⊆ E (H ) and if i j ∈ E (G ), then j i ∈ E (H ). In other words, G is completed with edges such that there is no directed edge anymore.

Chapter 2. Graphs

in both graphs: There is an edge (i , k) ( j , l ) in G 1 ⊗ G 2 if and only if i j is an edge of G 1 and kl an edge of G 2 (the same applies to directed edges). Mathematically, it translates to:

V (G 1 ⊗ G 2 ) =V (G 1 ) × V (G 2 ) E (G 1 ⊗ G 2 ) = (i , k)( j , l ) : i j ∈ E (G 1 ) , kl ∈ E (G 2 ) .
An example of Kronecker product is shown on Figure 2.2(b). The weights of such a product should be a function of the weights of both edges of G 1 and G 2 . Therefore, we do not give a definition of this function that depends on the application at hand.

Finally, the third graph product we consider is the strong product G 1 G 2 of G 1 and G 2 with yet again the Cartesian product V (G 1 ) × V (G 2 ) as its vertex set. The edge set is however the union of the edge sets of the two previous graph products. Mathematically, it translates to:

V (G 1 G 2 ) =V (G 1 ) × V (G 2 ) E (G 1 G 2 ) =E (G 1 G 2 ) E (G 1 ⊗ G 2 ) .
An example of strong product is shown on Figure 2.2(c). The definition of the weights depends again on the application at hand. We remark on Figure 2.2 that the symbol (resp. ⊗, ) actually refers to the edge pattern of the Cartesian product (resp. Kronecker product, strong product) of two path graphs.

Matrix Representations

From now on, we will assume a numbering of the vertices from 1 to N such that when we write "i ∈ V " i is an integer. We introduce now several matrix representations of graphs. In all these representations, an element of the matrix corresponds to an edge between two vertices. The numbering of the vertices are then the indices of these matrices.

Adjacency Matrix

The most basic matrix representation of a graph G = (V, E , w) is the adjacency matrix A (G ) defined as:

A (G ) i j =    w G (i j ) if i j ∈ E 0 otherwise.
By extension, we denote w (G ) i j the entries of the matrix A (G ) . If the graph is symmetric, the matrix

A (G ) is symmetric since w G (i j ) = w G ( j i ).
Using the adjacency matrix, the (out)-degree and in-degree of the vertex i ∈ V are:

d (G ) i = j ∈V w (G ) i j = (A1) i d (G ) i ,in = j ∈V w (G ) j i = 1 * A i ,
with 1 the all one column vector. The degree matrix is then defined as:

D (G ) =    d (G ) 1 . . . d (G ) N    .

Matrix Representations

Its counterpart for the in-degrees is the in-degree matrix D (G ) in whose diagonal elements are the in-degrees.

Graph Products

In [Brouwer and Haemers, 2012, Section 1.4.6], the authors show that the adjacency matrix of the Cartesian product G 1 G 2 is: .2) This addition corresponds to the disjoint union of the edges in G 1 and G 2 presented in section 2.2. Note that using this notation, the indexing of the vertices of G 1 G 2 are determined by those of the individual graphs. For example, if G 2 has N 2 vertices and i ∈ V (G 1 ), j ∈ V (G 2 ), then the vertex (i , j ) has index N 2 × (i -1) + j , with vertex indices starting at 1. The adjacency matrix of the Kronecker product G 1 ⊗ G 2 is [Brouwer and Haemers, 2012, Section 1.4.7]:

A (G 1 G 2 ) = A (G 1 ) ⊗ I + I ⊗ A (G 2 ) . ( 2
A (G 1 G 2 ) = A (G 1 ) ⊗ A (G 2 ) .
Finally, the adjacency matrix of the strong product is the sum of the two previous adjacency matrices.

Random Walk Matrix

When a graph has real non-negative weights, it is possible to define a random walk on the vertices of the graph using its weights. This is a random process whose realizations are paths on the graph. Generating a realization works as follows: Given a vertex of the path, the next vertex is iteratively chosen at random with probabilities proportional to the weights of the out-going edges of the current vertex. In other words, the probability that j is chosen after i is w i j /d i . This probability is well defined since we consider only connected graphs, such that vertex degrees are all positive.

The random walk matrix P (G ) is then the matrix such that P (G ) i j is the probability for j to be the next vertex of the path with i the current vertex. Let A (G ) be the adjacency matrix and D (G ) the out-degree matrix. The random walk matrix is then: .3) As a probability transition matrix, we have then:

P (G ) = D (G ) -1 A (G ) . ( 2
j P (G ) i j = 1 d (G ) i j w (G ) i j = 1.
The non-negative property is essential for a well-defined random walk matrix, and the graph does not need to be symmetric. Remark also that P is not necessarily symmetric even if G .

Laplacian Matrices for Symmetric Graphs

Another matrix representation of graphs that is of particular interest in this dissertation is the Laplacian matrix. This matrix is defined as the difference between the degree matrix and the adjacency matrix:

L (G ) = D (G ) -A (G ) ,
such that:

L (G ) i j =    d (G ) i if 2 i = j -w (G ) i j otherwise.
The normalized Laplacian matrix is a variation of the Laplacian matrix where weights are normalized by the vertex degrees:

L (G ) = (D (G ) ) -1/2 L (G ) (D (G ) ) -1/2 = I -(D (G ) ) -1/2 A (G ) (D (G ) ) -1/2 ,
such that:

L (G ) i j =      1 if i = j - w (G ) i j d (G ) i d (G ) j otherwise.
Obviously, L is not defined for graphs having a vertex degree equal to zero (i.e. a vertex connected to no other vertex), and may be complex if a degree is negative. This gives a first justification to the restriction in this dissertation to connected and non-negatively weighted graphs.

Graph Products First, using Equation 2.2, the degree matrix of a graph product is easy to compute: .4) with the Cartesian product on the left, and the Kronecker product on the right. It follows then: .5) Unfortunately, a similar formula cannot be devised for the Kronecker product (and then not for the strong product). Finally, the normalized Laplacian matrix does not have a simple writing with respect to the normalized Laplacian matrix of the two multiplied graphs mostly because inverting the matrix in Equation 2.4 is not a linear operation.

D (G 1 G 2 ) = D (G 1 ) ⊗ I + I ⊗ D (G 2 ) D (G 1 ⊗ G 2 ) = D (G 1 ) ⊗ D (G 2 ) , ( 2
L (G 1 G 2 ) = L (G 1 ) ⊗ I + I ⊗ L (G 2 ) . ( 2

Properties of Matrix Representations

We focus here on the three most important matrix representations of graphs, namely the adjacency matrix, the Laplacian matrix and the normalized Laplacian matrix. The results of this section are summarized in Table 2.1.

Adjacency Matrix

For directed graphs, even with non-negative weights, the adjacency matrix may not be diagonalizable. For example, the graph with 2 vertices connected by a single directed edge has the

w i j ≥ 0 A |γ i | ≤ d max L 0 ≤ λ i ≤ ρ G ≤ 2d max L 0 ≤ θ i ≤ 2
following adjacency matrix:

A = 0 1 0 0 .
This matrix is the Jordan block J 2 (0), which is not diagonalizable. Therefore, to decompose this matrix, only the Jordan decomposition and the Schur decomposition, including the block diagonalization computed from the Schur decomposition, can be used. We remark also that even with real weights, the adjacency matrix may have complex valued eigenvalues (see for example the directed cyclic graph in section 2. 5.2, page 24).

On the other hand, if the graph is symmetric, the adjacency matrix is also symmetric. Therefore, A is unitarily diagonalizable and has real eigenvalues: There exists a set of orthonormal eigenvectors u 0 , . . . , u N -1 associated to the eigenvalues 4

γ 0 ≥ • • • ≥ γ N -1 . Let U = [u 0 • • • u N -1 ]
and Γ = diag(γ 0 , . . . , γ N -1 ). We have then: .6) Negatively weighted graphs can have negative eigenvalues. For example, the graph with 2 vertices and an edge with weight -1 has the eigenvalue -1:

A = U ΓU * . ( 2
-1 A 1 1 = 0 -1 -1 0 1 1 = -1 1 1 .
Moreover, even with only non-negative weights the eigenvalues are not necessarily nonnegative. For example, -2 is an eigenvalue of the adjacency matrix of the path graph of length 3 (pictured on the left):

A 1 -2 1 = 0 1 0 1 0 1 0 1 0 1 -2 1 = -2 2 -2 = -2 1 -2 1 .

Eigenvalues Bounds

In [Brouwer and Haemers, 2012, Propositions 3.1.1 and 3.1.2], the authors show that if the graph is symmetric and non-negatively weighted, we have:

|γ i | ≤ d max ,
with d max = max i d i . Also, γ 0 = d max if and only if the graph is regular. Finally, the graph is bipartite if and only if γ N -1 = -γ 0 [Brouwer and Haemers, 2012, Proposition 3.4.1] such that the bounds -d max , d max are met only for regular bipartite graphs.

Graph Product

We can compute the eigenvalues and eigenvectors of A (G 1 G 2 ) with the eigenvalues and eigenvectors of A (G 1 ) and A (G 2 ) using Equation 2.2 and Equation 1.1:

u (G 1 G 2 ) N 2 i + j (N 2 k + l ) = u (G 1 ) i (k)u (G 2 ) j (l ) γ (G 1 G 2 ) N 2 i + j = γ (G 1 ) i + γ (G 2 ) j , with N 2 the number of vertices of G 2 .
The same is possible for the Kronecker product G 1 ⊗ G 2 , and the strong product G 1 G 2 :

u (G 1 ⊗ G 2 ) N 2 i + j (N 2 k + l ) = u (G 1 ) i (k)u (G 2 ) j (l ) γ (G 1 ⊗ G 2 ) N 2 i + j = γ (G 1 ) i γ (G 2 ) j u (G 1 G 2 ) N 2 i + j (N 2 k + l ) = u (G 1 ) i (k)u (G 2 ) j (l ) γ (G 1 G 2 ) N 2 i + j = γ (G 1 ) i + γ (G 2 ) j + γ (G 1 ) i γ (G 2 ) j .

Laplacian Matrix

Even if the Laplacian matrix is well defined for any kind of graph, its mathematical properties are tractable only for symmetric graphs with non-negative weights. In this context, L is symmetric real, and therefore unitarily diagonalizable with real eigenvalues. We denote χ 0 , . . . , χ N -1 its orthonormal eigenvectors and

λ 0 ≤ • • • ≤ λ N -1 its eigenvalues. Let χ = [χ 0 • • • χ N -1
] and Λ = diag(λ 0 , . . . , λ N -1 ). We have then: .7) Note that the eigenvalues of L are ordered increasingly contrary to those of A that are ordered decreasingly. This accounts for different properties of these representations and will be justified on the cyclic graph in section 2.5.2.

L = χΛχ * . ( 2
Using the relation x * Lx = i j w i j (x(i ) -x( j )) 2 , and w i j ≥ 0, we have ∀x, x * Lx ≥ 0. Therefore L is semi-definite positive: ∀i , λ i ≥ 0. Moreover, [Brouwer and Haemers, 2012, Proposition 1.3.7] shows that if G has k connected components then for all l < k λ l = 0 and λ k > 0. Since we supposed the graphs to be connected (i.e. k = 1), we have then λ 0 = 0 and λ 1 > 0. Additionally, in [Das, 2011], the author prove several upper bounds on λ N -1 . One in particular is of interest in this dissertation (Corollary 3.3): .8) with

λ N -1 ≤ max i 2d i (d i + d i ) = ρ G , ( 2
d i = 1 d i j w i j d j .
Equality holds in Equation 2.8 if and only if G is bipartite regular. Finally, given the maximum degree d max , we have

d i ≤ d max and d i ≤ 1 d i j w i j d max = d max , such that ρ G ≤ 2d max , with equality if G is regular.

Graph Product

We can compute the eigenvalues and eigenvectors of L (G 1 G 2 ) with the eigenvalues and eigenvectors of L (G 1 ) and L (G 2 ) using Equation 2.5 and Equation 1.1:

χ (G 1 G 2 ) N 2 i + j (N 2 k + l ) = χ (G 1 ) i (k)χ (G 2 ) j (l ) λ (G 1 G 2 ) N 2 i + j = λ (G 1 ) i + λ (G 2 ) j ,
with N 2 the number of vertices of G 2 . As seen in section 2. 3.3, there is no equivalent formula for the Kronecker product and the strong product.

Normalized Laplacian Matrix

As seen in section 2.3, the normalized Laplacian matrix may be ill-defined for directed or negatively weighted graphs. We study here its properties for symmetric non-negatively weighted graphs. For these graphs, the normalized Laplacian matrix is symmetric real, and therefore unitarily diagonalizable with real eigenvalues. We denote ψ 0 , . . . , ψ N -1 its eigenvectors and

θ 0 ≤ θ 1 ≤ • • • ≤ θ N -1 its eigenvalues. Let Ψ = [ψ 0 • • • ψ N -1
] and Θ = diag(θ 0 , . . . , θ N -1 ). We have then: .9) In [Chung, 1996, Lemma 1.6], the author shows that θ 0 = 0 and θ N -1 ≤ 2. Also if the graph has k connected components, then θ k-1 = 0 < θ k . This is proved using the corresponding eigenpairs of the Laplacian matrix:

L = ΨΘΨ * . ( 2
Lχ l = 0 ⇔ D -1/2 LD -1/2 D 1/2 χ l = 0 ⇔ L D 1/2 χ l = 0,
such that for all l < k, D 1/2 χ l is an eigenvector of L associated to the eigenvalue 0, giving a relation between χ l and ψ l when l < k. This proves that the eigenvalue 0 has the same multiplicity for L and L . Assuming a connected graph, we have then 0 = θ 0 < θ 1 . Finally, using [Chung, 1996, Lemma 1.7] ψ N -1 = 2 if and only the graph is bipartite.

Note that the correspondence between χ l and ψ l shown is only valid when the eigenvalue is 0. Showing a correspondence between other eigenvectors is far more difficult. Fortunately, such a correspondence is not useful, and the differences between the properties of the eigenvectors of L and L give rise to two alternative approaches when defining a Fourier transform on graph.

Chapter 2. Graphs

Example Graphs

In this section, we introduce several (unweighted) graphs whose structure is deterministic as opposed to those of the next section that depend on probability distributions. The graphs of this section will appear throughout this dissertation. Important results are summarized on Table 2.2. 2 cos

ω i 2 -2 cos ω i 1 -cos ω i Param. ω i = 2π i +1 2N +2 c i (k) = cos(kω i ), ω i = 2π i 2N c i (k) = cos(kω i ), ω i = 2π i 2N -2 T n Eigvect. 1, . . . , ζ k i ζ l j , . . . , ζ N -1 i ζ N -1 j T Eigval. 2 cos ω i + 2 cos ω j 2 -2 cos ω i + 2 -2 cos ω j 1 -cos ω i + 1 -cos ω j Param. ζ i = e ıω i and ω i = 2π i n G n
Eigvect. . . . , s i (k + 1)s j (l + 1), . . .

T

Eigval. 2 cos ω

i + 2 cos ω j 2 -2 cos ω i + 2 -2 cos ω j 1 -cos ω i + 1 -cos ω j Param. s i (k) = sin(kω i ), ω i = 2π i +1 2N +2 c i (k) = cos(kω i ), ω i = 2π i 2N c i (k) = cos(kω i ), ω i = 2π i 2N -2

Complete Graph

We denote K N the complete graph with N vertices having any pair of vertices connected by an edge. Notice that a complete graph is necessarily symmetric. [Brouwer and Haemers, 2012, Section 1.4.1] shows that A has the eigenvalues N -1 (associated to the all one eigenvector) and -1 with multiplicity N -1 (associated to eigenvectors δ 1δ i for i = 1). This reference shows also that L has the eigenvalues 0 (associated to the all one eigenvector) and N with multiplicity N -1 (associated to eigenvectors δ 1 -δ i for i = 1). Using the fact that L = (N -1)L , we obtain then the eigenvalues of L : 0 and N N -1 with multiplicity N -1 (same eigenvectors).

Cyclic Graph: Periodic Time Series

We denote C N the symmetric cyclic graph with N vertices and edges i i + 1 plus N 1 . Its directed counterpart is denoted C N and its directed edges are i i + 1 plus N 1 .

For symmetric cyclic graphs, the eigenvectors of A, L and L are the same, and equal to [Brouwer and Haemers, 2012, Section 1.4.3]. Note that we used the angular frequency ω i = 2πν i on Table 2.2 to save space.

u i = χ i = ψ i = 1, ζ 1 i , . . . , ζ N -1 i T , with ζ i = e ı2πν i and ν i = i N
Their eigenvalues are γ i = 2 cos(2πν i ) (matrix A), λ i = 2 -2 cos(2πν i ) (matrix L) and θ i = 1 -cos(2πν i ) (matrix L ). The adjacency matrix A of C N has the same eigenvectors, but its eigenvalues verify γ i = ζ i (remember that the matrices L and normalized Laplacian matrix are not defined for directed graphs).

We remark here that the eigenvalue ordering is crucial to observe the close relation between γ i , λ i and θ i . Indeed, each of them is written using ν i , and i is always associated to the same eigenvector u i = χ i = ψ i . Finally, the reader familiar with discrete signal processing may have recognized here the Discrete Fourier Transform (DFT) on N samples (that will be detailed in section 3.1): The eigenvectors are exactly its Fourier modes. The cyclic graph can be used to model the periodic time series with each vertex representing a time instant and each edge connecting one time instant to the next (and the previous if the graph is symmetric), thus representing a sampling period.

Path Graph

We denote P N the symmetric path graph with N vertices. This graph is essentially the cyclic graph C N with an edge removed. Its edges are then limited to i i + 1 (without N 1 ). Its directed counterpart is denoted P N with its directed edges limited to i i + 1 .

For P N , and contrary to the cyclic graph, [Brouwer and Haemers, 2012, Section 1.4.4] shows that the eigenvectors of A and L differ. The eigenvalues of A are γ i = 2 cos(2πν i ) with

ν i = i +1 2N +2 , associated to the eigenvectors u i = sin(1 × 2πν i ), . . . , sin(N × 2πν i ) T .
Using the results of the same reference, it can be shown that the eigenvalues of L are

λ i = 2 -2 cos(2πν i ) with ν i = i 2N associated to χ i = 1 + cos(2πν i ), . . . , cos (N -1)2πν i + cos N × 2πν i T .
Finally, the normalized Laplacian matrix L has the eigenvalues [Chung, 1996, Example 1.4].

θ i = 1 -cos(2πν i ) with ν i = i 2N -2
In the case of the directed graph P N , this graph is not strongly connected. In particular, no vertex is connected to vertex 1, and vertex N is not connected to any other vertex. The usefulness of such a graph is then limited. Also, the adjacency matrix is the Jordan block J N (0), therefore it is not diagonalizable and its only eigenvalue is 0.

The example of the path graph shows that by only removing one edge of a graph, namely the cyclic graph, the spectral decomposition of the matrix representations can be quite different than those of the original graph.

Torus: Image Pixel Graph

We denote T n the torus graph with N = n × n vertices arranged in a square grid layout with horizontally and vertically adjacent vertices connected to each other by edges, plus edges between the top and bottom rows and the left and right columns. We denote T n its directed counterpart. More specifically, if (i , j ) is a vertex of T n , then (i , j ) (i + 1, j ) and (i , j ) (i , j + 1) are its outgoing directed edges.

Using the fact that the torus graph T n is actually the Cartesian product C n C n of the two cyclic graphs C n (see Figure 2.2(a)), we can compute the eigenvalues and eigenvectors of both the adjacency matrix A (T n ) and Laplacian matrix L (T n ) using the relations:

u (T n ) ni + j (nk + l ) = u (C n ) i (k)u (C n ) j (l ) γ (T n ) ni + j = γ (C n ) i + γ (C n ) j
.

Chapter 2. Graphs

for the adjacency matrix, and:

χ (T n ) ni + j (nk + l ) = χ (C n ) i (k)χ (C n ) j (l ) λ (T n ) ni + j = λ (C n ) i + λ (C n ) j
, for the Laplacian matrix. For the normalized Laplacian matrix L (T n ) , using the regularity of the graph (∀i , d i = 4), we have

L (T n ) = 1 4 L (T n ) since each entry of L (T n ) is normalized by d i d j -1/2 = 1/4 to obtain L (T n ) . It follows that Ψ (T n ) = χ (T n ) and Θ (T n ) = 1 4 Λ (T n ) .
More intuitively, the eigenvectors of T n are combinations of horizontal and vertical modes (χ

(C n ) i and χ (C n ) j
), with the eigenvalues being the sum of the associated eigenvalues (λ i and λ j of the graph C n ). This example of the torus graph is fundamental as it encompasses the case of images where pixels are arranged in a grid structure. The reader familiar with signal processing on images will recognize here the 2D Fourier transform. However, eigenvalues do not represent a couple of frequencies, but rather a sum of frequencies. Therefore, part of the frequency information is lost here: for example, one cannot distinguish the frequencies ν i +ν j and ν j + ν i . This will be illustrated in section 3.5.3. The same computations can be applied to T n using C n .

Grid

We denote G n the grid graph with N = n × n vertices arranged in a square grid layout with horizontally and vertically adjacent vertices connected to each other by edges. We denote G n its directed counterpart.

Similarly to the torus graph, the eigenvectors and eigenvalues of G n and G n can be computed using those of the path graphs P n and P n since G n = P n P n and G n = P n P n .

Note that, similarly to the directed path graph, the directed grid graph is not strongly connected. The relations above show also the limits of the spectral decomposition of the adjacency matrix since limits of the path graph spectral properties are carried over to the grid graph. As a consequence the usefulness of the directed grid graph is very limited.

Example Graph Models

In addition to the deterministic graphs presented in the previous section, we define several stochastic graph models. These graph models are denoted with an underlined symbol for the model and without it for a realization of said model. For example, X is a model and X a realization of X . Realizations of these graph models are implemented in our toolbox presented in chapter 7.

Random Vertex Sampling on a Plane

We denote U N the graph model with N vertices randomly and uniformly chosen in a unit 2D square. In other words, each vertex i of a realization U N has two coordinates x i and y i uniformly chosen in [0, 1]. Choosing edges and giving weights can then be achieved using several approaches that we now explicit. First, we can use a Gaussian kernel of the Euclidean distance d(i , j ) between vertices i and j . We denote U G N a realization of U N with this weighting process such that A reads: .10) Working with all possible edges can be computationally intensive, and we generally prefer sparser graphs. The rationale behind working with less edges lies in the fact that vertices that are far from each other (i.e. with a small weight) have little influence over each other. We now show two methods to remove edges and sparsify a graph. We denote first U G,th N a realization with thresholded weights. The thresholding works by removing edges between vertices whose distance is greater than a threshold σ 2 , i.e. the weight is smaller than exp(-σ 2 2 /2σ 2 1 ). The adjacency matrix A reads then: .11) The second method consists in selecting for each vertex the edges to its k nearest neighbors (knn). We denote U G,knn N a realization with this edge selection scheme. Let KNN(i ) be the set of the k nearest neighbors of vertex i . The adjacency matrix A reads then:

Example Graph Models

w i j = e - d(i , j ) 2 2σ 2 1 , with σ 1 > 0. ( 2
w i j =      e - d(i , j ) 2 2σ 2 1 if d(i , j ) ≤ σ 2 0 otherwise. ( 2
w i j =      e - d(i , j ) 2 2σ 2 1 if i ∈ KNN( j ) or j ∈ KNN(i ) 0 otherwise.
(2.12)

Note that the set of neighbors of i in U G,knn N may not be reduced to its KNN(i ) closest vertices, but contains them. Indeed, i may be within the k closest vertices of a vertex i without j belonging to KNN(i ).

A simpler weighting scheme is to use unit weights and include the information of distance in the edge selection procedure. We denote U 1,th N and U 1,knn N the resulting graphs. However, this type of weights are prone to noise: For example if the distance between two vertices is close to σ 2 , then noise will have an influence on the appearance or not of an edge, whereas for U G,th N the edge has already a low weight such that its appearance (or not) has a low impact on Chapter 2. Graphs the overall structure.

Erdős-Rényi Graph Model

The Erdős-Rényi stochastic graph model has been introduced by Paul Erdős and Alfréd Rényi in [START_REF] Erdős | On Random Graphs I[END_REF]. We denote R N ,p the Erdős-Rényi graph model whose realization R N ,p has N vertices and each edge appears independently and uniformly with probability p.

Probabilistic spectral analysis of those graphs shows interesting properties of this model. In [Graham and Radcliffe, 2011, Theorem 7], the authors show that when pN log(N ), then the non-zero eigenvalues of L are close to 1 with high probability. In [Erdős et al., 2013, Section 2], the authors show that the eigenvalues γ l of A follow a semicircle law, centered around 0 with radius 2 N p(1 -p), except for γ 0 of higher value. Moreover, the eigenvectors of A are all delocalized on the graph. Examples of realizations of Erdős-Rényi graphs with N = 30 vertices and three probability p of edge appearance. Vertices are placed in 2D according to their connectivity using a spring model (all pair of vertices repulse each, and connected vertices attract each other with a force depending on the weight of the corresponding edge, such that the resulting stable state of the system gives the placement in 2D of the vertices).

Barabási-Albert Graph Model

The Erdős-Rényi graph model is useful to model uniform random graphs where vertices are statistically indistinguishable from one another. However, many real world graphs, and especially social networks, show a heavy tail (i.e. power law) distribution on the vertex degrees.

To model this behavior, the Barabási-Albert graph model B N ,m has been introduced in [START_REF] Albert | Statistical mechanics of complex networks[END_REF]. A realization of this model is generated using the preferential attachment algorithm.

The algorithm works as follows: Given an arbitrary graph with m vertices, at each step of the algorithm, a new vertex is added and attached to m vertices with probability proportional to the vertex degrees. This process is iterated until the graph has a total of N vertices. Naturally, the vertices with higher degrees will be favored. The resulting graph shows then a power-law distribution of the degrees.

Watts-Strogatz Graph Model

Finally, the Watts-Strogatz graph model S N ,d ,β is another model of graphs that aims at modeling real world networks [START_REF] Watts | Collective dynamics of 'small-world' networks[END_REF]]. Similarly to the Barabási-Albert graph model, generating a realization of this model is done iteratively. The procedure is as follows:

In the beginning, the N vertices are placed on a ring and connected to their d 2 neighbors on the left, and their d 2 neighbors on the right. Then vertices are considered sequentially: At each step i , all edges i j are rewired with probability β. If the rewiring accurs, the new edge i k is chosen such that it is not already an edge and it is not a loop (i = k). The resulting graph is a small-world network in the sense that the diameter of the graph is small (the maximum length of a shortest path between any pair of vertices is small).

Real Graphs

MOSAR Graph

During the course of the EU project Mastering hOSpital Antimicrobial Resistance (MOSAR) (i.e. from February 2007 to July 2012), the project produced a dataset of human contacts within an hospital. These contacts are snapshots of interactions at a distance of less than one meter every 30 seconds during several months. This represents a huge amount of snapshots that are advantageously represented using a dynamic contact graph. More specifically, we represent people using vertices and draw an edge between two people if a contact occurs.

On Figure 2.7(a), we show the five care units of the hospital with blue ellipses. These are denser regions of the graph with more contacts due to the spatial proximity of the people. shows all contacts between people within an hour. Finally, (c) shows all contacts within a day. Each edge is the support of the time series of the presence/absence of a contact. An edge appears on the figures when a contact occurs at least once during the time frame. The weighting scheme is not depicted. Ellipses on (a) correspond to group of people assigned to the same care units in the hospital.

Vertices in between care units are often staff that are not assigned to a care unit and move from one service to another.

Raw Data

Before actually using this dataset, we needed first to perform a preprocessing step to be able to work on clean data. First, the contacts were taken using RFID badges worn by people within the hospital. These badges recorded badges in their vicinity within 30 second periods. Several comments are in order here. First, the badges could not be turned off. For those worn by patients, this is not an issue. However, for the staff who are not present at all time in the hospital this leads to meaningless contacts at hours where those staff were absent. Fortunately, the dataset comprises an agenda of the staff with shifts. Second, during the course of the study, the batteries regularly died and three periods of three days were dedicated to changing these batteries. Third, noise with the radio channel leads to false contacts of lower radio power.

Therefore, the raw data needs to be cross-referenced with several data: First the assignation badge / person, and second the overall agenda of the staff. The overall procedure is the following:

1. Importation of raw data (contacts, agenda, badge assignation) 2. Deletion of contacts with lower power 3. Deletion of contacts outside the target period (actual time span of the study without the battery replacement periods) 4. Association of contacts to people (rather than badges in the raw data) 5. Deletion of staff to staff contacts with an absent staff (cross-referencing with agenda) 6. Same with staff to patient contacts 7. Merge consecutive contacts This last step is actually to overcome the hardware limitation of the badges (256 consecutive contacts for a contact entry recorded by the badge).

The original preprocessed data with this procedure showed discrepancies, and needed a complete rework. During this PhD a good amount of time have been dedicated in the beginning to first understand the pitfalls in the raw data and then rework this preprocessing using a relational database (PostgreSQL) with which it is easier to perform cross-referencing.

Real Graphs

The resulting dataset is much cleaner. It has been used in [START_REF] Martinet | Dynamic Contact Network Analysis in Hospital Wards[END_REF], Martinet, 2015, Obadia et al., 2015, Obadia, 2015], and in this dissertation.

The second part of this dataset comprises bacterial carriage information, with bacterial strains measured from weekly nose swabs. Unfortunately, the quantity of information within this dataset is not enough to test new methods of signal processing on graphs. Indeed, with twice as many patient as strains, a lot of strains were actually not present for very long in the hospital, and only a few strains spread widely. In turn, the spread data linked to these strains are essentially binary data: present v.s. not present. Within the context of signal processing, binary data are quite hard to represent and as such we will not make use of them here and focus on the contact data alone. This has been the subject of a Gretsi 2013 communication [5], an abstract [3] and a poster [1] where the objective was to study the correlation between graph wavelet atoms of [START_REF] Hammond | Wavelets on graphs via spectral graph theory[END_REF] and the diffusion of a bacterial strain.

Comments on the Resulting Graph

We denote M (t ) the cleaned up MOSAR contact graph at snapshot t . Since M (t ) is sparse and not connected, we will use an aggregation over a long enough period of time. For example, the aggregated graph can be M (t 1 , t 2 ) = t 2 t =t 1 M (t ) following the definition of Equation 2. 1. See Figure 2.7 for an illustration of a snapshot and two aggregations. Whereas the weights of M (t ) are all 1, those of M (t 1 , t 2 ) can be greater than one if at least one edge is shared by two snapshots. Therefore, the aggregating function is crucial to define meaningful weights.

Also, since the graph may not be connected after aggregation, we assume that it is split into connected subgraphs and the computations are performed on each of those subgraphs, or at least on its connected component having the maximum number of vertices.

Finally, in this dissertation, as said before, a given person is represented through a vertex of M . The corresponding vertex is always placed at the same position on the figures. Note that depending on the selected connected component, we do not need to draw all vertices, and except for Figure 2.7, only the necessary vertices will be drawn. In the representation of M , the position of each vertex is computed with all snapshots aggregated together. A spring model (see Figure 2.4) has then been used to compute the positions such that close vertices in the representation correspond to vertices with more contacts.

Molène Graph

We will be using a second real graph in this dissertation. The underlying data are weather reports in Brittany, France, during the month of January 2014. The advantage of this dataset over the previous one is that it has been published in open data5 thus allowing for reproducible results. The name of the dataset comes from the island of Molène at the far end of Brittany. We denote W the Molène graph. Vertices of W are the ground stations. Such a graph with geolocalized vertices is a good candidate for a realization of the graph model of section 2.6. 1. As described in that section, there are several possible weighting schemes for the edges. In this dissertation, we will be using the weights constructed from a Gaussian kernel of the geographical distance between ground stations. We will see in chapter 3 that these weights have good properties to define a graph Fourier transform. The weights read then:

w i j = e - d(i , j ) 2 2σ 2 1 ,
where d(i , j ) is the geographical distance between stations i and j , and σ 1 > 0. Same as in section 2. 6.1, we consider two edge selection scheme to get sparser graphs. We denote W th(σ 2 ) the thresholding scheme graph and W knn(k) the knn scheme. The effect of σ 2 is shown on Figure 2.8 (we will see the influence of knn in chapter 3).

In addition to this graph of ground stations, the dataset comprises hourly weather readings for those stations, among which temperatures, wind characteristics and rain are of interest in this dissertation. This represents 744 snapshots over the month of January 2014. This data will be studied in section 6.4.2.

Summary

We introduced in this chapter the classical graph structure and its associated definitions. In particular, we described the common matrix representations of graphs, namely the adjacency matrix, the Laplacian matrix, and the normalized Laplacian matrix. These matrices are essential to the definition of signal processing on graphs described in chapter 3. We also introduced the spectral decomposition of those matrices when the graph is symmetric:

A = U ΓU * L = χΛχ * L = ΨΘΨ * ,
with the following eigenpairs:

Au l = γ l u l Lχ l = λ l χ l L ψ l = θ l ψ l .
Finally, we reviewed classical graphs and graph models and their properties with respect to the spectral decomposition above. As the name suggests, the field of graph signal processing focuses on understanding and studying graph signals. More precisely, a graph signal is defined as values supported by the vertices of a graph. The graph structure shows then relations between these values, just like the values of a temporal signal are relating to each other through temporal relations.

Graph Fourier Transforms

Among many tools, the classical signal processing toolbox comprises several decompositions of signals into elementary signals. The most basic decomposition being the Fourier transform that decomposes a signal into oscillating components: the Fourier modes.

In the context of graph signal processing, the question that naturally arises is then the question of a decomposition of a graph signal into oscillating components, or at least into components ranging from slowly varying to quickly varying on the graph structure: the graph Fourier modes of a graph Fourier transform.

In this chapter, we recall the classical definitions of Fourier transform for temporal signals and present several definitions of graph Fourier transform from the literature. In the course we present several definitions and results related to the graph Fourier transform, ranging from a definition of graph frequency to a definition of convolution and its associated convolution theorem. Finally, we represent several Fourier modes of several graphs to give intuitions.

Chapter 3. Graph Fourier Transforms

Background: Continuous and Discrete Signals

Continuous Domain

We denote x a continuous 1D signal such that x(t ) is the value of the signal x at time t . Depending on the properties of the signal x, several forms of the Fourier transform exist: Definition 1 (Fourier Transform). The Fourier transform of the signal x(t ) and its inverse read:

x(ν) = F x (ν) = ∞ -∞ x(t )e -ı2πνt d t (3.1) x(t ) = F -1 x (t ) = ∞ -∞ x(ν)e ı2πνt d ν, (3.2)
with ν ∈ R the frequency.

If the signal x is periodic of period T , the Fourier transform has a simpler form. Indeed, we have:

x(ν) = ∞ -∞ x(s -T )e -ı2πνs d s = ∞ -∞ x(t )e -ı2πν(t +T ) d t = e -ı2πνT ∞ -∞ x(t )e -ı2πνt d t = e -ı2πνT x(ν), such that x(ν) is non zero only for νT integer, i.e. ν = k T , k ∈ Z.
Definition 2 (Fourier Transform of a T -periodic Signal).

x

(ν k ) = F T {x}(ν k ) = ∞ -∞ x(t )e -ı2πν k t d t (3.3) x(t ) = F -1 T {x}(t ) = ∞ k=-∞ x(ν k )e ı2πν k t , ( 3.4
)

with ν k = k T , k ∈ Z.
Equation 3.1 and Equation 3.3 are actually the projection of x on the Fourier basis e ν (t ) = e ı2πνt ν∈R using the canonical scalar product of the space of 1D functions, i.e.:

〈x, y〉 = ∞ -∞ x(t )y * (t )d t ,
such that x(ν) = 〈x, e ν 〉, and we have the following theorem:

Theorem 1 (Plancherel Theorem). The Fourier transform is a unitary operator. Given the two signals x and y the following equality holds: .5) The following corollary holds and states the energy preservation of the Fourier transform:

〈x, y〉 = 〈 x, y〉. ( 3
Corollary 1 (Parseval Identity).

x 2 = x 2 , with x 2 2 = 〈x, x〉.
We call operator a function whose input is a signal and output another signal. Let H be such an operator. If H is linear, i.e. if it verifies the superposition principle H {x +αy} = H {x}+αH {y}, then it can be written using its kernel h(t , s): .6) Moreover, if the kernel h is a function of the time difference ts, we obtain a convolutive operator H :

H {x}(t ) = ∞ -∞ x(s)h(t , s)d s. ( 3
H {x}(t ) = ∞ -∞ x(s)h(t -s)d s, ( 3.7)
where H {x} is the convolution of x and h. Indeed, the convolution of two signals x and y is defined as: .8) Performing the Fourier transform on the convolution of two signals yields the following theorem:

(x * y)(t ) = ∞ -∞ x(τ)y(t -τ)d τ. ( 3
Theorem 2 (Convolution Theorem).

x * y(ν) = x(ν) y(ν).

(

For a convolutive operator H , Theorem 2 interprets h as the frequency response of the operator H , while h is its impulse response (i.e. h = H {δ} where δ is the Dirac signal centered on t = 0). This framework of 1D signal processing can be extended to 2D signals, with s(x, y) a 2D signal: .11) We will use this 2D Fourier transform later to illustrate key differences with the graph Fourier transform.

Definition 3 (2D Fourier Transform). s(ν x , ν y ) = s(x, y)e -ı2πν x x e -ı2πν y y d xd y (3.10) s(x, y) = s(ν x , ν y )e ı2πν x x e ı2πν y y d ν x d ν y . ( 3

Discrete Domain

Let x(t ) be a continuous signal sampled with a time rate of τ s . This produces a discrete time series x[n] = x(nτ s ) or time series in short. Note that samples of time series are denoted with brackets [n] whereas those of continuous signals are denoted with parentheses (t ). Many of the tools of signal processing on continuous signals are carried over to the discrete framework by transforming the continuous sum into a discrete sum . We explain the subtleties and particularities of the discrete framework along the way.

The Fourier basis for time series is composed of the functions e ν [n] = e ı2π νn where ν = ν ν s is a reduced frequency. The direct adaptation of the continuous case reads then:

x( ν) = ∞ n=-∞ x[n]e -ı2π νn .
Due to the sampling process, the function x( ν) is periodic of period 1:

x( ν + 1) = ∞ n=-∞ x[n]e -ı2π( ν+1)n = ∞ n=-∞ x[n]e -ı2π νn = x( ν),
such that its study on the interval [0, 1] is sufficient. The corresponding Fourier transform and its inverse are then called the Discrete Time Fourier Transform: Definition 4 (Discrete Time Fourier Transform).

x( ν) = DTFT{x}( ν) = ∞ n=-∞ x[n]e -ı2π νn
(3.12) .13) x( .15)

x[n] = DTFT -1 { x}[n] = 1 0 x( ν)e ı2π νn d ν, ( 3
for ν ∈ [0, 1]. Finally, if the time series x[n] is periodic of period N , we have x( ν) is zero if ν N is not
ν k ) = DFT{x}( ν k ) = 1 N N -1 n=0 x[n]e -ı2π ν k n (3.14) x[n] = DFT -1 { x}[n] = 1 N N -1 k=0 x( ν k )e ı2π ν k n , ( 3
for ν k = k N .
The 1 N factors in Equation 3.14 and Equation 3.13 are necessary to ensure that DFT -1 DFT{x} = x and to verify Theorem 1, i.e. the operator DFT is unitary. As previously stated in section 2.5, we remark that the Fourier modes of the DFT are identical to the eigenvectors of the cyclic graph, with the eigenvalues of the graph as a function of the associated frequencies. This remark will be important for the introduction of the field of signal processing on graphs.

Finally, we look at time series operators. In general, the definition of a linear operator is similar to that of continuous signal. However, for N -periodic time series, we can simplify the boundaries without loss of generality. Let H be a linear operator and x an N -periodic time series. We have then:

H {x}[n] = N -1 k=0 x[k]h[n, k],
(3.16) .17) where

H {x}[n] = N -1 k=0 x[k]h[n -k], ( 3
h[-i ] = h[N -i ].
The matrix H is then circulant2 , i.e. all its rows are obtained by a shift of the elements of the row above it:

H =         h[0] h[N-1] h[1] h[1] h[0] . . . . . . . . . . . . . . . h[0] h[N-1] h[N-1] h[1] h[0]         .
As we can see, the vector h is the first column of H , coherent with it being the impulse response of H . Also, the convolution theorem is still valid and h is the frequency response of H .

Signal Processing on Graphs: Notations

In this section, we introduce a few notations that are common to all graph Fourier transform.

First of all, a graph signal is a function x : V → R or C. We assume that vertices are ordered such that x(i ) = x i with [x 1 , . . . , x N ] T a vector representing the graph signal x. Several approaches exist to define an equivalence of the Fourier transform. In all of them, the graph Fourier transform is a linear operator and as such can be represented by a matrix. Let F be the matrix of the graph Fourier transform, and F -1 its inverse:

Definition 6 (Graph Fourier Transform). The graph Fourier transform of the signal x is denoted

x. x and x verify then the following relations:

x = F x (3.18) x = F -1 x (3.19)
The columns of the matrix F -1 correspond then to the Fourier modes. Indeed, using Equation 3.19, we can show:

x i = N -1 l =0 F -1 i l x(l ).
Indices on the Fourier modes range from 0 to N -1 to stress the fact that the first Fourier mode actually represents the lowest frequency (0 in the temporal case). Notice also that we denote the l th Fourier component of a signal x(l ) as opposed to the value x i of x on vertex i to emphasize the fact that the elements of x are Fourier components. We associate each Fourier mode to a graph frequency ν l ∈ R. The exact definitions of the Fourier modes and their associated graph frequencies depends on the chosen approach and will be made explicit in the next sections.

A linear operator on graph signals can be represented by a matrix H such that the matrix multiplication H x is the output of H on the graph signal x. We then define H = F H F -1 . This matrix is the representation of H in the Fourier domain such that H x = H x. In other words, applying H to x is equivalent to applying H to x:

F -1 H x = F -1 H F x = H x.
We follow the convention of [START_REF] Shuman | The Emerging Field of Signal Processing on Graphs: Extending High-Dimensional Data Analysis to Networks and Other Irregular Domains[END_REF] to define the convolutive operator on graph signals as the linear operators verifying the convolution theorem:

H x(l ) = H x (l ) := h(l ) x(l ) , (3.20)
with l the Fourier mode index and h a graph signal. In this context, H x is the generalized convolution of h and x. It follows from the equation above that the matrix of a convolutive operator H must be diagonal: .21) In general, and contrary to time series, regardless of the vertex indices, the matrix H is not circulant, property which is specific to convolutive operator on periodic time series. This will be crucial to interpret results in chapter 5 and chapter 6.

H =    h(0) . . . h(N -1)    . ( 3
By enforcing the convolution theorem, this definition of convolution for graph signals is very interesting. Indeed, when the same graph frequency is associated to two Fourier modes we observe an interesting property. Let F l , F l +1 be these Fourier modes associated to the same graph frequency ν l . We will see that in all graph Fourier transforms we present, any linear combination of F l and F l +1 is also a Fourier mode associated to the frequency ν l . This is similar to the classical Fourier transform where cos(2πνt ) and sin(2πνt ) are of frequency ν and any linear combination of them is equal to A cos(2πνt + ϕ) and is also of frequency ν. We are therefore free to use alternative Fourier modes associated to the graph frequency ν l , for example3 : (F l + F l +1 )/2 and (F l -F l +1 )/2. For temporal signals, we can use either the basis {e ν } ν or the basis {c ν } ν ∪ {s ν } ν with c ν (t ) = cos(2πνt ) and s ν (t ) = sin(2πνt ). The relation between the two bases is obtained using Euler's formula:

e ı x = cos(x) + ı sin(x).
Unfortunately, the convolution theorem is not valid anymore on the basis of cos and sin:

H {c ν } = h(ν) + h(-ν) 2 c ν - h(ν) -h(-ν) 2ı s ν H {s ν } = h(ν) -h(-ν) 2ı c ν + h(ν) + h(-ν) 2 s ν
where H is a convolutive operator, and

H {e ν } = h(ν)e ν .
Similarly, when a graph frequency is associated to multiple Fourier modes, the convolution theorem (i.e. the matrix H is diagonal for a convolutive operator H ) is only valid for the chosen Fourier basis. A different Fourier basis leads to a different convolution theorem. Equivalently, a different Fourier basis leads to a different family of convolutive operators. Nevertheless, when a convolutive operator H is such that h(l ) = h(l + 1) for the two Fourier modes F l and F l +1 associated to the graph frequencies

ν l = ν l +1
, then this operator verifies the convolution theorem no matter the chosen Fourier basis since:

H (αF l + βF l +1 ) = αH F l + βH F l +1 = h(l )(αF l + βF l +1 ),
i.e. any linear combination of F l and F l +1 is an eigenvector of H associated to the eigenvalue h(l ). We can see the cyclic graph C N as a uniform sampling of a time interval, i.e. each vertex is a time instant with regularly spaced time instant in time, and two successive time instant are linked by an edge. First, let f be a T -periodic temporal signal. The three point stencil approximation of the continuous Laplacian ∆ f (t ) = f (t ) is then given by:

Laplacian Based Fourier Transform

∆ f (t ) ≈ f (t -h) + f (t + h) -2 f (t ) h 2 .
We can actually see the right-hand side of this equation as the Laplacian matrix applied to a graph signal. Let C N be the cyclic graph with N = T h and edges weighted by 1 h 2 . Using the periodicity of f , and the construction of C N we can show a :

(L f ) i = 2 f i -f i -1 + f i +1 h 2 , with f i = f (i h). It follows that: ∆ ≈ -L.
The same formula can be derived using the Euclidean plane, the torus graph T n and a five point stencil:

∆ f (x, y) ≈ f (x -h, y) + f (x + h, y) + f (x, y -h) + f (x, y + h) -4 f (x, y) h 2 (L f ) i , j = 4 f i , j -f i -1, j + f i +1, j + f i , j -1 + f i , j +1 h 2 .
Moreover, classical signal processing states that the continuous Fourier modes are the eigenvectors of the continuous Laplacian ∆ associated to squared frequencies:

∆e ν (t ) = -(2πν) 2 e ν (t ).
a We suppressed the modulo N on the indices for clarity.

Observation 2 [START_REF] Belkin | Towards a theoretical foundation for Laplacian-based manifold methods[END_REF] The authors of this communication prove that if a finite subset of a Euclidean space a is randomly sampled with uniform distribution, we can construct a graph with the N sampled points and weighted using a Gaussian kernel of the Euclidean distance (as in section 2. 6.1). The authors define then the quantity L N f (x) where L N is the point cloud Laplace operator and f a function on the space:

L N f (x) = f (x) 1 N j e - d 2 (x,x j ) 2σ 2 1 + 1 N j f (x j )e - d 2 (x,x j ) 2σ 2 1 ,
with x i the coordinates of the vertex i and x an arbitrary point of the space. It follows that

L N f (x i ) only depends on f = [ f (x 1 ), . . . , f (x N )]
T and the Laplacian matrix of the graph

with L N f (x i ) = (L N f ) i .
Under mild assumptions b , L N f converges to -∆ f as the number of sampled vertices N goes to infinity, where ∆ is the Laplace-Beltrami operator of the Euclidean space. Essentially, the Laplacian matrix converges to the Laplace-Beltrami operator. This formally justifies the observation above and the graph Fourier transform below.

a The exact hypothesis on the space is much weaker but not useful here. b such as the parameter σ 1 of the Gaussian kernel depends on N , and the Laplacian matrix L is rescaled by a factor depending on N and the dimension of the space

Definitions

Based on these observations, the graph Fourier transform is defined as the projection on the eigenvectors of the matrix L. Following the properties seen in section 2.4, namely L is semi-definite positive, we obtain the following definition: .23) such that x = F x and x = F -1 x:

Definition 7 (Laplacian-based Graph Fourier Transform). F =   χ 0 * . . . χ N -1 *   = χ * (3.22) F -1 = F * = [ χ 0 • • • χ N -1] = χ, ( 3
x(l ) = N i =1
x i χ * l (i )

x i = N -1 l =0 x(l )χ l (i )
Theorem 3 (Plancherel Theorem). The graph Fourier transform is unitary, and:

〈x, y〉 = 〈 x, y〉.

Proof. We have F = χ * , and χ unitary.

We also define the graph frequency using the analogy above [7]:

Definition 8 (Graph Frequency). If λ l is the eigenvalue associated to the eigenvector χ l , then the associted graph frequency is:

ν l = 1 2π λ l .
Using Definition 8, we obtain a perfect matching with the continuous Laplacian (up to the sign):

Lχ l = (2πν l ) 2 χ l .
Finally, similarly to the DFT where reduced frequencies are used to match the frequency interval [0, 1], we can rescale the graph frequencies into reduced graph frequencies using the upper bound ρ G on the eigenvalues of L:

Definition 9 (Reduced Graph Frequency). If λ l is the eigenvalue associated to the eigenvector χ l , then the associated reduced graph frequency is:

ν l = 1 2 λ l ρ G .
There are several comments to make on this definition. First, the factor 1 2 is the consequence of not being able to choose between positive and negative graph frequencies. Indeed, in the DFT, the reduced frequencies ν k lie in the interval [0, 1) such that the complex phases 2π ν k lie in the interval [0, 2π). Negative frequencies are then those in the interval [π, 2π), which is the same as the phases of [-π, 0), hence the name of negative frequencies. For time series we know that frequencies are coupled by opposite frequencies. For example, e ν + e -ν gives the signal 2 cos(2πνt ). Unfortunately, in the graph setting, the sign of the graph frequencies is unknown since we have access only to the squared graph frequencies through the eigenvalues of L. In this dissertation, we arbitrarily choose positive signs to graph frequencies and leave the study of a possible better choice of sign or even the meaning of a negative graph frequency for future works.

The second comment is on the normalization by the upper bound ρ G on the eigenvalues of L in Definition 9. As seen in section 2. 4.2, this upper bound is met only for bipartite regular graphs. The main example of such a graph is the cyclic graph C 2N with an even number of vertices. The corollary is that the bound is not met for C 2N +1 . This is coherent with the DFT where the reduced frequency ν N /2 = 1 2 exists if and only if N is even.

Discussion

First of all, the Laplacian matrix-based approach can also be used with the normalized Laplacian matrix instead of the Laplacian matrix. This leads to: .24) The associated reduced graph frequency are then: .25) With respect to Observation 1, and the fact that the graphs therein are regular, the normalized Laplacian matrix verifies L = h 2 2 L, with h the time between two samples (see Observation 1 for details) or L = h 2 4 L for the torus graph. Therefore, the first observation can intuitively justify working with either L or L. [Belkin and Niyogi, 2008, Section 5] provides also a convergence result of both the Laplacian matrix and the normalized Laplacian matrix in the context of sampled Euclidean space with an arbitrary sampling distribution, such that the second observation holds also. More precisely, given P the probability density function of the sampling process on the space (remember that in Observation 2, this sampling process was uniform), we have: .26) where B and C depend on the smoothness of P and f in the space (if either is smooth B and C are close to zero, e.g. a uniform probability leads to P constant and then B (x) = 0 = C (x)). In particular, we observe that the scaling factor P (x) is absent in the limit form corresponding to the normalized Laplacian matrix L . In other words, the normalized Laplacian matrix approximate the Laplace-Beltrami operator with an additive term, contrary to the Laplacian matrix that approximate it with a multiplicative term (and an additive term). Therefore, the normalized Laplacian matrix should be favored as soon as the sampling is not uniform or unknown.

F = Ψ * . ( 3
ν l = 1 2 θ l 2 . ( 3
L N f (x) ----→ N →∞ -P (x)∆ f (x) + B (x) L N f (x) ----→ N →∞ -∆ f (x) +C (x), ( 3
The second remark is the intrinsic limitation of the Laplacian matrix and the normalized Laplacian matrix: At the moment of writing this dissertation, strong results on these matrices are only available for non-negatively weighted symmetric graphs. Transportation networks (e.g. with one-way roads), or gene regulatory networks (e.g. genes producing proteins that in turn regulate positively or negatively the expression of those genes) are therefore not directly supported by this definition of graph Fourier transform.

The third remark is about Definition 8. In the absence of a better choice, all graph frequencies are chosen to be non-negative: Choosing between positive and negative graph frequencies for each of the eigenvalues is an open question that still needs to be elucidated. To give an intuition on the complexity of the task, we have to go back to the DFT where Fourier modes are almost all coupled by opposite frequencies. In the graph setting this coupling occurs only for very specific graphs such as the cyclic graph. Most graphs have Laplacian matrices with eigenvalues of multiplicity one.

In the rest of this dissertation, and unless otherwise stated, we use the graph Fourier transform defined using the Laplacian matrix such that F = χ * (Definition 7).

Adjacency Matrix Based Graph Fourier Transform

This approach has a different premise than the one based on the Laplacian matrix or the normalized Laplacian matrix. For the sake of completeness, we describe it and we will try to compare the two approaches as often as possible.

Adjacency Matrix Based Graph Fourier Transform Observation 3 [Sandryhaila and Moura, 2013]

The adjacency matrix of the directed cyclic graph is exactly the time shift operator of periodic time series with (Ax) i = x i -1 .

In turn, the Fourier modes e ν are eigenvectors of the time shift operator for periodic time series associated to the complex numbers e -ı2πν as eigenvalues.

Definition

The idea of the authors of [START_REF] Sandryhaila | Discrete Signal Processing on Graphs[END_REF] and [Sandryhaila and Moura, 2014a] is then to call A the graph shift and consider its eigenvectors as Fourier modes. This gives the following definition of graph Fourier transform:

Definition 10 (Adjacency Matrix Based Graph Fourier Transform). Given U the matrix of eigenvectors of A, we have:

F = U -1
(3.27)

F -1 = U (3.28)
Contrary to Definition 7, we cannot give a more explicit formula than x = F x. Indeed, if the graph is directed, then F is not unitary, i.e. F -1 = F * . On the other hand, if the graph is symmetric, then A is unitary, and F is unitary. We have therefore

F -1 = F * .
In [START_REF] Sandryhaila | Discrete Signal Processing on Graphs[END_REF] the authors suggest using the Jordan decomposition to compute the spectral decomposition of the graph shift. We saw however that this decomposition is not numerically stable. We will therefore use in this dissertation the block diagonalization using the Schur decomposition 4 . For eigenvalues with multiplicity 1, the eigenvectors are exactly the same between the two decompositions. For other eigenvalues, the important definition is that of stable subspace. Let γ l be an eigenvalue of A with multiplicity r > 1. Then the set of associated eigenvectors and generalized eigenvectors forms the subspace:

S l = l +r k=l α k u k : ∀k ∈ {l , . . . , l + r }, α k ∈ C ,
assuming that ∀k ∈ {l , . . . , l + r }, u k is an eigenvector or a generalized eigenvector associated with the eigenvalue γ k = γ l . S l is stable means that ∀x ∈ S l , Ax ∈ S l . In other words, S l is the subspace associated to the eigenvalue γ l . The Jordan decomposition and the block diagonalization using the Schur decomposition have the same stable subspaces associated to each eigenvalue. Moreover, for a given stable subspace, the associated set of eigenvectors or generalized eigenvectors is not unique in either decomposition such that any difference between these decompositions is not relevant.

Discussion

First of all, the meaning of the eigenvalues of the adjacency matrix A is unclear. Indeed, in the temporal setting, the frequencies are linked to the phase of the eigenvalues e ı2πν of the time shift. In general, the eigenvalues of the adjacency matrix are not unimodular complex exponential of the form e ı2πν , such that Observation 3 reaches its limits when linking frequencies and eigenvalues. Some interpretation is presented in [Sandryhaila and Moura, 2014a], where the authors show that smaller (resp. higher) real-valued eigenvalues correspond to bigger (resp. smaller) variation on the graph, where the considered metric is an adaptation of total variation to graph signal: TV G (x) = x -Ax 1 . Unfortunately, this does not solve the question of the complex-valued eigenvalues. Therefore, the question of graph frequencies remains.

The second drawback of the approach is the matrix F = U -1 itself. This matrix is in general neither unitary, nor proportional to a unitary matrix. Therefore, this definition of graph Fourier transform does not verify Plancherel Theorem or Parseval's identity.

Despite these drawbacks, the clear advantage of this approach is to allow the study of graph signal on directed graphs, with possibly negative weights, or even complex-valued weights since no structure is required on A to define this graph Fourier transform.

Examples of Fourier Modes

In this section, we look at the Fourier modes of a few graphs and graph models. The example studied here should give intuitions as to the meaning of the graph Fourier transform and interpret the weights on a graph. Second, as proven in section 2.5, the eigenvectors of all three graph Fourier transforms are sinewaves, and the eigenvalues are all of multiplicity 2, except for the first (and last if N is even). A close observation of modes 3 and 4 on Figure 3.1 (i.e. second and third columns) reveals that they are merely sine functions with a phase difference of π 2 , i.e. sine and cosine functions. More formally, the stable subspaces are the same for corresponding eigenvalues, the difference being the choice of eigenvectors to describe these subspaces.

Cyclic Graph: Periodic Time Series

Another interesting property of the cyclic graph and its periodic time series counterpart is the comparison between the Fourier transform and the Fourier series whose expression can be summarized using (N -periodic time series):

x[k] =          a 0 2 + ∞ n=1 a n cos 2πn N k + b n sin 2πn N k ∞ n=0 c n e ı 2πn N k
, with: The second and third columns correspond to 2 modes with equal graph frequencies. Finally, the fourth column corresponds to a high graph frequency.

a n = 2 N N -1 k=0 x[k] cos 2πn N k , b n = 2 N N -1 k=0 x[k] sin 2πn N k , c n = 1 N N -1 k=0 x[k]e -ı 2πn N k . L (a) λ 0 = -0.00 (b) λ 3 = 0.17 (c) λ 4 = 0.17 -0.3 0.0 0.3 (d) λ 9 = 1.00 L (e) θ 0 = -0.00 (f ) θ 3 = 0.09 (g) θ 4 = 0.09 -0.3 0.0 0.3 (h) θ 9 = 0.50 A (i) γ 0 = 2.00 (j) γ 3 = 1.83 (k) γ 4 = 1.83 -0.3 0.0 0.3 (l) γ 9 = 1.00
Both decomposition are performed with respect to an orthogonal family of time series, namely the cos 2πn N k n ∪ sin 2πn N k n and {e n/N (k)} n . This is the same for the graph Fourier modes.

But what is more interesting, is the impact the choice of the Fourier basis has on the convolution theorem. As seen in section 3.2, the convolution theorem for temporal signals is specific to the exponential Fourier basis. Since the matrices we manipulate are all real and unitarily diagonalizable, computing the graph Fourier transform through the eigenvectors necessarily yields real Fourier modes, and as such the basis of sine and cosine instead of the exponential basis. As such, and this is important, the convolution theorem for the cyclic graph and for N -periodic time series do not agree.

Cyclic Graph: The Importance of the Weights

The goal of this section is to study the influence of the weights on the graph Fourier transform. To this end, we consider periodic temporal signals, and sample the time line at irregular time instant. Figure 3.2 depicts this sampling with samples closer to each other at the beginning and the end of a period compared to samples in the middle of the period. We model this case study with a cyclic graph connecting the samples (the last sample is connected to the first through an edge to model the periodicity of the signal). In Figure 3.2 and Figure 3.3, the resulting cycle has 14 vertices. Note that on those figures, the cycle is unwrapped so as to better visualize the graph signal. The weights are then set using different schemes:

• C 14 denotes the cyclic graph with unit weights, • C (d) 14 denotes the cyclic graph with weights equal to the time difference between vertices, • C (G) 14 denotes the cyclic graph with weights equal to a Gaussian kernel of the time difference (σ 1 = 1).

• K (G) 14 denotes the complete graph with weights equal to a Gaussian kernel of the time difference between vertices (σ 1 = 1).

Our purpose in this framework is to observe clean sinewaves despite the non-uniform sampling. First of all, on C 14 , the sinewaves are stretched out where vertices are more distant, and compressed where vertices are closer. This is expected since weights are unitary. The second graph C (d) 14 shows that using the distance as weight is even worse, with compression / dilation particularly obvious for the eigenvectors of L and L (the eigenvectors of A seems to better behave although the mode u 2 is odd looking with two positive parts and one negative). Third, the graph C (G) 14 using the weighting scheme with a Gaussian kernel yields results closer to sinewaves especially for L and L. Finally, the graph K (G) 14 corresponding to the case study of [START_REF] Belkin | Towards a theoretical foundation for Laplacian-based manifold methods[END_REF] shows the best of all results for L as expected from a non-uniform sampling. Note however, that the sinewaves are not perfectly recovered, most certainly because of the additive term in Equation 3.26 (page 42).

C (G) 14 K (G)
This toy example shows first that weights based on similarities rather than distances give better results (C (d) 14 vs. C (G) 14 ). The second point is following the convergence result of [START_REF] Belkin | Towards a theoretical foundation for Laplacian-based manifold methods[END_REF], a Gaussian kernel of a distance as weights associated to the graph Fourier transform based on normalized Laplacian matrix gives the best results (graph K (G) 14 ).

Torus Graph: Image Pixel Graph

Vertices of a torus graph can be interpreted as the pixels of a square image. Indeed, the Fourier modes associated with an image are often taken as a combination of horizontal and vertical Fourier modes, just as in Definition 3. This section aims then at identifying key differences between the 2D Fourier transform and the graph Fourier transform.

The main difference lie in the Fourier modes. We can see on the fourth row of Figure 3.4 that the Fourier modes defined as a superposition of vertical and horizontal modes are not the Fourier modes we obtain by either of the three graph Fourier transforms. To understand this, suppose we work on the Laplacian matrix (first row). We know that the graph signals of the fourth row are eigenvectors of L since eigenvalues are the sum of eigenvalues of C n . Unfortunately, we already know that most of the eigenvalues of C n are of multiplicity 2. Given two eigenvalues ξ = ζ, the sum ξ + ζ can be obtained using the eigenvalues of C n with at least 8 different manners:

T n = C n C n , i.
(λ l , λ k ), (λ l , λ k+1 ), (λ l +1 , λ k ), (λ l +1 , λ k+1 ), (λ k , λ l ), (λ k , λ l +1 ), (λ k+1 , λ l ), (λ k+1 , λ l +1 ), with ξ = λ l = λ l +1 and ζ = λ k = λ k+1 .
In terms of stable subspaces, this means that the eigenvalue ξ + ζ is associated to at least 8 distinct eigenvectors. The diagonalization of L is therefore not unique and leads to strange yet correct eigenvectors as those of the first row of Figure 3. 4. The same justification works for the normalized Laplacian matrix and the adjacency matrix.

The non uniqueness of the choice of eigenvectors for the torus graph is the consequence of two properties. First, the cyclic graph C n already shows non uniqueness in the eigenvectors given an eigenvalue (see Table 2.2, with i and Ni , i.e. cos(ω i ) = cos(ω N -i )). Second, the torus graph is the Cartesian product of two identical graphs, such that the commutativity of the sum leads to at least twice the multiplicity of the eigenvalues.

In general, the graph products seems to be a very nice approach for computing Fourier modes efficiently since it only requires to compute the Fourier modes of the prime graphs to obtain the Fourier modes of their product. Such an approach is used in [Sandryhaila and Moura, 2014b]. However, its usefulness is limited to graphs which are known to be products of graphs. Indeed, most graphs are already prime graphs such that they do not benefit from the lower complexity of the spectral decomposition of a graph product. Moreover, if a real graph is a graph product, noise on the observation of the graph can easily make the observed graph a prime graph. Noise can here be the loss of an edge, or if weighted, a noisy observation of its weight, and these two types of noise lead to breaking the structure of a graph product.

Molène Graph: Sparsifying Graphs

The graph Fourier transforms we studied in the previous sections were those associated to regular structures that do not result from observations. We now focus on a real graph, and more specifically on the impact of removing edges of a graph. For that purpose, we use the Molène graph W . We defined this graph as the geographical graph of weather stations, with weights equal to a Gaussian kernel of the geographical distance between vertices. Naturally, we would like to remove edges with a limited impact on the Fourier modes. For that purpose we use two different schemes.

The first scheme is a thresholding scheme as presented in section 2.5: given a threshold σ 2 , all edges at distance greater than this threshold are removed from the graph. Note that this is equivalent to edges of weight less than exp(-σ 2 2 /2σ 2 1 ) being removed. Figure 3.6 shows then the impact of σ 2 on three Fourier modes compared to the original graph W of Figure 3.5 (left column). We can observe that a threshold of 100 kilometers has close to no impact on the spectral decomposition of the matrices, with a density (proportion of edges compared to the graph with all possible edges) of the graph reducing from 1 to 0.542. Increasing σ 2 to 50 km further reduces the density to 0.168 and we can observe some impact on the spectral decomposition.

The theory of matrix perturbation is able to give bounds on the error between the adjacency matrix A (t h) of the thresholded graph and the adjacency matrix A of the complete graph using the relation:

A = A (t h) + M , with = exp(-σ 2 2 /2σ 2 1
) and M a matrix with non-negative entries in [0, 1].

The second scheme is the k nearest neighbors scheme where each vertex selects the edges to its k nearest neighbors w.r.t. the geographical distance. The selected edges define the edges of the knn graph. Note that vertices may have more than k neighbors since a vertex i may be among the k closest vertex of j without the converse being true. When the graph is a uniform sampling of manifold, this scheme should give similar results to the thresholding scheme. However, when the graph degrees vary greatly among vertices, the two schemes lead to different graphs: A vertex with small degree in the complete graph will have its k closest neighbors relatively far away compared to a vertex with high degree. In the case of the Molène W th(100 km)

W th(50 km) graph, we can see that the knn scheme and the thresholding scheme lead to similar impacts on the Fourier modes.

MOSAR Graph: Aggregating Dynamic Graphs

This graph is more complex to study due to its nature of dynamical contact graph. Therefore, to apply the framework of graph signal processing presented above two determining steps needs to be performed beforehand:

1. Aggregation: Signal processing over graphs does not cover the dynamic graph case yet. Aggregation over a time period is therefore necessary. 2. Weighting: Given the quantities above, one needs to compute relevant weights depending on the aggregation process above.

We use in this section the same aggregation period as in Figure 2.7(c): over a day. This allows for a bigger connected component, while being able to have several graphs to compare their Fourier modes. A daily aggregation allows also to keep causalities in the graph: contacts occur at precise instants and are ordered in time such that aggregation looses this ordering. Note that from one day to another the size and the nodes of the biggest connected component may vary since the personnel in the hospital may not be the same, and some days may be less busy than others (see Figure 3.7). For completeness, we also show in Figure 3.8 a daily versus weekly aggregation with respect to the Fourier modes.

In Figure 3.7 and Figure 3.8, weights are set according to the total contact duration. Those graphs are denoted M (sum) (t , t + δt ). This quantity is a similarity rather than a distance since longer contacts correspond to closer people. Another weighting scheme is to consider the number of contacts between two people rather than the total duration 5 . The resulting graph is denoted M (nb) (t , t + δt ) (see Figure 3.9). We could devise alternative weighting schemes 5 A contact occurs at a time instant t and ends at time t + δ t . The total contact duration corresponds to the sum of δ t while the number of contacts corresponds to the number of time instant t . based on the particular content of this dataset, however the non-reproductibility of the results (the data is not public) and the specific nature of the dataset limit the interest of the dataset in the context of developing new methods for the signal processing on graphs toolbox.

M (sum) (t 1 , t 1 + 1 day) M (sum) (t 2 , t 2 + 1 day) -0.3 0
Overall, this social graph shows an example of a graph whose structure is not the image of a Euclidean space (or at least none that is obvious). The observation we make on this graph will then be applicable to many, if not all, social graphs since they share one important property of having power law distribution on the degrees. The key message here lies in the support of the Fourier modes. As we can see on Figure 3.7 through Figure 3.9, many Fourier modes are highly localized around specific vertices, even for modes with a low index that according to the intuition should be slowly varying. Also, a comparison of the three graph Fourier transforms shows that the adjacency matrix-based one has almost all its Fourier modes highly localized whereas the two based on the Laplacian matrix and the normalized Laplacian matrix have some Fourier modes much more evenly spread on the graph.

M (sum) (t 1 , t 1 + 1 day) M (nb) (t 1 , t 1 + 1 day) -0.3 0

Summary

We introduced in this chapter basic definitions of signal processing on graphs, with graph signals being defined as values supported by vertices of a graph:

x : V → R or C.
We then presented several definitions of the graph Fourier transform F based on either the adjacency matrix, the Laplacian matrix, or the normalized Laplacian matrix. The graph Fourier transforms are summarized into the following three relations:

F = U * F = χ * F = Ψ * .
We also presented the convolution operation x * y = F -1 x. y , where . is the componentwise product of vectors). Finally, we looked at a few graph Fourier modes on several graphs to give intuitions on what they are, and on how to give weights to a graph such that these Fourier modes agree with the intuition. In this chapter, we are interested in studying the structural properties of graphs. For that purpose, we take the path of leveraging the power of the signal processing toolbox for time series by embedding the structural properties of a graph into a time series. In other words, we are seeking a mapping from graphs to time series. We base our approach on a dual mapping proposed in [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] and improve upon it. In section 4.2, we present our approach based on Graph-Based Semi-Supervised Learning (GSSL) and achieving better smoothness in the resulting time series. In section 4.3, we interpret GSSL as a filter on graph signals, thus giving a frequency interpretation to the approach.

Graph Analysis

Background

Graph to Time Series Duality

In [START_REF] Campanharo | Duality between Time Series and Networks[END_REF], the authors introduced a stochastic approach to the problem of mapping time series to graphs and its inverse problem of mapping graphs to time series. The premises of these dual questions are twofold. First, we want to characterize time series through the structure of their dual graph using tools of graph theory. Second, we want to characterize the structure of graphs through the spectral properties of their dual time series.

In this dissertation, we are interested in the latter: studying the structure of a graph through the spectral properties of its dual time series.

We begin though with the presentation of the mapping from time series to graphs since this transformation gives insights on how to perform the inverse transformation. The first step is to define what will be the vertices of the graph. To do so, the statistics on the amplitudes of the time series, i.e. the set A = {x[n]} n , are used using a quantization of A into Q quantiles. In other words, A is partitioned into the intervals [a 0 , a 1 ], (a 1 , a 2 ], . . . , (a Q-1 , a Q ] such that there are the same number of samples x[n] into each quantile (interval). We now have the vertex set V = {1, . . . ,Q}. One thing is missing from the graph: the edges (and their weights). We remark first that the choice of quantiles in the first step actually leads to a mapping of each sample x[n] to a vertex i n ∈ V . We recognize in {i 0 , i 1 , . . . , } a path on the graph. The idea of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] is to consider this path as a random walk on the graph. An edge represents therefore the probability to move from one vertex to another. The random walk matrix P is therefore defined by the transition statistics of i n ∈ V . Assuming a constant vertex degree of 1, i.e. D = I , and using P = D -1 A, we obtain the adjacency matrix of the graph through P = A. In other words, all directed edges of the graph appear with a weight equal to the transition probability from the source vertex to the target vertex. Note that loops may appear in the resulting graph since i n = i n+1 if x[n] and x[n + 1] belong to the same quantile. Note also that the graph may not be symmetric since the statistics of the transition from one quantile to another may not be the same than their converse.

We can see the vertices of the graph as the image of the first order statistics of the time series, and the directed edges as the image of the second order statistics. Of course, for the graph to correctly represent the time series through these two statistics, one would need some hypothesis on the time series, for example some kind of stationarity in the sense that there should be no dependence in n when transitioning from an amplitude to another. [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] does not explicit the exact hypotheses, and our goal in this dissertation is to use the inverse mapping from graphs to time series. As such, we will not go further on this topic. [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. Given the input time series (a), the amplitudes are split into quantiles 1 through 5. The statistics are then computed to obtain the graph (b) (the loops are omitted). Finally, using a mapping from vertices to amplitudes, a time series (c) can be generated through a random walk on the graph.

The goal of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] was to build both the direct mapping and its inverse. The inverse is built using a symmetric approach. First, a random walk is performed on the graph. We obtain then the path i 0 , i 1 , . . . , i K of length K . The missing ingredient is then the mapping from vertices to amplitudes. This mapping is done in two steps. First, the vertices are sorted such that a(i ) is the new index of vertex i . We will come back to this central step in the next paragraph. The second step is as follows. The amplitudes are split into Q intervals [ j -1, j ].

The interval [ j -1, j ] is assigned to the j th vertex according to the ordering of the first step, i.e. the vertex i verifying a(i ) = j . These intervals sum to the interval [0,Q]. The n th sample of the time series is then chosen uniformly at random within the interval [a(i n ) -1, a(i n )]. These regular quantiles are chosen arbitrarily by the authors because we lack a priori knowledge on them.

We now focus on the first step we skipped above: the ordering of the vertices. Again, without a priori knowledge on the ordering, the authors assume smoothness in the resulting time series. If we consider the directed path graph P N , the ordering is quite obvious: One end of the path graph is the first vertex, the other end is the last vertex. Now, if we consider the directed cyclic graph C N , this is not obvious anymore (see Figure 4.3(a)). Indeed, if we order the vertices according to the directed edges, there will be a strong discontinuity in the resulting time series when the random walk goes from the last vertex of the ordering to the first. When the graph is more complex, with vertices with many neighbors, a good ordering ensuring smoothness is not at all obvious.

The idea of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] is to optimize the ordering a(i ) by minimizing the following cost function:

Cost(a) = i j ∈E w i j a(i ) -a( j ) . ( 4.1)
The rationale is simple: Close vertices should be mapped to close intervals of amplitudes, hence close in the ordering. The authors use a simulated annealing method to get an approximate solution of this optimization problem. Unfortunately, the exact parameters, or the code for that matter, were not published. We tried to implement this approach as best we could using the little information we had 1 . Two examples of resulting amplitude mappings and time series are shown on Figure 4.2 with the torus graph T 10 and the Molène graph W .

Drawbacks of [Campanharo et al., 2011]

As introduced, our goal is to map a graph to a time series and use the tools dedicated to the study of time series to understand and characterize the structural properties of the graph.

Unfortunately, the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] shows several drawbacks that we would like to address:

• Complexity: The optimization approach using Equation 4.1 is an optimization problem over permutations. As a consequence it is not solvable in polynomial time 2 , and the authors use an approximation algorithm, namely the simulated annealing, to get an approximated solution. This is costly, and computing the amplitudes of the graph with only 100 vertices and 200 edges of Figure 4.2(a) took almost 2 hours on an Intel Ivy Bridge core i7 laptop. • Regular Quantiles: The authors of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] propose a mapping that reduces to equally spaced quantiles. Unfortunately, the quantization step of the transformation from time series to graphs does not produce in general regularly spaced quantiles. As such, it is impossible to recover the input time series from the graph. T 10 and the Molène graph W using [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. The quantities 1 and 2 quantify the smoothness of the amplitudes on the graph and will be used to compare to our approach. This is explained in section 4.2. The color scale of (a) and (c) represent for each vertex i the value of a(i ).

• Non-uniqueness: For some graphs, the optimal solution may not be unique, or smooth. This is especially the case for the cyclic graph C N (see Figure 4.3). • Non-Optimal Smoothness: We will show that we can achieve better smoothness with our approach, especially for C N .

Our idea in section 4.2 is to relax the constraint of the mapping a being a permutation of the vertex indices. The only constraint we give is a(i ) ∈ [-1, 1]. This gives an implicit ordering of the vertices by increasing values of a(i ). With this, we see that the intervals [a(i ) -1, a(i )] are not suitable anymore. Without loss of generality, suppose the vertex indices verifies a(i ) ≤ a(i + 1). We then use the intervals a(i -1)+a(i ) 2 , a(i )+a(i +1) 2 for vertex i ∈ {2, . . . ,Q -1}, and a(1) + a(2)-a (1) 2 , a(1)+a (2) 2 for vertex 1, and

a(Q-1)+a(Q) 2 , a(Q) + a(Q)-a(Q-1) 2
for vertex Q. Amplitudes a(1) and a(Q) are therefore at the center of the intervals of vertices 1 and Q.

We now have all the ingredients, except for the mapping a, to generate a time series from a random walk. For a, our goal is to learn the mapping using a graph-based approach we now present.

Graph Based Semi-Supervised Learning

Machine learning is a vast field that can be summarized into the question of learning from data. We can subdivide machine learning approaches into three categories:

• Supervised Learning: We are given data samples, with each sample associated to a desired output. The task is then to generalize this association to any input. • Unsupervised Learning: We are given data samples, without any knowledge on the desired output. The task is to infer meaningful output from the structure. • Semi-Supervised Learning: We are given data samples, with a few of them associated to the desired output. The task at hand is to learn how to map the data to a meaningful output from both the structure of the data samples and the few desired output. 4. 1. Color scale indicate the value of a(i ) of the mappings for each vertex i . The 1 quantity corresponds to the cost function, and is equal in all three cases. The first two solutions are hand made while the third is obtained through the simulated annealing method of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF].

For example, when the output is the class of the input, we are in the presence of a classification problem.

In [START_REF] Avrachenkov | Generalized Optimization Framework for Graph-based Semi-supervised Learning[END_REF], the authors introduce a formalism generalizing several semisupervised learning techniques applied to classification on graphs and call it Graph-Based Semi-Supervised Learning. More precisely, the task is to find the class of all vertices given the graph (the structure) and the class of a few vertices (the desired output). To that end, their approach is divided into two steps. First, for each class, the level of belongings of a vertex to a class is computed. Second, the classification step consists in classifying a vertex into the class it belongs the most to. We now formalize this approach.

First, we are given a few vertices already classified into C classes. We denote x (c) the 0/1 vector with x (c) i = 1 if vertex i is known to be in class c, and x (c) i = 0 otherwise. Next, the level of belongings of vertex i to class c is denoted y (c) i . For each class c, we have then an input vector x (c) and an output vector y (c) . The classification step is given by:

class(i ) = argmax c y (c) i .
Definition 11 introduces the optimization function giving the level of belongings of the vertices w.r.t. each class, as presented in [START_REF] Avrachenkov | Generalized Optimization Framework for Graph-based Semi-supervised Learning[END_REF].

Definition 11 (Graph-Based Semi-Supervised Learning (GSSL)). We call GSSL the algorithm solving the following optimization problem:

argmin y i j w i j y i d 1-σ i - y j d 1-σ j 2 penalty +µ i d 2σ-1 i y i -x i 2 fitting term . (4.2)
with µ a regularization parameter, and σ a unifying parameter.

The novelty in [START_REF] Avrachenkov | Generalized Optimization Framework for Graph-based Semi-supervised Learning[END_REF] is the introduction of the unifying parameter σ. Its role is to encompass three known techniques of GSSL:

• Laplacian (σ = 1): corresponds to the following Laplacian-based optimization problem:

min y i j w i j y i -y j 2 + µ i d i y i -x i 2 .
• Normalized Laplacian (σ = 1 2 ): corresponds to the following normalized Laplacianbased optimization problem:

min y i j w i j y i d i - y j d j 2 + µ i y i -x i 2 .
• Page-Rank (σ = 0): corresponds to the following Page-Rank-based optimization problem:

min y i j w i j y i d i - y j d j 2 + µ i 1 d i y i -x i 2 .
The parameter µ balances between the fitting term that is minimal when y equals x, and the penalty term that is minimal when y is smooth on the graph. Small values of µ favor regularity on the graph, and high values favor fitness to the input. Finally, the authors prove that this minimization problem has a simple closed-form solution:

Property 1. The solution to the optimization problem of Definition 11 is:

GSSL{x} = µ 2 + µ I + 2 2 + µ D -σ AD σ-1 -1
x.

We see already that compared to the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF], we have a simple closed-form solution to an optimization problem. We would like to leverage this strength to define smoother amplitude mapping for the problem of mapping graphs to time series.

The Graph-Based Semi-Supervised Learning approach has been successfully used to classify peer to peer traffic in [START_REF] Avrachenkov | On the Choice of Kernel and Labelled Data in Semi-supervised Learning Methods[END_REF].

GSSL Rewriting

We now reformulate Definition 11 with a matrix notation that will be easier to deal with:

GSSL σ,α {x} = argmin y D σ-1 2 (y -x) 2 2 fitting term +α D σ-1 2 y 2 L penalty , (4.3)
where y L = y * L y is the normalized Laplacian matrix semi-norm 3 . While the advantage of the original formulation is to directly show the minimized quantity at the scale of an edge (first sum), our formulation is simpler and more symmetrical with respect to the degree matrix D σ-1 2 . The equivalence can be shown with µ = 2 α using:

GSSL σ,α {x} = argmin y D σ-1 2 (y -x) 2 + α D σ-1 2 y 2 L 3 .
L is not a norm since ψ 0 L = 0 and ψ 0 = 0 since it is an eigenvector of L .

Graph to Time Series Duality Using

GSSL = argmin y D σ-1 2 (y -x) * D σ-1 2 (y -x) + α D σ-1 2 y * L D σ-1 2 y = argmin y i y i d 1 2 -σ i - x i d 1 2 -σ i 2 + α D σ-1 2 y 2 -α D σ-1 2 y * D -1 2 AD -1 2 D σ-1 2 y = argmin y i d 2σ-1 i y i -x i 2 + α i d σ-1 2 i y i 2 -αy * D σ-1 AD σ-1 y = argmin y i d 2σ-1 i y i -x i 2 + α i y i d 1-σ i 2 d i - i j w i j y * i d 1-σ i y * j d 1-σ j = argmin y i d 2σ-1 i y i -x i 2 + α 2 2 i j w i j y i d 1-σ i 2 -2 i j w i j y * i d 1-σ i y j d 1-σ j = argmin y µ i d 2σ-1 i y i -x i 2 + i j w i j y i d 1-σ i - y j d 1-σ j 2
Also, we can recognize in the formalism of Equation 4.3 a Tikhonov regularization problem. This gives then an alternative proof of Property 1, which can be rewritten using α instead of µ:

Property 2. The solution to the optimization problem of Equation 4.3 is:

GSSL σ,α {x} = I + αD -(σ-1 2 ) L D σ-1 2 -1 x.

Graph to Time Series Duality Using GSSL

The approach we explain now has been published in [4] and [START_REF] Girault | Semi-Supervised Learning for Graph to Signal Mapping: a Graph Signal Wiener Filter Interpretation[END_REF] together with a poster [START_REF] Girault | Semi-Supervised Learning for Graph to Signal Mapping: a Graph Signal Wiener Filter Interpretation[END_REF].

In section 4. 1.2, we claimed that we can achieve better smoothness than the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. First of all, optimal smoothness for a simple graph as the directed cyclic graph is quite obvious: We would like to have the graph mapped to a sinusoidal time series. For more complex graphs, achieving smoothness in the time series is not a goal that can be as easily characterized. We remark in particular, that the random walk approach leads to singularities in the resulting time series, as shown on Figure 4.2.

As such, instead of characterizing the smoothness of the resulting time series, we focus on the mapping a from vertices to amplitudes and we would like to obtain a mapping as regular as possible on the graph. In other words, close vertices on the graph should be mapped to close amplitudes. First of all, as suggested in section 4. 1.2, our approach maps vertices to the interval [-1, 1] instead of [1, N ]. To remove the influence of this difference, we rescale a given mapping a before making any comparison:

a norm (i ) = 1 + (N -1) a(i ) -a min a max -a min ,
with a min the smallest amplitude and a max the highest. Then, we see that the premise of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] is exactly the same, as illustrated in Equation 4. 1. We then use their cost function as the first assessment of regularity, and denote it 1 since it uses the modulus of the difference |a(i ) -a( j )|. The second function we use is based on the squared modulus |a(i ) -a( j )| 2 , and denoted 2 : .5) The difference between the two quantities is that 1 will allow a few edges with a high value of w i j |a(i ) -a( j )|, while 2 will favor global regularity. We can already see on Figure 4.3 that 2 favors the output of Figure 4.3(b) that globally have less variation of amplitudes from one vertex to another, while 1 does not differentiate between the three solutions of Figure 4.3.

1 (a) = i j ∈E w i j a norm (i ) -a norm ( j ) (4.4) 2 (a) = i j ∈E w i j a norm (i ) -a norm ( j ) 2 . ( 4

Achieving Smoothness using GSSL with 2 Classes

Instead of the cost function of Equation 4.1, we leverage the strength of GSSL and its closedform solution to build a cost function that is much simpler to optimize. We remind that the problem we try to solve is to associate amplitudes to vertices such that a random walk on the graph is associated to a smooth time series.

We remark first that amplitudes have always a minimum and a maximum. Our method intuitively associate these extrema to the two farthest vertices in the graph, in the sense of the shortest path. These two vertices constitute then the centroids of two classes: the small amplitudes and the high amplitudes. We denote x (1) and x (2) the delta vectors centered on the centroids. These correspond to the desired output in the context of Semi-Supervised Learning. The strength of GSSL is to give more than the class of a vertex: It gives the level of belongings y (c) of each vertex to each class c:

∀c, y (c) = I + αD -(σ-1 2 ) L D σ-1 2 -1 x (c) .
In the context of classification, the class of a vertex i is given by the index of its component of highest magnitude. Here we directly use the values y (1) i and y (2) i . Note that y (c) i ≥ 0 since ∀i , ∀c, x (c) i ≥ 0 and

I + αD -(σ-1 2 ) L D σ-1 2
-1 is a positive matrix for α > 0.

Without loss of generality, we assume that class 1 is the class associated to the minimum amplitude, and class 2 to the maximum amplitude. We propose to define the amplitude of vertex i as:

a(i ) =
y (2) i -y (1) i y (2) i + y (1) i

. ( 4.6)
The idea is to measure how much a vertex is in between the centroids of classes 1 and 2.

The numerator ensures that if a vertex i does not belong to class 1, then y (2) i = 1 and y (1) i = 0 such that a(i ) = 1 is maximal. On the other hand, if i does not belong to class 2, we have a(i ) = -1 minimal. The amplitudes span therefore the interval [-1, 1]. The rescaling by y (2) i + y (1) i -1 helps vertices away from both centroids to better span the interval [-1, 1] instead of aggregating around the zero amplitude.

On the same example than Figure 4.2(c), our approach gives the results of Figure 4.4. On this example, according to both 1 and 2 , our solution gives better smoothness at the cost of smaller extreme amplitudes. Note that amplitudes can always be rescaled similarly to the definition of 1 and 2 .

Discussion

Directed graphs To handle directed graphs, we use the same approach than [START_REF] Avrachenkov | Generalized Optimization Framework for Graph-based Semi-supervised Learning[END_REF] and symmetrize the graph using A sym = A+A * 2 prior to applying the GSSL method. Note that computing a random walk does not suffer from the same limitation and can be performed on the original directed graph. Moreover, the symmetrization A sym is equivalent to considering that closeness of two vertices is the mean of how close they are to one another. This should be a good input for the computation of the amplitudes.

Weights

Similarly to [START_REF] Campanharo | Duality between Time Series and Networks[END_REF], only non-negative weights can be considered. This is due to the random walk matrix being only defined for non-negative weights.

Complexity Computing the matrix GSSL involves a matrix inversion that can be costly. A simpler approach is to consider Property 2 as the following system of linear equations:

(1 + αL ) D σ-1 2 y = D σ-1 2 x .
For example, using MATLAB, we can use the function mldivide (or the operator \) such that:

D σ-1 2 y = (1 + αL )\(D σ-1 2
x), and no matrix is inverted.

Non Linearity Equation 4.6 uses the output of GSSL, but its amplitudes A are not the output of a linear operator because of the division. This non-linearity is however not an issue in practice and the results are interesting, as we show in the next section.

Dual Mapping

The approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] allows for an inverse mapping that, given a long enough time series, will recover the probability transition matrix of the random walk. This is made possible through the invertible mapping from vertices to amplitudes.

In our case, we do not force the invertibility property of the mapping. The drawback is then that we may obtain several vertices mapped to the same amplitudes (especially in the case of graphs showing high regularity such as the cyclic graph or the torus graph). However, our initial goal has been to study the structural properties of a graph through the properties of a time series, and as such the inverse transformation is not essential.

Example Graphs and Dual Time Series

In this section we apply our approach to several graphs and study the influence of the parameters α and σ. First, we start with the cyclic graph. Using the directed graph C N or the symmetric graph C N does not change the amplitude mapping of our method (due to the symmetrization of the adjacency matrix). We use first the graph C N to remove the influence of the randomness in the random walk and visually assess the regularity of the mapping on the time series (vertices of C N have a single outgoing directed edge). The results are shown on Figure 4.5 using C 30 . Note that the influence of the parameter σ is not studied for C N since the degree matrix associated to the symetrized adjacency matrix verifies D sym = 0.5I , and as such GSSL σ,α does not depend on σ. We observe on the same figure that the resulting dual time series ranges from a square signal (α small) to a triangle signal (α = 100). Carefully choosing the value of α (e.g. α = 10) yields a dual time series that is smooth in the sense that it is close to a sinusoid. Quantitatively, 1 does not discriminate between the GSSL solutions and that of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. On the other hand, it is clear that the triangle signal has the best 2 and the result from the simulated annealing approach is not optimal with respect to 2 .

In Figure 4.6, we performed the same generation of dual time series with the symmetric cyclic graph C 30 . This shows that even if the amplitudes are smooth on the graph, the resulting dual time series is not as smooth as a sinusoid. This is a consequence of the random walk that being able to get back to the vertex it just left, thus creating a singularity in the generated time series. This is unavoidable and especially clear for symmetric graphs.

The next graph is the torus graph T 10 . This graph introduces a regular structure that is inherently two dimensional. Similarly to C N and C N , this graph is regular such that σ plays no role in the resulting amplitudes or the generated time series. The results are shown on Figure 4.7. These time series get smoother when α increases, as expected, except for the unavoidable singularities. With respect to 1 and 2 , our approach has a clear advantage as soon as α is high enough, i.e. the regularity term in Equation 4.3 is sufficiently weighted. Note that if α is too high, the amplitudes will spread to a very small interval, and as such a compromise should be found between high values of α and enough richness in the amplitudes.

The next graph, whose resulting amplitudes are shown on Figure 4.8, is the Molène graph W . Similarly to the torus graph, this graph results from a two dimensional space. However, it is an irregular sampling of this space rather than a regular sampling. Noticeably, this space has boundaries (mostly the coast lines) where degrees are smaller (fewer neighbors). All this leads to a non regular graph, i.e. with vertices of different degrees. Therefore, the unifying parameter σ influences the results. This can be observed with α = 100: The set of vertices mapped to amplitudes close to 0 (white vertices) are closer to the south east when σ is higher. In other words, the Page-Rank method (σ = 0) maps vertices to amplitudes more evenly. The influence of α is similar to previous graphs. Compared to the output of the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] shown on Figure 4.2(c) (page 58), the resulting time series is much smoother. This is illustrated with smaller 1 and smaller 2 if α is not too high. Overall, σ = 1 and α = 100 seems to produce the best results with respect to 1 and 2 .

Finally, the mapping between vertices of the MOSAR graph M (sum) (t , t + 1 week) and amplitudes are shown on Figure 4.9. First, as suggested before, the placement of the vertices may not be representative of the overall structure of the graph for the given period since the structure of the contacts evolves over time. This is particularly visible on Figure 4.9 where 4.8 for the MOSAR graph M (sum) (t , t + 1 week). Finding the farthest vertices is done with respect to the pseudo-distances 1 -e -w i j . This ensures a better scaling of distances since the contacts total duration can vary greatly from one edge to another. This figure does not include the values of 1 and 2 since we do not have a reference with the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] due to the computational complexity of the approach and our results are not completely satisfactory (see section 4.2.4).

the vertex associated to class 2 is lost within a group of vertices in class 1 (see for example the amplitude of high magnitude when α = 500). Overall the results show that the chosen centroids are far away from any other vertex such that amplitudes on the graph are quite small. Indeed, for all these vertices lying in between centroids, we have y (1) i and y (2) i close to 0.

Limits of the GSSL Approach with 2 Classes

As we can see on Figure 4.5 through Figure 4.9, our approach has some limitations. First of all, similarly to [START_REF] Campanharo | Duality between Time Series and Networks[END_REF], only static graphs can be studied. Some extension to dynamic graphs has been introduced in [Singh Mor, 2013] where we studied the impact of adding vertices one after the other and updating the amplitudes on the resulting time series. Phase changes in the resulting time series show then major structural changes in the graph such as a path connecting into a cycle. We did not investigate further for lack of time.

The second limitation can be observed on Figure 4.9. As written in section 2.7.1, the MOSAR graph shows five main communities that correspond to care units in the hospital. Our method could not recover all these care units, and only one care unit stands out with respect to the amplitudes (the positives ones). From this it appears that choosing differently the two centroids could lead to a better mapping. For example, if the outlying vertices are not considered during the shortest path computation, we could place the centroids in the middle of actual classes of vertices. Another extension of our approach have been studied in [Fromherz, 2014] where we tried to extend the 2 classes approach to N classes (one per vertex) and took the output of the GSSL operator as feature vectors for each vertex. A hierarchical clustering approach can then be used to perform community detection (see [Tremblay and Borgnat, 2014] for a comprehensive work where feature vectors are based on graph wavelet atoms and community detection is done through hierarchical clustering). Moreover, using different values of α allows for different communities, and eventually multiscale communities.

Beyond the Classification Approach

In this section, we leave the problem of finding a mapping from vertices to amplitudes and we study GSSL as a filter on graph signals.

Graph Signal Operator Interpretation: Wiener Filter

The most remarkable feature of GSSL lies in its nature of graph operator. Indeed, as we already remarked in the background section, GSSL is an operator having as input a graph signal x and a graph signal y as output. We show in this section that GSSL is actually a generalized convolution with a very specific frequency response. First, we recall the matrix representation of the operator:

GSSL σ,α = I + αD -(σ-1 2 ) L D σ-1 2 -1 .
We then simplify it into:

GSSL σ,α = D -(σ-1 2 ) I + αL -1 D σ-1 2 .
Let

H α = D σ-1 2 GSSL σ,α D -(σ- 1
2 ) , we have:

H α = I + αL -1 .
And the frequency response of H α , with F = Ψ (the graph Fourier transform based on the normalized Laplacian matrix L ), is:

H α = 1 I + αΘ .
Since Θ is diagonal, H α is diagonal, and using Equation 3.21 (page 38), H α is a convolutive operator. Using the reduced graph frequency ν l = 1 2 θ l /2 associated to the graph Fourier transform F = Ψ, we have:

h α (l ) = 1 1 + αθ l = 1 1 + 8α ν l 2 . ( 4.7)
The reader familiar with stochastic signal processing will have recognized here the frequency response of an optimal Wiener filter (in the sense of minimizing the Mean Square Error (MSE)). To see this, we consider a stationary temporal signal x observed with an additive uncorrelated stationary noise n: s = x + n. The optimal filter taking s as input and minimizing the mean squared sum of errors E |x -y| 2 between its output y and x is the Wiener filter g with frequency response equal to [Papoulis, 1991]:

g (ν) = 1 1 + S n (ν) S x (ν)
, with S n (ν) = E n(ν) n * (ν) the power spectrum density of the noise n and S x (ν) = E x(ν) x * (ν) the power spectrum density of the signal. The fraction S x (ν) S n (ν) is here the Signal to Noise Ratio (SNR):

g (ν) = 1 1 + 1 SNR(ν)
.

We see that if the SNR equals 1 8αν 2 , then the frequency response of the corresponding Wiener filter has the same form than Equation 4. 7. This is the case for n a white noise of power spectrum density S n (ν) = 8α and for a signal x of power spectrum density S x (ν) = 1 ν 2 . We will give more details on this in chapter 6, including Wiener filters and a formalism for Wiener filters in the case of graph signals.

Back to the GSSL filter, when σ = 1 2 , we have H α = GSSL 1 2 ,α such that GSSL 1 2 ,α is a convolutive operator. For other values of σ, none of the graph Fourier transform we defined in chapter 3 allows to interpret GSSL σ,α as a convolutive operator. However, we have the following spectral decomposition of the GSSL matrix:

GSSL σ,α D -(σ-1 2 ) ψ l = 1 1 + αθ l D -(σ-1 2 ) ψ l .
Therefore, we can define a class of graph Fourier transform defined as

F (σ) = D -(σ-1 2 )
Ψ, with the same graph frequency as the normalized Laplacian matrix 4 . Using this graph Fourier transform F (σ) , GSSL σ,α is a convolutive operator. This however comes at the cost of F (σ) not being unitary:

F (σ) (F (σ) ) * = D -2σ+1 .

Application

An immediate application of this interpretation of GSSL as a graph filter is to apply it to a non-smooth graph signal. For example, given the task we started with of mapping vertices to amplitudes, suppose we have a mapping that is not smooth enough. This can be for example the mapping of Figure 4.3(c) produced by the simulated annealing approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. Note first, that compared to the approach of section 4.2, the input of GSSL is not anymore two 0/1 graph signals, but a single graph signal x.

From Equation 4.7, we know that GSSL is a low pass filter, and we identified it to a Wiener filter. Figure 4.10 shows the frequency responses, as written in Equation 4.7, for the three graphs we were able to perform the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. Also, since we used the reduced graph frequencies in Equation 4.7 all the frequency responses range from ν 0 = 0 to ν N -1 ≤ 0.5. Therefore, from one graph to another, the only differences with respect to the frequency response lie in the localization of the reduced graph frequencies (red dots on Figure 4.10).

We now apply GSSL to the mapping produced by the simulated annealing approach. The first graph, is the cyclic graph C 30 on Figure 4.11. While 1 is not modified5 by GSSL, 2 is clearly reduced when the filter is applied. Also, there seems to be an optimal value of α for which 2 is minimal. Overall, GSSL improves then the smoothness of the mapping of vertices to amplitudes.

Similar results can be observed on the torus graph T 10 on Figure 4.12. While 1 only slightly increases with α, 2 slightly decreases, especially for α = 10. It can be observed on Figure 4.10(b) that the input of GSSL is actually already a lowband signal with high frequency components of very low magnitude, such that the input of GSSL is only slightly altered by the operator.

On the other hand, we see on Figure 4.10(c) that the mapping on the Molène graph has some noticeable high frequency components. As a consequence, the mapping produced by GSSL filtering of the mapping of Figure 4.2(c) is much better as illustrated on Figure 4.13. Indeed, 2 is always smaller, and for σ = 1, 1 is always smaller. The best results are achieved for σ = 1 and α = 10 since the case α = 100 shows very flat amplitudes. However, compared to Figure 4.8, 1 and 2 are still higher such that our original approach based on two classes should be favored.

Finally, we can apply our interpretation in terms of filter to the amplitude mapping we proposed in Equation 4. 6. Indeed, we can rewrite this equation into:

a(i ) = GSSL σ,α x (2) -x (1) (i )
GSSL σ,α x (2) + x (1) (i ) , using the linearity of the GSSL operator. In the equation above, the numerator is then a filtering of the graph signal with 1 on the vertex associated to the maximum amplitude, -1 to the vertex of minimum amplitude, and 0 elsewhere. GSSL smoothes out this signal. On the other hand the denominator is the operator GSSL applied to the graph signal with 1 on the vertices of minimum and maximum amplitude, and 0 elsewhere: This denominator boost 4.7 on three different graphs (scale on the left). The top curve is for α = 0.1, the middle one for α = 10, and the bottom one for α = 100. (Blue) Moduli of the Fourier components of the graph signal obtained using the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF], with respect to the graph Fourier transform based on the normalized Laplacian matrix L (scale on the right). Vertical dotted lines identify to reduced graph frequencies of the graphs. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 5.1.1 Time Shift . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 5.1.2 2D Translation (Images) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 5.1.3 Graph Shift and Algebraic Signal Processing Theory . . . . . . . . . . . . . . 77 5.1.4 Generalized Translations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77 5.2 Preliminary Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 5.3 Graph Translation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 5.3.1 .4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89 We now leave the setting of graph analysis of chapter 4 to the setting of operators on graph signals. More precisely, we look at the definition of an operator on graph signals that would be the equivalent of the classical time shift operator. Finding a good candidate for this equivalence would allow the study of numerous equivalences between graph signals and temporal signals where the time shift is at the very core of time: Time shifting is the time passing by. The time shift operator led to the study of time invariant operators, to Time-Frequency decompositions, to wavelets, and to the formalism of stationary signals to cite only a few. We begin by recalling the standard time shift operator, and list several proposed equivalents operators for graph signals. We build on this background to propose our own operator enforcing a key property: isometry.

Graph Translation

Background

Time Shift

The time shift T is an operator acting on a signal x(t ) by shifting the time reference by one time unit:

T {x}(t ) = x(t -1). (5.1) This operator is linear (T {αx + y} = αT {x}+T {y}), and a convolution of its input with the Dirac mass located at t = 1:

T {x} = x * δ 1 , since: (x * δ 1 ) (t ) = x(u)δ 1 (t -u)d u = x(t -1).
Using the convolution theorem in Theorem 2 (page 35), the action of the time shift on a Fourier mode is simple and only depends on the Fourier mode: .2) Using the linearity of the time shift, the Fourier transform of a shifted signal T {x} is a complex phase shift of the Fourier transform x of x. Note also that, due to the convolution property, the Fourier modes are eigenvectors of the time shift.

T {e ν } = e -ı2πν e ν . ( 5
As remarked in section 3.2, the convolution theorem is valid using the basis {e ν } ν and gives Equation 5.2. Using the Fourier basis {c ν } ν ∪ {s ν } ν the convolution theorem is not valid. For the time shift, this is illustrated by:

T {c ν } = cos(2πν)c ν -sin(2πν)s ν T {s ν } = -sin(2πν)c ν + cos(2πν)s ν
All the properties we described in this section apply also to the shifting operator T τ , shifting a signal by τ time units instead of just one:

T τ {x}(t ) = x(t -τ).
Note that shifting by τ time units is actually equivalent to shifting τ times when τ is integer. Finally, the time shift operator for a time series x[n] is given by:

T k {x}[n] = x[n -k].

2D Translation (Images)

Shifting the pixels of an image can be done using two axes: horizontally or vertically. Combining these two shifts yields a shift in any possible direction. We denote s(x, y) a 2D signal, i.e. the sample at coordinates (x, y) of the signal s. We denote e ν x ,ν y (x, y) = e ı2πν x x e ı2πν y y the 2D Fourier mode of horizontal frequency ν x and vertical frequency ν y . Definition 3 gives then:

e ν x ,ν y (ν x , ν y ) =    1 if ν x = ν x and ν y = ν y 0 otherwise.
We denote then T x and T y the shift operators according to the horizontal and vertical axes:

T τ x x {s}(x, y) = s(x -τ x , y) T τ y y {s}(x, y) = s(x, y -τ y ), such that T τ x
x {e ν x ,ν y } = e ı2πν x τ x e ν x ,ν y T τ y y {e ν x ,ν y } = e ı2πν y τ y e ν x ,ν y . These 2D shift operators verify then the convolution theorem, and they are isometric. We remark also that these two operators commute T x T y {s} = T y T x {s} since they apply indepen-

Background

dently to each axis. It follows also from the commutative property that the operator T τ x x T τ y y corresponds actually to an infinite number of paths in the 2D space. Indeed, shift in x and in y can be alternated (composed) such as at the end a total of τ x and τ y shifts have been performed. This is different from the time shift where there is only one axis to shift on.

Graph Shift and Algebraic Signal Processing Theory

In [START_REF] Pueschel | Algebraic Signal Processing Theory[END_REF], the authors introduce the so called Algebraic Signal Processing framework. One of the result of this framework is to define a Fourier transform from a given translation operator. For example, given the matrix representation of the time shift operator for periodic time series, one can compute the DFT matrix using the diagonalization of this matrix (see section 2.5).

This leads the authors of [START_REF] Sandryhaila | Discrete Signal Processing on Graphs[END_REF] to seek the equivalent of a time shift in the framework of Graph Signal Processing. As seen in section 3.4, they use the operator whose matrix representation is the adjacency matrix A of a graph, and call it the Graph Shift. Note that this approach is only valid with the graph Fourier transform defined using the adjacency matrix otherwise the graph shift is not a convolutive operator.

The graph shift agrees with the time shift only on the directed cyclic graph C N . On the symmetric cyclic graph C N , the graph shift diffuses its input according to the edges. Therefore, the value at vertex i is diffused to vertices i -1 and i + 1. A similar property can be verified on both the directed and symmetric torus graphs T n and T n , where values are diffused both horizontally and vertically at once. We will see also that the graph shift does not in general preserve the energy of the signal, except for the directed cyclic graph C N .

Generalized Translations

Another approach developed in [START_REF] Shuman | The Emerging Field of Signal Processing on Graphs: Extending High-Dimensional Data Analysis to Networks and Other Irregular Domains[END_REF], Shuman et al., 2015] stems from the observation that the translation T τ of a signal x is actually the convolution of x with the Dirac mass δ τ . They generalize this property shown by the time shift operator to the graph setting using the generalized translation T i :

(T i {x}) j = (x * δ i ) j = N -1 l =0 x(l ) χ l * i χ l j , (5.3) since δ i (l ) = (F δ i ) (l ) = χ * δ i (l ) = χ * l i = χ i l * = χ l i * = χ l *
i , with χ l the l th Fourier mode. Note, that the last equality of Equation 5.3 uses the Laplacian matrix approach of graph Fourier transform. A similar formula can be obtained using the normalized Laplacian matrix.

From our discussion on the convolution theorem of section 3.2 (page 38), we know that the exponential basis and the sine/cosine basis are not equivalent. Also, the graph Fourier transform based on the Laplacian matrix gives a real Fourier basis, hence not the exponential one on the cyclic graph C N . Finally, since the time shift does not verify T (ν) = T (-ν), then the time shift is not a convolution with respect to the real Fourier basis given by the graph Fourier transform on L. Therefore, the generalized translation and the time shift do not agree on the cyclic graph C N .

Finally, we will see that the generalized translations are not isometric operators.

Preliminary Results

In this chapter, we are interested in the definition of a translation operator for graph signals.

Compared to the literature, out aim is to have an isometric translation operator that we shall use in chapter 6 to define stationary graph signals.

Such an operator should first be a function H : C N → C N , i.e. an input signal is mapped to an output signal. In order to make sense and be rightfully called a translation operator it should share key properties with the time shift, or at least some of them. We summarized several of them in the following list:

(i) Linear: H {αx + y} = αH {x} + H {y} (superposition principle) (ii) Convolutive: ∀ν k = ν l , H kl = 0 (iii) Energy preserving (isometry): E x = E H {x} (iv) Localization preserving: ∀i , ∃k i : H {δ i } = δ k i (a delta is mapped to another delta)
For item (iii), we will be using the accepted convention of graph signal energy [Shuman et al., 2013b, Sandryhaila andMoura, 2014a]:

E x = x * x = i ∈V |x i | 2 .
( 5.4) Whether this norm is adapted to the graph framework is still an open question. For example, if the vertices are a regular sampling of a Euclidean space, then each of them represent the signal at the position of the vertex. However, if we look at the definition of energy in such a space, we have (in 1D):

E x = b a |x(t )| 2 d t , if the support of x is [a, b] (i.e.
x is non zero only on [a, b]). This can be approximated using a Riemann integral. Given a = t 0 < t 1 < • • • < t T = b, we have the following approximation of the energy:

T -1 i =0 x t i +1 + t i 2 2 (t i +1 -t i ) .
We see that using a regular sampling of time, we obtain a constant value for t i +1 -t i , and the discrete definition of energy is a good approximation of the energy of the continuous signal. However, graphs are irregular sampling, similar to a non-constant t i +1 -t i . To have a better energy we should then consider how big is "t i +1 -t i " for a given vertex, and weight |x i |2 by this quantity. This will be the subject of future works. In the absence of a better norm, we will use the one of Equation 5.4 in this dissertation.

For item (iv), the localization preservation is actually a more general property involving the spread of a signal: Intuitively, if a signal is localized on a given vertex and its energy spreads slightly around it, then its translated signal spreads slightly around another vertex. Quantifying the spread of a signal is however the subject of an ongoing research and the definition has not yet converged [START_REF] Agaskar | A Spectral Graph Uncertainty Principle[END_REF], Pasdeloup et al., 2015, Tsitsvero and Barbarossa, 2015] 1 . We then use the characterization above that only ensures that a delta signal is mapped to another delta signal.

We derive now a few preliminary results on the general form a translation operator H should have. In particular, we show that we need to drop at least one of the properties above to properly define a translation operator on graphs.

Property 3. Let H be a linear (item (i)) and isometric (item (iii)) operator. Then H is unitary.

Proof. We show here that H * H = I . First, we have:

E H x = x * H * H x = x * x = E x .
Using the delta signal δ i we show that the diagonal of H * H is one:

E δ i = 1 = δ * i H * H δ i = H * H i i .
Then, using the signal δ i + δ j , we have:

E δ i +δ j = 2 = (δ i + δ i ) * H * H (δ i + δ i ) = 2 + H * H i j + H * H j i ,
such that (H * H ) i j + (H * H ) j i = 0. Finally, the signal δ i + ıδ j allows us to conclude:

E δ i +ıδ j = 2 = (δ i + ıδ i ) * H * H (δ i + ıδ i ) = 2 + ı H * H i j -ı H * H j i , such that (H * H ) i j -(H * H ) j i = 0.
[ Horn and Johnson, 2013, Theorem 2.1.4] proves also the equivalence. with Ω the matrix of a convolutive real operator.

Proof. We have H = diag( h(0), . . . , h(N -1)) (H convolutive) such that:

H * H = I ⇔ H * H = I since F is unitary, and H = F H F * ⇔ ∀l , h(l ) = e ıω l ⇔ H = e ı Ω with Ω = diag(ω 0 , . . . , ω N -1 ) ⇔ H = e ıΩ

Isometry and Localization Invariance

We now study whether a linear (item (i)), isometric (item (iii)) and localization preserving (item (iv)) operator on graph signals can still be a between the two vertices u and u 0 . The exact definition of this distance is subject to ongoing research. Nevertheless, if d(u, u 0 ) = 0 implies u = u 0 (natural assumption if the graph is connected), then the only graph signal with a zero spread are the delta signals. As such, an operator preserving the graph spread maps a delta signal to another delta signal. Another approach to define the graph spread is to use the 0 -norm of the signal. This leads to a similar result on the image of a delta signal, up to rescaling of the image, which does not hinder our results.

convolutive operator (item (ii)). Let H be such an operator. We have them (item (i) and item (iv)):

H = [H δ 1 • • • H δ N ] = δ k 1 • • • δ k N H * =    δ * k 1 . . . δ * k N    .
It follows from the isometry (item (iii)) that

H * H = I such that δ * k i δ k j = 0 if i = j .
H is then a permutation matrix. Such a matrix is diagonalizable with very specific eigenvalues and eigenvectors. Assume first that the permutation has a single cycle 2 . Assuming also that the vertices are indexed according to this cycle (i.e. vertex i is mapped to i + 1 though the operator H ), then H is exactly the time shift operator for N -periodic time series. As a consequence, the eigenvectors of H are the Fourier modes of the DFT. Since the graph Fourier modes are not the DFT Fourier modes (except for the cyclic graph, for example), and the eigenvectors of a convolution are the graph Fourier modes, then H is not a convolution.

When H is a permutation matrix with more than one cycle, the analysis we developed can be applied independently to all cycles. In conclusion, all four properties of the time shift can be preserved if and only if the graph Fourier transform is essentially a DFT, or multiple independent DFTs 3 . The graph Fourier transform F being rarely equal to the DFT, there is no operator on graph signals verifying all four properties. An equivalent translation operator for graph signals will therefore not verify one of the properties, and we make the choice in the following to enforce only the first three: linearity, convolution, and isometry (with respect to the chosen energy definition).

Graph Translation

Now that we have set the basic elements of the theory, we can introduce our own equivalence of the time shift, study its properties and compare it to the literature. The choice we develop here is to preserve item (i) through item (iii) without trying to enforce the localization preservation of our operator.

Definition

In the rest of this chapter, we assume that the graph Fourier transform is unitary. Let N = diag(ν 0 , . . . , ν N -1 ) be the diagonal matrix of the graph frequencies. At this point, these graph frequencies are generic. They can be for example reduced graph frequency. We use the notation N for this matrix only for the following definition.

Definition 12. The graph translation on the graph G with respect to the unitary graph Fourier transform F and diagonal graph frequencies matrix N is the operator T τ

G verifying:

T τ G = e -ı2πτN

Graph Translation

In particular, when τ = 1 and using the Laplacian matrix based graph Fourier transform with the associated reduced graph frequency ν l = 1 2 λ l /ρ G (see section 3.3, page 39), we obtain the following graph translation operator [7]: .5) Similarly for the normalized Laplacian matrix based graph Fourier transform we obtain:

T τ G = e -ıπτ L ρ G . ( 5
T τ G = e -ıπτ L 2 .
( 5.6) In the rest of this dissertation, T G is the graph translation defined using L (Equation 5.5), and T G is the graph translation defined using L (Equation 5.6).

Properties

We now go through the important properties the graph translation verifies. First, as seen in section 5.2, the graph translation is an isometric convolutive operator. Also, as we did not enforce it, the graph translation does not preserve the signal localization in general.

Another very important property is that the graph translation is a complex operator. This is quite different from other operators from the literature, but necessary as we saw in section 5. 2. The time shift does not suffer from this because of the particular structure of the DFT matrix. This raises several questions. First, dealing with complex signals leads to the definition of analytic signals, and whether we could define the analytic graph signal corresponding to a real graph signal. The second question lies in both the notion of convolution and the graph Fourier transform. We are using a real graph Fourier transform, whereas the Fourier transform that interprets the time shift as a convolutive operator is complex (i.e. the Fourier basis is the basis of complex exponential, and the basis of sine and cosine does not work for interpretation purposes). Therefore, could we define a complex graph Fourier transform such that the graph translation is a real convolutive operator? This is an open question.

Back to the classical setting of graph Fourier transform and graph translation, the consequence of the convolution and isometric properties is the preservation of the Power Spectrum Density (PSD), i.e. | T τ G x| 2 = | x| 2 , ∀τ, since only the complex phase of x is modified by the graph translation and not the amplitude.

Another consequence of this property lies in the energy of the input-output difference:

x -T τ G x 2 2 = N -1 l =0 2 1 -cos(2πτ ν l ) | x(l )| 2 .
(5.7)

Since the reduced graph frequency ν l increases with l , we can show, for small values of τ, that the energy above increases with l when applied to the Fourier mode χ l . A similar result can be shown with T G . This follows the intuition where a small shift of a slowly varying signal (low frequency) will be very similar to the original signal whereas a small shift of a quickly varying signal (high frequency) can be quite different.

Similarly to the time shift, we can define a mathematical group on the graph translations whose law is the compositional one. In other words, T

τ 1 G T τ 2 G = T τ 1 +τ 2 G is itself a graph translation, just like T τ 1 T τ 2 = T τ 1 +τ 2 is itself a time shift.
One property that is obvious for the time shift is its cyclic nature: T N x = x (for a N -periodic time series). For the graph translation, this property is not obvious anymore. Let τ be such that ∀x, T τ G x = x. Using the graph Fourier transform, we have then:

∀l , e -ı2πτν l = 1, which leads to:

∀l , ∃k l ∈ Z, τν l = k l .
First, we assume that all graph frequencies are rational: ∀l , ν l = p l q l . Then the least common multiplier τ = lcm(q 0 , . . . , q N -1 ) verifies the equation above. This case can be extended to ∀l , ν l = ζ p l q l with ζ an irrational number independent of l , in which case τ = ζ lcm(q 0 , . . . , q N -1 ) verifies the equation above. On the other hand, if any ratio of graph frequencies is irrational, which can easily happen when dealing with a real graph, then the only solution to the equation above is τ = 0. In other words, the graph translation is not cyclic in general.

Additionally to not being cyclic, T τ G x does not converge when τ goes to infinity (except for DC signals, i.e. signals invariant by graph translation). This last property agrees with the time shift which is cyclic, hence not converging when iterated. We will see that this is not the case for either operators from the literature because of their non-isometric nature.

Comparison to Other Translation Operators on Graphs

In section 5.1, we described the two operators on graphs from the literature that share similar properties with the time shift, namely the generalized translation and the graph shift. Results are summarized in Table 5. 1. First, the generalized translation T i is defined as the convolution by a delta centered on vertex i . As such, it is obviously a convolutive operator. Unfortunately, it is not isometric since δ i does not verify δ i = 1. Finally, in [Shuman et al., 2015, Example 2], the authors show a localization property of the generalized translations, but not a localization preservation. More precisely, given the heat kernel g κ having exponentially decreasing power spectrum:

g κ (l ) = e -κ(2πν l ) 2 , ( 5.8)
bounds are given on the spread of the translated signal T i g κ . In Equation 5.8, κ controls the band of the signal: Bigger values correspond to low-band graph signals with only low frequencies, and smaller values to wide-band graph signals with a lot of frequencies. The bound is quite complex but reduces to an energy decay from i to its neighbors that is steeper for smaller values of κ. Actually, we can see T i g κ as the convolution of δ i by g κ , i.e. as the filtering of δ i with a low pass filter, thus showing that the steepness of the decay depends on κ in a natural way: small κ correspond to an all-pass filter, in which case T i g κ ≈ δ i .

On the other hand, in general T i δ j = δ j * δ i can be quite delocalized on the graph, and is not equal to δ i + j , or any other δ k for that matter, as one would have expected with the time shift [START_REF] Shuman | Vertexfrequency analysis on graphs[END_REF]. This also shows that the composition of two generalized translations is not a generalized translation.

The next operator is the graph shift defined as the adjacency matrix A. With respect to the graph Fourier transform defined using A, the graph shift is a convolutive operator.

In [Sandryhaila and Moura, 2014a], the authors introduce the normalization of A by the eigenvalue of highest magnitude: A norm = 1 γ max A. While this normalization ensures that the graph shift does not alter the eigenvector associated to the eigenvalue of highest magnitude, all other Fourier modes have a vanishing output on the graph shift since γ l γ max ≤ 1. Therefore the graph shift is not isometric. Finally, item (iv) is preserved if and only if each vertex has degree 1 and only one neighbor. Assuming a connected graph, this is only verified for the directed cyclic graph.

Table 5.1. Properties of several translation operators on graphs. The time shift operator T on a 1D regular space is separated from the three graph signal operators to highlight the difference in the type of space they apply to.

Convolution

Mathematical Isometry Localization Group Invariance

T (time shift) {T τ } τ T G {T τ G } τ A {A τ } τ †
T i † if and only if A is non singular, i.e. A -1 is well defined.

Illustration on Graph Signals

We introduced the graph translation operator as a complex phase shifting operator in the graph Fourier domain. As such the comparison with the time shift with respect to the Fourier domain is immediate. Its behavior in the vertex domain is however not given by its spectral properties. In this section, we illustrate on several toy signal examples the response of the graph translation and the responses of the two operators from the litterature that were derived from the properties of the time shift: the graph shift and the generalized translation.

Toy Signals

To that end, we propose to use three toy graph signals. First, the delta signal δ i that is classically used to identify the impulse responses of convolutive operator. All the operators being linear, the reponses to a linear combination of delta signals can be derived from each of the impulse response. One key property of a delta signal is its Fourier transform. Indeed, delta signals are delocalized in the Fourier domain, and as such comprise many Fourier modes, all being translated at once by the operators we wish to study. It is worth noticing however that some graph can show localization, to some extent, in both the vertex and the Fourier domains. We refer the interested reader to [START_REF] Pasdeloup | Toward An Uncertainty Principle For Weighted Graphs[END_REF] where the authors study an uncertainty principle applied to graph signals.

The second signal is the heat kernel g κ . This signal is a low band signal with only Fourier modes of low frequencies. Provided κ is large enough, this signal is localized in the Fourier domain, and as such the exact opposite of the previous delta signal. In Figure 5.1 through Figure 5.4, we chose to use the graph Fourier transform based on Laplacian matrix and the associated graph frequencies to define g κ .

The third signal is even further delocalized in the vertex domain and equal to the constant signal 1. Depending on the graph Fourier transform, this constant signal can be highly localized in the Fourier domain (with the Laplacian matrix approach, 1 is a DC component),

X T 5 X AX /γ 0 T G X T G X δ 1
T 5 G X T 10 G X T 20 G X T 50 G X T 100 G X δ 1 1 1 1 1 0.0 0.5 1

Graph Translation

X T 35 X AX /γ 0 T G X T G X δ 1
T 5 G X T 10 G X T 20 G X T 50 G X T 100 G X δ 1 1 1 1 1 0.0 0.3 0.6
T 5 G X T 10 G X T 20 G X T 50 G X T 100 G X δ 1 1 1 1 1 1 0.0 0.3 0.6
T 5 G X T 10 G X T 20 G X T 50 G X T 100 G X δ 1 1 1 1 1 0.0 0.3 0.6 1 g 1e -4
0.0 or more or less delocalized depending on the graph. Noticeably, for regular graphs, using any of the three graph Fourier transforms does not change the localization in the Fourier domain of the constant graph signal 1. Indeed, we already know that it is a DC graph signal for the Laplacian matrix-based approach. Since d L = L = d I -A with d the degree for regular graphs, then 1 is also an eigenvector of L and A, and it is actually associated to the lowest frequency in both cases, therefore 1 is a DC graph signal in all three graph Fourier transform.

Description of the Figures

In Figure Each figure is divided into three blocks. The first top block shows one application of the operators to each of the three signals. The second middle block shows ten applications of the operators. And the third bottom block shows the graph translation T G iterated 5, 10, 20, 50 and 100 times. Note that for the generalized translation, iterating the operator does not make sense [START_REF] Shuman | Vertexfrequency analysis on graphs[END_REF], and we choose another one instead of iterating the first one ten times.

Before actually commenting on those figures, we need a warning about Figure 5. 4. The vertices of this graph have been placed in the 2D plane using GraphViz (see section 2.7.1). It follows that this representation may not be accurate. For example, the proximity relations of vertices of a d > 2 dimensional space can not be accurately rendered in 2D. This is also the case for social graphs such as the MOSAR graph. Therefore, on Figure 5.4, it is more meaningful to look at the spread of a signal rather than at where it is localized.

Comments on the Responses

We can now comment the figures. First, the response of the generalized translation is as expected. Indeed, the heat kernels are localized on the target vertices whereas the responses on the delta signal are not particularly interesting. Notice that the constant signal stays constant through T i since we used a generalized translation built from the graph Fourier transform based on L. We see on the MOSAR graph the illustration of the previous remark where the localized heat kernel is not cleanly localized in the 2D space. Nevertheless, guarantees exist on the localization such that what we observe on the figures is the result of a non accurate representation rather than a misbehaving operator.

Next, the graph shift. Without rescaling the operator, the energy of its output is not contained within reasonable bounds. The rescaling of the output in the figures reflects the use of A norm instead of A to keep the energy of the output from exploding. Another drawback of this operator that we observe on the figure is its convergence, as opposed to the cyclic time shift and the non-converging graph translation. When iterated on a signal, the output of A norm converges to the DC component of the signal. The DC component is meant here with respect to the adjacency matrix-based graph Fourier transform. As seen with the Molène

Discussion

graph and the MOSAR graph, this DC component can be highly localized. For example, the delocalization of 1 is not preserved at all, and the delta signal vanishes.

The graph translation on the other hand shows first evidence of a diffusion operator, in the sense of an isotropic translation operator. In particular, the one step translation of the delta signal using any of the two graph translation diffuses the energy of the signal to the neighboring vertices. After a few steps, the energy is spread across all vertices, corresponding then to a wide spread diffusion in the graph. The second remark concerns the difference between the response to a delta signal and to the heat kernel. While the first response is noticeably different than the original delta signal, for the second, input and output are much more alike. This observation is related to the fact that variations of the signal through graph translation are higher for higher frequencies, as proved in Equation 5. 7. The next observation is on the third blocks of the figures, and more particularly on the response to the heat kernel where we can observe a continuum of low band signals. By enforcing the preservation of the power spectrum (i.e. through isometry and convolution), we allowed the graph translation to preserve the bandwidth of the signal through multiple applications of the operator. Last but not least, the main difference between the two graph translations lie in the responses to the signal 1. While T G keeps this signal unaltered (it is a DC component), T G alters it. Indeed, 1 is not a DC component w.r.t. the graph Fourier transform based on L . Therefore, using one or the other graph translation depends, in particular, on the meaning of a DC component.

Discussion

We built the graph translation as an equivalent of the time shift for graph signals. As such, the graph translation can be identified as an operator in a one dimensional space: Translating is performed in only one (very abstract) direction. This contrasts with the two translations T x and T y of a 2D space that can translate in two directions. The question that naturally arises concerning this operator is whether it is possible to devise an alternative definition that would consider the graph as a higher dimensional space. For example, the torus graph, the Molène graph, or the graph built from a random sampling of a 2D plane are definitely graph representations of a 2D space. With these graphs, the goal would be to define two graph translations: one for the x axis, and one for the y axis. Unfortunately, this is an open question since we have access to only one set of frequencies (and Fourier modes) instead of the two sets of frequencies (and Fourier modes) for the 2D Fourier transform of Definition 3 (page 35).

Also, as proved in section 3.5.1 for the cyclic graph and section 2. 5.4 for the torus graph, we cannot get the Fourier bases {e ν } ν or {e ν x ,ν y } ν x ,ν y with the diagonalization of A, L or L . However, the 1D and 2D shift operators only verify the convolution theorem w.r.t. the exponential bases. Therefore, given the graph Fourier transform of the cyclic graph and the torus graph, no convolutive operator is the time shift or the 2D shift. It follows that the definition of graph translation can not recover the classical 1D and 2D shift operators on the cyclic graph and the torus graph using the definition of generalized convolution of [START_REF] Shuman | The Emerging Field of Signal Processing on Graphs: Extending High-Dimensional Data Analysis to Networks and Other Irregular Domains[END_REF].

The prerequisite on the graph translation were a convolutive and isometric operator. In the course of this chapter, we used the classical definitions of convolution and energy for graph signals, while showing their limits. As showed, this leads to the impossibility to recover classical shift operators. There are however several paths to explore to improve the current results. First, the definition of energy, which is at the core of the isometric property, can be improved to better take into account the irregularity of the graph structure. Second, since we are dealing with a complex operator, the question of a definition of analytic graph signals, negative frequencies, and complex Fourier modes are of interest to better cope with the complex nature of the graph translation. Third, defining multiple set of frequencies and Fourier modes, similarly to the case of graph products but for general graphs, would lead to a better representation of a graph signal. These are a few open questions worth investigating.

Summary

In this chapter, we introduced a new translation operator defined by similarity to the time shift. This is the graph translation. Depending on the graph Fourier transform chosen, this leads to two matrix representations of the operator:

T τ G = e -ıπτ L ρ G T τ G = e -ıπτ L 2 ,
with the Laplacian matrix-based approach on the left and the normalized Laplacian matrix-based approach on the right. The important properties of these operators are the convolutive nature and the isometry that lead to the power spectrum preservation.

Finally, we showed on examples that these operators are to some extent comparable to diffusion operators, and they preserve the slowly varying nature of low band graph signals through iterations of the graph translation.

Stationarity

Background

Temporal Signal Statistics

We start by introducing the basic tools from statistics that we will be using afterwards. Let x(t ) be a stochastic signal. Its statistics are described by its first order probability distribution:

F x (x, t ) = P[x(t ) ≤ x],
and its order k distributions:

F (k) x (x 1 , . . . , x k ; t 1 , . . . , t k ) = P[x(t i ) ≤ x i , ∀i ∈ {1, . . . , k}].
Differentiating according to all x i yields the probability density functions:

f (k) x = ∂ k F (k)
x ∂x 1 . . . ∂x k .

We will use several statistical properties to study stochastic signals, starting with the first order moment or mean:

η x (t ) = E[x(t )] = R x f x (x, t )d x,
and the autocorrelation:

R x (t 1 , t 2 ) = E[x(t 1 )x * (t 2 )] = R 2 x 1 x * 2 f x (x 1 , x 2 ; t 1 , t 2 )d x 1 d x 2 .
The average power (second order moment) of the signal is given by the autocorrelation through

E[|x(t )| 2 ] = R x (t , t ).
The autocovariance is built on the two quantities above to remove the influence of the mean:

C ov x (t 1 , t 2 ) = R x (t 1 , t 2 ) -η x (t 1 )η * x (t 2 ).
When the mean of the signal is identically zero, we have then

R x (t 1 , t 2 ) = C ov x (t 1 , t 2 ).
From the autocovariance we obtain the variance (and the standard deviation) through:

σ 2 x (t ) = C ov x (t , t ).
Finally, normalizing the autocovariance by the standard deviations yields then the correlation coefficient:

C or x (t 1 , t 2 ) = C ov x (t 1 , t 2 ) σ x (t 1 )σ x (t 2 ) .
More details on these quantities and the study of stochastic signals can be found in [Papoulis, 1991].

Stationary Temporal Signals

Apart from the quantities above, one important concept to study stochastic signals is that of stationary temporal signal. Such a signal is usually characterized by the "independence of [its] statistical properties relative to an absolute time" [Flandrin, 1999, p. 36]. This leads to the following definition of Strict Sense Stationary (SSS) signal: Definition 13 (SSS signal). A stochastic signal x is SSS if and only if it verifies for any τ:

T τ {x} d = x,
with T the time shift operator verifying T τ {x}(t ) = x(tτ).

In Definition 13, the sign d = means the equality in distribution of the two signals, i.e. x d = y if and only if ∀k, F (k) x = F (k) y . In practice, verifying that a signal is SSS is impractical since it involves verifying the equality of the distributions of all orders. We usually are satisfied with a weaker form of stationary signal where only its first and second order statistics are to be preserved through a change of the absolute time. These are the Wide (or Weak) Sense Stationary (WSS) signals:

Definition 14 (WSS signal). A stochastic signal x is WSS if and only if it verifies for any τ:

η T τ {x} = η x i.e. ∀t , E[x(t -τ)] = E[x(t )] (6.1) R T τ {x} = R x i.e. ∀t 1 , t 2 , E[x(t 1 -τ)x * (t 2 -τ)] = E[x(t 1 )x * (t 2 )] (6.2)
We remark first that using the autocorrelation or the autocovariance in Equation 6.2 does not change the definition since using Equation 6.1 we have η T τ {x} η * T τ {x} = η x η * x . Our second remark is that Equation 6.2 has the following equivalent expression:

E[x(t + τ)x * (t )] = R x (τ), (6.3)
with R x (τ) the autocorrelation function. In other words, the autocorrelation R x (t 1 , t 2 ) of the WSS signal x depends only on the time difference τ = t 1 -t 2 . These two equations are equivalent for temporal signals, but as we will see later in this chapter, this will not be the case for graph signals.

We can now study the frequency components of a WSS signal using the Wiener-Khintchin theorem:

Theorem 4 (Wiener-Khintchin Theorem). Let x(t ) be a WSS signal. There exists F x (ν) a monotone function on R called the power spectrum distribution such that (Fourier-Stieljes integral):

R x (τ) = R e ı2πτν d F x (ν).
Under mild assumptions, we can differentiate F x (ν) to obtain the power spectrum density S x (ν) and the following Fourier transform of the autocorrelation function: .5) Another way of looking at the spectral components of stochastic signals stems from the Cramér harmonic decomposition of a signal x(t ): .6) with the equality interpreted as a Mean Square limit (MS limit). The Fourier transform is then (MS limit):

R x (τ) = R S x (ν)e ı2πτν d ν (6.4) S x (ν) = R R x (τ)e -ı2πντ d τ. ( 6
x(t ) = R x(ν)e ı2πνt d ν, ( 6
x(ν) = R x(t )e -ı2πνt d t , ( 6.7)
From this, we can compute the autocorrelation of the Fourier transform (without assuming yet that x(t ) is WSS):

E x(ν 1 ) x * (ν 2 ) = R 2 R x (t 1 , t 2 )e -ı2π(ν 1 t 1 -ν 2 t 2 ) d t 1 d t 2 .
If we now assume that x is WSS, we have R x (t 1 , t 2 ) = R x (τ), hence:

E x(ν 1 ) x * (ν 2 ) = R 2 R x (τ)e -ı2π(ν 1 (τ+t 2 )-ν 2 t 2 ) d τd t 2 = S x (ν 1 ) R e -ı2π(ν 1 -ν 2 )t 2 d t 2 = S x (ν 1 )δ(ν 1 -ν 2 ).
We can also show the following for the mean of the Fourier transform:

E [ x(ν)] = R η x e -ı2πνt d t = η x δ(ν).
The converse is also true, leading to the following spectral characterization of WSS signals [Papoulis, 1991, page 418]:

Theorem 5 (WSS Spectral Characterization). A stochastic signal x(t ) is WSS if and only if: .9) WSS signals are therefore characterized by a double orthogonality principle [Flandrin, 1999]:

E [ x(ν)] = η x δ(ν) (6.8) E x(ν 1 ) x * (ν 2 ) = S x (ν 1 )δ(ν 1 -ν 2 ). ( 6
• Fourier Modes: The Fourier basis {e ν } ν is an orthogonal basis of functions: R e ı2πνt e -ı2πν t d t = 0 if ν = ν .

• Spectral Components: The spectral components of a stochastic signal are orthogonal in a statistical sense:

E x(ν) x * (ν ) = 0 if ν = ν .
We refer the interested reader to [Papoulis, 1991] and [Flandrin, 1999] for details.

Discrete and Periodic Time Series

Let x[n] be an N -periodic time series. We represent x through the column vector x[0], . . . , x[N -1] T . The mean can be represented by the vector

η x = E[x] = η x [0], . . . , η x [N -1]
T , and the autocorrelation of x is actually an autocorrelation matrix

R x = E[xx * ] such that R x nm = E x[n]x * [m] .
The definitions of stationary signals do not change when dealing with discrete time series:

x[n] is SSS if it is equal in distribution to its shifted time series x[n -k], for all k. x[n] is WSS if η x = η1 and R x nm = R x [n -m].
The WSS property implies a specific structure on the autocorrelation matrix: R x is circulant (similarly to circular convolutive operators of page 37):

R x =        R x [0] R x [1] . . . . . . R x [N -1] R x [N -1] R x [
[0] R x [1] R x [1] . . . . . . R x [N -1] R x [0]        .
This remark is important since our definition of WSS graph signals will not carry this property to the graph setting. We can also define the autocorrelation matrix of the spectral components of x, i.e. of x: S x = E[ x x * ]. For a WSS signal x, the matrix S x is then diagonal since Theorem 5 gives S x [ ν k , ν l ] = 0 if k = l . The spectral characterization of WSS periodic time series is then given by: Theorem 6 (WSS Time Series Spectral Characterization). x is a WSS time series if and only if:

E[ x( ν l )] = η x δ( ν l ) S x = diag (S x ( ν 0 ), . . . , S x ( ν N -1
)) .

Mean Square Estimation

We saw in chapter 4 that the Graph-Based Semi-Supervised Learning (GSSL) method can actually be interpreted as a Wiener filter in the setting of graph signal. We recall in this section the classical Wiener filters.

This type of filtering stems from the problem of performing estimation by minimizing the MSE in the context of WSS signals. Formally, having an observation s(t ) of a signal x(t ) on the time interval [a, b], we want to estimate the signal x(t ) using the estimator x(t ). The goal is now to find the estimator x(t ) minimizing the MSE:

E (x(t ) -x(t )) 2 .
For that purpose, we can use a linear estimator (filter):

x(t ) = b a s(u)h(t , u)d u.
The optimal choice of h minimizing the MSE is achieved when xx and s are orthogonal [Papoulis, 1991, page 481]: .10) which is equivalent to: .11) We assume now that [a, b] = R, and we would like to estimate the signal x from a noisy observation s = x + n. We suppose also that x and n (and then s) are WSS signals. To estimate x, we suppose that x is obtained using a convolutive noncausal operator: .12) The noncausal property stems from the fact that x(t ) depends on s(u) for u > t (in the future compared to the time instant t ). This is the simplest case for solving the problem of recovering x from s. From Equation 6.10, we know that x is optimal when xx is orthogonal to s. Setting α = tτ in Equation 6.11, we obtain:

∀α ∈ [a, b], E x(t ) - b a s(u)h(t , u)d u s(α) = 0, ( 6
∀α ∈ [a, b], R xs (t , α) = b a R s (u, α)h(t , u)d u. ( 6
x(t ) = R s(u)h(t -u)d u. ( 6
R xs (t , t -τ) = R R s (u -t + τ)h(t -u)d u = R R s (u)h(τ -u)d u = R xs (τ).
Using the convolution theorem in Theorem 2, we obtain the expression of h:

h(ν) = S xs (ν) S s (ν)
The convolution by h in Equation 6.12 is then a noncausal Wiener filter. Finally, if we also suppose that x and n are orthogonal (non correlated), h as a simpler form. We have R xs = R x , and R s = R x + R n , which leads to:

h(ν) = S x (ν) S x (ν) + S n (ν) = 1 1 + S n (ν) S x (ν) = 1 1 + 1 SNR(ν)
with SNR(ν) the Signal to Noise Ratio of s. In other words, knowing the power spectrum density of both the signal to recover and the noise, the noncausal Wiener filter is entirely defined and it gives an optimal estimator of the signal in the sense of the MSE.

Stationary Graph Signal

Now that we have introduced statistical properties of stochastic signals, we can present the framework of stationary graph signals that we proposed by analogy to the classical case. In this chapter, we always suppose the graph Fourier transform to be unitary, i.e. F -1 = F * .

Graph Signal Statistics

We will use similar quantities to describe the statistics of a graph signal than for a temporal signal. Remember first that a graph signal x is represented by a column vector, and as such, an operation between two graph signals has to be understood as an operation between two vectors. The distribution function of order k of x is:

F (k) x (i 1 , . . . , i k ; x 1 , . . . , x k ) = P[x i j ≤ x j , ∀ j ∈ {1, . . . , k}],
and the density function function of order k of x is:

f (k) x (i 1 , . . . , i k ; x 1 , . . . , x k ) = ∂ k F (k)
x ∂x 1 . . . ∂x k .

We have then the mean and the autocorrelation:

η x = E [x] R x = E xx * .
Applying the graph Fourier transform on the autocorrelation matrix yields the spectral autocorrelation matrix:

S x = F R x F -1 = E x x * .
In particular, S x is the autocorrelation matrix of the spectral components of x. We will call it the spectral autocorrelation matrix of x in the rest of this chapter. Remark that, contrary to R x , the matrix S x depends on the structure of the graph through the Fourier matrix F . In other, words, given a graph signal, its spectral autocorrelation matrix changes depending on the graph structure chosen, or the graph Fourier transform. This will be crucial to understand why a graph signal is stationary or not.

Definitions

We now leverage the unitarity of both the graph Fourier transform and the graph translation to define a tractable definition of stationary graph signals.

∀k, F (k) x = F (k) T τ G x .
Similarly to the classical case, verifying the SSS property is impractical. We therefore propose the WSS definition: Definition 16 (WSS Graph Signal). x is a Wide (or Weak) Sense Stationary (WSS) graph signal if and only if its mean and its autocorrelation matrix are invariant through graph translation, i.e. we have for any τ: .14) We can now study the spectral characterization of WSS graph signals. Let x be a WSS graph signal. We have then for the mean of the spectral components:

η T τ G x = η x i.e. ∀i , E[(T τ G x) i ] = E[x i ] (6.13) R T τ G x = R x i.e. ∀i , j , E[(T τ G x) i (T τ G x) * j ] = E[x i x * j ] ( 6
∀τ, E[ x(l )] = E T τ G x l = e -ı2πτ ν l E[ x(l )].
Therefore, E[ x(l )] = 0 if ν l = 0. Using both the graph Fourier transform based on the Laplacian matrix or the normalized Laplacian matrix, we have ν l = 0 if and only if l = 0 since we supposed that graphs are connected. We obtain then:

E[ x] = E[ x(0)]δ 0 .
Next, we study the spectral correlations of x. First, we have:

S T G x = F R T G x F * = F E (T G x)(T G x) * F * = F T G F * F E[ x x * ]F * F T * G F * = T G S x T G * ,
since T G := F T G F * . Then, using Equation 6.14, we obtain S T G x = S x , such that:

S x = T G S x T G * .
Since the graph translation is a convolutive operator, we have T G diagonal, and the elements of the right hand side matrix above are easy to compute. We obtain:

(S x ) kl = e ı2π( ν k -ν l ) (S x ) kl .

Using the same argument on the reduced graph frequencies as for the mean, we obtain:

∀ ν k = ν l , (S x ) kl = 0.
We can now state the spectral characterization of WSS graph signals:

Theorem 7 (WSS Spectral Characterization). A stochastic graph signal x is WSS if and only if:

E[ x] = E[ x(0)]δ 0 (6.15) ∀ ν k = ν l , (S x ) kl = E[ x(k) x(l ) * ] = 0 (6.16)
If all reduced graph frequencies are distinct, Equation 6.16 is actually equivalent to S x diagonal, which is exactly the same characterization as for time series (Theorem 6).

Discussion

There are several comments to make on these definitions and on the spectral characterization of Theorem 7. First, we consider the SSS definition (Definition 15). In the classical case, we have x n d = x m . For graph signals, we do not have the direct equivalent of x i d = x j for two vertices i and j . This is due to the graph translation not being an operator translating energy from one vertex to another, contrary to the time shift that translates energy from one vertex (time instant) to another. Also, compared to temporal signals where temporal samples are regularly spaced, vertices of a graph can show very different neighborhoods from one vertex to another due to the irregularity of the structure. As such, if two vertices show the same distribution but with different contexts, it is natural to consider that the signal is not equivalent on both vertices.

This remark leads to the definition of WSS graph signals (Definition 16). We begin with Equation 6.13 and the mean of the graph signal. As shown in Theorem 7, the mean of a WSS graph signal is reduced to a DC graph signal, and as such it depends on the graph Fourier transform chosen. For the Laplacian approach we have (η x ) i = (η x ) j similarly to temporal signals (χ 0 ∝ 1). On the other hand, for the normalized Laplacian approach the mean depends on the vertices such that d

-1/2 i η x i = d -1/2 j η x j (ψ 0 ∝ [ d 1 , . . . , d N ] T ).
Contrary to WSS temporal signals, WSS graph signals may not show a constant mean accross the graph. Depending, on the task at hand, we will choose the graph Fourier transform that best fits.

The second remark about Definition 16 is on Equation 6.14 and the autocorrelation matrix. While in the classical case, a WSS signal is associated to an autocorrelation function R x (τ), in the graph framework, there is no such function. This results from the fact that R x is not circulant. Intuitively, if it were to be circulant, the difference ij between two vertices i and j would be meaningful. In the temporal setting, the difference between two time instant is the time difference, i.e. a measure of how far (in time) the two samples are. In the graph setting, ij is not in any way a measure of the distance between vertices i and j . This also shows that using an analogy of Equation 6.3 to define WSS graph signals, i.e. through the relation:

∀i , j , E[(x) i + j x * i ] = R x ( j ),
does not make sense since i + j may not be a vertex, and even if it is, j is not a measure of how far vertices i + j and i are.

Our next remark is on the invariance operator chosen to define stationarity. We used here the graph translation for that purpose. We could actually use a different operator and still have the same spectral characterization of WSS graph signals. Let H be a convolutive and isometric operator with 1 ∈ Sp(H ) (i.e. there exists an x = 0 such that x is invariant through H : H x = x). Using Property 4 (and the assumption of a unitary graph Fourier transform), we have H = exp(ıΩ), with Ω = diag(ω 0 , . . . , ω N -1 ). Using the same technique as for the graph translation, we obtain the following spectral characterization of WSS graph signals with respect to H rather than T :

E[ x(l )] = E[ x(l 0 )]δ l 0 (l ) ∀ω k = ω l , (S x ) kl = E[ x(k) x(l ) * ] = 0,
with l 0 such that e ıω l 0 = 1. Except for ω l 0 = 0, all other ω l can be freely chosen while obtaining the exact same spectral characterization than the graph translation. The choice of the ω l we made when defining the graph translation stemmed from a parallel with the time shift, and as we saw in chapter 5, the resulting operator has a behavior we can justify.

We focus now on the complex nature of the graph translation. Indeed, as we saw in chapter 5, this operator is complex, and applied to a real graph signal, it gives a complex graph signal. However, we observe that even with a complex valued graph signal T G x, its autocorrelation matrix R T G x is real since R T G x = F * S T G x F and both F and S T G x are real matrices. In particular, this shows that the statistics of the real and imaginary parts of

T G x are linked through R T G x = E[(T G x)(T G x) * ]
. Therefore, T G x is not an arbitrary complex valued graph signal, and a WSS complex valued graph signal is not guaranteed to be the translated of a real graph signal. Furthermore, we can actually characterize WSS graph signals using the complex nature of the autocorrelation matrix. We haven proven above that a WSS real graph signal verifies R T G x real. We now look at the opposite to prove the equivalence, i.e. whether a real graph signal x with R T G x real is a WSS graph signal. Let x be a real stochastic graph signal whose translated has a real autocorrelation matrix R T G x . We obtain then for the elements of the spectral autocorrelation matrix:

S T G x kl = e ı2π( ν l -ν k ) (S x ) kl
But since R T G x is real, S T G x kl is real. Finally (S x ) kl is real, and we have:

∀ ν l = ν k , (S x ) kl = 0,
which is exactly Equation 6. 16. If we suppose also that E[ x] = E[ x(0)]δ 0 (Equation 6.15), we obtain the WSS spectral characterization of Theorem 7. We have then another characterization of WSS graph signals, this time in the vertex domain:

Theorem 8 (Vertex Domain WSS Characterization). A stochastic graph signal x is WSS if and only if:

E[x] ∝ F 0 R T G x is real,
where F 0 is the first Fourier mode (χ 0 for the Laplacian matrix-based approach and ψ 0 for the normalized Laplacian matrix-based approach).

Finally, we look at SSS graph signals. While a real graph signal does not hinder the power of the WSS definition, this is not the case for the SSS definition. Let x be a real SSS graph signal. Definition 15 gives the following relation between probability density functions:

∀k, f (k) x = f (k) T G x .
We claim (see Appendix A for a proof) then that we have:

∀l ∈ {0, . . . , N -1}, ∀ x ∈ R, f (1) x (l ; x) = f (1) T G x (l ; x).

Using T G x (l ) = e -ı2π ν l x(l ) such that x(l ) = e ı2π ν l T G x (l ), we finally have:

∀l ∈ {0, . . . , N -1}, ∀ x ∈ R, f (1) x (l ; x) = f (1) x (l ; e ı2π ν l x).

Since x is real, x is real and f (1) x (l ; x) is non zero only when x is real. However, when ν l ∈ (0, 1 2 ), e ı2π ν l x is not real, such that f x (l , x) is non zero only if ν l = 0 or ν l =1 2 . In other words, using the Laplacian matrix-based approach, the graph signal x = aχ l , with a a random real number, is not a SSS signal if 1 0 < l < N -1. This is quite different from temporal signals where the signal x(t ) = cos(2πνt + ϕ), where ϕ follows a uniform law in [0, 2π] and ϕ does not depend on t , is a SSS signal without enforcing ν = 0.

Finally, we can extend the theory of MSE to stochastic graph signals. In particular, if s = x + n is a noisy observation of the graph signal x. The optimal solution in the sense of minimizing the MSE is then achieved for xx and s orthogonal. Following the definition of convolution for graph signals, we can show that the optimal filter verifies:

h(l ) = S xs (l ) S s (l ) .
We call this filter a Wiener filter on graph. This shows that the name of Wiener filter used in section 4.3.1 is coherent with the theory of stationary graph signals developed in this chapter.

Sources of Non-Stationarity

First of all, similarly to temporal signals, stationarity of graph signals is characterized by a double orthogonality of both the basis of Fourier modes and the spectral components of a signal. This can be observed with WSS graph signals using Theorem 7:

• Fourier Modes: The Fourier basis {F l } l is an orthogonal basis of functions since F is unitary.

• Spectral Components: The spectral components of a stochastic signal are orthogonal in a statistical sense:

E x(l ) x * (l ) = 0 if ν l = ν l .
This shows that leaving the cocoon of WSS graph signals can be done either by relaxing the orthogonal spectral components thus allowing spectral correlations, or by keeping the double orthogonality but changing the basis of decomposition [Flandrin, 1999, page 38].

This remark is the basis of interpretation of the non-stationarity of graph signals. Indeed, if we first look at the basis of decomposition of a graph signal, we remark that the graph Fourier transform depends on the graph structure. Indeed, if we change the structure to something completely different, the graph Fourier transform also changes completely. In this context, while the autocorrelation matrix R x of a graph signal x does not depend on the graph structure, the spectral autocorrelation matrix S x does and it is highly unlikely that S x is diagonal regardless of the structure. We will experiment with this in the next section.

There are however graph signals with a spectral autocorrelation matrix that is independent of the structure: the Gaussian white noise. We define them as graph signals with a flat PSD S = σ 2 I . By unitarity of the graph Fourier transform, we have R = σ 2 I and no matter the unitary graph Fourier transform chosen, S = σ 2 I . This is actually quite natural: Since R = σ 2 I , there is no correlation between vertices such that the structure plays no role in explaining the missing correlations.

Note that one question arises from the discussion on the influence of the graph structure: Can we compute the graph that describes best as a WSS graph signal a given set of realizations of a stochastic graph signal? In other words, can we change the structure such that the data are realizations of a stochastic WSS graph signal? Up until now, the graph structure is given beforehand by the relation between vertices, but we can imagine the structure being given by the data itself, or refined by the data. Such an approach is different from Principal Component Analysis (PCA) (Karhunen-Loève Transform) that focuses on building an orthogonal basis to analyze the data and from which the underlying graph structure is not obvious. Here the goal is to define the edges and their weights such that the graph Fourier transform they define is such that the data corresponds to realizations of a WSS graph signal. This will be the subject of future works.

Coming back to the double orthogonality principle above, a poorly chosen graph structure is not the only explanation to a non-stationary graph signal. Suppose we have a WSS graph signal x. Suppose also that we observe this graph signal with one perturbed reading on a given vertex i . For example, the perturbation can be a higher variance: .17) with r ∼ N (0, 1). This does not introduce any perturbation on the mean. On the other hand, we have R y = R x + E i i , with E i i the matrix with zeros everywhere except for a 1 at row i and column i . For the spectral components, we obtain then:

y = x + rδ i , ( 6
S y = S x + F E i i F * = S x + (F δ i )(F δ i ) * .
However, since δ i is localized on a single vertex, we do not expect F δ i to be localized on a single graph frequency, even if some graphs show very high localization of Fourier modes (see [START_REF] Agaskar | A Spectral Graph Uncertainty Principle[END_REF]). Therefore, S y is not diagonal, and y is not WSS.

Experimental Studies

Now that we have introduced the definitions and characterizations of stationary graph signals, we experiment with these definitions. This section is divided into two subsections. The first one deals with synthetic data, both graphs and graph signals. The second subsection has been the subject of a communication [6], and introduces a methodology to study series of graph signals.

Toy Signals

We are interested in this section in illustrating stationary graph signals, but also their nonstationary counterpart. Indeed, as we saw in the previous section, non-stationarity pertains to both the data studied and the graph used to compute the graph Fourier transform, and we now illustrate this on examples. We first show how to generate data following a given mean and autocorrelation matrix, and then implement this with several stochastic graph signals and graphs.

Stochastic Signal Synthesis

Before actually experimenting with graphs and toy graph signals, we explain how to generate realizations of Gaussian graph signals. The Gaussian hypothesis is necessary to completely describe the graph signal statistics using only its mean and its autocorrelation matrix. The only constraint we have is that the matrices C ov x = R xη x η * x and C ov x = S xη x η * x need to be semi-definite positive to be valid autocovariance and spectral autocovariance matrices. This classical result on autocovariance matrices is proved for a signal x using the following inequality:

∀y, y * C ov x y = E ((x -η x ) * y) * ((x -η x ) * y) = E |(x -η x ) * y| 2 ≥ 0.
The same inequality can be derived for C ov x .

The semi-definite positive property allows to use the Cholesky decomposition on these matrices. First, we consider the Cholesky decomposition of the matrix C ov x = LL * , with L a lower triangular matrix given by the Cholesky decomposition (see section 1.2.4, page 12, for details). To generate a realization of x, we use the graph signal y such that the y i are independent and identically distributed according to the standard Gaussian distribution

Experimental Studies

N (0, 1). In other words, η y = 0 and C ov y = R y = I . Note that y is a Gaussian white noise graph signal. We have then η x + Ly = x since:

E[µ x + Ly] = η x + LE[y] = η x E[(η x + Ly)(η x + Ly) * ] = η x η * x + LL * = η x η * x +C ov x = η x η * x + R x -η x η * x = R x ,
and each of the x i is Gaussian with x i ∼ N (µ x ) i , (C ov x ) i i . The same approach can be used with C ov x to generate realizations of x instead of x. This generation procedure does not assume that x is stationary, only that it is Gaussian.

Stationary Graph Signals

We (arbitrarily) choose two stochastic graph signals x 1 and x 2 to illustrate key properties of WSS graph signals. They are defined using their PSD: We chose a non-zero mean to illustrate the influence of the zero frequency on the autocorrelation matrix R. Indeed, all three graph Fourier transform presented in chapter 3 are unitary graph Fourier transform, and as such can be used to define the graph translation and stationarity. Figure 6.2 through Figure 6.5 show then the same stochastic graph signal on four different graphs (the cyclic graph, the torus graph, the Molène graph and the MOSAR graph). These pages are split into a total of six columns. The first division into 2 columns corresponds to the two stochastic graph signals x 1 and x 2 . Note that since they are defined in the spectral domain through Equation 6.18, their exact definition in the vertex domain depends on the graph Fourier transform. Each of these transforms corresponds then to a subcolumn, hence the total of 6 columns. Finally, we arbitrarily chose the following definition of reduced graph frequency for the graph Fourier transform based on the adjacency matrix:

S (1) x (l ) = 1 1 + ν l 2 S (2) x (l ) = 1 1 + | ν l -0.2| 2 , ( 6
ν l = 1 4 1 - γ l γ max ,
such that they correspond to the interval [0, 0.5]. This choice is indeed arbitrary, but it is not essential to the interpretation of the figures as long as the results from the graph Fourier transforms are not compared to one another. (first) and the signal x i + a(δ 1 + δ 10 ) with a ∼ N (0, 1) (second). Finally, the two last rows show the impact of a "wrong" graph (Erdős-Rényi graph with p = 0.15) on the spectral autocorrelation matrix, with an illustration of the structure on the last row.

x 1 x 2 L L A L L A Realizations -1 0 1 -1 0 1 -1 0 1 -1 0 1 -1 0
x The first three rows of each figure show three realizations of the stochastic graph signals described above. Their empirical autocorrelation matrix R x (forth row) and their spectral autocorrelation matrix S x (fifth row) are computed on 10 000 realizations. As expected from these WSS graph signals, the matrix S x is diagonal.

1 x 2 L L A L L A Realizations -1 0 1 -1 0 1 -1 0 1 -1 0 1 -1 0
On the other hand the matrix R x is more interesting. Indeed, on the cyclic graph C 30 and the torus graph T 10 , the autocorrelation matrix is circulant. This is actually expected. First, we observe that the PSD is a function of the reduced graph frequencies, i.e. if two Fourier modes are associated to the same frequency ν l = ν l +1 , then their PSD is equal: S (i ) x (l ) = S (i ) x (l + 1). We are therefore free to choose any Fourier basis within a particular stable subspace and keep the same matrix S x (see discussion of section 3.2, page 38, on the link between convolutive operators and graph Fourier transform). In turn, this observation allows us to use the DFT to show that S x corresponds to a circulant matrix R x (see section 6.1.3).

Also, the circulant structure in the torus graph stems from its nature of Cartesian product of two cyclic graphs. On the other hand, the matrix R x is not at all circulant on the Molène graph and the MOSAR graph (Figure 6.4 and Figure 6.5) as expected from real irregular graphs where the difference ij between two vertices is meaningless. In this case, the underlying structure of the autocorrelation matrix R x is not obvious, and it is only captured by the diagonal spectral autocorrelation matrix S x .

Non-Stationary Graph Signals

The next graph signals we study are non-stationary graph signals. For that purpose, we start with the two graph signals x 1 and x 2 , and introduce two correlations in the vertex domain.

The first correlation is actually a variance increase on vertex 1 (chosen arbitrarily). This is the graph signal of Equation 6.17 (page 101). In the temporal setting, such an increase of variance y(t ) = x(t ) + aδ τ (t ) at a single time τ on a WSS signal x clearly breaks the stationarity since the shifted signal is then not statistically the same: T y)(t ) = x(t ) + aδ τ-1 (t ). In the graph setting, this is also the case, as discussed in section 6.2.3, and we can observe a non-diagonal spectral autocorrelation matrix S x on the figures (sixth row).

The second correlation we consider is a correlation between two arbitrarily chosen vertices: 1 and 10. For that purpose, we generate realizations of the signal y = x + a(δ 1 + δ 10 ), with a and x orthogonal, and a following the standard Gaussian distribution N (0, 1). It follows that:

η y = η x R y = E[(x + a(δ 1 + δ 10 ))(x + a(δ 1 + δ 10 )) * ] = R x + (δ 1 + δ 10 )(δ 1 + δ 10 ) * .
Similarly to the previous correlation, in the temporal setting this would break the stationarity of x, and we observe the non-diagonal spectral autocorrelation matrices that are characteristic to non-stationary graph signals on the figures (seventh row).

Influence of the Graph Structure

Finally, the last illustration of non-stationarity is on the influence of the structure. Indeed, as discussed in section 6.2.3, stationarity is defined with respect to a given graph Fourier transform. If we change the graph structure, the graph Fourier transform changes. In this section we change the structure of the graphs of Figure 6.2 through Figure 6.5 to something completely different and look at the resulting spectral autocorrelation matrix of the WSS graph signals x 1 and x 2 .

These "wrong" graph structures are an Erdős-Rényi graph (with p = 0.15) for the cyclic graph, a random vertex sampling of the plane for the torus graph (random coordinates in 2D are assigned to each vertex of the torus and a Gaussian kernel of the Euclidean distance is used to weigh edges), and finally, arbitrary cyclic graphs for the Molène and MOSAR graphs. The last row of the figures show realizations of x 1 and x 2 in the context of the "wrong" graph. We observe on the matrices above them that the spectral autocorrelation matrices are not diagonal with a lot of off-diagonal elements spread everywhere in the matrix. This illustrates that stationarity is indeed dependent on the underlying graph, and non-stationarity may as well be the image of a poorly chosen graph, i.e. a structure that does not explain the correlations of the graph signal.

Weather Dataset

We introduced in section 2.7.2 a dataset comprising geographical locations of weather stations together with weather reports. We have been using extensively the graph of the weather stations that we call the Molène graph W . In this section, we are interested in the weather reports in the context of W . Among many readings, three are of interest to us since they span the largest number of stations and snapshots: the temperature t, the wind speed w and the rain r. The graph shown until now has been the graph associated to the temperature, where vertices are weather stations with temperature readings for all 744 snapshots (one month of hourly readings). Each reading is associated to a different graph where the vertices are the weather stations with the maximal number of readings (a few stations are missing readings of either the temperature, the wind, or the rain). The graphs are shown on Figure 6.6. This section has been the subject of a communication [6].

Pre-Processing

We now have N time series of length 744 corresponding to one month of hourly readings. Our goal in this section is to consider these N time series as 744 graph signals, i.e. one graph signal per snapshot. Furthermore, we want to consider these 744 graph signals as 744 realizations of the same stochastic graph signal.

As illustrated on Figure 6.6, the data we want to study have a significant trend (slightly decreasing temperature over the month, slightly decreasing wind speed) and is cyclostationary, with daily cycles for example. As a consequence, assuming the ergodicity of the data collected to estimate statistical quantities such as the mean and the autocorrelation matrix is not reasonable.

To solve this problem, we propose to use the Empirical Mode Decomposition (EMD) [START_REF] Huang | The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[END_REF]]. Essentially, this decomposition interprets a given time series with several Intrinsic Mode Functions (IMFs). These IMFs are time series ranging from low frequency (last IMF) to high frequency (first IMF) and summing up into the original time series. In this section, we are interested in the first IMF, on which an ergodicity assumption is reasonable since its high frequency nature removes the trend and cyclostationarity influences. In other words, we consider the temporal samples of the first IMF as realizations of the same random variable, and these random variables (one per vertex) constitute a stochastic graph signal. Eventually, we study the stationarity of this graph signal through its empirical mean and autocorrelation.

For a given vertex i, we denote IMF (k) t i the k th IMF of the temperature data, IMF (k) w i that of the wind speed, and IMF (k) r i that of the rain. Figure 6.6 shows the original temporal data of the three circled vertices (second row) together with their first IMF (third row) and second IMF (bottom row). Figure 6.7, Figure 6.8, and Figure 6.9 show several snapshots of the original data with its temporal mean (first row), together with IMF (2) and its mean (second row), and IMF (1) and its mean (third row). Weather readings over time of three ground stations of (a) the temperature, (b) the wind speed, and (c) the rain. The first row shows the weather stations considered on the Molène graph. Note that the vertex set changes between readings since not all stations have readings over the whole period. The second row shows the original data, the third row shows IMF (2) , and the fourth row IMF (1) .

At this point, there is one important remark to make about the rain data. Indeed, the EMD fails to decompose some of the time series because they show saturation phenomena: There are too many snapshots without precipitation in Brittany. This is illustrated on Figure 6.6(c) where we observe sinusoid segments in IMF (1) r and in IMF (2) r . Unfortunately, when no precipitation occurs, we would like a first IMF that is zero for all vertices, rather than sinusoids having no real meaning from the data. These data show the limit of our method based on the EMD to obtain a good ergodicity hypothesis. In the rest of this chapter, we are not considering the rain readings for this reason.

The question we have now is whether the IMF (1) time series of the stations can actually The IMF (1) (first row), the IMF (1) with a cyclic graph as underlying graph structure (second row), the PSD model based on the Laplacian (third row), and the PSD model based on the normalized Laplacian (fourth row). See Figure 6.11 for these models. Note that for all those matrices, the moduli of their elements is shown. Finally, the color scales on R and S are different to better account for the differences in values of these matrices. This does not hinder the interpretation of the figure since R and S are two different objects that cannot be compared. be seen as 744 realizations of a WSS graph signal. A related question is whether the graph structure we have chosen can interpret these correlations as the signature of a WSS graph signal.

Qualitative Study of Stationarity

First, we look at the mean of the data on Figure 6.7 and Figure 6.8. We see here that the original data does not have a mean related to the vertex degrees. For example, on Figure 6.7, the station on an island at the south of Brittany is far from other stations (smaller degree), yet it has the same mean temperature than the group of stations (higher degree) directly to its east. In other words, the mean of the original graph signals is not a DC component since it is neither flat (Laplacian based graph Fourier transform) nor proportional to d i (normalized Laplacian based graph Fourier transform). On the other hand, the mean is much closer to zero with IMF (2) and IMF (1) , where by construction, only the last IMF has a non zero mean (equal to the mean of the global trend of the decomposed time series), with lower values for the first IMF.

As a consequence, we study the stationarity of IMF (1) only. Note that, we could remove the mean from the original data, but then, the trends and the cyclostationarities remain and still prevent us from using an ergodicity hypothesis.

Now that Equation 6.15 has been verified (the mean of the considered stochastic graph signal is a DC component), we look at Equation 6.16 (and whether the spectral autocorrelation matrix is diagonal). To this end, we compute the autocorrelation matrices R and the spectral autocorrelation matrices S for each data and for both the Laplacian based and normalized Laplacian based graph Fourier transforms. To better appreciate stationarity, we also compute the correlation coefficient matrices of both the graph signals and their graph Fourier transforms. We denote them C or and C or . Note that the notation with an hat is actually an abuse of notation, and we cannot use the graph Fourier transform F in its usual way to obtain C or from C or . Indeed, the correlation coefficient is obtained with the autocovariance matrix C ov and its diagonal Σ through:

C or = Σ -1/2 C ovΣ -1/2
which is not a linear operation, and since C ov = F C ovF * , then C or = F C or F * . Figure 6.7 and Figure 6.8 are read as follows: The matrix S gives on its diagonal the PSD, while the matrix C or gives the correlation coefficients between spectral components. The same applies for R giving the variance on its diagonal (remember that the mean is 0 such that R = C ov) and C or giving the correlation coefficients between vertices.

We see on these matrices that both IMF (1) of the temperature and the wind data are not realizations of a WSS graph signal since there are non zero spectral correlations. Nevertheless, there are several observations to make. First, we observe a PSD in 1/ν α on the temperature (diagonal of S for IMF (1) ). We will look at this in the next section. Second, even if there are spectral correlations, they are limited. In particular, we observe for the temperature small amplitudes on the correlation coefficient matrix C or , and very few correlations for the wind data. Furthermore, we observe that two highly correlated spectral components (l = 7 and l = 8) in Figure 6.8 (circled in blue) are actually between two Fourier modes whose support are overlapping. These modes are shown on Figure 6.10. Finally, to highlight the fact that the graph explains many correlations in the vertex domain, we look at the same data on a different graph: an arbitrary cyclic graph. This graph breaks the structural meaning of the Molène graph by connecting the vertices in an arbitrary order. On Figure 6.7 and Figure 6.8 (second row of matrices) the spectral autocorrelation matrices with respect to the new cyclic graph show two things. First, and this is particularly clear for the wind, there are more correlations between spectral components. Second, the PSD on the spectral autocorrelation matrix is flatter (almost constant) such that the information lies in the spectral correlations rather than in the PSD, as opposed to when we use the geographical graph of the stations. Our choice of graph is therefore more relevant than the cyclic graph.

Models

To close this section on weather data, we look at the PSD of the IMF (1) . We observed on the temperature data what looks like a 1/ν α function. On the wind data, the maximum of the PSD is shifted to a non zero frequency. We then model these PSDs with two functions β/ 1 + γ ν l α and β/ 1 + γ| ν l -δ| α . We then fit2 the PSD through linear regression to find the parameters α, β and γ (and δ for the wind). Results are shown on Figure 6.11. The fit is quite good for the temperature, but less perfect on the wind. This may be due to the irregularities of the structure (and of the graph frequencies), or the limited number of samples to compute S, or the mere fact that we are dealing with real data. This example is for illustrative purposes and future works will refine the method. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 7.2 Basic Elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 7.2.1 Presentation of the Toolbox . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 7.2.2 Graphs and Graph Signals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 7.2.3 Operators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 7.3 Visualization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 7.3.1 Layout . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 7.3.2 Plotting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 7.3.3 High Level Visuals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 7.4 Interface with L A T E X . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120 During the course of this PhD, we developed an extensive MATLAB toolbox1 dedicated to handling and studying both graphs and graph signals. Most noticeably, this toolbox comprises functions to quickly show graph signals. Additionally, we defined several L A T E X commands to display directly in L A T E X/TikZ a graph signal. All figures of this dissertation have been generated using our toolbox and these L A T E X commands. In this short chapter, we present their strengths.

Matlab Toolbox

State of the Art

When this PhD started, only one toolbox related to signal processing on graphs was available: the Spectral Graph Wavelet Transform (SGWT) toolbox of [START_REF] Hammond | Wavelets on graphs via spectral graph theory[END_REF]. This toolbox is dedicated to efficiently computing the SGWT. Inspiration for our toolbox came from it. In contrast, our goal has been to generate a wide range of graphs and efficiently work with graph signals without worrying about how to display both a graph and one or more graph signals.

The second toolbox that can be found as of today is the Graph Signal Processing Toolbox presented in [START_REF] Perraudin | GSPBOX: A toolbox for signal processing on graphs[END_REF]. This toolbox has similar goals than ours, except for the graph plotting part that is less extensive.

Basic Elements

Our toolbox is split into two classes of functions, the elementary functions presented in this section, and those dedicated to plotting in the next. Elementary functions comprise toolbox initialization, graph creation, signal handling, and several operators.

Presentation of the Toolbox

The toolbox is initialized the first time using grasp_install. This function downloads the third party libraries necessary for the toolbox to run. Using this approach, we do not need to bundle the third party libraries with the toolbox and software licensing is not an issue. To start the toolbox, the function grasp_start sets MATLAB path with third party libraries and the various directories of the toolbox. The toolbox has an example script at the root directory named grasp_translation_page.

This script generated the content of the page http://perso.ens-lyon.fr/benjamin.girault/pages/ graph_translation-en.html that shows the operators of chapter 5 (the generalized translations, the graph shift and the graph translation operators). This should be a starting point for anyone wanting to use our toolbox.

Graphs and Graph Signals

We have implemented all graph structures presented in this dissertation, including generators of the graph models of section 2.6. A few other graphs are also available. The graph representation of our toolbox and that of the GSP toolbox are similar with the adjacency matrix being at the core of the representation. Combining the two toolboxes to have access to more graphs is possible with just a few MATLAB instructions. Writing this function is part of the future features of our toolbox.

Additionally, we provide functions to manipulate graphs such as building the Gaussian kernel of the distance, thresholding and KNN edge selection schemes. Other matrices such as the degree matrix, the Laplacian matrix and the normalized Laplacian matrix, and their diagonalization is directly available. We also implemented the definitions of graph frequencies of this dissertation.

Finally, several graph signals used in this dissertation are available, from the simple deltas and heat kernels to the realization of a stochastic Gaussian graph signal using its mean η and its autocorrelation matrix R.

Operators

The operators used in this dissertation are available, and comprise the convolutive operator of section 3.2, the graph Fourier transforms and their inverse of chapter 3, the simulated annealing approach discussed in section 4. 1.1, the GSSL closed form solution of section 4. 1.4, the complete graph to time series mapping and its inverse mapping of chapter 4, the generalized translation of section 5. 1.4, and our graph translation of section 5.3.1. We are using the same natural convention of graph signal and operator than the two toolboxes of section 7.1, i.e. a column vector for the former and a matrix for the latter, and as such, we can use the operators of one toolbox with the other.

Visualization

The strength of our toolbox lies in this drawing capabilities. We explain now the visuals it can generate.

Layout

First of all, not all graphs have a natural layout of the vertices. For example, giving 2D coordinates to vertices of a social network is a very complex task. When no layout is provided with the graph structure, we use an interface with GRAPHVIZ, and delegate the task of finding a good layout to this tool. By default, a spring algorithm is used, i.e. edges are replaced with springs attracting the two end vertices and all non connected vertices repulse each other. Several other algorithm from GRAPHVIZ can be chosen.

Plotting

We now reach the core of the toolbox: the function show_graph. This function takes the layout of the graph provided, plots the vertices at the given coordinates and draws edges, with arrows for directed edges, between vertices. If a graph signal is provided, the color scale on the vertices reflects the interval spanned by the values (or any other interval provided). A slower option can associate values to edges, but it has not been optimized yet and does not scale well with the number of edges. Additional options comprise a background image, custom labels, highlighted vertices, or a custom vertex size.

For matrices, we implemented the function show_matrix that allows for irregular row heights and column widths contrary to the MATLAB built-in function imshow.

High Level Visuals

There is one function that is very handy to quickly study the Fourier modes of a given graph:

show_fouriermodes. This function displays either all modes, or a selected number of them in a panel of figures. Several options of show_graph can be given to show_fouriermodes that forward them to show_graph, such as the interval of values for the color scale or the size of the vertices.

Another function generates an animated GIF image file from several graph signals on a graph: grasp_generate_gif. This function is used on http://perso.ens-lyon.fr/benjamin. girault/pages/graph_translation-en.html to show iterations of the graph translation and the graph shift.

When writing show_graph, we tried to optimize as much as possible the number of MATLAB low-level graphics calls, such that updating a figure is fast. This is essential for grasp_generate_gif, but also for the GUI in animate_translations_gui (see Figure 7.1).

This GUI iterates an operator on a graph signal in real time and draws it on screen. Optimization is done on both the iteration part, and the drawing step such that the number of iteration per second is high. The speed still depends on the number of vertices, such that a high number of vertices slows down the process. Similar techniques can be applied to have an animation of a different set of graph signals.

Interface with L A T E X

We wrote in the introduction of this chapter that the figures of this dissertation are all generated by L A T E X. To do so, we rely on two parts. First, in MATLAB, the toolbox has two functions grasp_exportcsv and grasp_exportcsv_signal to generate CSV files describing the graph for the former and a graph signal for the latter. These CSV files are then the input of three L A T E X/TikZ commands:\flatgraph to plot a graph (the most common figure in this dissertation), \stemgraph to plot a graph in a slanted plane with vertical bars on vertices to describe a graph signal (see section 3.5.2, page 45), and \backgroundflatgraph to display a graph with a background image (such as the Molène graph). Note that these L A T E X commands can be used with CSV files generated by a different software than MATLAB.

The color scales used are defined such that printing in black and white should still show the necessary information. Finally, matrices are plotted with \tikzmatrix. These commands may be slow on huge graphs, but they only need to be compiled once with the use of the TikZ library externalize.

Summary

This rather short chapter explained the strengths of the plotting capabilities of the MAT-LAB toolbox we developed during the course of this PhD. Using L A T E X/TikZ commands, we were also able to use L A T E X to generate figures without the need to go back and forth between L A T E X and MATLAB during the writing of a communication.

(PSD) preserving operator. We then showed on several examples that this operator acts as a diffusion-like operator with a complex output.

With this graph translation, we have been able in chapter 6 to propose a tractable definition of stationary graph signals through a statistical invariance of a stochastic graph signal. Using the Wide (or Weak) Sense Stationary (WSS) definition, we showed that the spectral characterization of WSS graph signal is the same than in the the temporal setting. We then looked at the possible interpretations of non-stationarity in data and showed that, in contrast to the temporal setting, the graph structure also have an influence over the stationarity and it should be correctly chosen. After several illustrations on toy graph signals, we proposed a method to study a dataset of structured temporal data: weather reports. We saw that the graph signals we looked at were not WSS, but a PSD model was able to capture most of the correlations in the weather readings.

Challenges

During this dissertation, we faced many challenges, solved some of them, but a lot of challenges remain. First, we saw in chapter 2 that the Laplacian matrix and the normalized Laplacian matrix are only defined for symmetric graphs. In [Chung, 2005], the author introduces a normalized Laplacian matrix for directed graphs with non-negative weights. This approach is introduced in the context of random walks, but would be very interesting to use in the context of signal processing on graphs since this normalized Laplacian matrix for directed graphs is still Hermitian, hence unitarily diagonalizable.

Then, we saw two challenges on the topic of Fourier modes. First, we have observed that some Fourier modes (w.r.t. either graph Fourier transform) can have a very limited support, i.e. can be localized. Whereas Fourier modes in the temporal setting have an infinite support that is both their strength (natural spectral characterization of stationary signals) and their drawback (need for Time-Frequency or Time-Scale decompositions), we would like a similar basis to study graph signals. This is linked to the mutual coherence of the delta and Fourier bases that should be as low as possible, i.e. the delta basis is highly localized and the Fourier basis is highly delocalized. Coherence appears in [START_REF] Shuman | Vertexfrequency analysis on graphs[END_REF] and should be a start to characterize a good Fourier basis. The second question on the Fourier modes is whether it is possible to use a graph Fourier transform on a dynamic graph. Taking the graph Fourier transform at each snapshot is a start, but the link between Fourier modes of two consecutive snapshots is not obvious.

When we studied the problem of mapping graphs to time series, we found many open questions. First of all, we did not look at the properties of the resulting time series. It should be interesting to link graph properties and properties of time series. Second, several extensions were investigated such as increasing the number of centroids to create the mapping of vertices to amplitudes, or a mapping from a dynamic graph to a time series. These two subjects raise many interesting questions: Given more than 2 classes, how do we set the amplitudes based on the level of belongings to each class? or, given as many classes as vertices, is it possible to perform multiscale community detection with good results? or can we characterize the dynamics of a dynamic graph through a time series? Finally, we used GSSL to obtain the amplitudes, but there are alternative approaches such as Multi Dimensional Scaling used in [START_REF] Hamon | From graphs to signals and back: Identification of graph structures using spectral analysis[END_REF]] that can give one amplitude per vertex following the structure of the graph.

In chapter 5, we also raised many open questions. First, the complex nature of the graph translation leads to the question of defining analytic signals for graphs, but also to the question of graph Fourier transform with complex modes instead of real ones, and finally, are there graph frequencies that should be negative and not positive as we arbitrarily chose? We discussed about the definition of the energy of a graph signal and showed that the usual definition is most certainly not the best since it does not take into account the irregularity of the graph structure. This raises also the question of the energy x of the graph Fourier transform x of a signal x since the graph frequencies are not regularly spaced, similarly to the vertices. Eventually, we will have to study Parseval's identity x = x w.r.t. the new definitions of energies. Furthermore, we could use the same approach than [START_REF] Sandryhaila | Discrete Signal Processing on Graphs[END_REF] and define the class of filters invariant through graph translation to have an analogy with the Linear Time Invariant operators of the temporal setting.

Last but not least, the definition of stationarity we gave relies on the definition of graph translation of chapter 5. As such most questions on stationarity are linked to questions on the graph translation operator. We also remark that in the simpler case of images, the question of a stationary image is still not completely solved in the literature, and requires, for example, the use quaternions, with the Quaternionic Fourier Transform, and heavy tools from mathematics. Going from a 1D shift to a 2D shift is far from obvious, hence in the case of graphs the problem is even more complex. As a consequence, we do not expect from our definition to solve this complex problem. On the other hand, several challenges are offered to us. First of all, we have used a very simple stationarity characterization in chapter 6. We are now looking at a better, statistically correct, stationarity characterization, and the surrogates approach of [START_REF] Borgnat | Stationarization via surrogates[END_REF] is a promising approach to properly build a stationarity test. The next challenge is to compute a graph that describes best a given stochastic graph signal as a WSS graph signal, as suggested in section 6.3. To do that, we can for example start with the given graph structure and refine it through optimization. The advantage of such an approach is then to interpret it as structure denoising (as opposed to graph signal denoising), or to compare the resulting structure to the original one and derive conclusions on the locality of the studied phenomena. We have hopefully set the first stone for future works.

Proof of the Claim Page 100 APPENDIX A

In section 6.2.3, we claimed that an SSS graph signal verifies: ∀l ∈ {0, . . . , N -1}, ∀ x ∈ R, f (1) x (l ; x) = f (1) T G x (l ; x).

To show this, we remark first that x = F x. Using [Papoulis, 1991, Eq. (8-8) and (9-9), p. 183], we obtain:

f (N )
x ( x 0 , . . . , x N -1 ) = f (N ) x (x 1 , . . . , x N ) det(F ) , (A.1) with x 1 , . . . , x N T = F -1 x 0 , . . . , x N -1 T , and det(F ) the determinant of the graph Fourier transform matrix. Note that we wrote f (N ) x (x 1 , . . . , x N ) for of f (N ) x (1, . . . , N ; x 1 , . . . , x N ). Using the SSS property on x, we have:

f (N ) x = f (N ) T G x ,
and using Equation A.1, we obtain: det(F ) f (N ) x ( x 0 , . . . , x N -1 ) = det(F ) f (N ) T G x ( x 0 , . . . , x N -1 ), since F T G x = T G x. We finally have:

f (N )
x ( x 0 , . . . , x N -1 ) = f (N ) T G x ( x 0 , . . . , x N -1 ), and using the marginals of f (N ) x , we obtain:

∀l ∈ {0, . . . , N -1}, ∀ x ∈ R, f (1) x (l ; x) = f (1) T G x (l ; x).

Graph Signal Processing

Fourier transform

Operator computing the frequency components of a signal defined as x(t )e x ∂x 1 ...∂x k . 92 probability distribution Function F x (x, t ) = P[x(t ) ≤ x] (first order) or F (k) x (x 1 , . . . , x k ; t 1 , . . . , t k ) = P[x(t i ) ≤ x i , ∀i ] (order k) for a stochastic process x. 91

Figure 1 . 1 .

 11 Figure 1.1.A random matrix A with entries uniformly drawn in[-1, 1]. These two figures illustrate two color scales that are used throughout this dissertation for both matrices and graphs. The first color scale is essential to illustrate amplitudes with darker tones for higher amplitudes and lighter tones for smaller amplitudes. Yet the sign of the entries can be identified with red (positive) and blue (negative) colors. When all elements of the matrix are positive, a better color scale is the second one with more tones from white to black.

1. 1 .

 1 Fundamental Definitions and Properties matrix A = [a 1 . . . a n ] are then orthonormal:

Figure 2 . 1 .

 21 Figure 2.1. Example graphs illustrating two graph properties.

Figure 2 . 2 .

 22 Figure 2.2. Examples of graph products. In each of the sub-figures, the graph in the center is the product of the top and left graphs.

T

  . . . , c i (k) + c i (k + 1) c j (l ) + c j (l + 1) , .. .

Figure 2 . 3 .

 23 Figure 2.3. The graphs shown here all have N = 50 vertices uniformly sampled in the 2D unit square. Edges are all weighted by a Gaussian kernel w i j = expd(i , j ) 2 /2σ 2 1 of the Euclidean distance d(i , j ). In (a) all edges are present. In (b) only edges between close enough vertices are kept. Finally, in (c) only k edges to the closest vertices are kept, for each vertex.

Figure 2 . 4 .

 24 Figure 2.4.Examples of realizations of Erdős-Rényi graphs with N = 30 vertices and three probability p of edge appearance. Vertices are placed in 2D according to their connectivity using a spring model (all pair of vertices repulse each, and connected vertices attract each other with a force depending on the weight of the corresponding edge, such that the resulting stable state of the system gives the placement in 2D of the vertices).

Figure 2 . 5 .

 25 Figure 2.5. Examples of realizations of Barabási-Albert graphs of 30 vertices with three values of m. Placement of vertices is performed similarly to Figure 2.4.

Figure 2 . 6 .

 26 Figure 2.6. Examples of realizations of Watts-Strogatz graph of 30 vertices with two mean degree d and two rewiring probability β. Placement of vertices is performed similarly to Figure 2.4.

Figure 2 . 7 .

 27 Figure 2.7. Three aggregations of the MOSAR graph. (a) shows a single snapshot of the dynamic graph M . (b)shows all contacts between people within an hour. Finally, (c) shows all contacts within a day. Each edge is the support of the time series of the presence/absence of a contact. An edge appears on the figures when a contact occurs at least once during the time frame. The weighting scheme is not depicted. Ellipses on (a) correspond to group of people assigned to the same care units in the hospital.

Figure 2 . 8 .

 28 Figure 2.8. Molène graph of ground weather stations with several edge selection schemes. (a) has all possible edges (complete graph). Edges of (b) (resp. (c)) are limited to ground stations at a maximal distance of σ 2 = 100km (resp. σ 2 = 50km). The resulting graphs get sparser as σ 2 decreases. The choice of σ 2 should be small enough to get sparse matrices and reduce the complexity of the computations while being high enough for the graph to stay connected and keep its structural properties.

CHAPTER 3 Contents 3 . 1

 331 Background: Continuous and Discrete Signals .. . . . . . . . . . . . . . . . . . . . . 34 3.1.1 Continuous Domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 3.1.2 Discrete Domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 3.2 Signal Processing on Graphs: Notations . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3.3 Laplacian Based Fourier Transform . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.3.1 Definitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 3.3.2 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3.4 Adjacency Matrix Based Graph Fourier Transform . . . . . . . . . . . . . . . . . . . . 42 3.4.1 Definition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.4.2 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.5 Examples of Fourier Modes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.5.1 Cyclic Graph: Periodic Time Series . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.5.2 Cyclic Graph: The Importance of the Weights . . . . . . . . . . . . . . . . . . 45 3.5.3 Torus Graph: Image Pixel Graph . . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.5.4 Molène Graph: Sparsifying Graphs . . . . . . . . . . . . . . . . . . . . . . . . . 48 3.5.5 MOSAR Graph: Aggregating Dynamic Graphs . . . . . . . . . . . . . . . . . . 50

  integer. We obtain the Discrete Fourier Transform (DFT) of N -periodic time series: Definition 5 (Discrete Fourier Transform (DFT)).

Observation 1 [

 1 Hammond et al., 2011, Section 3.1]

Figure 3 .

 3 Figure 3.1 shows several Fourier modes of the cyclic graph C 30 associated to the three graph Fourier transforms we presented. These Fourier modes show key properties of the graph Fourier transforms. First of all, there is no difference between the two first rows (graph Fourier transform based on the Laplacian matrix or on the normalized Laplacian matrix) since L = 2L . The factor 2 is observed solely on the eigenvalues.

Figure 3 . 1 .

 31 Figure 3.1. Several Fourier modes of the cyclic graph C 30 with N = 30 vertices. Each row corresponds to a given graph Fourier transform, with the first row (a)-(d) using the Laplacian matrix L, the second row (e)-(h) using the normalized Laplacian matrix L , and the third row (i)-(l) using the adjacency matrix A. The corresponding eigenvalue is shown in the legend of the subfigures. The first column is the DC Fourier mode (lowest frequency).The second and third columns correspond to 2 modes with equal graph frequencies. Finally, the fourth column corresponds to a high graph frequency.

2 Figure 3 . 2 .

 232 Figure 3.2. Comparison of three Fourier modes of the cyclic graph using unit weights on the left and vertex distance in the 1D line on the right. Vertices are closer to each other at the beginning and the end of the cycle and farther from each other in the middle.

14 λ 2 Figure 3 . 3 .

 233 Figure 3.3. Comparison of three Fourier modes of the cyclic graph using a Gaussian kernel of the distance between vertices as weights on the left and the same weights on a complete graph on the right.

---Figure 3 . 4 .

 34 Figure3.4. Six Fourier modes of the torus graph T 10 associated to the same eigenvalue. Each row corresponds to Fourier modes obtained using (i) the Laplacian matrix, (ii) the normalized Laplacian matrix, (iii) the adjacency matrix, and (iv) the fact that T 10 is a Cartesian product of two cyclic graph C 10 .

γ 2 Figure 3 . 5 .

 235 Figure 3.5. Three Fourier modes of the complete Molène graph and the knn Molène graph. The first row uses the Laplacian matrix, the second the normalized Laplacian matrix and the third the adjacency matrix.

γ 2 Figure 3 . 6 .

 236 Figure 3.6. Three Fourier modes of the Molène graph thresholded at 100 kilometers (W th(100 km) ) and 50 kilometers (W th(50 km) ). The first row uses the Laplacian matrix, the second the normalized Laplacian matrix and the third the adjacency matrix.

γ 2 Figure 3 . 7 .

 237 Figure 3.7. Variations of three Fourier modes of the MOSAR graph between two aggregations over two arbitrary days of the study. The weights on the edges are the total contact duration between two vertices (people).

MFigure 3 . 8 .

 38 Figure 3.8. Variations of three Fourier modes of the MOSAR graph between two aggregations over a day and a week of the study. The weights on the edges are the total contact duration between two vertices (people).

γ 2 Figure 3 . 9 .

 239 Figure 3.9. Variations of three Fourier modes of the MOSAR graph between two aggregations over the same day of the study. The weights on the edges are the total contact duration between two vertices (people) and the number of individual contacts.

  Figure 4.1(a) shows a time series whose amplitudes are subdivided into 5 quantiles. On this example, the quantiles have the same size because each sample is uniformly drawn in[-1, 1]. In general this is not the case. The vertices obtained for the example are shown on Figure4.1(b).

Figure 4 . 1 .

 41 Figure 4.1. Illustration of the dual mapping presented in[START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. Given the input time series (a), the amplitudes are split into quantiles 1 through 5. The statistics are then computed to obtain the graph (b) (the loops are omitted). Finally, using a mapping from vertices to amplitudes, a time series (c) can be generated through a random walk on the graph.

  W : generated time series.

Figure 4 . 2 .

 42 Figure 4.2. Examples of resulting mapping of vertices to amplitudes and generated time series for the torus graph

Figure 4 . 3 .

 43 Figure 4.3. Illustration of three optimal mappings from vertices to amplitudes minimizing the cost function of Equation4.1. Color scale indicate the value of a(i ) of the mappings for each vertex i . The 1 quantity corresponds to the cost function, and is equal in all three cases. The first two solutions are hand made while the third is obtained through the simulated annealing method of[START_REF] Campanharo | Duality between Time Series and Networks[END_REF].

Figure 4 . 4 .

 44 Figure 4.4. (a) Application of our approach to mapping vertices to amplitudes using GSSL, with parameters α = 200 and σ = 0.5. The values of 1 and 2 are written in bold to highlight better results compared to Figure 4.2(c). The color scale illustrate the values of a(i ) according to Equation 4.6. The two centroids of the classes are circled in blue to the south east (class 1) and red to the north west (class 2). (b) A corresponding time series is shown on the right.

Figure 4 . 5 .

 45 Figure 4.5. Amplitudes computed using GSSL and Equation 4.6 on the directed cyclic graph C 30 . Several values of the regularization parameter α are used. The unifying parameter σ plays no role here since the graph is regular. Examples of realizations of the dual time series are shown on the bottom.

Figure 4 . 6 .

 46 Figure 4.6. Same as Figure 4.5 for the symmetric cyclic graph C 30 .

  α = 0.1.

  α = 10.

  α = 100.

Figure 4 . 7 .

 47 Figure 4.7. Same as Figure 4.5 for the torus graph T 10 .

Figure 4 . 8 .

 48 Figure 4.8. Same as Figure 4.5 for the Molène graph W . Note that this time the unifying parameter σ influences the results since the graph is not regular. The figure is organized as follows. First row: σ = 0 (Page-Rank), second row: σ = 0.5 (normalized Laplacian matrix), third row: σ = 1 (Laplacian matrix).

Figure 4 . 9 .

 49 Figure 4.9. Same as Figure4.8 for the MOSAR graph M (sum) (t , t + 1 week). Finding the farthest vertices is done with respect to the pseudo-distances 1 -e -w i j . This ensures a better scaling of distances since the contacts total duration can vary greatly from one edge to another. This figure does not include the values of 1 and 2 since we do not have a reference with the approach of[START_REF] Campanharo | Duality between Time Series and Networks[END_REF] due to the computational complexity of the approach and our results are not completely satisfactory (see section4.2.4).

Figure 4 .

 4 Figure 4.10. (Red) Frequency response of the GSSL operator according to Equation4.7 on three different graphs (scale on the left). The top curve is for α = 0.1, the middle one for α = 10, and the bottom one for α = 100. (Blue) Moduli of the Fourier components of the graph signal obtained using the approach of[START_REF] Campanharo | Duality between Time Series and Networks[END_REF], with respect to the graph Fourier transform based on the normalized Laplacian matrix L (scale on the right). Vertical dotted lines identify to reduced graph frequencies of the graphs.

Figure 4 .

 4 Figure 4.11. Amplitudes of Figure 4.3(c) taken as a graph signal on the cyclic graph C 30 and given as input to GSSL. Three different values of α are considered. Similarly to Figure 4.5, σ plays no role since the graph is regular. (Top) new amplitudes. (Bottom) dual time series.

  α = 0.1.

  α = 10.

  α = 100.

Figure 4 .

 4 Figure 4.12. Same as Figure 4.11 for the torus graph T 10 (input on Figure 4.2(a)).

Figure 4 . 13 .

 413 Figure 4.13. Same as Figure 4.11 for the Molène graph W (input on Figure 4.2(c)). The unifying parameter σ influences the results since the graph is not regular. First row: σ = 0 (Page-Rank). Second row: σ = 0.5 (normalized Laplacian matrix). Third row: σ = 1 (Laplacian matrix).

Property 4 .

 4 Let F be a unitary graph Fourier transform. Let H be a convolutive (item (ii)) and isometric (item (iii)) operator. H is then of the form: H = exp(ıΩ),

Figure 5 . 1 .

 51 Figure 5.1. Translations of the cyclic graph C 30 . Four operators are applied to three different signals. The operators are the generalized translation T i , the graph shift A (with (γ 0 = 2), and the graph translations T G (graph Fourier transform based on L) and T G (graph Fourier transform based on L ). The three signals are the delta signal δ 1 localized on vertex 1, the heat kernel g 3 , and the constant signal 1. The figure is split into three blocks: (top) one iteration, (middle) ten iterations (or a different i for the generalized translation), (bottom) several iterations of the graph translation T G .

Figure 5 . 2 .

 52 Figure 5.2. Same as Figure 5.1 for the torus graph T 10 (γ 0 = 4).

Figure 5 . 3 .

 53 Figure 5.3. Same as Figure 5.1 for the Molène graph W (γ 0 = 2.96).

Figure 5 . 4 .

 54 Figure 5.4. Same as Figure 5.1 for the MOSAR graph M .

  .1 through Figure5.4, the three signals we described are translated using a generalized translation based on L, the graph shift, the graph translation based on L, and the graph translation based on L . The modulus is shown. The reason behind the use of the modulus is the complex nature of the graph translation, and the phases of the translated signals are not particularly interesting and as such are omitted from figures that are already quite dense.

6 Figure 6 . 1 .

 661 Figure 6.1. PSD of the two graph signals considered to illustrate stationarity and non-stationarity. The low band graph signal x 1 is represented in red, and the bandpass x 2 in blue.

Figure 6 . 2 .

 62 Figure 6.2. Stochastic signal analysis on the cyclic graph C 30 .The figure is vertically split in two, with x 1 on the left, and x 2 on the right, and then split into three subcolumns, one per graph Fourier transform. The three first rows comprise realizations of the graph signal, followed by the autocorrelation matrix R and the spectral autocorrelation matrix S. The next two rows show the spectral correlation of the signal x i + aδ 1 with a ∼ N (0, 1) (first) and the signal x i + a(δ 1 + δ 10 ) with a ∼ N (0, 1) (second). Finally, the two last rows show the impact of a "wrong" graph (Erdős-Rényi graph with p = 0.15) on the spectral autocorrelation matrix, with an illustration of the structure on the last row.

Figure 6 . 3 .

 63 Figure 6.3. Same as Figure 6.2 for the torus graph T 10 , and a random sampling of the plane as "wrong" graph (all pair of vertices are connected such that we have not shown the edges to simplify the last row).

Figure 6 . 4 .

 64 Figure 6.4. Same as Figure 6.2 for the Molène graph W , and an arbitrary cyclic graph as "wrong" graph.

Figure 6 . 5 .

 65 Figure 6.5. Same as Figure 6.2 for the MOSAR graph M , and an arbitrary cyclic graph as "wrong" graph. To reduce the number of elements, each element of the matrices shown is the mean of a 2x2 block of the original matrix.

  Figure 6.6.Weather readings over time of three ground stations of (a) the temperature, (b) the wind speed, and (c) the rain. The first row shows the weather stations considered on the Molène graph. Note that the vertex set changes between readings since not all stations have readings over the whole period. The second row shows the original data, the third row shows IMF(2) , and the fourth row IMF(1) .

Figure 6 . 7 .

 67 Figure 6.7. Study of the stationarity of the temperature readings of the Molène dataset. The top three rows show several snapshots of the raw data (first row), the second IMFs (second row), and the first IMF (third row). The next half of the figure shows the autocorrelation matrices R and correlation coefficient matrices C or (first two columns), spectral autocorrelation matrix S and correlation coefficient matrix C or for both the Laplacian based (third and fourth columns) and the normalized Laplacian based (last two columns) graph Fourier transforms. Four stochastic graph signals are considered for these matrices: The IMF(1) (first row), the IMF(1) with a cyclic graph as underlying graph structure (second row), the PSD model based on the Laplacian (third row), and the PSD model based on the normalized Laplacian (fourth row). See Figure6.11 for these models. Note that for all those matrices, the moduli of their elements is shown. Finally, the color scales on R and S are different to better account for the differences in values of these matrices. This does not hinder the interpretation of the figure since R and S are two different objects that cannot be compared.

Figure 6 . 8 .Figure 6 . 9 .

 6869 Figure 6.8. Same as in Figure 6.7 with the wind readings.

8 Figure 6 . 10 .

 8610 Figure 6.10. Molène graph with the vertices of the wind data. The graph signals shown are two Fourier modes χ l of the Laplacian matrix based graph Fourier transform and corresponding to the circled spectral correlations of Figure 6.8. On the north part of the region, the Fourier modes are very similar, while the variations on the south coast line are slightly different.

CHAPTER 7 Contents 7 . 1

 771 State of the Art .

Figure 7 . 1 .

 71 Figure 7.1. GUI from our MATLAB toolbox showing the iterations of an operator. Here, the graph translation operator has been iterated on a heat kernel of a torus graph.

  the standard deviation. 92 Mean Square limit (MS limit)y are equal x in the MS limit sense if E |x -y| 2 = 0. 93 mean (first order moment) E[x(t )] = R x f x (x, t )d x.92 power spectrum density Fourier transform of the autocorrelation function. 93 power spectrum distribution Integral of the power spectrum density. 93 probability density function Differential of the probability distribution: f (k)x = ∂F(k)

Part I A Journey from Linear Algebra to Signal Processing on Graphs

  Contents 1.1 Fundamental Definitions and Properties . . . . . . . . . . . . . . . . . . . . . . . . . 7 1.2 Matrix Decompositions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 1.2.1 QR Decomposition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 1.2.2 Jordan Decomposition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 1.2.3 Schur Decomposition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 1.2.4 Cholesky Decomposition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12

	CHAPTER 1
	Matrices, Vectors and Matrix
	Decompositions

  2

Chapter 1. Matrices, Vectors and Matrix Decompositions U

  • • • u n ] with arbitrary orthonormal vectors u i . Since the columns of 1 are orthonormal, it follows that:

Table 2 . 1 .

 21 Properties of several matrix representations of graphs. In the following table, the largest magnitude of the vertex degrees is d max = max i |d i |. The symbols and represent a verified and a non-verified properties. On the other hand, the symbol ? represents a missing result.

				diagonalizable		eigenvalues
					real non-negative	bounds
	directed graph	w i j ∈ R i j ≥ 0 w	A L L A L L	? ? ? ?	? ? ? ?	? ? ? ?	? ? ? ?
	symmetric graph	i j ∈ R w	A L L	?	?	3 ?	? ? ?

Table 2 .2. Eigenpairs

 2 of matrix representations of several symmetric graphs. We use the angular frequencies ω = 2πν instead of the classical frequency ν to keep the width of the table to a minimum. Note also that whereas the graphical representation of a cyclic graph is clear, that of a torus graph is actually indistinguishable from the representation of a grid graph since edges going back from one end to the other (horizontally or vertically) overlay the other edges. Nevertheless, they are indeed present in the graph.

		A	L	L
	C N	Eigvect.

  Contents4.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 4.1.1 Graph to Time Series Duality . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 4.1.2 Drawbacks of [Campanharo et al., 2011] . . . . . . . . . . . . . . . . . . . . . 57 4.1.3 Graph Based Semi-Supervised Learning . . . . . . . . . . . . . . . . . . . . . 58 4.1.4 GSSL Rewriting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 4.2 Graph to Time Series Duality Using GSSL . . . . . . . . . . . . . . . . . . . . . . . . . 61 4.2.1 Achieving Smoothness using GSSL with 2 Classes . . . . . . . . . . . . . . . . 62 4.2.2 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 4.2.3 Example Graphs and Dual Time Series . . . . . . . . . . . . . . . . . . . . . . 64 4.2.4 Limits of the GSSL Approach with 2 Classes . . . . . . . . . . . . . . . . . . . 68 4.3 Beyond the Classification Approach .. . . . . . . . . . . . . . . . . . . . . . . . . . . 68

	CHAPTER 4

4.3.1 Graph Signal Operator Interpretation: Wiener Filter . . . . . . . . . . . . . . 68 4.3.2 Application . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70

  Contents6.1 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 6.1.1 Temporal Signal Statistics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 6.1.2 Stationary Temporal Signals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 6.1.3 Discrete and Periodic Time Series . . . . . . . . . . . . . . . . . . . . . . . . . 94 6.1.4 Mean Square Estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95 6.2 Stationary Graph Signal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 6.2.1 Graph Signal Statistics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 6.2.2 Definitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 6.2.3 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 6.3 Sources of Non-Stationarity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 6.4 Experimental Studies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 6.4.1 Toy Signals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 6.4.2 Weather Dataset . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 Until this chapter, we were interested in studying deterministic graph signals, i.e. graph signals whose samples are certain, not subject to uncertainty. We now aim at the study of stochastic graph signals subject to probabilities. We denote a stochastic signal, being temporal or graphical, using a bold character x. We propose in this chapter a definition of stationary graph signals by analogy to the classical case. In the temporal setting, stationarity is a fundamental definition giving a powerful framework to study stochastic signals. In particular, understanding of non-stationarity is crucial to study local phenomena. After summarizing classical results and introducing our definition of stationary graph signals, we explain how non-stationarity can be explained and illustrate our results on both toy graph signals and real-world examples.

	CHAPTER 6

  We saw that temporal stochastic signals are stationary if and only if they are statistically invariant through the time shift. We present in this section a novel definition of stationary graph signals we proposed in[2]. For that purpose, we define stationary graph signals as statistically invariant through the graph translation operator.

	Definition 15 (SSS Graph Signal). x is a Strict Sense Stationary (SSS) graph signal if and only if
	its probability distribution functions of any order are invariant through graph translation, i.e.
	we have for any τ:

  -ı2πνt d t , or x[n]e -ı2πνn for a discrete time series. 33 frequency response For a convolutive operator H , function h(ν) such that H {e ν } = h(ν)e ν . 35 impulse response For a convolutive operator H , function h(t ) such that H {δ} = h. 35 Operator computing a signal from its frequency components and equal to x(ν)e ı2πνt d ν or x(ν)e ı2πνn d ν for a discrete time series. 34 kernel For a linear operator H , function h(t , s) such that H {x}(t ) = x(s)h(t , s)d s. 35 Generalized convolution of a graph signal with a delta signal localized on a vertex. 77 generalized convolution Operation between two graph signals defined as x y(l ) = x(l ) y(l ). 38 Adjacency matrix seen as the equivalent of the time shift operator for . 43 graph frequency Real number ν associated to a Fourier mode that should reflect how fast is varying the Fourier mode. 37 Function mapping the vertex set V to scalar values. 37 Operator computing the frequency components of a graph signal using the graph Fourier matrix: x = F x. 35 heat kernel Signal g κ such that g κ (l ) = e -κν 2 l . 82 Rescaling of all graph frequencies such that lie in the interval [0, 2π]. 41 stationary graph signal Stochastic graph signal which is statistically invariant through graph translation (SSS) or whose first two statistical moments are invariant through graph translation (WSS). 97 (t 2 )] = R 2 x 1 x * 2 f x (x 1 , x 2 ; t 1 , t 2 )d x 1 d x 2 . 92 autocovariance C x (t 1 , t 2 ) = R x (t 1 , t 2 )η x (t 1 )η * x (t 2 ). 92 average power (second order moment) E[|x(t )| 2 ] = R |x| 2 f x (x, t )d x. 92

	Statistics
	graph signal
	graph Fourier transform
	inverse Fourier transform reduced graph frequency
	operator
	Operator transforming a signal into a new signal. 35 Statistics
	power spectrum autocorrelation
	Squared modulus of the Fourier transform. 82 E[x(t 1 )x
	reduced frequency
	Frequencies ν = ν τ Graph Signal Processing
	generalized translation
	graph translation
	Operator acting on graph signal whose action is defined such as to ensure similar
	behavior than the time shift. 80
	graph shift

s for a discrete time series with sampling period τ s . 35 signal energy Integral of the squared modulus of the signal |x(t )| 2 d t , or |x[n]| 2 for a discrete time series. 34 stationary signal Stochastic signal which is statistically invariant through time shift (SSS) or whose first two statistical moments are invariant through time shift (WSS). 91, 92 time shift

Operator T acting on a signal x(t ) such that T {x}(t ) = x(t -1). 75 *

p A is a polynomial of degree N , and as such has exactly N roots in C: the eigenvalues.

Remember that G has no loop, and as such ∀i , w i i = 0.

Given a graph non-negatively weighted, inverting the sign of all weights is equivalent to inverting the sign of all entries of L, hence inverting the sign of the eigenvalues of L. Therefore, there are symmetric graphs whose Laplacian matrix has negative eigenvalues.

Numbering of the eigenvalues and eigenvectors starts at 0 instead of 1 because of their interpretation in terms of Fourier modes and frequencies that we explicit in chapter 3. In short, the index 0 represents the lowest possible frequency on the graph, same as the Direct Current (Frequency 0) (DC) component of a temporal signal is of frequency 0.

In http://data.gouv.fr.

In the context of N -periodic time series, linear convolutive operators where h[n, k] = h[nk] depend on the sample difference nk (without the modulo) are not interesting in this dissertation since they do not verify the convolution theorem.

This is a special kind of Toeplitz matrix. Toeplitz matrices are defined as matrices with constant diagonals and correspond to linear convolutive operator.

These two eigenvectors are orthonormal if F l and F l +1 are.

Function bdschur of MATLAB.

Our implementation is part of the toolbox presented in chapter7.

When the graph is unweighted this is the Minimum Linear Arrangement Problem which is a NP-hard problem[START_REF] Arora | A new rounding procedure for the assignment problem with applications to dense graph arrangement problems[END_REF]. Therefore, the weighted version is at least as difficult, such that it is also NP-hard.

Note that F(1) = D -1 2 Ψ = χ, i.e. F(1) is the not graph Fourier transform associated to the Laplacian matrix, even though σ = 1 corresponds to the Laplacian-based GSSL approach.

Intuitively, in the case of C 30 , 1 is equal to twice a max -a min (a is increasing on one half of the period and decreasing on the other), divided by a max -a min , hence constant.

The spread of a graph signal is often defined as ∆

G = min u 0 u d(u, u 0 ) 2 |X u | 2 , with d(u, u 0 ) a given distance

A cycle of the permutation can be thought of as a connected component of the graph whose edges are the mapping i → k i of the permutation.

In[Grady and Polimeni, 2010, Section 5.1], a shift-invariant graph is a graph whose graph Fourier transform based on the Laplacian matrix is the DFT. They show also that a graph is shift invariant if and only if there exists a permutation of the vertex indices such that its Laplacian matrix is circulant. We can show the same equivalence on the adjacency matrix. Therefore, using the adjacency matrix-based or the Laplacian matrix-based graph Fourier transform, only the shift invariant graphs can verify the four time shift properties we presented.

For a graph which is not bipartite regular, we have λ N -1 < ρ G , and ν N -1 < 1

, such that the condition can be extended to 0 < l ≤ N -1.

Due to the very small number of frequencies, many other model may fit very well. We are satisfied with this one as it is common in classical signal processing.

Available here: https://gforge.inria.fr/projects/grasp/.

Acknowledgements

Summary

In this chapter, we introduced basic definitions of matrices and vectors that are essential to this dissertation. The most important elements of this chapter are:

• An eigenvector x of A associated to the eigenvalue λ verifies Ax = λx.

• An Hermitian matrix A is unitarily diagonalizable, and we have A = U ΛU * with Λ diagonal and U unitary.

Part II

New Methods for Signal Processing on Graphs

Summary

In this chapter, we used Graph-Based Semi-Supervised Learning (GSSL) to map vertices of a graph to amplitudes of a time series. We saw that this approach achieves better smoothness than the previous approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF]. In the course, we interpreted GSSL as a graph filter with similarities with respect to its frequency response to a Wiener filter (Equation 4.7):

We then use these models to generate realizations of the corresponding WSS graph signals (with zero mean, like the first IMF). Doing so, we model the data with a stationary graph signal. The third and fourth rows of Figure 6.7 and Figure 6.8 show then the results using 744 realizations (same number as in the data). We observe in particular that these models are close to the real data in terms of autocorrelation matrix and spectral autocorrelation matrix. Indeed, even if the correlation coefficient matrix C or is not perfectly recovered, we can see several (spatial) correlations explained by the models. In particular, when looking at the correlation coefficient matrix of the wind, several high correlations between vertices are recovered.

Finally, we notice that the variance on the vertices (shown on the diagonal of R) is not explained by these models. One of the reason could be the fact that the sensors are not exactly the same, and as such do not show the same variance on the readings. As such, non-stationarity can be the image of a different sensor. Future works definitely include the detection and correction of abnormal readings with respect to the graph structure. Wind. L: (1.70,1.09,16.24,0.135), L : (1.86, 1.07, 14.48, 0.164).

Summary

In this chapter, we introduced two definitions of Strict Sense Stationary (SSS) and Wide (or Weak) Sense Stationary (WSS) graph signals using our graph translation operator. Using the WSS definition we showed that WSS graph signals are characterized in the spectral domain with a mean equal to a DC component, and uncorrelated spectral components, thus having the same spectral characterization as WSS temporal signals.

We showed also that the non-stationarity stems from the graph structure or the graph signal itself, or both. This has been evident on synthetic graph signals. Finally, we studied a real dataset of weather reports and showed that they contained stationarity in some sense.

Conclusion

Summary

In this dissertation, after quickly reminding in chapter 1 the essential results from linear algebra on matrices and vectors and in particular matrix decompositions, we reviewed important results from spectral graph theory in chapter 2, among which the diagonalizations of the adjacency matrix, the Laplacian matrix and the normalized Laplacian matrix. We presented several regular graphs on which we computed the diagonalizations, several graph models, and finally two real graphs derived from two distinct datasets.

Chapter 3 takes advantage of the matrix representations and decompositions of graphs to introduce three graph Fourier transforms and the analogies to the classical case that justify them. We also present the convolution of two graph signals and its associated convolution theorem. A notion of graph frequencies is introduced for latter use in chapter 5. We close this chapter with Fourier modes of several graphs: the weighted and unweighted cyclic graphs corresponding to the periodic time series and showing the influence of the weights on the resulting Fourier modes, the torus graph corresponding to the pixel graph of an image and showing the impact of multiple Fourier modes with the same frequency, the Molène graph corresponding to an irregular sampling of a 2D space and showing the influence of edge removal schemes, and finally the MOSAR graph corresponding to a dynamic graph and showing the impact of aggregating snapshots into a static graph.

Part II focuses then on four contributions to the field of signal processing on graphs. In chapter 4, we developed a graph to time series mapping embedding the structural properties of a graph into a temporal signal. We started with the approach of [START_REF] Campanharo | Duality between Time Series and Networks[END_REF] and improved upon it by increasing the regularity of the vertex to amplitude mapping. We did so by leveraging the strength of the closed form solution to the Graph-Based Semi-Supervised Learning (GSSL) method. We obtained smoother results on several graphs. We then observed that this closed form solution is the application of a convolutive operator to a graph signal. Using the results of chapter 6, we showed that this filter is a Wiener filter for graph signals, hence a smoothing operator. With this understanding of GSSL, we applied it to a few nonsmooth vertex to amplitude mapping and observed better smoothness in the filtered mapping.

In chapter 5, we first observed that the operators from the literature that have been defined by similarity to the time shift operator for temporal signal, lack a crucial property to obtain a tractable definition of stationarity: the isometry. We introduced then a novel convolutive operator, the graph translation enforcing the isometry, hence being a Power Spectrum Density Linear Algebra

List of Figures Glossaries

-norm

Norm of a vector defined as

Linear Algebra block diagonal matrix

A matrix that can be written as a diagonal of submatrices. 8

canonical scalar product

Operation between two vectors defined as 〈x, y〉 = x y * = i x i y * i . 8 canonical basis Set of vectors {e i } i verifying e i (i ) = 1 and ∀ j = i , e i ( j ) = 0 and that constitute a basis of the vector space C n . 8

characteristic polynomial

Determinant of the matrix x I -A: p A (x) = det(x I -A). 9

Cholesky decomposition

Computing the matrix decomposition A = LL * for a Hermitian matrix A with L a lower triangular matrix. 12 definite positive matrix A matrix verifying ∀x, x , x * Ax > 0. 9

diagonalizable matrix

A matrix A that can be written as A = V DV -1 with D a diagonal matrix. 9

diagonal matrix

A matrix whose only non-zero elements lie on its diagonal. 8 eigenpair An eigenvalue λ and an eigenvector x verifying Ax = λx. 9 eigenvalue A real or complex number λ verifying Ax = λx with x an associated eigenvector. 9 eigenvector A non-zero vector x verifying Ax = λx with λ its associated eigenvalue. 9 generalized eigenvector A non-zero vector x verifying Ax = λx + y with λ its associated eigenvalue and y another generalized eigenvector associated to the eigenvalue λ. 11

Hermitian matrix

A matrix equal to its transpose conjugate matrix. 8

identity matrix

A diagonal matrix with only ones on its diagonal. 8

invertible matrix

A matrix A such that there exists a matrix B verifying AB = I . B is denoted A -1 . 8

Jordan chain

A set of generalized eigenvector {x i } i associated to an eigenvalue λ of a matrix A and such that Ax i = λx i + x i -1 . 11

Jordan block

A matrix whose diagonal is constant and its upper diagonal equal to 1. 10

Jordan decomposition

Computing a similar matrix which is block diagonal matrix with Jordan block. 10

Kronecker matrix product

Matrix operation between two matrices A and B such that A ⊗ B can be written as a block matrix of block identical to B and with the block i j equal to A i j B . 10 multiplicity Multiplicity of eigenvalue as the zero of the characteristic polynomial. 9

orthogonal vectors

Two vectors u and v verifying 〈u, v〉 = i u i v * i = 0, or a set of vectors orthogonal to each other. 8

orthonormal vectors

Orthogonal vectors of unitary norm. 8

QR decomposition

Computing the decomposition A = QR with Q a unitary matrix and R an upper triangular matrix. 10

Rayleigh quotient

Quotient x * Ax x * x for an Hermitian matrix A and a vector x, denoted R A (x). 9

Schur decomposition

Computing a triangular unitarily similar matrix T of A: A = U T U * . 11

semi-definite positive matrix

A matrix verifying ∀x, x , x * Ax ≥ 0. 9

similar matrices

Two matrices A and B such that there exists an invertible matrix V verifying A = V BV -1 . 9

spectral decomposition

The set of all eigenpairs of a given matrix. 9

spectrum

The set of all eigenvalues of a given matrix, denoted Sp. 9

symmetric matrix

A matrix equal to its transpose matrix. 8

transpose conjugate operator ( * )

Matrix operator returning the complex conjugate of the transpose matrix. 8

transpose operator ( T )

Matrix operator inverting lines and columns of a matrix such that (A T ) i j = A j i . 8

triangular matrix

A matrix whose only non-zero elements lie on its diagonal or its upper right half. 8

unitarily diagonalizable matrix

A matrix A that can be written as A = U DU * with D a diagonal matrix and U a unitary matrix. 9

Graphs unitarily similar matrices

Two matrices A and B such that there exists a unitary matrix U verifying A = U BU * . 9

unitary matrix

A matrix U verifying UU * = I . 8

Graphs adjacency matrix

Square matrix whose entries are the weights of the edges between vertices, or zero if there is no such edge. 18

Barabási-Albert graph

Graph model whose realization is constructed iteratively beginning with a graph with d vertices and adding a vertex at each step. This vertex is connected to d vertices with probability proportional to the degree. 28

bipartite graph

Graph whose vertex set can be partitioned into two subsets such that all edges lie between those sets. 16

Cartesian product of graphs

Graph whose vertex set is the Cartesian product of the vertex sets, and an edge i , k j , l appears if either i = j and kl is an edge or k = l and i j is an edge. 17

complete graph

Graph whose edges link all pair of vertices. 24

connected component

A connected subgraph. 16

connected graph

Graph having a path between any two vertices. 16

cyclic graph

Graph whose edges link the vertices one after the other with the last vertex connected to the first to form a loop. 21, 24

directed edge

directed (or non symmetric) graph

A graph having at least an oriented edge, i.e. an edge i j such that j i is not an edge. 16

(non-oriented) edge

Relation between two vertices in a graph. v

edge weight

A (real) value associated to an edge. 16

Erdős-Rényi graph

Graph model whose realization has edges drawn uniformly and independently with a fixed probability. 28

Graphs graph

A set of vertices with edges between them. 15

graph transformation

Function of one or several graphs whose output is a new graph obtained by union of the vertices and edges with weights function of the input weights. 17

grid graph

Graph whose vertices are arranged in a square grid with edges between horizontally or vertically adjacent vertices. 24

in-degree matrix

Diagonal matrix with the in-degree of vertex i on its i th element. 18

Kronecker product of graphs

Graph whose vertex set is the Cartesian product of the vertex sets, and an edge i , k j , l appears if i j and kl are edges in the original graphs. 17

Laplacian matrix

Matrix representation of a graph defined as the difference between the degree matrix and the adjacency matrix. 19 loop An edge linking a vertex to itself. 15

normalized Laplacian matrix

Matrix representation of a graph defined as D -1/2 LD -1/2 with D the degree matrix and L the Laplacian matrix. 20

(out-)degree matrix

Diagonal matrix with the degree of vertex i on its i th element. 18

path

A succession of edges such that the target of an edge is the source of the next. 16

path graph

Graph whose edges link the vertices one after the other, thus constituting a single path. 25

path length

The number of edges in a path.