Subgraph-based filterbanks for graph signals

I Introduction

Graphs are a modeling tool suitable to many applications involving networks, may they be social, neuronal, or driven from computer science, molecular biology … Data on these graphs may be defined as a scalar (or vector) on each of its nodes, forming a so-called graph signal . In a sense, a graph signal is the extension of the 1-D discrete classical signal (where the signal is defined on the circular graph, each node having exactly two neighbors) to any arbitrary discrete topology where each node may have an arbitrary number of neighbors. Temperature measured by a sensor network, age of the individuals in a social network, Internet traffic in a router network, etc. are all examples of such graph signals.

Adapting classical signal processing tools to signals defined on graphs has raised significant interests in the last few years . For instance, the graph Fourier transform, the fundamental building block of signal processing, has several possible definitions, either based on the diagonalisation of one of the Laplacian matrices , or based on Jordan’s decomposition of the adjacency matrix . Building upon this graph Fourier transform, authors have defined different sampling and interpolation procedures , windowed Fourier transform , graph empirical mode decomposition , different wavelet transforms, including spectral graph wavelets , diffusion wavelets , and wavelets defined via filterbanks . Among the applications of graph signal processing, one may cite works on fMRI data , on multiscale community detection , image compression , etc. In fact, graph signal processing tools are general enough to deal with many types of irregular data .

Graph filterbanks using downsampling have been initially defined for bipartite graphs because: i) Bipartite graphs, by definition, contain two sets of nodes that are natural candidates for the sampling operations; ii) Downsampling followed by upsampling (which forces to zero the signal on one of the two sets of nodes) can be exactly written as a filter in the graph Fourier space. This enables to write exact anti-aliasing equations for the low-pass and high-pass filters to cancel the spectral folding phenomenon due to sampling . However, for arbitrary graphs, one needs to decompose the graph in a (non-unique) sum of bipartite graphs , and analyze each of them separately. Another solution for arbitrary graph is based on downsampling according to the polarity of the graph Fourier mode of highest frequency .

We propose a significantly different way of defining filterbanks. Instead of trying to find an exact equivalent of both the decimation operator, hereafter $\bm{(\downarrow)}$, and a filtering operator $\bm{C}$, we directly define a decimated filtering operator $\bm{L}=\bm{(\downarrow)}\bm{C}$; following here the notations of , where $\bm{L}$ is not to be confused with the Laplacian operator of the graph, noted $\bm{\mathcal{L}}$. Consider the 1-D straight-line graph where each node has two neighbors, and a partition of this graph in subgraphs of pairs of adjacent nodes. The classical Haar low-pass (resp. high-pass) channel samples one node per subgraph and defines on it the local average (resp. difference) of the signal. By analogy, we consider a partition of the graph in connected subgraphs, not necessarily of same size. Creating one “supernode” per subgraph, the low-pass channel (resp. high-pass channels) defines on it the local, i.e., over the subgraph, average (resp. differences) of the signal. The coarsened graphs on which those downsampled signals are defined, are then derived from the connectivity between the subgraphs: two supernodes are linked if there are edges between the associated subgraphs.

With this approach, we design a critically-sampled, compact-support biorthogonal filterbank that is valid for any partition in connected subgraphs. Depending on the application at hand, one has the choice on how to detect such partitions. For compression and denoising, an adequate way is to use a partition in communities, i.e. groups of nodes more connected with themselves than with the rest of the network . This community structure is indeed linked to the low frequencies of graph signals . For hierarchical clustering trees, multiresolution bases on graphs have been explored in . As a difference here, we not only take into account a hierarchical clustering in groups, but also the local intra-cluster topology in each group when defining the analysis atoms.

Section II recalls the definition of the graph Fourier transform we use. In Section III, after detailing the difficulties to extend classical filterbanks to graph signals, we discuss the state-of-the-art of graph filterbanks. The main contribution is in Section IV, first discussed as an analogy to the Haar filterbank, before presenting fully the proposed graph filterbank design. Section V proposes how to obtain a relevant partition in connected subgraph. Section VI shows applications, in compression and denoising. We conclude in Section VII.

II The Graph Fourier Transform

Let $\mathcal{G}=(\mathcal{V},\mathcal{E},\mathbf{A})$ be a undirected weighted graph with $\mathcal{V}$ the set of nodes, $\mathcal{E}$ the set of edges, and $\mathbf{A}$ the weighted adjacency matrix such that $\mathbf{A}_{ij}=\mathbf{A}_{ji}\geq 0$ is the weight of the edge between nodes $i$ and $j$. Let $N$ be the total number of nodes. Let us define the graph’s Laplacian matrix $\bm{\mathcal{L}}=\mathbf{D}-\mathbf{A}$ where $\mathbf{D}$ is a diagonal matrix with $\mathbf{D}_{ii}=\mathbf{d}_{i}=\sum_{j\neq i}\mathbf{A}_{ij}$ the strength of node $i$. $\bm{\mathcal{L}}$ is real symmetric, therefore diagonalizable: its spectrum is composed of $\left(\lambda_{l}\right)_{l=1\dots N}$ its set of eigenvalues that we sort: $0=\lambda_{1}\leq\lambda_{2}\leq\lambda_{3}\leq\dots\leq\lambda_{N}$; and of $\mathbf{\bm{Q}}$ the matrix of its normalized eigenvectors: $\bm{Q}=\left(\bm{q}_{1}|\bm{q}_{2}|\dots|\bm{q}_{N}\right)$. Considering only connected graphs, the multiplicity of eigenvalue $\lambda_{1}=0$ is 1 . By analogy to the continuous Laplacian operator whose eigenfunctions are the continuous Fourier modes and eigenvalues their squared frequencies, $\bm{Q}$ is considered as the matrix of the graph’s Fourier modes, and $\left(\sqrt{\lambda_{l}}\right)_{l=1\dots N}$ its set of associated “frequencies” . For instance, the graph Fourier transform $\bm{\hat{x}}$ of a signal $\bm{x}$ defined on the nodes of the graph reads: $\bm{\hat{x}}=\bm{Q}^{\top}\bm{x}$.

