ASWT-SGNN: Adaptive Spectral Wavelet Transform-based
Self-Supervised Graph Neural Network

Following the success of CNNs in computer vision and natural language processing, researchers have sought to extend CNNs to the graph domain. The key challenge lies in defining the convolution operator for graphs. Graph convolutional neural networks can be broadly categorized into two main approaches: spectral and spatial. Spectral methods employ the graph Fourier transform to transfer signals from the spatial domain to the spectral domain, where convolution operations are performed. Spectral GNN (Bruna et al. 2013) is the first attempt to implement CNNs on graphs. ChebyNet (Defferrard, Bresson, and Vandergheynst 2016) introduce a parameterization method using Chebyshev polynomials for spectral filters, which enables fast localization of spectral filters. GCN (Kipf and Welling 2016) proposes a simplified version of ChebyNet, which achieves success in graph semi-supervised learning tasks. However, these spectral methods face challenges with generalization, as they are limited by fixed graph sizes, and larger filter sizes result in increased computational and memory costs. Spatial methods draw inspiration from weighted summation in CNNs to extract spatial features from topological graphs. MoNet (Monti et al. 2017) provides a general framework for designing spatial methods by utilizing the weighted average of functions defined within a node’s neighborhood as a spatial convolution operator. GAT (Velickovic et al. 2017) proposes a self-attentive mechanism to learn the weighting function. However, these methods employ unlearnable aggregators, and the localization of the convolution operation remains uncertain.

In graph analysis, contrastive learning was initially introduced by DGI (Velickovic et al. 2019) and InfoGraph (Sun et al. 2019), drawing inspiration from maximizing local-global mutual information. Building upon this, MVGRL (Hassani and Khasahmadi 2020) incorporate node diffusion into the graph comparison framework. GCA (Zhu et al. 2021) learns node representations by considering other nodes as negative samples, while BGRL (Thakoor et al. 2021) proposes a no-negative-sample model. CCA-SSG (Zhang et al. 2021) optimizes feature-level objectives in addition to instance-level differences. GRADE (Wang et al. 2022b) investigates fairness differences in comparative learning and proposes a novel approach to graph enhancement. Several surveys (Xie et al. 2022; Ding et al. 2022; Kumar, Rawat, and Chauhan 2022) summarize recent advancements in graph contrastive learning. Despite these methods’ notable achievements, most rely on GCN and its variants as the base models, inheriting the limitations of GCN. These limitations restrict the performance of these graph contrastive learning methods in tasks that require preserving fine-grained node features.

The wavelet transform exhibits favorable structural properties by utilizing finite length and attenuation basis functions. This approach effectively localizes signals within both the spatial and spectral domains. Additionally, it is noteworthy that the basis and its inverse in the wavelet transform often showcase sparsity, contributing to its utility and efficiency. To construct the wavelet transform on graphs, Hammond et al. (Hammond, Vandergheynst, and Gribonval 2011) propose a method that approximates the wavelet using Chebyshev polynomials. This approach effectively avoids the need for eigen-decomposition of the Laplace matrix. Building upon this, GWNN (Xu et al. 2019) redefines graph convolution based on graph wavelets, resulting in high efficiency and sparsity. M-GWCN (Behmanesh et al. 2022) applies the multi-scale graph wavelet transform to learn representations of multimodal data. Collectively, these works showcase the value of graph wavelets in signal processing on graphs. However, these methods primarily employ wavelet transforms in supervised or semi-supervised tasks and heavily rely on labeled data.

Let $\mathcal{G}=(\mathcal{V},\mathcal{E})$ denotes a graph, where $\mathcal{V}=\{v_{1},\cdots,v_{N}\}$ represents the set of nodes and $\mathcal{E}\subseteq\mathcal{V}\times\mathcal{V}$ represents the set of edges. $\mathcal{G}$ is associated with a feature matrix $\bm{X}\in\mathbb{R}^{N\times p}$, $\bm{X}=[x_{1},\dots,x_{n}]$ and an adjacency matrix $\bm{A}\in\mathbb{R}^{N\times N}$, where $\bm{A}_{ij}=1$  if $(v_{i},v_{j})\in\mathcal{E}$ and $\bm{A}_{ij}=0$  otherwise. We define the Laplacian matrix of the graph as $\bm{L}=\bm{D}-\bm{A}$, where $\bm{D}=Diag(d_{1},\dots,d_{N})$, $d_{i}={\textstyle\sum_{j}\bm{A}_{ij}}$. The symmetric normalized Laplacian is defined as $\bm{L}_{sym}=\bm{D}^{-\frac{1}{2}}\bm{L}\bm{D}^{-\frac{1}{2}}=\bm{U}\bm{\Lambda}\bm{U}^{\top}$. $\bm{U}=\left[{\bm{u}}_{1},\ldots,{\bm{u}}_{N}\right]$, where $\bm{\bm{u}}_{i}\in\mathbb{R}^{N}$ denotes the $i$-th eigenvector of $\bm{L}_{sym}$ and $\bm{\Lambda}=Diag(\lambda_{1},\dots,\lambda_{N})$ is the corresponding eigenvalue matrix.

