Stationary signal processing on graphs

Nathanaël Perraudin; Pierre Vandergheynst

Source-linked AI summary

Stationary signal processing on graphs

Nathanaël Perraudin, Pierre Vandergheynst

arXiv:1601.02522v5Jan 11, 2016cs.DSstat.APstat.ML

TL;DR

The paper addresses how to model statistical regularity for signals on irregular domains where classical translation-based stationarity is unsuitable. It defines graph stationarity through localization, develops scalable PSD estimation and Wiener-type procedures, and applies them to noisy and partially observed signals. The framework supports denoising, regression, and imputation, while graph structure can remain unrecoverable from stationary signals with uninformative covariance.

Problem
Classical stationarity requires a regular structure, motivating a graph-based formulation for signals on irregular domains without an underlying regular continuum.
Method
The paper defines graph stationarity through a graph localization operator, generalizes Welch PSD estimation, and uses the PSD in Wiener-based graph regularization.
Results
The framework provides scalable PSD estimation and Wiener-type procedures for denoising, regression, and missing-data imputation on graph signals.
Takeaways & Limitations
Graph stationarity offers a model for exploiting relationships encoded by graph edges in signal analysis and estimation.
Takeaways & Limitations
Stationary signals with uncorrelated unit-variance samples have covariance Σx = I, which loses all graph information despite stationarity.

Abstract

from arXiv · show

Graphs are a central tool in machine learning and information processing as they allow to conveniently capture the structure of complex datasets. In this context, it is of high importance to develop flexible models of signals defined over graphs or networks. In this paper, we generalize the traditional concept of wide sense stationarity to signals defined over the vertices of arbitrary weighted undirected graphs. We show that stationarity is expressed through the graph localization operator reminiscent of translation. We prove that stationary graph signals are characterized by a well-defined Power Spectral Density that can be efficiently estimated even for large graphs. We leverage this new concept to derive Wiener-type estimation procedures of noisy and partially observed signals and illustrate the performance of this new model for denoising and regression.

1 Introduction

The paper generalizes stationarity from regular domains to signals on weighted undirected graphs, using graph structure to model relationships that classical methods may miss. It develops graph-based PSD estimation and Wiener-style optimization for prediction, denoising, regression, and missing-data imputation.

Classical stationarity relies on statistical regularity under translations, which requires a regular structure such as time.
The paper extends stationarity to irregular graph domains without assuming an underlying regular continuum.Signals are scalar values on vertices connected by weighted undirected edges.
The paper motivates graph stationarity with examples where values vary smoothly along graph edges, including survey responses among socially connected users.This structure can support predicting answers for users who never completed the survey.
Graph localization generalizes translation by adapting a spectral kernel to each node and the graph topology.On a ring graph, localization reduces to traditional cyclic translation.
The framework generalizes Welch PSD estimation, providing a scalable and robust estimator that improves covariance estimation when signals are limited.
A Wiener-based regularization term depends on noise level and signal PSD, replacing the traditional Dirichlet prior for graph signal optimization.The resulting scheme is presented as a MAP estimator and is intended for denoising, regression, and semi-supervised learning.

2 Background theory

The paper develops graph-based signal-processing foundations by defining signals, Laplacian spectral representations, and localized graph filters. It then connects localization to generalized translation and extends wide-sense stationarity concepts from periodic signals to graph vertices.

Graph and signal representation: A graph consists of vertices connected by weighted undirected edges, and a graph signal assigns one scalar value to each vertex.
Graph and signal representation: The combinatorial Laplacian L = D − W is symmetric positive semidefinite, with orthonormal eigenvectors and non-negative eigenvalues that generalize squared frequencies.For connected graphs, the zero eigenvalue is unique; its multiplicity equals the number of connected components.
Spectral representation and filtering: The graph Fourier transform uses Laplacian eigenvectors as a basis, with ˆf = U ∗f and inverse f = U ˆf; on cyclic graphs it matches the discrete Fourier transform.
Spectral representation and filtering: Graph filtering extends convolution by multiplying graph-Fourier coefficients pointwise through a spectral kernel g(λℓ).Direct computation generally requires O(N^3) complexity, motivating approximations for large graphs.
Localization and translation: Because irregular graphs lack an intuitive shift, the localization operator applies a spectral kernel at a node while adapting its shape to graph topology.For a K-order polynomial kernel, its support lies within a radius-K ball around the node; in irregular graphs, localized filters can vary across vertices.
Stationarity: Classical wide-sense stationarity requires translation-invariant first two moments, whereas graph stationarity uses graph vertices as the underlying structure.For temporal WSS signals, autocorrelation depends only on time difference and is linked to the Power Spectral Density through the Wiener-Khintchine theorem.