III State of the art

In classical setting, the decimation $\bm{(\downarrow 2)}$ operator by 2 is paramount. It keeps one every two nodes and follows what we call the “one every two nodes paradigm”, as seen on Fig. 1a). The classical design of a filterbank is to find a set of operators, e.g. a low-pass filter $\bm{C}$ and a high-pass filter $\bm{D}$, that combine well with decimation such that perfect recovery is possible from the decimated low- and high-pass filtered signals .

The usual Haar filterbank will be used as a leading example to expose the main issues encountered when attempting to extend filterbanks to graph signals. Consider the 1-D discrete signal $\bm{x}$ of size $N$. Let us recall that, at the first level of the classical Haar filterbank, $\bm{x}$ is decomposed into : - its approximation $\bm{x}_{1}$ of size $N/2$: $\bm{x}_{1}=\bm{(\downarrow 2)}\bm{C}\bm{x},$ where $\bm{C}$ is the sliding average operator, here in matrix form:

- its detail $\bm{x}_{2}$ of size $N/2$: $\bm{x}_{2}=\bm{(\downarrow 2)}\bm{D}\bm{x},$ where $\bm{D}$ is the sliding difference operator:

This Haar filterbank is orthogonal and critically sampled. Our objective is to generalize it to signals on arbitrary graphs.

III-B Adapting filterbanks to graph signals

For graph signals, a key difficulty is the design of a suitable decimation operator, and it comes in two separate problems:

how to wire together the nodes that are kept, so as to create the downsampled graph?

On a straight line or a regular grid, issue i) is solved by the one every two nodes paradigm and issue ii) does not exist: the structure after downsampling is exactly the same as the original (straight line or 2D grid); see Figs. 1 a) and b).

To tackle these issues, the following works propose a way to adapt the one every two nodes paradigm to arbitrary graphs.

Narang and Ortega consider first signals defined on bipartite graphs (i.e. two-colourable graphs). In this particular case, one may still downsample the graph by naturally keeping one every two nodes, as shown in Fig. 1c). For non-bipartite graphs, they develop a preprocessing step of the structure where the graph is decomposed in a sum of bipartite subgraphs, on which the filterbank is successively applied. Other methods to create bipartite graphs have been proposed, by oversampling , or with maximum spanning tree . For issue ii), the downsampled structure has edges between two nodes if they have at least a common neighbor in the initial graph. As seen in Fig. 1c), a downsampled bipartite graph is not necessarily bipartite and the preprocessing of the structure is mandatory to iterate the filterbank design. Still, an interesting property of this design is the specific behavior of bipartite graphs Laplacian’s eigenvalues, that enables the authors to write specific anti-aliasing equations for the design of the filters $\bm{C}$ and $\bm{D}$ (see Eq. (13) of ).

III-B2 A downsampling based on the Laplacian’s last eigenvector

Shuman et al. focus on the eigenvector associated to the Laplacian’s largest eigenvalue. They create two sets of nodes depending on this eigenvector’s sign, as illustrated in Fig. 1d). According to the Fourier interpretation of the Laplacian’s eigenbasis, the last eigenvector corresponds to the “highest frequency” of a graph signal. This idea, inspired by graph coloring studies , generalizes the fact that, for structured grids and bipartite graphs, the sign of this eigenvector does alternate every two nodes. To tackle issue ii), the authors in rely on the Kron reduction of the initial graph to obtain new graphs, and post-process them to remove links from otherwise very dense downsampled graphs, implying some degree of arbitrary choices.

In summary, there are many choices (and some stochasticity) in the pre- or post-processing steps to obtain suitable decimated graphs on which the filterbank can be cascaded.

IV Filterbanks on connected subgraphs

The idea we explore lets go of the “one every two nodes” paradigm, and concentrates on graph coarsening: given a partition in connected subgraphs of the initial graph, the approximation and detail(s) will be obtained on each “supernode” that represent each connected subgraph. Hence we will not attempt to define separately analogies of downsampling $\bm{(\downarrow 2)}$ and filtering $\bm{C}$ and $\bm{D}$. Instead, we directly define analogies to graph signals of the decimated sliding average operator $\bm{L}=\bm{(\downarrow 2)}\bm{C}$ and of the decimated sliding difference operator $\bm{B}=\bm{(\downarrow 2)}\bm{D}$. They read for the Haar filterbanks:

In Section IV-A, we first take a close look at the effect of these two Haar operators $\bm{L}$ and $\bm{B}$ on the input signal $\bm{x}$, to give insight in the fundamental analogy that is further formalized in Sections IV-B to IV-D. In Section IV-E, we detail the analysis atoms created by the proposed filterbanks.

