Reduction

The pygsp.reduction module implements functionalities for the reduction of graphs’ vertex set while keeping the graph structure.

tree_multiresolution(G, Nlevel[, …])

Compute a multiresolution of trees

graph_multiresolution(G, levels[, sparsify, …])

Compute a pyramid of graphs (by Kron reduction).

kron_reduction(G, ind)

Compute the Kron reduction.

pyramid_analysis(Gs, f, **kwargs)

Compute the graph pyramid transform coefficients.

pyramid_synthesis(Gs, cap, pe[, order])

Synthesize a signal from its pyramid coefficients.

interpolate(G, f_subsampled, keep_inds[, …])

Interpolate a graph signal.

graph_sparsify(M, epsilon[, maxiter, seed])

Sparsify a graph (with Spielman-Srivastava).

pygsp.reduction.graph_multiresolution(G, levels, sparsify=True, sparsify_eps=None, downsampling_method='largest_eigenvector', reduction_method='kron', compute_full_eigen=False, reg_eps=0.005)[source]

Compute a pyramid of graphs (by Kron reduction).

‘graph_multiresolution(G,levels)’ computes a multiresolution of graph by repeatedly downsampling and performing graph reduction. The default downsampling method is the largest eigenvector method based on the polarity of the components of the eigenvector associated with the largest graph Laplacian eigenvalue. The default graph reduction method is Kron reduction followed by a graph sparsification step. param is a structure of optional parameters.

Parameters
GGraph structure

The graph to reduce.

levelsint

Number of level of decomposition

lambdfloat

Stability parameter. It adds self loop to the graph to give the algorithm some stability (default = 0.025). [UNUSED?!]

sparsifybool

To perform a spectral sparsification step immediately after the graph reduction (default is True).

sparsify_epsfloat

Parameter epsilon used in the spectral sparsification (default is min(10/sqrt(G.N),.3)).

downsampling_method: string

The graph downsampling method (default is ‘largest_eigenvector’).

reduction_methodstring

The graph reduction method (default is ‘kron’)

compute_full_eigenbool

To also compute the graph Laplacian eigenvalues and eigenvectors for every graph in the multiresolution sequence (default is False).

reg_epsfloat

The regularized graph Laplacian is \(\bar{L}=L+\epsilon I\). A smaller epsilon may lead to better regularization, but will also require a higher order Chebyshev approximation. (default is 0.005)

Returns
Gslist

A list of graph layers.

Examples

>>> from pygsp import reduction
>>> levels = 5
>>> G = graphs.Sensor(N=512)
>>> G.compute_fourier_basis()
>>> Gs = reduction.graph_multiresolution(G, levels, sparsify=False)
>>> for idx in range(levels):
...     fig, ax = Gs[idx].plot(title='Reduction level: {}'.format(idx))
pygsp.reduction.graph_sparsify(M, epsilon, maxiter=10, seed=None)[source]

Sparsify a graph (with Spielman-Srivastava).

Parameters
MGraph or sparse matrix

Graph structure or a Laplacian matrix

epsilonfloat

Sparsification parameter, which must be between 1/sqrt(N) and 1.

maxiterint, optional

Maximum number of iterations.

seed{None, int, RandomState, Generator}, optional

Seed for the random number generator (for reproducible sparsification).

Returns
MnewGraph or sparse matrix

New graph structure or sparse matrix

References

See [SS11], [Rud99] and [RV07]. for more informations

Examples

>>> from pygsp import reduction
>>> from matplotlib import pyplot as plt
>>> G = graphs.Sensor(100, k=20, distributed=True, seed=1)
>>> Gs = reduction.graph_sparsify(G, epsilon=0.4, seed=1)
>>> fig, axes = plt.subplots(1, 2)
>>> _ = G.plot(ax=axes[0], title='original')
>>> Gs.coords = G.coords
>>> _ = Gs.plot(ax=axes[1], title='sparsified')
../_images/reduction-1.png
pygsp.reduction.interpolate(G, f_subsampled, keep_inds, order=100, reg_eps=0.005, **kwargs)[source]

Interpolate a graph signal.

Parameters
GGraph
f_subsampledndarray

A graph signal on the graph G.

keep_indsndarray

List of indices on which the signal is sampled.

orderint

Degree of the Chebyshev approximation (default = 100).

reg_epsfloat

The regularized graph Laplacian is $bar{L}=L+epsilon I$. A smaller epsilon may lead to better regularization, but will also require a higher order Chebyshev approximation.

Returns
f_interpolatedndarray

Interpolated graph signal on the full vertex set of G.

References

See [Pes09]

pygsp.reduction.kron_reduction(G, ind)[source]

Compute the Kron reduction.

This function perform the Kron reduction of the weight matrix in the graph G, with boundary nodes labeled by ind. This function will create a new graph with a weight matrix Wnew that contain only boundary nodes and is computed as the Schur complement of the original matrix with respect to the selected indices.

Parameters
GGraph or sparse matrix

Graph structure or weight matrix

indlist

indices of the nodes to keep

Returns
GnewGraph or sparse matrix

New graph structure or weight matrix

References

See [DB13]

pygsp.reduction.pyramid_analysis(Gs, f, **kwargs)[source]

Compute the graph pyramid transform coefficients.

Parameters
Gslist of graphs

A multiresolution sequence of graph structures.

fndarray

Graph signal to analyze.

h_filterslist

A list of filter that will be used for the analysis and sythesis operator. If only one filter is given, it will be used for all levels. Default is h(x) = 1 / (2x+1)

Returns
candarray

Coarse approximation at each level

pendarray

Prediction error at each level

h_filterslist

Graph spectral filters applied

References

See [SFV13] and [Pes09].

pygsp.reduction.pyramid_synthesis(Gs, cap, pe, order=30, **kwargs)[source]

Synthesize a signal from its pyramid coefficients.

Parameters
GsArray of Graphs

A multiresolution sequence of graph structures.

capndarray

Coarsest approximation of the original signal.

pendarray

Prediction error at each level.

use_exactbool

To use exact graph spectral filtering instead of the Chebyshev approximation.

orderint

Degree of the Chebyshev approximation (default=30).

least_squaresbool

To use the least squares synthesis (default=False).

h_filtersndarray

The filters used in the analysis operator. These are required for least squares synthesis, but not for the direct synthesis method.

use_landweberbool

To use the Landweber iteration approximation in the least squares synthesis.

reg_epsfloat

Interpolation parameter.

landweber_itsint

Number of iterations in the Landweber approximation for least squares synthesis.

landweber_taufloat

Parameter for the Landweber iteration.

Returns
reconstructionndarray

The reconstructed signal.

candarray

Coarse approximations at each level

pygsp.reduction.tree_multiresolution(G, Nlevel, reduction_method='resistance_distance', compute_full_eigen=False, root=None)[source]

Compute a multiresolution of trees

Parameters
GGraph

Graph structure of a tree.

NlevelNumber of times to downsample and coarsen the tree
rootint

The index of the root of the tree. (default = 1)

reduction_methodstr

The graph reduction method (default = ‘resistance_distance’)

compute_full_eigenbool

To also compute the graph Laplacian eigenvalues for every tree in the sequence

Returns
Gsndarray

Ndarray, with each element containing a graph structure represent a reduced tree.

subsampled_vertex_indicesndarray

Indices of the vertices of the previous tree that are kept for the subsequent tree.