For a discrete graph signal $\bm{X}$, its Fourier transform has the following form (Shuman et al. 2013):

where $\bm{U}_{\ell}^{*}$ denotes the conjugate transposition of the eigenvalues $\lambda_{\ell}$ corresponding to the eigenvectors. Since the complex number is not involved in the scope of this work, it can be regarded as a common transposition, that is, $\bm{U}^{*}=\bm{U}^{T}$.

Graph Fourier transform provides a means to define the graph convolution operator using the convolution theorem. The filter signal, denoted as $f_{o}$, can be mathematically expressed as:

where $g(\bm{\Lambda})$ denotes the filter function on the eigenvalues, and $\bm{U}g(\bm{\Lambda})\bm{U}^{T}$ is the graph filtering matrix.

The implementation of graph convolution using the Fourier transform lacks spatial localization despite its ability to achieve localization in the spectral domain. Additionally, its localization is significantly influenced by computational efficiency. We employ the spatial and spectral concentration metric proposed by Tsitsvero et al. (Tsitsvero, Barbarossa, and Di Lorenzo 2016) to establish a more comprehensive and precise notion of localization.

where the square of the Euclidean paradigm $\left\|\bm{x}\right\|_{2}^{2}$ denotes the total energy of the signal $\bm{x}$. $\bm{B}$ is the diagonal matrix representing the node restriction operator. Given a node subset $\mathcal{S}\subseteq\mathcal{V}$, $\bm{B}={Diag}\left\{\mathbb{I}(i\in\mathcal{S})\right\}$, where $\mathbb{I}(\cdot)$ is an indicator function. $\bm{C}=\bm{U}\Sigma_{\mathcal{F}}\bm{U}^{-1},$ is the band-limiting operator. Given a matrix $\bm{U}$ and a subset of frequency indices $\mathcal{F}\subseteq\mathcal{V}^{*}$, where $\mathcal{V}^{*}=\{1,\cdots,N\}$ denotes the set of all frequency indices. $\Sigma_{\mathcal{F}}$ is a diagonal matrix defined as $\Sigma_{\mathcal{F}}={Diag}\left\{\mathbb{I}(i\in\mathcal{F})\right\}$.

More specifically, considering a pair of vertex sets $\mathcal{S}$ and frequency sets $\mathcal{F}$, where $a^{2}$ and $b^{2}$ represent the percentage of energy contained within the sets $\mathcal{S}$ and $\mathcal{F}$ respectively, our objective is to determine the balance between $a$ and $b$ and identify signals that can achieve all feasible pairs. The resulting uncertainty principle is formulated and presented in the following theorem (Tsitsvero, Barbarossa, and Di Lorenzo 2016).

where $\lambda_{max}(\bm{BC})$ is the maximum eigenvalue of $\bm{BC}$.

Considering the uncertainty principle, we retain the entire frequency band if we follow the graph localization method in Eq. (2). However, real graph signals show a non-uniform distribution in the frequency domain, which suggests that $\bm{C}$ can be chosen more efficiently. Moreover, a deeper network leads to a wider propagation of the signal in the graph domain, which limits the width of the frequency bands and leads to the dominance of low-frequency signals (smoothing signals), which produce over-smoothing.

As the graph Fourier transform, the graph wavelet transform also necessitates a set of suitable bases to map the graph signal to the spectral domain. In this case, we denote the wavelet operator as $\bm{\Psi}_{s}=\bm{U}g_{s}(\bm{\Lambda})\bm{U}^{\top}$, where $s$ is the scaling parameter. The wavelet transform breaks down a function $g$ into a linear combination of basis functions localized in spatial and spectral. This paper employs the Heat Kernel wavelet as the low-pass filter (denoted as $g^{l}$). In contrast, the Mexican-Hat wavelet is the band-pass filter (denoted as $g^{b}$). These filters are defined as follows:

In this work, we integrate filter functions to achieve the combined effect of a low-pass filter and a wide band-pass filter. Figure 1 exemplifies how the combination of low-pass and band-pass filters can generate a more intricate wavelet filter capable of capturing signal components at various scales. This design enables us to obtain a richer representation of the signal, encompassing its frequency information more nuancedly.