Let us rephrase the Haar filterbank from the proposed new point of view of operators on connected subgraphs.

Consider the 1-D classical signal $\bm{x}$ of even size $N$ defined on the straight line graph $\mathcal{G}$, of size $N$, where each node has two neighbors. We consider the partition $\bm{c}$ of this graph in $K=N/2$ connected subgraphs $\left\{\mathcal{G}^{(k)}\right\}_{k\in\{1,K\}}$ connecting neighbors two-by-two: we call it the Haar partition, and it reads (when coded as a vector):

where $\bm{c}(i)$ is the label of node $i$’s subgraph.

IV-A2 Interpret operators 𝑳𝑳\bm{L} and 𝑩𝑩\bm{B} in terms of local Fourier modes

Consider subgraph $\mathcal{G}^{(k)}$ and $\bm{x}^{(k)}$ the restriction of $\bm{x}$ to this subgraph. Define $\mathcal{G}^{(k)}$’s local adjacency matrix $\bm{A^{(k)}}$ :

Its Laplacian matrix is diagonalisable with two local Fourier modes: $\bm{q}_{1}^{(k)\top}=\frac{1}{\sqrt{2}}\left(1,\quad 1\right)$ of associated eigenvalue $\lambda_{1}^{(k)}=0$ and $\bm{q}_{2}^{(k)\top}=\frac{1}{\sqrt{2}}\left(-1,~{}~{}1\right)$ of associated eigenvalue $\lambda_{2}^{(k)}=2$.

The actual effect of the operation $\bm{x}_{1}=\bm{L}\bm{x}$ in Haar filterbank is to assign to each subgraph $\mathcal{G}^{(k)}$ the first local Fourier component of $\bm{x}^{(k)}$ :

Similarly, the actual effect of the operation $\bm{x}_{2}=\bm{B}\bm{x}$ may be rewritten as :

In other words $[\bm{x}_{1}(k)\quad\bm{x}_{2}(k)]^{\top}$ is the local (reduced to $\mathcal{G}^{(k)}$) Fourier transform of $\bm{x}^{(k)}$.

IV-A3 Analogy for graph signals

Consider a graph $\mathcal{G}$ and a partition $\bm{c}$ of this graph in $K$ connected subgraphs $\left\{\mathcal{G}^{(k)}\right\}_{k\in\{1,K\}}$. Consider one of these subgraphs $\mathcal{G}^{(k)}$ of size $N_{k}$. Consider $\hat{\bm{x}}^{(k)}$ the local graph Fourier transform of $\bm{x}^{(k)}$, the graph signal reduced to $\mathcal{G}^{(k)}$. To this end, we diagonalize $\mathcal{G}^{(k)}$’s local Laplacian matrix to find its $N_{k}$ eigenvectors (a.k.a. local Fourier modes) sorted w.r.t. their eigenvalues and compute the successive inner products. We propose the following fundamental analogy: the first coefficient of $\hat{\bm{x}}^{(k)}$ will contribute to the approximation $\bm{x}_{1}$ of the signal, and the following coefficients to its successive details $\bm{x}_{2},\ldots,\bm{x}_{N_{k}}$.

IV-A4 Graph support of the decimated components

For 1-D Haar filterbanks, the two components are defined on straight-line graphs of size $N/2$. For arbitrary graphs, all subgraphs have not necessarily the same size. This implies that only subgraphs of size at least $l$ contribute to $\bm{x}_{l}$. For instance, a three-node subgraph ${\mathcal{G}^{(k)}}$ will have three eigenvectors for its local Laplacian and its first (resp. second, third) eigenvector will contribute to $\bm{x}_{1}$ (resp. $\bm{x}_{2}$, $\bm{x}_{3}$).

The proposed analogy naturally defines the graphs on which are defined the downsampled signals $\bm{x}_{l}$. Let us introduce a supernode $k$ standing for each subgraph $\mathcal{G}^{(k)}$. For $l=1$, the approximation signal $\bm{x}_{1}$ lies naturally on a graph of adjacency matrix $\bm{A}_{1}$ where $\bm{A}_{1}(k,k^{\prime})$ is the sum of the weights of the edges connecting subgraph $k$ to subgraph $k^{\prime}$ in the original graph. Then, for the subsequent graph of adjacency matrix $\bm{A}_{l}$ on which is defined the detail signal $\bm{x}_{l}$, only supernodes standing for subgraphs of size at least $l$ are needed and $\bm{A}_{l}$ is defined the same way by summing the edges between involved subgraphs. We formalize this analogy in the following.

IV-B Formalization of the operators necessary to the design

To help the assimilation of the definitions introduced here, the reader may in parallel look at Appendix -A, where a trivial concrete example is fully detailed.

Applying $\bm{C}^{(k)\top}$ to a graph signal $\bm{x}$, one obtains the graph signal reduced to subgraph $\mathcal{G}^{(k)}$, i.e. $\bm{x}^{(k)}$. Conversely, applying $\bm{C}^{(k)}$ expands a signal defined on $\mathcal{G}^{(k)}$ to a graph signal on $\mathcal{G}$ with zero-padding.

As a consequence, the intra-subgraph adjacency matrix $\bm{A}_{int}$, that is, the adjacency matrix of the graph that contains only the intra-subgraph edges and none of the inter-subgraph edges, may be written as:

The complement to obtain the full adjacency matrix $\bm{A}$ is called the inter-subgraph adjacency matrix and is defined as $\bm{A}_{ext}=\bm{A}-\bm{A}_{int}$; it keeps only the links connecting subgraphs together.

IV-B2 Subgraph Laplacian operators

On each $\mathcal{G}^{(k)}$, let us define $\bm{\mathcal{L}_{int}}^{(k)}$ the local Laplacian matrix, computed from $\bm{A}_{int}^{(k)}$. It is diagonalisable:

with $\bm{\Lambda}^{(k)}$ the diagonal matrix of sorted eigenvalues ($\lambda_{1}^{(k)}$ is the smallest):

and $\bm{Q}^{(k)}$ the basis of local Fourier modes:

and $\bm{P}^{{(k)\top}}=\left(\bm{Q}^{(k)}\right)^{-1}$ with:

We normalize the $\bm{q}_{i}^{(k)}$ with the $L_{p}$ norm:

Note that in the specific case where $p=2$ for this normalization, we end up with $\bm{P}^{(k)}=\bm{Q}^{(k)}$. We will see that in this case, the filterbank simply codes for an orthogonal transform. We discuss the choice of $p$ (usually 1 or 2) in Section IV-E.

For each $\bm{q}_{i}^{(k)}$ of size $N_{k}$ defined on the local subgraph $\mathcal{G}^{(k)}$, let $\bar{\bm{q}}_{i}^{(k)}$ be its zero-padded extension to the whole global graph:

Similarly $\bar{\bm{p}}_{i}^{(k)}$ stands for the zero-padded extension of $\bm{p}_{i}^{(k)}$.

IV-B3 Analysis, synthesis and group operators

This means that operator $\bm{\Theta}_{l}$ groups together all local Fourier modes associated to the $l$-th eigenvalue of all subgraphs containing at least $l$ nodes.

This means that $\bm{\Omega}_{l}$ groups together indicator functions of subgraphs containing at least $l$ nodes.

IV-B4 On the operators’ uniqueness

Operators as we have defined them are not unique if no further rules are enforced. In fact, for each eigenvector $\bm{q}_{i}^{(k)}$, its opposite $-\bm{q}_{i}^{(k)}$ is also an eigenvector. Moreover, in the case of eigenvalue multiplicity, associated eigenvectors are not unique. To enforce uniqueness, any set of deterministic rules to extract eigenvectors will work. To solve the orientation issue, one may for instance decide to set the first non-zero coefficient of all vectors to be positive. For eigenvalues with multiplicity, we discuss a possible set of rules in Appendix -B that guarantees uniqueness.

IV-C The filterbank design

each of them defined on a graph whose adjacency matrix reads:

By this formula, the convention is that there is no self-loop, i.e., for $k=k^{\prime}$, $\bm{A}_{l}(k,k)=0$. Also, adjacency matrices $\bm{A}_{l}$ for $l\geq 2$ are only needed if one decides to cascade the filterbank on detail channels; here, the cascade will only be done on $\bm{A}_{1}$ (see Fig.3 and Section IV-D).

IV-C2 A remark on the storage of structural information

An important question is whether the total amount of stored information pre- and post-analysis equal or not? In terms of signal information only (i.e., discarding the structural information), the total amount of stored information is equal on both sides of the analysis block as each of the downsampled signals $\bm{x_{l}}$ is of size $|\mathcal{I}_{l}|$ and $\sum_{l}|\mathcal{I}_{l}|=N$.

On the other hand, in terms of structural information (i.e. the information of the adjacency matrices), one needs to keep both the structural information pre- and post-analysis. Indeed, the post-analysis structural information is not enough to reconstruct the original graph (unlike the signal $\bm{x}$ who can be perfectly reconstructed from its approximation and details, as we will see in Section IV-C3). The amount of stored structural information therefore increases after analysis and this is –at least for now– an irreducible storage price to pay. This is a common downfall of all graph filterbanks yet proposed, e.g. . Finding ways to critically sample both the structure and the signal defined on it is part of our ongoing research and is not in the scope of this article.

In the following, all graph structures are stored in the form $\bm{A}=\bm{A}_{int}+\bm{A}_{ext}$, as this does not increase the amount of information (the number of links) but still subtly encodes the connected subgraph structure: indeed the partition $\bm{c}$ can be exactly recovered by detecting the connected components of $\bm{A}_{int}$. In Narang et al. , authors also need to keep an information equivalent to $\bm{c}$: the bipartite graph decomposition of $\bm{A}$ and, for each bipartite graph, the information of the two sets of nodes. It is also the case in Shuman et al. ’s work, where authors need to keep the downsampling vector $\bm{m}$.

IV-C3 Synthesis block

$\bar{\bm{p}}_{l}^{\mathcal{I}_{l}(j)}\bar{\bm{q}}_{l}^{\mathcal{I}_{l}(j)\top}$ is a matrix of size $N\times N$, with non zero coefficients only for indices in subgraph $\mathcal{I}_{l}(j)$. In this non-zero block, it equals $\bm{p}_{l}^{\mathcal{I}_{l}(j)}\bm{q}_{l}^{\mathcal{I}_{l}(j)\top}$. One finally obtains:

the Identity, as $\bm{P}^{(k)\top}=\left(\bm{Q}^{(k)}\right)^{-1}$, which ends the proof. ∎

