Reduction¶
The pygsp.reduction
module implements functionalities for the reduction
of graphs’ vertex set while keeping the graph structure.
|
Compute a multiresolution of trees |
|
Compute a pyramid of graphs (by Kron reduction). |
|
Compute the Kron reduction. |
|
Compute the graph pyramid transform coefficients. |
|
Synthesize a signal from its pyramid coefficients. |
|
Interpolate a graph signal. |
|
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')
- 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
- 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.