3 Stationarity of graph signals

The paper generalizes wide-sense stationarity from regular domains to signals on weighted undirected graphs by replacing translation with graph localization. GWSS signals have graph-structured covariance, diagonal graph spectral covariance, and a PSD, while filtering preserves stationarity.

Localization gives graph stationarity a scale-sensitive interpretation because localized kernels encode how correlation depends on graph distance.Polynomial kernels of degree K are exactly localized within graph-radius K, making low-degree kernels especially useful.
GWSS signals have constant first moments and covariance obtained by localizing a graph kernel over the graph.The localization operator adapts the kernel to local graph structure, producing correlations between vertices and uncorrelated spectral components.
A graph signal is GWSS exactly when its covariance matrix is jointly diagonalizable with the graph Laplacian.This characterization links the covariance structure directly to the Laplacian eigenvectors and graph frequencies.
The graph spectral covariance matrix is U*ΣxU; for GWSS signals it is diagonal, defining the graph power spectral density.The PSD entries are given by the diagonal spectral covariance values γx(λℓ).
On a regular cyclic graph, the localization operator becomes traditional translation, recovering classical cyclic-stationarity results while extending stationarity to irregular domains.The framework also establishes equivalence with Girault’s isometric-translation definition, though repeated eigenvalues make the proposed definition more restrictive.
Graph filtering preserves GWSS, changes the mean to mxg(0), and reshapes the PSD according to the filter.Filtering white noise can therefore generate stationary graph signals with a prescribed PSD, including when the white noise is non-Gaussian.

4 Estimation of the signal PSD

The paper estimates graph-signal PSDs by averaging normalized squared coefficients of graph-windowed Fourier transforms, avoiding Laplacian diagonalization. The estimator extends Welch-style spectral averaging to irregular graphs and scales with graph edges.

Method: The method generalizes Welch’s PSD estimation by replacing time-domain translated windows with spectral kernels shifted by multiples of τ.The shifted kernels are localized at every graph node through the windowed graph Fourier transform.
Method: Averaging squared transform coefficients over vertices estimates the PSD, with normalization by each graph window’s Frobenius norm to account for irregular spectra.The resulting discrete estimates are interpolated into a continuous PSD estimator.
Method: A single signal realization can be used, while additional realizations may be averaged to reduce estimator variance.The procedure designs a filterbank, computes the estimator, approximates filter norms with random white signals, then interpolates the resulting spectral samples.
Bias–variance trade-off: The estimator trades bias against variance: sharper spectral localization reduces bias, whereas smaller vertex-domain windows can reduce variance when correlations are localized.The balance depends on the kernel localization and the spatial extent of signal correlations.
Experiments and complexity: On a 20′000-node graph with K = 1 realization, the algorithm successively estimates the PSD, and its computation time scales linearly with the number of edges.The experiment uses m = 30 filters; the complexity is O(Oc(K + K2)M|E|).
Error analysis: Approximation error comes from estimator bias, finite-sample estimation, and polynomial filtering, with bias controlled by σ and finite-sample error improving as 1/K1.Polynomial-filtering error is usually negligible for rough approximations where σ ≫ λmax/N but may become large otherwise.

5 Graph Wiener filters and optimization framework

The section extends Wiener estimation to graph signals and introduces Wiener optimization for general linear measurement operators. The framework incorporates graph-frequency signal and noise structure, admits optimality results, and supports scalable iterative solutions.