We show in Fig. 2 a schematic representation of the analysis and synthesis blocks of the proposed graph filterbanks.

IV-C4 Critical sampling and biorthogonality

as shown by Lemma 1. When we choose $p=2$ for the normalization of Eq. (16), this filterbank is orthogonal.

IV-D The analysis cascade

Fig. 4 shows the analysis cascade on a toy signal defined on a simple graph of size 14. Let us take a close look at subgraph number 5 of the original graph: it contains 2 nodes and the signal on each of its node is of same absolute value but of opposite signs. As expected, $\bm{x}_{1}^{(1)}(5)$, the approximation signal on the corresponding supernode at level $(1)$ is null (average of the 2 original values); and $\bm{x}_{2}^{(1)}(5)$, the first detail signal is large: it is the difference between the 2 original values. Moreover, as this subgraph contains only 2 nodes, its local Laplacian does not have a third eigenvector: this subgraph does not participate to the second and third detail signals and its associated supernode does not appear in $(\bm{x}_{3},\bm{A}_{3})^{(1)}$ nor $(\bm{x}_{4},\bm{A}_{4})^{(1)}$.

The analysis cascade decomposes the original signal in 3 detail signals (of sizes 5, 3 and 1) at level $(1)$, 2 detail signals (of sizes 2 and 1) at level $(2)$, 1 detail signal and one approximation signal (both of size 1) at level $(3)$. From these 7 downsampled signals of total size 14, one may perfectly reconstruct the original signal, using the synthesis operators defined in Section IV-C3.

IV-E Atoms of analysis and the choice of normalization

To study the effect of the filterbank, one may look at the dictionary of analysis atoms, and of recovery atoms. With additional assumptions, it could be wavelets. However, we will avoid using the term wavelet and prefer the more general term of “atom”, as they do not necessarily have wavelet’s properties: they are for instance not related to translation on the graph.

To each output of the analysis cascade (approximations and details at all levels) is associated an analysis atom. Approximation analysis atoms (“scaling function-like”) are associated to the first channel of each level: at level $(j)$, they are the columns of $\bm{\Theta}_{1}^{(j)}$, upsampled back to the original graph’s size:

Detail analysis atoms (“wavelet-like”) are obtained from all but the first channel at each level. More precisely at level $(j)$, the detail analysis atoms associated to channel $l\neq 1$ are the columns of $\bm{\Theta}_{l}^{(j)}$, upsampled back to the original graph’s size:

The choice of the $L_{p}$ norm in Eq. (16) is now dictated by the desired properties of the atoms. A first possible choice is $p=1$, i.e. normalization in $L_{1}$, as it is the only normalization that ensures that the detail analysis atoms $\bm{\Psi}_{l}^{(j)}$ have zero mean – a desirable feature to have atoms as close as possible to a wavelet interpretation. Another possibility would be to normalize in $L_{2}$ (as for the Haar filterbank). In this case, detail atoms do not have zero mean in general, however the energy of the modes is constant. In Section VI, the normalization will be application-dependent. The default normalization is with $L_{1}$.

In Fig. 5, we show the 14 analysis atoms corresponding to the analysis cascade of Fig. 4: one approximation atom (from the approximation channel at the last level of the cascade); and 13 detail atoms that represent the other channels.

A property of the atoms is that their support is always compact: each is defined and non-zero only on one subgraph. On the other hand, in the global Fourier domain, the atoms are exactly localized only if the decomposition in subgraphs corresponds exactly to different connected components of the whole graph (in this case, the global Fourier matrix is the concatenation of all local Fourier matrices). If not, the further away is the graph from this disconnected model, the less localized are the atoms in the global Fourier domain.

V Detecting a partition of connected subgraphs

The proposed filterbank explicitly integrates the graph structure in connected subgraphs. A central question arises: how does one choose a particular partition $\bm{c}$ of the graph in connected subgraphs? The partition choice has a strong influence on what the filterbank achieves, as shown in Fig. 6 where we compare the effect of downsampling for two different partitions on a toy graph. The practitioner has the choice among a wide variety of options to find such a partition: he or she could follow graph partitioning techniques of or , or use graph nodal domains – either very high frequency ones as in or others — or any other solution… While the proposed filterbank is well-defined for any of these partitions, the final decision regarding the partitioning algorithm will depend on what the user wants the filterbanks to achieve.

In the following, we show applications for compression and denoising. We seek to typically transform the original signal into a sparser one after analysis. For that, we look for partitions that separate the graph into groups of nodes more connected to themselves than with the rest of the graph: they are known as communities. Indeed, as in image or video compression, we suppose that low-frequencies contain the useful information of the signal. Approximating a community of nodes, each one with its signal value, by a supernode on which is the average over the community is a way to keep such low-frequencies.

Literature is abundant on community detection (see the survey ). To detect non-overlapping communities, we use the greedy Louvain method . It maximizes (approximatively) over all the possible partitions $\mathbf{c}$, the so-called modularity (see ), a well-known objective function that measures the quality of a partition in communities $\mathbf{c}$, defined as:

where $d_{i}=\sum_{j}\mathbf{A}_{ij}$ and $2m=\sum_{i}d_{i}$. The Louvain method iteratively repeats two main phases, starting from an initial situation where each node is in its own community: 1) Select a node and group it with its adjacent node that causes the largest increase of modularity; do this sequentially with all other nodes, until no individual move can improve the modularity; 2) Aggregate each community in a “supernode” and build a new adjacency matrix of this “supernode” graph. Phase 1 is then applied to this new graph, and so on and so forth. The algorithm stops when phase 1 is not able to increase the modularity anymore. We modify this algorithm and have two different implementations:

The SC (Small Communities) implementation. It consists in performing phase 1 only once: this implementation ensures that the partition separates the graph in small communities (typically smaller than 10 nodes).

The LC (Large Communities) implementation. When performing the usual algorithm, a stopping criterion is added: the algorithm is stopped (if not already stopped thanks to the first criterion) before the size of the largest community becomes larger than a given threshold $\tau$. In fact, iterating both phases, communities become gradually larger; and recall that our proposal relies on the diagonalisation of the local Laplacian matrices, which has a cubic computation cost. In order to control computation time, we do not allow communities larger than the threshold, hereafter $\tau=1000$ nodes.

For comparison, in Section VI-B2, we will show some results obtained with another famous multiscale community detection algorithm, called Infomap . With this algorithm also, one may define analog SC and LC implementations.

Note on stochasticity: The Louvain and the Infomap algorithms are stochastic: they do not necessarily output the same partition at every run on the same data. This implies that the output of the analysis cascade of our filterbank may differ from one realisation to another (this is also the case of other methods such as the filterbanks based on bipartite graphs). Stochasticity is not an issue as synthesis operators are built according to the solutions found during the analysis: reconstruction is always perfect. For the results in Table I and Figures 10, 11 and 14, we show the median computed over 10 realisations.

V-B Choice of adjacency matrix

When performing community detection, one may choose to use only the original adjacency matrix $\bm{A}$ as it is, or incorporate some information about the graph signal $\bm{x}$ to follow more closely its evolution. We explore two choices:

CoSub, short for Connected Subgraphs, is based on simply applying the Louvain algorithm on the adjacency matrix $\bm{A}$;

EdAwCoSub, short for Edge AwareWe use the term edge aware in relation to usage in the Signal processing community; the reader can think of it as “signal-adapted” if preferred. Connected Subgraphs, takes the signal $\bm{x}$ into account for subgraph partitioning and modifies the adjacency matrix into

where $\sigma_{x}=\mbox{std}(\{|\bm{x}(i)-\bm{x}(j)|\}_{i\sim j})$ ($i\sim j$ means $i$ neighbor to $j$ in $\mathcal{G}$). This choice of $\sigma_{x}$ is classical in the clustering literature . The Louvain algorithm is then applied on $\bm{A_{x}}$.

The obtained partition enables us to write $\bm{A}=\bm{A}_{int}+\bm{A}_{ext}$ in both cases. Edge-awareness may also be implemented, as in , by adapting image segmentation methods to graph signals; such an advanced comparison between edge-awareness methods is left for future work. In Section VI, we compare the 4 implementations of the proposed filterbank: CoSub SC and LC, EdAwCoSub SC and LC; to methods from the literature.

V-C Complexity of the algorithm

At a given level of the analysis cascade, computing the analysis atoms requires: i) to run the partitioning algorithm: the Louvain algorithm has a linear complexity $O(N)$ ; ii) the diagonalisation of the Laplacian associated to $\bm{A}_{int}$, i.e. a block diagonal matrix containing as many blocks as there are detected communities. Given that the diagonalization of a matrix of size $N$ costs $O(N^{3})$, the diagonalization of a block diagonal matrix containing $K$ blocks of same size thus costs $O(N^{3}/K^{2})$. Overall, at each level of the cascade, computing the analysis atoms costs $O(N+N^{3}/K^{2})$. Typically, if $K=N/\alpha$ with $\alpha$ an average small number of nodes per community (see end of Sec. VI-B for typical values of $\alpha$), the complexity turns out to be $O((\alpha^{2}+1)N)$. Cascading the analysis on all levels thus costs $O(\alpha^{2}N\log{N})$. This is to compare to the global graph Fourier analysis that costs $O(N^{3})$.

VI Applications

All the reported examples are computed using a developed Matlab toolbox that is available for downloadURL: http://perso.ens-lyon.fr/pierre.borgnat/Codes/CoSubFBtoolbox.zip. The comparisons with methods from the literature use the implementations from the original authors, when they are available.

Fig. 7 shows successive approximated signals $\{(\bm{x_{1}},\bm{A_{1}})^{(l)}\}_{l=1:5}$ of a smooth signal defined on the Minnesota traffic graph , using the CoSub SC implementation. Notice how the last level’s approximated signal, even if small in size (12 nodes) still captures the original signal’s information remarkably well.

VI-A2 A small image

