An -Channel Critically Sampled Filter Bank for Graph Signals
We investigate an -channel critically sampled filter bank for graph signals where each of the filters is supported on a different subband of the graph Laplacian spectrum. We partition the graph vertices such that the set comprises a uniqueness set for signals supported on the subband. For analysis, the graph signal is filtered on each subband and downsampled on the corresponding set of vertices. However, the classical synthesis filters are replaced with interpolation operators, circumventing the issue of how to design a downsampling pattern and graph spectral filters to ensure perfect reconstruction for signals that do not reside on bipartite graphs. The resulting transform is critically sampled and graph signals are perfectly reconstructable from their analysis coefficients. We empirically explore the joint vertex-frequency localization of the dictionary atoms and sparsity of the analysis coefficients, as well as the ability of the proposed transform to compress piecewise-smooth graph signals.
AN -CHANNEL CRITICALLY SAMPLED FILTER BANK FOR GRAPH SIGNALS
|Yan Jin and David I Shuman††thanks: This research has been funded in part by a grant to Macalester College from the Howard Hughes Medical Institute through the Precollege and Undergraduate Science Education Program. ††thanks: The authors would like to thank Federico Poloni for providing a proof to Proposition 2, and Andrew Bernoff for helpful discussions about the matrix partitioning problem discussed in Section 3. ††thanks: MATLAB code for all numerical experiments in this paper is available at http://www.macalester.edu/~dshuman1/publications.html. It leverages the open access GSPBox , into which it will soon be integrated.|
|Department of Mathematics, Statistics, and Computer Science|
|Macalester College, St. Paul, MN 55105|
|Email: firstname.lastname@example.org, email@example.com|
Index Terms— Graph signal processing, filter bank, sampling, interpolation, wavelet, compression
In graph signal processing , transforms and filter banks can help exploit structure in the data, in order, for example, to compress a graph signal, remove noise, or fill in missing information. Broad classes of recently proposed transforms include graph Fourier transforms, vertex domain designs such as [3, 4], top-down approaches such as [5, 6, 7], diffusion-based designs such as [8, 9], spectral domain designs such as [10, 11], windowed graph Fourier transforms , and generalized filter banks, the last of which we focus on in this paper. For further introduction to dictionary designs for graph signals, see .
The extension of the classical two channel critically sampled filter bank to the graph setting was first proposed in . Fig. 1 shows the analysis and synthesis banks, where and are graph spectral filters , and the lowpass and highpass bands are downsampled on complementary sets of vertices. For a general weighted, undirected graph, it is not straightforward how to design the downsampling and the four graph spectral filters to ensure perfect reconstruction. One approach is to separate the graph into a union of subgraphs, each of which has some regular structure. For example, [14, 15] show that the normalized graph Laplacian eigenvectors of bipartite graphs have a spectral folding property that make it possible to design analysis and synthesis filters to guarantee perfect reconstruction. They take advantage of this property by decomposing the graph into bipartite graphs and constructing a multichannel, separable filter bank, while  adds vertices and edges to the original graph to form an approximating bipartite graph. References [17, 18] generalize this spectral folding property to -block cyclic graphs, and leverage it to construct -channel graph filter banks. Another class of regular structured graphs is shift invariant graphs [19, Chapter 5.1]. These graphs have a circulant graph Laplacian and their eigenvectors are the columns of the discrete Fourier transform matrix. Any graph can be written as the sum of circulant graphs, and [20, 21, 22] take advantage of this fact in designing critically sampled graph filter banks with perfect reconstruction. Another approach is to use architectures other than the critically sampled filter bank, such as lifting transforms [23, 24] or pyramid transforms .
Our approach in this paper is to replace the synthesis filters with interpolation operators on each subband of the graph spectrum. While this idea was initially suggested independently in , we investigate it in more detail here. Our construction leverages the recent flurry of work in sampling and reconstruction of graph signals -. The key property we use is that any signal whose graph Fourier transform has exactly non-zero coefficients can be perfectly recovered from samples of that signal on appropriately selected vertices (see, e.g., [26, Theorem 1] [35, Proposition 1]).
2 -Channel Critically Sampled Filter Bank
We consider graph signals residing on a weighted, undirected graph , where is the set of vertices, is the set of edges, and is the weighted adjacency matrix. Throughout, we take to be the combinatorial graph Laplacian , where is the diagonal matrix of vertex degrees. However, our theory and proposed transform also apply to the normalized graph Laplacian , or any other Hermitian operator. We can diagonalize the graph Laplacian as , where is the diagonal matrix of eigenvalues of , and the columns of are the associated eigenvectors of . The graph Fourier transform of a signal is , and filters a graph signal by . We use the notation to denote the submatrix formed by taking the columns of associated with the Laplacian eigenvalues indexed by . Similarly, we use the notation to denote the submatrix formed by taking the rows of associated with the vertices indexed by the set .
We start by constructing an ideal filter bank of graph spectral filters, where each filter is equal to 1 on a subset of the spectrum, and 0 elsewhere. We choose the filters so that for each , for exactly one . Fig. 2 shows an example of such an ideal filter bank. Equivalently, we are forming a partition of and setting
The next step, which we discuss in detail in the next section, is to partition the vertex set into subsets such that forms a uniqueness set for .
Definition 1 (Uniqueness set ).
Let be a subspace of . Then a subset of the vertices is a uniqueness set for if and only if for all , implies . That is, if two signals in have the same values on the vertices in the uniqueness set , then they must be the same signal.
The following equivalent characterization of a uniqueness set is often useful.
The set of vertices is a uniqueness set for if and only if the matrix whose columns are is nonsingular, where each is a Kronecker delta centered on a vertex not included in .
The channel of the analysis bank filters the graph signal by an ideal filter on subband , and downsamples the result onto the vertices in . For synthesis, we can interpolate from the samples on to . If there is no error in the coefficients, then the reconstruction is perfect, because is a uniqueness set for . Fig. 3 shows the architecture of the proposed -channel critically sampled filter bank with interpolation on the synthesis side.
3 Partitioning the Graph Into Uniqueness Sets for Different Frequency Bands
In this section, we show how to partition the set of vertices into uniqueness sets for different subbands of the graph Laplacian eigenvectors. We start with the easier case of and then examine the general case.
First we show that if a set of vertices is a uniqueness set for a set of signals contained in a band of spectral frequencies, then the complement set of vertices is a uniqueness set for the set of signals with no energy in that band of spectral frequencies.
On a graph with vertices, let denote a subset of the graph Laplacian eigenvalue indices, and let . Then is a uniqueness set for if and only if is a uniqueness set for .
This fact follows from either the CS decomposition [36, Equation (32)]) or the nullity theorem [37, Theorem 2.1]. We also provide a standalone proof below that only requires that the space spanned by the first columns of is orthogonal to the space spanned by the last columns, not that is an orthogonal matrix.
Proof of Proposition 1.
We assume without loss of generality that . Suppose first that the set is a uniqueness set for , but is not a uniqueness set for . Then by Lemma 1, the matrix
has full rank, and the matrix
is singular, implying
equation (1) implies that there must exist a vector such that
Yet, implies is orthogonal to , and, similarly, implies is orthogonal to In matrix notation, we have , so has a non-trivial null space, and thus the square matrix is not full rank and is not a uniqueness set for , a contradiction. We conclude that if is full rank, then must be full rank and is a uniqueness set for , completing the proof of sufficiency. Necessity follows from the same argument, with the roles of and interchanged. ∎
The issue with using the methods of Proposition 1 for the case of is that while the submatrix is nonsingular, it is not necessarily orthogonal, and so we cannot proceed with an inductive argument. The following proposition and corollary circumvent this issue by only using the nonsingularity of the original matrix. The proof of the following proposition is due to Federico Poloni , and we later discovered the same method in , [41, Theorem 3.3].
Let be an nonsingular matrix, and be a partition of . Then there exists another partition of with for all such that the square submatrices are all nonsingular.
First consider the case , and let . Then by the generalized Laplace expansion ,
where the sign of the permutation determined by and is equal to 1 or -1. Since , one of the terms in the summation of (2) must be nonzero, ensuring a choice of such that the submatrices and are nonsingular. We can choose as the desired partition. For , by induction, we have
where again is always 1 or -1 and one of the terms in the summation in (3) must be nonzero, yielding the desired partition. ∎
For any partition of the graph Laplacian eigenvalue indices into subsets, there exists a partition of the graph vertices into subsets such that for every , and is a uniqueness set for .
Corollary 1 ensures the existence of the desired partition, and the proof of Proposition 2 suggests that we can find it inductively. However, given a partition of the columns of into two sets and , Proposition 2 does not provide a constructive method to partition the rows of into two sets and such that the submatrices and are nonsingular. This problem is studied in the more general framework of matroid theory in , which gives an algorithm to find the desired row partition into two sets. We summarize this method in Algorithm 1, which takes in a partition of the spectral indices and constructs the partition of the vertices. In Fig. 4, we show two examples of the resulting partitions.
While Algorithm 1 always finds a partition into uniqueness sets, such a partition is usually not unique, and the initial choices of in each loop play a significant role in the final partition. In the numerical experiments in the next section, we use the greedy algorithm in [31, Algorithm 1] to find an initial choice for , and use row reduction after permuting the complement of to the top to find an initial choice for . We discuss this issue further in Section 5.2.
4 Illustrative Examples
We can represent the proposed filter bank as an dictionary matrix that maps graph signals to their filter bank analysis coefficients.
4.1 Joint vertex-frequency localization of atoms and example analysis coefficients
We start by empirically showing that the dictionary atoms (the columns of ) are jointly localized in the vertex and graph spectral domains. On the Stanford bunny graph  with 2503 vertices, we partition the spectrum into five bands, and show the resulting partition into uniqueness sets in Fig. 6(c). The first row of Fig. 5 shows five different example atoms whose energy is concentrated on different spectral bands. We see that these atoms are generally localized in the vertex domain as well. Some atoms such as the one shown at wavelet scale 3 are more spread in the vertex domain, possibly as a result of using ideal filters in the filter bank, an issue we revisit in Section 5.1.2. The second row of Fig. 5 shows the spectral content of all atoms in each band, with the average for each represented by a thick black line. As expected, the energies of the atoms are localized in the spectral domain as well. While the transform is not orthogonal, each atom is orthogonal to all atoms concentrated on other spectral bands. One possible extension is to find a fast algorithm to orthogonalize the atoms within each band. Note also that the wavelet atoms at all scales have mean zero, as they have no energy at eigenvalue zero.
Next we apply the proposed transform to a piecewise-smooth graph signal that is shown in the vertex domain in Fig. 6(a), and in the graph spectral domain in Fig. 6(b). The full set of analysis coefficients is shown in Fig. 6(d), and these are separated by band in the third row of Fig. 5. We see that with the exception of the lowpass channel, the coefficients are clustered around the two main discontinuities (around the midsection and tail of the bunny). The bottom row of Fig. 5 shows the interpolation of these coefficients onto the corresponding spectral bands, and if we sum these reconstructions together, we recover exactly the original signal in Fig. 6(a).
4.2 Sparse approximation
Next, we compress a piecewise-smooth graph signal via the sparse coding optimization
where is a predefined sparsity level. After first normalizing the atoms of various critically-sampled dictionaries, we use the greedy orthogonal matching pursuit (OMP) algorithm [45, 46] to approximately solve (4). For the -channel filter bank, the partition into uniqueness sets is shown in Fig. 4, and the filter bank is shown in Fig. 2. We show the reconstruction errors in Fig. 7(d).
5 Ongoing Work
We now briefly outline a number of extensions that we are currently investigating for the proposed transform.
5.1 Computational approximations to improve scalability
In the numerical examples in the previous sections, we have computed a full eigendecomposition of the graph Laplacian and subsequently used it for all three of the sampling, filtering, and interpolation operations; however, such an eigendecomposition does not scale well with the size of the graph as it requires operations with naive methods. We briefly mention a few potential computational approximations that would improve the scalability of the proposed filter bank.
Reference  has a nice review of the computational complexities of the various routines for identifying uniqueness sets. The random sampling method proposed in  does not require the full eigendecomposition of the graph Laplacian and therefore scales significantly better with the size of the graph. We are currently examining whether there is a way to extend the approach of  to non-uniformly randomly sample in a manner that leads with high probability to uniqueness sets for higher bands of the graph Laplacian spectrum.
One possibility for efficiently designing the filter bank without full knowledge of the spectrum is to use the tight warped filter bank design of . These can either be adapted to the width of the spectrum using only an estimate of the maximum eigenvalue, or to an approximation of the eigenvalue distribution based on the spectrum slicing method of [47, Section 3.3]. In either case, the filters can then be applied using the polynomial approximation method discussed in [10, 48]. An advantage of the approximation is that the value of the filtered signal at a given vertex is a weighted average of the signal values in a neighborhood around that vertex . However, the interpolation step needs to be adjusted accordingly.
If are the analysis coefficients of the branch, the standard least squares reconstruction is to solve , and let and ; however, this requires a full eigendcomposition of to get . Our ongoing work includes developing an efficient optimization algorithm that directly uses a function of rather than to stably reconstruct signals supported on a specific spectral band from its samples.
5.2 Reconstruction robustness to noisy or thresholded filter bank coefficients
For any partition into uniqueness sets, the proposed transform perfectly recovers the signal from the analysis coefficients. However, the choice of the downsampling patterns and interpolation method can dramatically affect the stability of the reconstruction if the analysis coefficients are corrupted by noise or are only partially available, or approximate filters are used to speed up the transform (see, e.g., [35, Section III.B] for more on reconstruction stability). Thus, we plan to investigate the role of the initial choice of the ’s in Algorithm 1 in improving the stability of the reconstruction.
5.3 Iterated filter bank
As with the classical wavelet construction, it is possible to iterate the filter bank on the output from the lowpass channel. This could be beneficial, for example, in the case that we want to visualize the graph signal at different resolutions on a sequence of coarser and coarser graphs. An interesting question is how iterating the filter bank with fewer channels at each step compares to a single filter bank with more channels and each filter supported on a smaller region of the spectrum.
5.4 Signals that are sparsely represented by the MCSFB transform
We plan to more formally characterize the relationships between the decay of the analysis coefficients under the proposed transform and different properties of graph signals, as well as the underlying graph structure. Globally smooth signals trivially lead to sparse analysis coefficients because the coefficients are only nonzero for first set of vertices in the partition. The same line of reasoning applies more generally to signals that are sparse in the graph spectral domain. More interesting is a theory for the coefficient decay for piecewise-smooth graph signals such as the ones shown in Fig. 6(a) and Fig. 7(a).
-  N. Perraudin, J. Paratte, D. I Shuman, V. Kalofolias, P. Vandergheynst, and D. K. Hammond, “GSPBOX: A toolbox for signal processing on graphs,” https://lts2.epfl.ch/gsp/.
-  D. I Shuman, S. K. Narang, P. Frossard, A. Ortega, and P. Vandergheynst, “The emerging field of signal processing on graphs: Extending high-dimensional data analysis to networks and other irregular domains,” IEEE Signal Process. Mag., vol. 30, no. 3, pp. 83–98, May 2013.
-  M. Crovella and E. Kolaczyk, “Graph wavelets for spatial traffic analysis,” in Proc. IEEE INFOCOM, Mar. 2003, vol. 3, pp. 1848–1857.
-  W. Wang and K. Ramchandran, “Random multiresolution representations for arbitrary sensor network graphs,” in Proc. IEEE Int. Conf. Acc., Speech, and Signal Process., May 2006, vol. 4, pp. 161–164.
-  A. D. Szlam, M. Maggioni, R. R. Coifman, and J. C. Bremer, Jr., “Diffusion-driven multiscale analysis on manifolds and graphs: top-down and bottom-up constructions,” in Proc. SPIE Wavelets, Aug. 2005, vol. 5914, pp. 445–455.
-  M. Gavish, B. Nadler, and R. R. Coifman, “Multiscale wavelets on trees, graphs and high dimensional data: Theory and applications to semi supervised learning,” in Proc. Int. Conf. Mach. Learn., Haifa, Israel, Jun. 2010, pp. 367–374.
-  J. Irion and N. Saito, “Hierarchical graph Laplacian eigen transforms,” JSIAM Letters, vol. 6, pp. 21–24, 2014.
-  R. R. Coifman and M. Maggioni, “Diffusion wavelets,” Appl. Comput. Harmon. Anal., vol. 21, no. 1, pp. 53–94, 2006.
-  M Maggioni, J. C. Bremer, R. R. Coifman, and A. D. Szlam, “Biorthogonal diffusion wavelets for multiscale representations on manifolds and graphs,” in Proc. SPIE Wavelet XI, Sep. 2005, vol. 5914.
-  D. K. Hammond, P. Vandergheynst, and R. Gribonval, “Wavelets on graphs via spectral graph theory,” Appl. Comput. Harmon. Anal., vol. 30, no. 2, pp. 129–150, Mar. 2011.
-  D. I Shuman, C. Wiesmeyr, N. Holighaus, and P. Vandergheynst, “Spectrum-adapted tight graph wavelet and vertex-frequency frames,” IEEE Trans. Signal Process., vol. 63, no. 16, pp. 4223–4235, Aug. 2015.
-  D. I Shuman, B. Ricaud, and P. Vandergheynst, “Vertex-frequency analysis on graphs,” Appl. Comput. Harmon. Anal., vol. 40, no. 2, pp. 260–291, Mar. 2016.
-  S. K. Narang and A. Ortega, “Local two-channel critically sampled filter-banks on graphs,” in Proc. Int. Conf. Image Process., Hong Kong, Sep. 2010, pp. 333–336.
-  S. K. Narang and A. Ortega, “Perfect reconstruction two-channel wavelet filter-banks for graph structured data,” IEEE. Trans. Signal Process., vol. 60, no. 6, pp. 2786–2799, Jun. 2012.
-  S. K. Narang and A. Ortega, “Compact support biorthogonal wavelet filterbanks for arbitrary undirected graphs,” IEEE Trans. Signal Process., vol. 61, no. 19, pp. 4673–4685, Oct. 2013.
-  A. Sakiyama and Y. Tanaka, “Oversampled graph Laplacian matrix for graph filter banks,” IEEE Trans. Signal Process., vol. 62, no. 24, pp. 6425–6437, Dec. 2014.
-  O. Teke and P. P. Vaidyanathan, “Graph filter banks with M-channels, maximal decimation, and perfect reconstruction,” in Proc. IEEE Int. Conf. Acc., Speech, and Signal Process., Mar. 2016, pp. 4089–4093.
-  O. Teke and P. P. Vaidyanathan, “Extending classical multirate signal processing theory to graphs – Part II: M-channel filter banks,” IEEE Trans. Signal Process., vol. 65, no. 2, pp. 423–437, Jan. 2017.
-  L. J. Grady and J. R. Polimeni, Discrete Calculus, Springer, 2010.
-  V. N. Ekambaram, G. C. Fanti, B. Ayazifar, and K. Ramchandran, “Circulant structures and graph signal processing,” in Proc. Int. Conf. Image Process., Melbourne, Australia, Sep. 2013.
-  V. N. Ekambaram, G. Fanti, B. Ayazifar, and K. Ramchandran, “Critically-sampled perfect-reconstruction spline-wavelet filterbanks for graph signals,” in Proc. IEEE Glob. Conf. Signal and Inform. Process., Dec. 2013, pp. 475–478.
-  M. S. Kotzagiannidis and P. L. Dragotti, “The graph FRI-framework- spline wavelet theory and sampling on circulant graphs,” in Proc. IEEE Int. Conf. Acc., Speech, and Signal Process., Shanghai, China, Mar. 2016.
-  M. Jansen, G. P. Nason, and B. W. Silverman, “Multiscale methods for data on graphs and irregular multidimensional situations,” J. R. Stat. Soc. Ser. B Stat. Methodol., vol. 71, no. 1, pp. 97–125, 2009.
-  S. K. Narang and A. Ortega, “Lifting based wavelet transforms on graphs,” in Proc. APSIPA ASC, Sapporo, Japan, Oct. 2009, pp. 441–444.
-  D. I Shuman, M. J. Faraji, and P. Vandergheynst, “A multiscale pyramid transform for graph signals,” IEEE Trans. Signal Process., vol. 64, no. 8, pp. 2119–2134, Apr. 2016.
-  S. Chen, R. Varma, A. Sandryhaila, and J. Kovačević, “Discrete signal processing on graphs: Sampling theory,” IEEE Trans. Signal Process., vol. 63, no. 24, pp. 6510–6523, Dec. 2015.
-  I. Pesenson, “Sampling in Paley-Wiener spaces on combinatorial graphs,” Trans. Amer. Math. Soc, vol. 360, no. 10, pp. 5603–5627, 2008.
-  S. K. Narang, A. Gadde, and A. Ortega, “Signal processing techniques for interpolation in graph structured data,” in Proc. IEEE Int. Conf. Acoust., Speech and Signal Process., May 2013, pp. 5445–5449.
-  A. Anis, A. Gadde, and A. Ortega, “Towards a sampling theorem for signals on arbitrary graphs,” in Proc. IEEE Int. Conf. Acoust., Speech and Signal Process., May 2014, pp. 3864–3868.
-  A. Gadde and A. Ortega, “A probabilistic interpretation of sampling theory of graph signals,” in Proc. IEEE Int. Conf. Acoust., Speech and Signal Process., Apr. 2015, pp. 3257–3261.
-  H. Shomorony and A. S. Avestimehr, “Sampling large data on graphs,” in Proc. IEEE Glob. Conf. Signal and Inform. Process., Dec. 2014, pp. 933–936.
-  G. Puy, N.Tremblay, R. Gribonval, and P. Vandergheynst, “Random sampling of bandlimited signals on graphs,” arXiv e-Prints, vol. abs/1511.05118, 2015.
-  S. Chen, A. Sandryhaila, J. M. F. Moura, and J. Kovačević, “Signal recovery on graphs: Variation minimization,” IEEE. Trans. Signal Process., vol. 63, no. 17, pp. 4609–4624, Sep. 2015.
-  M. Tsitsvero, S. Barbarossa, and P. Di Lorenzo, “Signals on graphs: Uncertainty principle and sampling,” IEEE Trans. Signal Process., vol. 64, no. 18, pp. 4845–4860, Sep. 2016.
-  A. Anis, A. Gadde, and A. Ortega, “Efficient sampling set selection for bandlimited graph signals using graph spectral proxies,” IEEE Trans. Signal Process., vol. 64, no. 14, pp. 3775–3789, Jul. 2016.
-  C.C. Paige and M. Wei, “History and generality of the CS decomposition,” Linear Algebra Appl., vol. 208-209, pp. 303–326, Sep. 1994.
-  G. Strang and T. Nguyen, “The interplay of ranks of submatrices,” SIAM Review, vol. 46, no. 4, pp. 637–646, 2004.
-  Wikipedia, https://en.wikipedia.org/wiki/Steinitz_exchange_lemma.
-  F. Poloni, “Partitioning an orthogonal matrix into full rank square submatrices,” MathOverflow, http://mathoverflow.net/q/243990.
-  C. Greene, “A multiple exchange property for bases,” Proc. AMS, vol. 39, no. 1, pp. 45–50, Jun. 1973.
-  C. Greene and T. L. Magnanti, “Some abstract pivot algorithms,” SIAM J. Appl. Math., vol. 29, no. 3, pp. 530–539, Nov. 1975.
-  Wikipedia, https://en.wikipedia.org/wiki/Laplace_expansion.
-  D. Gleich, “The MatlabBGL Matlab library,” http://www.cs.purdue.edu/homes/dgleich/packages/matlab_bgl/index.html.
-  Stanford University Computer Graphics Laboratory, “The Stanford 3D Scanning Repository,” http://graphics.stanford.edu/data/3Dscanrep/.
-  J. A. Tropp, “Greed is good: Algorithmic results for sparse approximation,” IEEE Trans. Inf. Theory, vol. 50, no. 10, pp. 2231–2242, Oct. 2004.
-  Michael Elad, Sparse and Redundant Representations, Springer, 2010.
-  B. N. Parlett, The Symmetric Eigenvalue Problem, SIAM, 1998.
-  D. I Shuman, P. Vandergheynst, and P. Frossard, “Chebyshev polynomial approximation for distributed signal processing,” in Proc. Int. Conf. Distr. Comput. Sensor Sys., Barcelona, Spain, Jun. 2011.