where $\theta=\{s_{0},s_{1},...,s_{l}\}$ is the set of scale parameters.

Using the graph wavelet transform instead of the graph Fourier transform in Eq. (2), we get the graph convolution as follows:

where $\bm{G}$ is a diagonal matrix, which acts as a filter, the scale set $\theta$ of wavelet coefficients is omitted for simplification.

The model formulation discussed in the preceding sections necessitates an eigen-decomposition of the Laplace operator of the input graph $\mathcal{G}$. However, this process presents a computational complexity of $\mathcal{O}(N^{3})$, rendering it infeasible for larger graphs. In order to overcome this constraint, we employ the polynomial approximation method previously introduced by Hammond et al. (Hammond, Vandergheynst, and Gribonval 2011). This method involves expressing the wavelet filter as $g_{\theta}(\lambda)\approx p_{\theta}(\lambda)=\gamma_{0}+\gamma_{1}\lambda+\cdots+\gamma_{m}\lambda^{m}$. This allows rewriting the wavelet operator as $\bm{\Psi_{\theta}}=\bm{U}p_{\theta}(\bm{\Lambda})\bm{U}^{\top}=p_{\theta}(\bm{L}_{sym})$.

While existing methods rely on Chebyshev polynomial approximations, we aim to optimize scales. Therefore, we utilize a least squares approximation that can parameterize the set of wavelet scales $\theta$, which can be expressed as follows:

Where $\bm{V}_{\bm{\Lambda}}\in\mathbb{R}^{N\times(m+1)}$ is the Vandermonde matrix of $\bm{\Lambda}$ from order $0$ to order $m$, $N$ is the number of eigenvalues.

Accurate calculation of the eigenvalues requires an expensive eigen-decomposition of the graph Laplace operator, which is not feasible in our case. As an alternative, we transform the set of eigenvalues into a sequence of linearly spaced points $\bm{\xi}=\{\xi_{i}\}_{i=1}^{K}$ within the interval $[0,2]$ on the spectral domain. However, in graphs that exhibit multiscale features, the eigenvalues do not follow a uniform distribution on the spectral domain. Instead, they display spectral gaps corresponding to different scales within the data. To address this non-uniform distribution, we incorporate the estimated spectral density $\bm{\omega}$ as the weight for each of the $K$ sample points $\bm{\xi}$ on the spectral domain (Fan et al. 2020). Consequently, we can compute the weighted least squares coefficients,

where the spectral density $\bm{\omega}=\{\omega_{j}\}_{j=1}^{K}$, and ${\omega}_{i}=\frac{1}{N}\sum_{j=1}^{N}\{\mathbb{I}(\lambda_{j}=\xi_{i})\}$.

The goal of spectral density estimation is to approximate the density function without expensive graph Laplacian eigen-decomposition. To achieve this, we determine the number of eigenvalues less than or equal to each $\xi_{i}$ in the set $\bm{\xi}$. It can be achieved by computing an approximation to the trace of the eigen-projection matrix $\bm{P}$ (Di Napoli, Polizzi, and Saad 2016). In practice, directly obtaining the projector $\bm{P}$ is often not feasible. However, it can be approximated efficiently using polynomials or rational functions of the Laplacian operator $\bm{L}$. In this approximation, we interpret $\bm{P}$ as a step function of $\bm{L}$, which can be expressed as follows:

Although it is impossible to compute $h(\lambda)$ exactly cheaply, it can be approximated using a finite sum of Jackson-Chebyshev polynomials, denoted as $\phi(\lambda)$. Please refer to Appendix AA for detailed information on the approximation method. In this form, it becomes possible to estimate the trace of $\bm{P}$ by an estimator developed by Hutchinson (Hutchinson 1989) and further improved more recently (Tang and Saad 2012). Hutchinson’s stochastic estimator relies solely on matrix-vector products to approximate the matrix trajectory. The key idea is to utilize Rademacher random variables, where each entry of a randomly generated vector $\bm{R}\in\mathbb{R}^{N}$ takes on the values $-1$ and $1$ with equal probability of $\frac{1}{2}$. Thus, an estimate of the trace $tr(\bm{P})$ can be obtained by generating $n_{r}$ samples of random vectors $\bm{R}_{k}$, $k=1,\cdots,n_{r}$ and computing the average over these samples. The estimator can be expressed as follows:

We obtain the approximation $\bm{\Omega}$ to the cumulative spectral density function.