Images can be studied as graph signals defined on the two-dimensional regular grid (each pixel is a node, and each node has four neighbors), and may therefore be analyzed by the proposed graph-based filterbank. Consider the small $32\times 32$ image of Fig. 8: it is a graph signal defined on a regular graph of size $1024$. Its graph Fourier transform is represented on the same figure. We analyze this image with the EdAwCoSub SC algorithm and the rest of Fig. 8 represents a selection of atoms of the filterbank, shown both in the node and the global graph Fourier domain, and partially reconstructed images from the projection of the original image on these atoms. For the interested reader, a dedicated PDF file in our Toolbox shows all 1024 atoms. Note that the support of the subgraphs are clearly impacted by edge-awareness. We see that each atom is compactly supported (in the node domain) and only (very) approximately localized in the global graph Fourier domain. Moreover, we see how, within a given level $(j)$, the mean frequency of $\bm{\Psi}_{l}^{(j)}$ increases as $l$ increases. The cause of the frequency delocalisation is that regular grids are not decomposable in a sum of almost disconnected subgraphs: the local Fourier modes on which we base our design are therefore far from localized in the global Fourier domain. In fact, regular grids are a typical structure for which our method (and the graph partition in communities) is not very appropriate. Nevertheless, we still show results on images for pedagogical purposes and in order to compare performance with other methods from the literature.

VI-B Graph signal reconstruction via non-linear approximation

One of the use of classical filterbanks is compression. The main idea relies on the fact that natural signals are approximately smooth at different scales of analysis, and have therefore a sparse representation on filterbanks’ atoms. One may thus transform the signal with the filterbanks, keep the low-pass coefficients and a fraction of high-pass coefficients while setting the others to zero, and still obtain a decent reconstruction of the original signal. In the following, we apply this non-linear approximation (NLA) scheme on images and on the Minnesota traffic graph.

Typical filterbank comparisons using NLA look at reconstruction results after three levels of analysis. In our case, as we do not know beforehand in how many subgraphs the partitioning algorithm will cut the graph, we cannot predict how many low-pass coefficients will be left after a given number of levels of the analysis cascade. Thus, for a comparison with other methods, and for a given compression ratiothe compression ratio is defined as the ratio of the size of the original data over the size of the compressed data, we will compute all non-linear approximations corresponding to all different levels of the cascade, and keep the level for which the reconstruction result is the best.

Consider for instance the benchmark image cameraman shown in Fig. 9 (left). Its size is $256\times 256$ ($N=65536$ is the size of the associated graph signal). Table I compares the reconstruction details after NLA of CoSub LC and SC, EdAwCoSub LC and SC, to the classical image filterbank CDF 9/7, the Graph Bior filterbank with Zero DC, graphBior(6,6) filters and Gain Compensation (GrBior); and the same Graph Bior filterbank but including edge-awareness (EdAwGrBior). Also, Fig. 10 recaps the PSNR of reconstruction for each of the three benchmark images of Fig. 9. We see that our method, without edge-awareness, does not perform as well as GrBior. This is due to the fact that our method is not best suitable to 2D grids as they do not have a natural structure in communities and, on the contrary, GrBior is best suitable to bipartite graphs, of which 2D grids are an example. On the other hand, when adding edge-awareness, the graph becomes more structured and we obtain results similar to EdAwGrBior.

Note on the LC implementation. We see here that the best reconstruction for the LC implementation is always obtained after only one level of analysis. In this case, the filterbanks may hardly be seen as a multiscale analysis, but rather as a graph windowed Fourier transform, where the window is simply an indicator function on each subgraph. Note that this window changes from one subgraph to another and it is hence different from the proposition of for windowed graph Fourier transform. To observe a multiscale analysis with the LC implementation, one needs to either decrease the threshold $\tau$ or increase the graph’s size.

VI-B2 Reconstruction performance for graph signals

The graph signal model underlying the NLA scheme is that graph signals should be smooth with respect to the topology on which they are defined. For instance, let us consider the graph signal of Fig. 7 (left): it is by construction smooth with respect to the underlying graph as it is the sum of the first five eigenvectors of its Laplacian matrix. We compare in Fig. 11 (left) the reconstruction performance for our filterbank implementations, to Shuman’s Laplacian pyramid filterbanks . This method was not originally written with edge-awareness, but one can simply consider $\bm{A}_{x}$ (as in Eq. (31)) instead of $\bm{A}$ to make it edge-aware and define what we call the edge-aware Laplacian pyramid method (EdAw Lap. Pyr.). Also, up to our knowledge, GrBior filterbanks have only been implemented for one-level cascades on arbitrary graphs, which explains why we do not consider them here. The full black line on the same Figure shows the performances for a random Gaussian signal of zero mean and variance 1, normalized to have the same energy as the smooth signal (all methods collapse on the same black line). As expected, random signals are dense on any analysis atoms, and reconstruction is comparatively poor. Moreover, we see that our proposed filterbanks really have an edge for signals who are smooth compared to the community structure of the underlying graph. On the right of Fig. 11 are represented the performances obtained with the Infomap algorithm rather than the Louvain algorithm (see Section V-B). In this particular case, performances with Infomap are better. Empirically, we find that using the Louvain algorithm or the Infomap algorithm yields in general similar results.

Note on the typical size of communities: In this Minnesota example (resp. cameraman example), the typical community size of the first level of the cascade is 2 (resp. 2) for CoSub SC, 5 (resp. 5) for EdAwCoSub SC, 80 (resp. 200) for CoSub LC, and 40 (resp. 100) for EdAwCoSub LC.

VI-C Application in denoising, on the Minnesota traffic graph