Graph Wiener filters: Graph Wiener filters replace classical frequencies with graph Laplacian eigenvalues and estimate signals through graph filters.The Wiener estimation scheme models a stationary signal generated by s(L), measurements filtered by h(L) with additive noise, and recovery through g(L).
Wiener optimization: The proposed Wiener optimization handles general linear operators H that are not assumed to be diagonalizable with the graph Laplacian.Because direct Wiener filtering is unavailable in this case, the method uses a regularized optimization problem incorporating the signal PSD and noise PSD.
Wiener optimization: The optimization addresses fixed-regularization limitations by weighting graph frequencies according to their signal-to-noise ratios.Low-SNR frequencies receive stronger penalization, while the weighting depends on the signal and noise PSDs rather than only a constant parameter.
Model classes: The framework generalizes smooth and band-limited graph-signal assumptions through frequency-dependent PSD choices.Smooth signals use PSD proportional to 1/λℓ, while band-limited signals have zero PSD above a selected eigenvalue.
Theoretical guarantees: For GWSS Gaussian signals with i.i.d. Gaussian noise, the optimization problem is a MAP estimator, and with white noise it yields the linear minimum mean square estimator.When H is jointly diagonalizable with L, the problem is optimal for the weighting and can be solved by a single filtering operation.
Algorithmic solution: A forward-backward splitting scheme solves the optimization efficiently when the Wiener operator is ill-conditioned or too costly to form explicitly.The method applies gradient and proximal steps iteratively, using the corresponding Wiener filter for the proximal operation.

6 Evidence of graph stationarity: illustration with USPS

The USPS experiments test whether data become stationary when represented on data-adaptive nearest-neighbor graphs rather than regular grids. Individual digit classes show especially strong stationarity on the nearest-neighbor graph.

Experimental setup: The study compares circularly shifted digits, original digits, and digit classes 3, 7, and 9 using grid and 20-nearest-neighbor graphs.Pixel means are removed before covariance-based stationarity analysis.
Stationarity results: Nearest-neighbor graphs adapt to the data, and individual digits are stationary on these graphs.The graph connects pixels using feature-vector distances and exponential edge weights.
Stationarity results: The stationarity measure increases as the covariance becomes more diagonal, with nearest-neighbor graphs outperforming grids for specific digit classes.Less universal data require a finer structure than a regular grid.
Digit-3 illustration: For digit 3, covariance and Laplacian eigenvectors are very similar and exhibit a 3-like shape.Using only 50 samples, the associated graph and PSD can generate digit-shaped samples by filtering i.i.d. Gaussian noise.
Graph construction: The experiment motivates nearest-neighbor graphs as proxies for covariance structure because strongly correlated variables receive strong edge weights.This links graph construction from signal realizations to the spectral content of the covariance.

7 Experiments

Experiments evaluate graph Wiener deconvolution, in-painting, and prediction on synthetic and real data. Wiener methods outperform classical graph regularizers across the reported synthetic and meteorological comparisons, including nonstationary temperature data.

Graph Wiener deconvolution: In synthetic graph deconvolution, Wiener filtering recovers signals from heat-kernel measurements with additive noise and outperforms Tikhonov and TV.The comparison attributes the advantage to the Wiener method’s more appropriate prior assumption.
Graph Wiener in-painting: Synthetic in-painting estimates an unknown PSD from 50 signals and compares Wiener optimization with Tikhonov, TV, and a covariance-based MAP estimator.The experiment evaluates relative in-painting error as the number of measurements changes.
Meteorological dataset: The meteorological experiments use hourly temperature and humidity observations from 744 time steps on a geographical nearest-neighbor graph.Station means are removed, and the observations provide realizations for covariance estimation and stationarity checks.
Prediction - Temperature: Wiener optimization always outperforms graph TV and Tikhonov for temperature prediction with 50% measurements and an initial average SNR of 13.4 dB.The temperature data are not really stationary on the distance graph, yet the method remains superior in the reported experiment.
Prediction - Humidity: Humidity prediction uses a 50% measurement mask and varying noise levels, producing conclusions similar to the temperature experiment.The testing framework otherwise remains identical.
Image datasets: USPS and ORL face experiments extend evaluation to image in-painting and denoising with PSD estimation from limited training images.The ORL images contain 10,304 pixels, making explicit covariance computation and Gaussian MAP estimation difficult.

8 Conclusion

The paper extends stationarity to graph signals and uses the resulting PSD to construct Wiener regularization for large-graph estimation problems.

The contribution extends the common concept of stationarity to signals defined on graphs.
The proposed PSD-based Wiener regularization framework addresses denoising, regression, and semi-supervised learning on large graphs.The conclusion presents these as applications enabled by efficient PSD estimation.
The work experimentally connects graph-signal stationarity with popular nearest neighbor graph construction.The authors identify clarifying this connection and inferring graph structure as future work.

A Convex models

This section formulates graph-based convex models for denoising, deconvolution, and inpainting, using smoothness or sparsity of graph gradients as priors.