Finally, the cumulative spectral density $\bm{\Omega}$ is differentiated to obtain an approximation of the spectral density $\bm{\omega}$, ie $\bm{\omega}=\frac{d}{d\bm{\xi}}\bm{\Omega}$. Using the estimated spectral density $\bm{\omega}$ as a weight for each sample point $\xi_{i}$ on the spectral domain, the weighted least squares coefficient $\bm{\gamma}_{\theta}$ can be calculated by substituting it into Eq. (10).

These coefficients are used to approximate the wavelet filter matrix $\bm{\Psi}_{\theta}$, which can be expressed as follows:

In this paper, we construct a graph wavelet multilayer convolutional network (ASWT-SGNN). Based on the above, we define the $l$-th layer of ASWT-SGNN as

where $\sigma$ is the activation function, $\bm{H}^{l^{\prime}}\in\mathbb{R}^{N\times p}$ is the results of graph convolution of the $l$-th layer, $\bm{W}^{l}\in\mathbb{R}^{p\times q}$ is the weight of the $l$-th layer, $p$ is the number of features in current layer, and $q$ is the number of features in next layer. $\bm{H}^{l^{\prime}}$ in Eq. (15) is described as follows:

where $1-\alpha$ represents the proportion of the original features $\bm{H}^{l}$ in $l$-th layer, and $0\leq\alpha\leq 1$. By incorporating this residual connection, we guarantee that regardless of the number of layers we add, the resulting representation will always retain a portion of the initial features. To enhance the encoder’s capability to aggregate both local and global information effectively, we define the diffusion operator $\bm{F}^{l}$ of $l$-th layer in Eq. (16) as

where $\beta$ stands for the ratio of the graph wavelet term and $0\leq\beta\leq 1$. ${\theta}$ is the learnable multi-scale parameters. $\bm{G}^{l}$ is the learnable diagonal filter matrix of $l$-th layer. $\tilde{\bm{A}}=\bm{A}+\bm{I}$ represents the adjacency matrix of the self-loop graph of $\mathcal{G}$, where $\bm{I}$ is the identity matrix.

Typical GCL methods involve generating augmented views and subsequently optimizing the congruence between their encoded representations. In this paper, we generate two augmented graphs, $\bm{z}$ and $\bm{o}$, by using feature augmentation. We randomly sample the mask vector $\bm{m_{f}}\in\{0,1\}^{B}$ to hide part of the dimensions in the node feature. Each element in mask $\bm{m_{f}}$ is sampled from Bernoulli distribution $Ber(1-f_{d})$, where the hyperparameter $f_{d}$ is the feature descent rate. Therefore, the augmented node feature $\hat{\bm{X}}$ calculated by the following formula:

We use the comparison objective for the node representation of the two graph augmentation obtained. For node $v_{a}$, the node representations from different graph augmentation $\bm{z}_{a}$ and $\bm{o}_{a}$ form a positive pair. In contrast, the node representations of other nodes in the two graph augmentation are considered negative pairs. Therefore, we define the paired objective of each positive pair $(\bm{z}_{a},\bm{o}_{a})$ as

where the critic $\theta(\bm{z},\bm{o})$ is defined as $sim(\bm{H_{z}},\bm{H_{o}})$, and $sim(\cdot,\cdot)$ refers to cosine similarity function. $\bm{H_{z}}$ is the graph embedding generated by graph augmentation through the proposed method ASWT-SGNN. The overall objective to be maximized is the average of all positive pairs,

The first step of ASWT-SGNN is wavelet coefficients approximation. Eq. (8) shows that the immediate solution requires eigen-decomposition of the Laplace matrix to obtain the eigenvalues and eigenvectors of the matrix. However, the complexity of the direct solution is very high. For example, the time complexity of the quick response (QR) algorithm is $\mathcal{O}(N^{3})$, and the space complexity is $\mathcal{O}(N^{2})$. Therefore, we use the least squares approximation to approximate the solution in this step. If the $m$-order Chebyshev polynomial approximation is used, the sum of the complexity of each item in the polynomial is calculated as $\mathcal{O}(m\times|E|)$.

The second step of ASWT-SGNN is to use the wavelet transform for graph convolution. Spectral CNN (Bruna et al. 2013) has high parameter complexity $\mathcal{O}(N\times p\times q)$. ChebyNet (Defferrard, Bresson, and Vandergheynst 2016) approximates the convolution kernel by the polynomial function of the diagonal matrix of the Laplace eigenvalue, reducing the parameter complexity to $\mathcal{O}(m\times p\times q)$, where $m$ is the order of the polynomial function. GCN (Kipf and Welling 2016) simplifies ChebyNet by setting m=1. In this paper, the feature transformation is performed first, and the parameter complexity is $\mathcal{O}(p\times q)$. Then the graph convolution is performed, and the parameter complexity is $\mathcal{O}(N)$.