Another application of filterbanks is denoising. We consider first a piece-wise constant graph signal (that has only two possible values: +1 and -1) defined on the Minnesota traffic graph, as shown in Fig. 12a; so as to compare our proposition with previously published methods. We corrupt this signal with an additive Gaussian noise of standard deviation $\sigma$. Fig. 12b shows such a corrupted signal with $\sigma=1/4$. We then attempt to recover the original image by computing the first level of the analysis cascade, and reconstructing the signal from all low-pass coefficients and thresholded high-pass coefficients having absolute value higher than $T=3\sigma$. In order for such a thresholding scheme to be justified for denoising, the energy of coefficients associated to noise should have a constant variance in all the details after analysis. For that, we use here a $L_{2}$ normalization of the local Fourier modes (rather than $L_{1}$) for this denoising experiment (see the discussion in Section IV-E).

Fig. 12 compares results obtained with EdAwGrBior and EdAw Lap. Pyr. filterbanks to our proposition, for $\sigma=1/4$. Fig. 14 (left) summarizes SNR results for different values of $\sigma$. These results may be compared to those obtained by Sakiyama et al. and summarized in Table 5 of . We study also the denoising on the smooth signal of Fig. 7. Results are shown in Fig. 13 for $\sigma=1/4$ and summarized in Fig. 14 (right) for different values of $\sigma$.

All four of our implementations outperform GrBior. Moreover, CoSub LC and the Laplacian pyramid obtain similar results; our method performing slightly better at high noise level. Edge-awareness helps in the case of the piece-wise constant signal and not so much for the smooth signal.

VII Conclusion

While previous methods are based on global filters defined in the global Fourier space of the graph, we defined local filters based on the local Fourier spaces of each connected subgraph. Thanks to this paradigm, a simple form of filterbanks is designed, that one could call Haar graph filterbank.

We first illustrated this for compression on images, mainly for pedagogical and state-of-the-art comparison purposes. In fact, without edge-awareness, our proposition is not really appropriate for such regular structures. Edge-awareness, on the other hand, by giving structure to the network raises performance to the state-of-the-art. The improvement over existing methods becomes truly apparent for irregular graphs for which a community structure exists. Existence of communities is a very common, if not universal, property of real-world graphs; and our filterbanks rely on this particular organization of complex networks. For such graphs, our proposition outperforms existing ones on non-linear approximation experiments and equals state-of-the art on denoising experiments.

Within this framework, future work will concentrate on extending the local filters to more sophisticated filters, and on finding ways to critically sample jointly the graph structure and the graph signal defined on it.

Consider the trivial graph composed of five nodes shown in Fig. 15. Three of them form a closed triangle. The other two are connected. The triangle and the pair are connected to each other with only one link. Its adjacency matrix $\bm{A}$ reads:

In this example, we consider the partition that separates the triangle (subgraph $\mathcal{G}^{(1)}$) from the pair (subgraph $\mathcal{G}^{(2)}$):

Therefore $\Gamma^{(1)}=(1,2,3)$ is the list of the nodes composing the triangle, and $\Gamma^{(2)}=(4,5)$ is the list of the nodes composing the pair. The subsampling operators associated to $\mathcal{G}^{(1)}$ and $\mathcal{G}^{(2)}$ read:

The intra- and inter- subgraph adjacency matrices read:

The local Laplacian operators $\bm{\mathcal{L}_{int}}^{(1)}$ and $\bm{\mathcal{L}_{int}}^{(2)}$ read:

In the following, we choose the $L_{1}$ normalisation for the $\bm{Q}^{(k)}$. One may diagonalize $\bm{\mathcal{L}_{int}}^{(1)}$ and obtain $\bm{Q}^{(1)}$ and $\bm{P}^{{(1)}}$:

as well as $\bm{\Lambda}^{(1)}=\mbox{diag}(0,3,3)$. One may also diagonalize $\bm{\mathcal{L}_{int}}^{(2)}$ and obtain $\bm{\Lambda}^{(2)}=\mbox{diag}(0,2)$ as well as $\bm{Q}^{(2)}$ and $\bm{P}^{{(2)}}$:

Here, the approximated graph’s adjacency matrix reads:

Let us add a second level of analysis where $\bm{c}^{(2)}=(1,1)$: we group together the two nodes of the approximated graph. The second-level approximated graph is thereby reduced to one node, and the analysis operators are:

Therefore, the 4 detail analysis atoms read:

The approximation analysis atom at level $(2)$ reads:

If using a $L_{2}$ normalization, $\bm{\Psi}_{2}$ would read:

-B A proposition for uniqueness of the operators

To enforce uniqueness of the graph Fourier basis in the case of eigenvalue $\lambda$ having multiplicity $m>1$, a possible rule can be set as follows. Consider the first vector of this eigenspace. All we know is its orthogonality with vectors of all other eigenspaces, i.e. $N-m$ vectors. We decide to arbitrarily force its last $m$ coefficients to zero, and then find the unique $N-m$ coefficients that respects orthogonality with other known vectors and proper normalization. Note that, if at least one of the vectors of the other eigenspaces have non-zero values only on these last $m$ coefficients, we then look for the set of $m$ coefficients closest possible to the last one such that uniqueness is guaranteed. For the second vector, it has to be orthogonal to the already decided $N-m+1$ vectors: we arbitrarily force its $m-1$ last coefficients to zero and find the unique set of its set of coefficients thanks to orthogonality. And so on and so forth up to the multiplicity $m$.

Note that, for practical implementations, classical functions for eigenvector computation (for instance eig or svd in Matlab) empirically output the same choice of eigenvectors when run on two exactly identical inputs, even when there are eigenvalues with multiplicity.

References