Graph inpainting models use xT Lx to penalize graph-gradient energy and encourage signals that vary smoothly along graph edges.The formulation includes a masking operator M and noise-dependent constraint parameter ϵ.
Known Gaussian noise allows the constraint parameter to be set as ϵ = n√#y, avoiding a greedy search over γ.Here n is the noise standard deviation and #y is the number of measurements.
The same prior framework is extended to combine deconvolution, inpainting, and denoising.The deconvolution formulation introduces a convolution kernel h.
USPS and ORL experiments illustrate graph construction, PSD behavior, and inpainting reconstructions using competing priors.In the ORL example, Wiener optimization reaches SNR 14.42 dB, compared with 13.53 dB for classic TV and 12.12 dB for Tikhonov.
An ℓ1 graph-gradient penalty favors a small number of major signal changes and is therefore suited to piecewise-smooth graph signals.The graph gradient is defined on edges, while the ℓ2 penalty corresponds to xT Lx.
The optimization problems are solved with proximal splitting methods.

B Proof of Theorem 3

The proof derives the posterior distribution for Gaussian graph signals under noisy linear measurements using Bayes’ law.

The proof assumes x is a Gaussian random multivariate signal and models the measurements conditionally on x.
For simplicity, s2(L) is assumed invertible, although the proof states that this assumption is not necessary.
Bayes’ law combines the signal and measurement probabilities to obtain the posterior expression.The resulting weighting is written as w(L) = σs−1(L).

C Proof of Theorem 5

The proof shows that the graph-domain optimization becomes separable in the graph Fourier domain and yields Wiener filtering when the regularization spectrum is chosen from the signal and noise PSDs.

Assuming zero mean, the proof rewrites the optimization in the graph Fourier domain using Parseval’s identity.
Because h(Λ) is diagonal, the optimization separates across graph eigenvalues λℓ.
The Fourier-domain solution is expressed as ¯x = h w2 + h2 ˆy after suppressing eigenvalue-index notation.
Substituting ˆy = hˆx + ˆwn yields an error expression involving signal and noise terms weighted by w2 and h2.
Minimizing expected error gives the optimal spectral choice z(λℓ) = n(λℓ) when z is constrained to be positive.The alternative infinite solutions are either infeasible or optimal only when s2 = 0.
Replacing w2(λℓ) with n(λℓ) s2(λℓ) identifies the solution as the Wiener filter associated with h(L) = H.

D Proof of Theorem 4

The proof derives a Wiener estimator for graph-wide-sense stationary signals under noisy linear measurements and shows its equivalence to the linear minimum mean square error estimator.

Measurement model: For GWSS x, the covariance is modeled as Σx = s2(L), with noisy measurements y having covariance Hs2(L)H∗ + σ2I.The measurement model uses i.i.d. noise with PSD σ2 and mean my = Hmx.
Wiener optimization: The Wiener optimization is formulated after centering the signal and measurements by their means.The proof introduces ˜x = x − mx and ˜y = y − mx before solving the optimization problem.
Wiener optimization: Assuming s2(L) and Hs2(L)H∗ + σ2I are invertible simplifies the derivation, although the proof states this assumption is unnecessary.The derivation uses the Woodbury, Sherman and Morrison formula to transform the solution.
Estimator equivalence: The resulting Wiener linear estimator is equivalent to the linear minimum mean square error estimator for y = Hx + wn.The proof explicitly identifies the two estimators as equivalent.
Estimator equivalence: The derived estimator is therefore also a linear minimum mean square estimator.This conclusion follows from the stated equivalence between the Wiener solution and the linear minimum mean square error estimator.

E Development of equation 21

The development relates pairwise squared distances among graph-signal samples to their covariance or Gram matrix through centering and matrix substitution.

Sample distances: The construction considers K samples of a random multivariate variable x defined on an N-vertex graph.The samples are denoted {x1, x2, . . . xK}.
Sample distances: The squared-distance matrix is defined from the sample coordinates and normalized by K.The proof denotes the entries by D[i, j].
Matrix transformation: Double centering transforms the distance matrix into a covariance-related Gram matrix.The displayed derivation includes the centering matrix J and the covariance matrix ¯Σx.
Matrix transformation: Substituting D[i, j] = ¯Σx[i, i] + ¯Σx[j, j] − 2¯Σx[i, j] yields the covariance-based form used in the equation development.This identity connects pairwise squared distances with diagonal and cross terms of ¯Σx.