This lemma shows that nodes with similar network neighborhoods have similar wavelet coefficients.

The theorem stated above indicates that the proposed method can map nodes with the same label to an area centered around the expectation in the embedding space. This holds true for any graph in which each node’s feature and neighborhood pattern are sampled from distributions that depend on the node label.

Datasets We evaluate the approach on eight benchmark datasets, which have been widely used in GCL methods. Specifically, citation datasets include Cora, CiteSeer and PubMed (Yang, Cohen, and Salakhudinov 2016), co-purchase and co-author datasets include Photo, Computers, CS and Physics (Suresh et al. 2021). Wikipedia dataset includes WikiCS (Mernyei and Cangea 2020).

The proposed method of approximating wavelet operators is applied to real-world graph data. We evaluate its performance using two distinct metrics: the similarity between the approximated wavelet filter $g_{\theta}$ and the precisely computed wavelet filter, and the Mean Absolute Error (MAE) between the approximated wavelet operator $\bm{\Psi}_{\theta}$ and the actual wavelet operator. Figure 2 presents a specific example showing the accurate and approximate multiscale filters, computed precisely and approximately, demonstrating a significant overlap between them. The MAEs for different scaled training sets are depicted in Figure 3. The experimental results indicate that the proposed wavelet operator approximation method reduces the approximation error in the high-density spectral domain, while avoiding the need for computationally expensive eigen-decomposition.

Table 1 displays the results of node classification accuracy. It is worth noting that our proposed ASWT-SGNN achieves state-of-the-art (SOTA) performance on six out of the eight graphical benchmarks. Specifically, compared to other self-supervised methods across the eight datasets, ASWT-SGNN outperforms the best self-supervised method, MaskGAE, by an average of 0.8%, and it outperforms the worst self-supervised method, GMI, by 3.2%. Furthermore, experiments are carried out within a semi-supervised setting, revealing that the proposed method consistently outperforms the best semi-supervised benchmark, GCNII, by an average margin of 3.1%. Additionally, it surpasses the best self-supervised benchmark, MaskGAE, by an average of 2.0%. These results further demonstrate the effectiveness of the proposed method ASWT-SGNN in node classification tasks.

The proposed method incorporates two pivotal hyperparameters: $\alpha$ and $\beta$. By adjusting these parameters, ASWT-SGNN can be simplified to its core forms: when $\alpha$ is set to 1 and $\beta$ is set to 0, ASWT-SGNN aligns with GCN; setting both $\alpha$ and $\beta$ to 1 makes ASWT-SGNN behave similarly to GWNN (Xu et al. 2019); and when $\alpha\neq 1$ and $\beta$ is set to 0, ASWT-SGNN exhibits similarities to GCNII (Chen et al. 2020). We comprehensively compare ASWT-SGNN and transformation models, including GCN, GWNN, and GCNII. Clear observations can be made from Figure 4: as the number of layers increases, the performance of GCN and GWNN significantly declines. In contrast, the performance of GCNII and ASWT-SGNN remains relatively stable even with more layers stacked. Notably, due to the incorporation of graph wavelet bases, ASWT-SGNN model outperforms GCNII in semi-supervised node classification tasks.

To comprehensively investigate the impact of parameters $\alpha$ and $\beta$ on model performance, we systematically vary their values from 0 to 1. The results are presented in Figure 5, illustrating several noticeable trends. When $\alpha$ is set to 0, indicating the exclusion of residual connections, accuracy is significantly decreased. On the other hand, larger values of $\alpha$ result in fixed node representations, leading to lower performance. $\beta$ represents the proportion of the graph wavelet base; if it is tiny, it may not effectively extract local information. Conversely, excessively large values of $\beta$ could result in the neglect of global information, ultimately reducing performance. These findings highlight the complex relationship between parameters $\alpha$ and $\beta$ in shaping the model’s performance. This delicate balance allows our model to synthesize local and global information, achieving optimal performance effectively.

Furthermore, we use t-SNE (Van der Maaten and Hinton 2008) to visualize the node embeddings. Figure 6 illustrates that compared to MA-GCL, ASWT-SGNN exhibits more pronounced gaps between different classes. This suggests that ASWL-SGNN captures more detailed class information and clarifies the boundaries between samples from different classes. Further extensive experimental analyses, including sparse, robustness, and uncertainty analyses, are exhaustively presented in Appendix BB.