Year: 2021

Abstract:
Recently, graph neural networks (GNNs) have shown prominent performance in graph representation learning by leveraging knowledge from both graph structure and node features. However, most of them have two major limitations. First, GNNs can learn higher-order structural information by stacking more layers but can not deal with large depth due to the over-smoothing issue. Second, it is not easy to apply these methods on large graphs due to the expensive computation cost and high memory usage. In this paper, we present node-adaptive feature smoothing (NAFS), a simple non-parametric method that constructs node representations without parameter learning. NAFS first extracts the features of each node with its neighbors of different hops by feature smoothing, and then adaptively combines the smoothed features. Besides, the constructed node representation can further be enhanced by the ensemble of smoothed features extracted via different smoothing strategies. We conduct experiments on four benchmark datasets on two different application scenarios: node clustering and link prediction. Remarkably, NAFS with feature ensemble outperforms the state-of-the-art GNNs on these tasks and mitigates the aforementioned two limitations of most learning-based GNN counterparts.

heading: 1 INTRODUCTION
text: In recent years, graph representation learning has been extensively applied in various application scenarios, such as node clustering, link prediction, node classification, and graph classification (Kipf & Welling, 2016b;a; Hamilton et al., 2017; Bo et al., 2020; Hettige et al., 2020; Wang et al., 2016; Wu et al., 2020; Abu-El-Haija et al., 2019). The goal of graph representation learning is to encode graph information to node embeddings. Traditional graph representation learning methods, such as DeepWalk (Perozzi et al., 2014), Node2vec (Grover & Leskovec, 2016), LINE (Tang et al., 2015), and ComE (Cavallari et al., 2017) merely focus on preserving graph structure information. GNN-based graph representation learning has attracted intensive interest by combining knowledge from both graph structure and node features. While most of these GNN-based methods are designed based on Graph AutoEncoder (GAE) and Variational Graph AutoEncoder (VGAE) (Kipf & Welling, 2016b), these methods share two major limitations:
Shallow Architecture. Previous work shows that although stacking multiple GNN layers in Graph Convolutional Network (GCN) (Kipf & Welling, 2016a) is capable of exploiting deep structural information, applying a large number of GNN layers might lead to indistinguishable node embeddings, i.e., the over-smoothing issue (Li et al., 2018). Therefore, most state-of-the-art GNNs resort to shallow architectures, which hinders the model from capturing long-range dependencies.
Low Scalability. GNN-based graph representation learning methods can not scale well to large graphs due to the expensive computation cost and high memory usage. Most existing GNNs need to repeatedly perform the computationally expensive and recursive feature smoothing, which involves the participation of the entire graph at each training epoch. Furthermore, most methods adopt the same training loss function as GAE, which introduces high memory usage by storing the dense-form adjacency matrix on GPU. For a graph of size 200 million, its dense-form adjacency matrix requires a space of roughly 150GB, exceeding the memory capacity of the current powerful GPU devices.
To tackle these issues, we propose a new graph representation learning method, which is embarrassingly simple: just smooth the node features and then combine the smoothed features in a node-adaptive manner. We name this method node-adaptive feature smoothing (NAFS), and its
goal is to construct better node embeddings that integrate the information from both graph structural information and node features. Based on the observation that different nodes have highly diverse “smoothing speed”, NAFS adaptively smooths each node feature and takes advantage of both loworder and high-order neighborhood information of each node. In addition, feature ensemble is also employed to combine the smoothed features extracted via different smoothing operators. Since NAFS is training-free, it significantly reduces the training cost and scales better to large graphs than most GNN-based graph representation learning methods.
This paper is not meant to diminish the current advancements in GNN-based graph representation learning approaches. Instead, we aim to introduce an easier way to obtain high-quality node embeddings and understand the source of performance gains of these approaches better. Feature smoothing could be a promising direction towards a more simple and effective integration of information from both graph structure and node features.
Our contributions are as follows: (1) New perspective. To the best of our knowledge, we are the first to explore the possibility that simple feature smoothing without any trainable parameters could even outperform state-of-the-art GNNs; this incredible finding opens up a new direction towards efficient and scalable graph representation learning. (2) Novel method. We propose NAFS, a nodeadaptive feature smoothing approach along with various feature ensemble strategies, to fully exploit knowledge from both the graph structure and node features. (3) State-of-the-art performance. We evaluate the effectiveness and efficiency of NAFS on real-world datasets across various graph-based tasks, including node clustering and link prediction. Empirical results demonstrate that NAFS performs comparably with or even outperforms the state-of-the-art GNNs, and achieves up to two orders of magnitude speedup. In particular, on PubMed, NAFS outperforms GAE (Kipf & Welling, 2016b) and AGE (Cui et al., 2020) by a margin of 9.0% and 3.8% in terms of NMI in node clustering, while achieving up to 65.4× and 88.6× training speedups, respectively.

heading: 2 PRELIMINARY
text: In this section, we first explain the notations and problem formulation. Then, we review current GNNs and GNN-based graph representation learning.

heading: 2.1 NOTATIONS AND PROBLEM FORMULATION.
text: In this paper, we consider an undirected graph G = (V ,E) with |V| = n nodes and |E| = m edges. Here we suppose that m ∝ n as it is the case in most real-world graphs. We denote by A the adjacency matrix of G. Each node can possibly have a feature vector of size f , which stacks up to an n× f feature matrix X. The degree matrix of A is denoted as D = diag (d1, d2, · · · , dn) ∈ Rn×n, where di = ∑ vj∈V Aij . We denote the final node embedding matrix as Z, and evaluate it in both the node clustering and the link prediction tasks. The node clustering task requires the model to partition the nodes into c disjoint groups G1, G2, · · · , Gc, where similar nodes should be in the same group. The target of the link prediction task is to predict whether an edge exists between given node pairs.

heading: 2.2 GRAPH CONVOLUTIONAL NETWORK.
text: Based on the assumption that locally connected nodes are likely to enjoy high similarity (McPherson et al., 2001), each node in most GNN models iteratively smooths the representations of its neighbors for better node embedding. Below is the formula of the l-th GCN layer (Kipf & Welling, 2016a):
X(l) = δ ( ÂX(l−1)Θ(l) ) , Â = D̃−1/2ÃD̃−1/2, Ã = A + In, (1)
where X(l) is the node embedding matrix at layer l, X(0) is the original feature matrix, Θ(l) are the trainable weights, and δ is the activation function. Â is the smoothing matrix that helps each node to smooth representations of neighboring nodes.
As shown in Eq. 1, each GCN layer contains two operations: feature aggregation (smoothing) and feature transformation. Figure 1 shows the framework of a two-layer GCN. The l-th layer in GCN firstly executes feature smoothing on the node embedding X(l−1). Then, the smoothed feature X̃(l−1) is transformed with trainable weights Θ(l) and activation function δ to generate new node embedding X(l). Note that GCN will degrade to MLP if feature smoothing is removed from each layer.

heading: 2.3 GNN-BASED GRAPH REPRESENTATION LEARNING.
text: GAE (Kipf & Welling, 2016b), the first and the most representative GNN-based graph embedding method, adopts an encoder to generate node embedding matrix Z with inputs Â and X. A simple inner product decoder is then used to reconstruct the adjacency matrix. The final training loss of GAE is the binary entropy loss between A′ and Ã, the reconstructed adjacency matrix and the original adjacency matrix with self loop added.
A′ = sigmoid(Z · ZT), L = ∑
1≤i,j≤n
−Ãi,j logA′i,j − (1− Ãi,j) log(1−A′i,j)). (2)
Motivated by GAE, lots of GNN-based graph representation learning methods are proposed recently. MGAE (Wang et al., 2017) presents a denoising marginalized autoencoder that reconstructs the node feature matrix X. ARGA (Pan et al., 2018) adopts the adversarial learning strategy, and its generated node embeddings are forced to match a prior distribution. DAEGC (Wang et al., 2019) exploits side information to generate node embeddings in a self-supervised way. AGC (Zhang et al., 2019) proposes an improved filter matrix to better filter out the high-frequency noise. AGE (Cui et al., 2020) further improves AGC by using the similarity of embedding rather than the adjacency matrix to consider original node feature information. Compared with GNN-based graph representation learning methods that rely on the trainable parameters to learn node embeddings, our NAFS is training-free and thus enjoys higher efficiency and scalability.

heading: 3 OBSERVATION AND INSIGHT
text: In this section, we make a quantitative analysis on the over-smoothing issue at the node level and then provide some insights when designing NAFS on graphs.

heading: 3.1 FEATURE SMOOTHING IN DECOUPLED GNNS
text: Recently, many works (Wu et al., 2019; Zhu & Koniusz, 2021; Chen et al., 2020; Zhang et al., 2021) propose to decouple the feature smoothing and feature transformation in each GCN layer for scalable node classification. Concretely, they execute the feature smoothing operation in advance, and the smoothed features are then fed into a simple MLP to generate the final predicted node labels. Under this framework, the predictive node classification accuracy of these methods is comparable with or even higher than the one of coupled GNNs, and these works claim that the true power of GNNs lies in feature smoothing rather than feature transformation.
We split the framework of these decoupled GNNs into two parts: feature smoothing and MLP training. Feature smoothing aims to combine the graph structural information and node features into better features for the subsequent MLP; while MLP training only takes in the smoothed feature and is specially trained for a given task. As stated by previous decoupled GNNs (Wu et al., 2019; Zhu & Koniusz, 2021), the true success of GNNs lies in feature smoothing rather than feature transformation. Correspondingly, we propose to remove feature transformation and preserve the key feature smoothing part alone for simple and scalable node representation.
There is another branch of GNNs that also decouple the feature smoothing and feature transformation. The most representative method of this category is APPNP (Klicpera et al., 2018). It first feeds the raw node features into an MLP to generate intermediate node embeddings; then the personalized PageRank based propagation operations are performed on the node embeddings to produce final
prediction results. However, compared with scalable decoupled GNNs mentioned in the previous paragraph, this branch of GNNs still have to recursively execute propagation operations in each training epoch, which makes it impossible to perform on large-scale graphs. In the remaining part of this paper, the terminology “decoupled GNNs” refers particularly to the scalable decoupled GNNs mentioned in the previous two paragraphs.

heading: 3.2 MEASURING SMOOTHING LEVEL
text: To capture deep graph structural information, a straightforward way is to simply stack multiple GNN layers. However, a large number of feature smoothing operations in a GNN model would lead to indistinguishable node embeddings, i.e., the over-smoothing issue (Li et al., 2018). Concretely, if we execute ÂX for infinite times, the node embeddings within the same connected component would reach a stationary state. When adopting Â = D̃r−1ÃD̃−r, Â∞ follows
Â∞i,j = (di + 1)
r(dj + 1) 1−r
2m+ n , (3)
which shows that the influence from node vi to vj is only determined by their degrees. Under the extreme condition that r = 0, all the nodes within one connected component have exactly the same representation, making it impossible to apply the node embeddings to subsequent tasks.
Here we introduce a new metric, “Over-smoothing Distance”, to measure each node’s smoothing level. A smaller value indicates that the node is closer to the stationary state, i.e., closer to over-smoothing.
Definition 3.1 (Over-smoothing Distance). The Over-smoothing Distance Di(k) parameterized by node i and smoothing step k is defined as
Di(k) = Dis([Â kX]i, [Â ∞X]i), (4)
where [ÂkX]i denotes the ith row of ÂkX, representing the representations of node vi after smoothing k times; [Â∞X]i denotes the ith row of Â∞X, representing the stationary state of node vi; Dis(·) is a distance function or a function positively relative with the difference, which can be implemented using Euclidean distance, the inverse of cosine similarity, etc.

heading: 3.3 DIVERSE SMOOTHING SPEED ACROSS NODES.
text: To examine the factors that affect Di(k), we divide nodes in the PubMed dataset into three different groups according to their degrees. In Figure 2, we show the trends of nodes in the first group and the second group where the group-averaged Di(k) changes with the number of smoothing step increases. The trend of Di(k) averaged over all the nodes is also provided, and the Euclidean distance is chosen as the distance function. Figure 2 shows that the Di(k) of nodes with degrees larger than 60 drops more quickly than nodes with degrees smaller than 3, which implies that nodes with high degrees approach the stationary state more rapidly than the nodes with low degrees.
The quantitative analysis empirically illustrates that the degree of each node plays an essential role in one’s optimal smoothing step. Intuitively, nodes with high degrees should have relatively small smoothing steps than nodes with low degrees. In addition, in Appendix A, we have made a detailed theoretical analysis about the graph sparsity, another factor that influences the smoothing speed.

heading: 3.4 DESIGN INSIGHTS OF NAFS
text: Though using feature smoothing operations inside the previous decoupled GNNs is scalable for large graph representation learning, it will lead to sub-optimal node representation. According to the observation in Sec. 3.3, it is sub-optimal to execute feature smoothing for all the nodes indiscriminately as previous decoupled GNNs do since nodes with different structural properties have diverse smoothing speeds. Therefore, node-adaptive feature smoothing (i.e., different nodes have different smoothing levels or mechanisms with equivalent effect) must be adopted to satisfy each node’s diverse needs of smoothing level. In our proposed method NAFS, we utilize the metric Di(k) defined in Def. 3.1 to help accomplish the aim of node-adaptive feature smoothing.

heading: 4 PROPOSED METHOD
text: In this section, we present NAFS, a training-free method for scalable graph representation learning. We first compute the smoothed features with the feature smoothing operation. Then the feature ensemble operation is used to combine the smoothed features generated by different smoothing strategies. The pseudo code of NAFS is provided in Appendix D.

heading: 4.1 NODE-ADAPTIVE FEATURE SMOOTHING
text: Figure 3 provides an overview of NAFS. When repeatedly executing X(l) = ÂX(l−1), the smoothed node embedding matrix X(l−1) contains deeper graph structural information with l increases. The multi-scale node embedding matrices {X(0),X(1), ...,X(K)} (K is the maximal smoothing step) are then combined into a single matrix X̂ such that both local and global neighborhood information are reserved for each node.
The analysis in Sec. 3.3 illustrates that the speed each node achieving its stationary state is highly diverse, which suggests that nodes should be treated individually. To this end, we define “Smoothing Weight” based on Di(k) introduced in Def. 3.1 for each node so that the smoothing operation can be performed in a node-adaptive manner. Definition 4.1 (Smoothing Weight). The Smoothing Weight wi(k) parameterized by node vi and smoothing step k is defined based on the softmax value of {Di(0), Di(1), · · · , Di(K)}:
wi(k) = e Di(k)/ K∑ l=0 eDi(l), (5)
where K is the maximal smoothing step.
To calculate Di(k) more efficiently, an alternative is to replace [Â∞X]i in Eq. 4 with Xi and implement Dis(·) as the cosine similarity. Larger Di(k) in this case means that node vi is farther from the stationary state and [ÂkX]i intuitively contains more relevant node information. Therefore, for node vi, the smoothed feature with larger Di(k) (i.e., larger wi(k)) should contribute more to the final node embedding. The smoothing weight can be formulated in the following matrix form: Definition 4.2 (Smoothing Weight Matrix). The Smoothing Weight Matrix W(k) parameterized by smoothing step k is defined as the diagonal matrix derived from η(k) ∈ Rn:
W(k) = Diag(η(k)), η(k)[i] = wi(k), ∀1 ≤ i ≤ n. (6)
Given the maximal smoothing step K, the multi-scale smoothed features {X(0),X(1), ...,X(K)}, and the corresponding Smoothing Weight Matrices {W(0),W(1), ...,W(K)} , the final smoothed
feature X̂ can be represented as: X̂ = K∑ k=0 W(k)ÂkX.

heading: 4.2 FEATURE ENSEMBLE
text: Different smoothing operators actually act as different knowledge extractors. For example, by setting r = 0.5, 1 and 0, Â = D̃r−1ÃD̃−r represents the symmetric normalized adjacency matrix
D̃−1/2ÃD̃−1/2 (Kipf & Welling, 2016a), the random walk transition probability matrix D̃−1Ã (Xu et al., 2018), and the reverse random walk transition probability matrix ÃD̃−1 (Zeng et al., 2019), respectively. These variants of Â captures and reserves different scales and dimensions of knowledge from both graph structures and node features.
To achieve the same effect, the feature ensemble operation has multiple knowledge extractors: we vary the value of r in the normalized adjacency matrix Âr = D̃r−1ÃD̃−r to easily acquire different knowledge extractors. These knowledge extractors are adopted inside the feature smoothing operation to generate different smoothed features. Concretely, the value of r controls the normalized weight of each edge. So different values of r generate different weight values for all the edges in the graph, which would increase the diversity of our smoothed features evidently. The ablation study in Sec. 6.5 shows that the increased diversity contributes a lot to the high performance of our generated node embeddings when applied to downstream tasks.
The detailed procedure of feature ensemble operation is as follows: Given {r1, r2, ..., rT }, we firstly perform the feature smoothing operation to generate corresponding smoothed features {X̂(1), X̂(2), ..., X̂(T )}. Then, we combine them as Z← ⊕i∈{1,2,...,T}X̂(i), where⊕ is the ensemble strategy, which can be implemented as concatenating, mean pooling, max pooling, etc.

heading: 5 ADVANTAGES OVER TRADITIONAL APPROACHES
text: NAFS generates node embeddings in a training-free manner, making it highly efficient and scalable. Moreover, the node-adaptive smoothing strategy enables it to capture deep structural information. In this subsection, we analyze the advantages of our NAFS over GAE and its variants.
Deep Structural Information. By assigning each node with personalized smoothing weights, NAFS can gather deep structural information without encountering the over-smoothing issue and keep the time and memory cost low. While for GAE and its variants, they either 1) have a coupled structure that cannot go deep due to low efficiency and high memory cost (e.g., GAE) or 2) are unable to adaptively capture structural information and encounter the over-smoothing issue when going deeper (e.g., AGE).
Efficiency. Compared with GAE and its variants, our proposed NAFS does not have any trainable parameters, giving it a significant advantage in efficiency. When generating node embeddings, NAFS only needs to execute the feature smoothing and feature ensemble operations, which has a time complexity of O(Kmf +Kn), where K denotes the maximal smoothing step. On the other hand, the asymptotic training time complexity of GAE and its variants is O(nf2 + n2f), which is one magnitude larger than NAFS. In Sec. 6.4, we further present efficiency comparison in detail on real-world datasets between NAFS and GAE.
Memory Cost. Our proposed NAFS enjoys not only high efficiency but also low memory cost. NAFS only requires to store the sparse adjacency matrix and the smoothed features, and thus the memory cost is O(m+nf), which grows linearly with graph size n in typical real-world graphs. For GAE and its variants, they need additional space to store the training parameters and especially the reconstructed dense adjacency matrix, making their memory cost O(n2 +m+ nf), which is also one magnitude larger than NAFS. Scalability comparison conducted on synthetic datasets between NAFS and GAE can be found in Appendix C.2.

heading: 6 EXPERIMENTAL RESULTS
text: In this section, we conduct extensive experiments to evaluate the proposed NAFS. We first introduce the considered baseline methods, used datasets, and experimental settings. Then, we demonstrate the advantages of NAFS from the following three perspectives: (1) end-to-end comparison with the state-of-the-art methods, (2) scalability and efficiency, and (3) effectiveness along with ablation studies.

heading: 6.1 DATASETS AND BASELINE METHODS
text: Four widely-used network datasets (i.e., Cora, Citeseer, PubMed, and Wiki) are used in our experiments. In Cora and Citeseer, the node features are binary word vectors; while in PubMed and Wiki, each node has a TF-IDF weighted word vector. We include the properties of these datasets in Appendix B.1.
For different downstream tasks, the considered baselines are as follows:
• Node clustering: GAE and VGAE (Kipf & Welling, 2016b), MGAE (Wang et al., 2017), ARGA and ARVGA (Pan et al., 2018), AGC (Zhang et al., 2019), DAEGC (Wang et al., 2019), and AGE (Cui et al., 2020).
• Link prediction: Spectral Clustering (SC) (Ng et al., 2002), DeepWalk (Perozzi et al., 2014), GAE and VGAE (Kipf & Welling, 2016b), ARGA and ARVGA (Pan et al., 2018), GALA (Park et al., 2019), and AGE (Cui et al., 2020).
For NAFS, we investigate four variants: NAFS-simple, NAFS-mean, NAFS-max, and NAFS-concat. NAFS-simple only contains the feature smoothing operation. Besides, it only uses the smoothed features at maximal smoothing step K with r = 0.5 and excludes the node-adaptive strategy. And for NAFS-mean, -max, and -concat, they all have the complete NAFS framework and only differ in the ensemble strategy adopted in feature ensemble operation.

heading: 6.2 EXPERIMENTAL SETTINGS
text: Node Clustering. For the node clustering task, we apply K-Means (Hartigan & Wong, 1979) to node embeddings to get the clustering results. Three widely-used metrics are used for evaluation: Accuracy (ACC), Normalized Mutual Information (NMI), and Adjusted Rand Index (ARI).
Link Prediction. For the link prediction task, 5% and 10% edges are randomly reserved for the validation set and the test set. Once the node embeddings have been generated, we follow the same procedure as GAE to reconstruct the adjacency matrix. Two metrics - Area Under Curve (AUC) and Average Precision (AP) are used in the evaluation of the link prediction task.
Hyperparameters. When generating node embeddings, we use the values of r in [0, 0.1, 0.2, 0.3, 0.4, 0.5] to get six different normalized adjacency matrix Â. The only exception is the link prediction task on the PubMed dataset, where we use values of r in [0.3, 0.4, 0.5] instead. The optimal value of
Figure 4: Efficiency comparison on citation networks.
Figure 5: Ablation study.
the maximal smoothing steps ranges from 1 to 70. Hyperparameters for all the baseline methods are tuned following the settings in their original paper. We run all the methods using 200 epochs and repeat 10 times on all the datasets, and report the mean value of each evaluation metric. The details of the experimental environment setting are in Appendix B.2.

heading: 6.3 END-TO-END COMPARISON
text: Link Prediction. Table 1 shows the performance of different methods on the link prediction task. On the three datasets, the proposed NAFS consistently achieves the best or the second-best performance compared with all the baseline methods. Remarkably, although NAFS-simple falls behind the other NAFS variants on Cora and Citeseer, it achieves state-of-the-art performance on PubMed, and outperforms the current state-of-the-art method - ARGA by 0.6% and 0.1% on AUC and AP, respectively.
Node Clustering. The node clustering results of each method are shown in Table 2. Among three ensemble strategies, NAFS-concat has the overall best performance across the three datasets, which also consistently outperforms the strongest baseline - AGE. For example, the best of NAFS variants exceeds AGE by 2.6%, 1.8%, 0.7%, and 3.5% on Cora, Citeseer, PubMed, and Wiki, respectively. Although NAFS-concat has the overall best performance, it falls behind NAFS-mean or NAFS-max on some metrics of some datasets; and it needs more CPU memory to store the node embeddings and longer time for inference.
It is quite interesting to find that the simplified version of NAFS, NAFS-simple, outperforms trainingbased GAE on all the four datasets, and it even outperforms the current SOTA method, AGE, on Citeseer. The competitive performance of NAFS-simple further illustrates that it is possible to achieve decent performance on graphs with only feature smoothing without any training parameters.

heading: 6.4 EFFICIENCY AND SCALABILITY ANALYSIS
text: To validate the efficiency and scalability of NAFS, we compare it with GAE and AGE, measuring their overall running time for generating node embeddings. The comparison is conducted on the three citation networks. For fairness, all methods only use CPUs for computation. Figure 4 showcases the results along with the speedup ratio of NAFS against AGE.
Results from Figure 4 show that NAFS is significantly faster than the considered methods on the three datasets. For example, NAFS is almost two magnitudes faster than AGE on the relatively large dataset - PubMed. The high efficiency of NAFS is attributed to no trainable parameters, while GAE and AGE are both training-based methods. Although AGE executes the feature smoothing step in advance, it adopts a time-consuming ranking policy to select positive and negative examples. This alteration in AGE brings both superior performance and longer running time compared with GAE. In addition, we also include more details in Appendix C.2 to demonstrate the excellent scalability of NAFS (See Figure 6).

heading: 6.5 ABLATION STUDY
text: To thoroughly investigate the proposed NAFS, ablation studies on the node clustering task are designed to analyze the effectiveness of feature smoothing and feature ensemble in NAFS. The experiments are conducted on the PubMed dataset, and Normalized Mutual Information (NMI) is used to measure the performance. To interpret its effectiveness better, we also visualize the node embedding in Appendix C.4.
Different Weighting Strategies. In NAFS, we use node-adaptive weight to average the smoothed features of different steps. Here we change the “Adaptive Weight” in our method to “Single Hop” (only the smoothed feature at the last smoothing step is reserved) and “Naive Average” (smoothed features at every smoothing step has equal weight) and evaluate their performance. Figure 5(a) shows the performance of the three different weighting strategies.
From Figure 5(a), the weighting strategy “Single Hop” performs the worst among the three, since it only makes use of the information at the last smoothing step, which could lead to the over-smoothing issue when the maximal smoothing step becomes large. On the other hand, the weighting strategy “Adaptive Weight” shows a better performance than “Naive Average”. Besides, when the number of maximal smoothing steps becomes large, the performance of “Naive Average” begins to drop, while “Adaptive Weight” does not. “Naive Average” assigns the same weight across all the different steps of smoothed features, which also leads to over-smoothing when it tries to exploit extremely deep structural information. Instead, by assigning adaptive weights, the “Adaptive Weight” strategy in NAFS could exploit such deep information and avoid the over-smoothing issue, thus improving the quality of the generated node embedding.
Different Ensemble Strategies. We use different ensemble strategies to obtain the final node embeddings in our proposed method - NAFS, which include “Mean”, “Max”, and “Concat”. To evaluate the impacts of these different ensemble strategies, we change the maximal smoothing step, K, from 15 to 40, and evaluate corresponding node clustering performance. For reference, the performance of NAFS without feature ensemble, “r=0.3 only” is also reported. The experimental results in Figure 5(b) illustrate that at most times, “r=0.3 only” performs worse than others, which shows the effectiveness of the three different ensemble strategies. However, if we limit the comparison within the ensemble strategies, the performance superiority is unclear. In Table 2 and 1 we can also have that different ensemble strategies perform diversely across different datasets and tasks. It indicates that these three ensemble strategies all have necessities in their own way.

heading: 7 CONCLUSION
text: This paper presents NAFS, a novel graph representation learning method. Unlike other GNNbased approaches, NAFS focuses on improving the feature smoothing operation in the GNN layer and generate node embeddings in a training-free manner. NAFS proposes Node-Adaptive Feature Smoothing to generate smoothed features with adaptivity to each node’s individual properties; it further employs feature ensemble to combine multiple smoothed features from diverse knowledge extractors effectively. Experiments results on typical tasks demonstrate that NAFS performs comparably with or even outperforms the state-of-the-art GNNs, and at the same time enjoys high efficiency and scalability where NAFS shows its absolute superiority.

heading: 8 REPRODUCIBILITY STATEMENT
text: The source code of NAFS can be found in Anonymous Github (https://anonymous.4open. science/r/NAFS). To ensure reproducibility, we have provided the overview of datasets and baselines in Section 6.1 and Appendix B.1. The settings for each task and baseline can be found in Section 6.2. Our experimental environment is presented in Appendix B.2, and please refer to “README.md” in the Github repository for more reproduction details.

heading: A THEORETICAL ANALYSIS
text: The essential kernel of NAFS is the Smoothing Weight, which determines the output results. We now analyze the factors affecting the value of Smoothing Weight. To simplify our analysis, we suppose r = 0 in the normalized adjacency matrix and apply Euclidean distance as the distance function in Definition 3.1. Thus we have
Di(k) = ||[ÂkX]i − [Â∞X]i||2, (7)
where || · ||2 symbols two-norm. Theorem A.1. For any node i in graph G, there always exists
Di(k) ≤ λk2 √√√√√ n∑j=1(d̃j ||Xj ||22) d̃i , (8)
where d̃i = di + 1, d̃j = dj + 1, ||Xj ||2 denotes the two-norm of the initial feature of node j, and 0 < λ2 < 1 denotes the second largest eigenvalue of the normalized adjacency matrix Â.
To prove Theorem A.1, we introduce the following lemma.
Lemma A.2.
|(eiÂk)j − (eiÂ∞)j | ≤
√ d̃j
d̃i λk2 , (9)
where ei denotes a one-hot row vector with its ith components as 1 and other components as 0, λ2 is the second largest eigenvalue of Â and d̃i denotes the degree of node i plus 1 (to include itself).
d̃i = di + 1, d̃j = dj + 1,
The proof of Lemma A.2 can be found in (Chung & Graham, 1997). Next we will prove Theorem A.1.
Proof. According to equation 7, we can have that
Di(k) = ||[ÂkX]i − [Â∞X]i||2 = ||(eiÂk − eiÂ∞)X||2
= √√√√ n∑ j=1 ((eiÂk)j − (eiÂ∞)j)2X2j
≤ √√√√√√ λ2k2 n∑ j=1 d̃j f∑ p=1 X2jp
d̃i = λk2 √√√√√ n∑j=1(d̃j ||Xj ||22) d̃i ,
(10)
where Xjp denotes the pth feature of node j.
Based on Lemma A.2 we then consider the smoothing distance for weighted averaged embedding features to the station state. ||[ÂkX]i − [Â∞X]i||2 (11) Therefore there holds Theorem A.1
Di(k) ≤ λk2 √√√√√ n∑j=1(d̃j ||Xj ||22) d̃i . (12)
We denote the constant n∑ j=1 (d̃j ||Xj ||22) as cdx because it is independent with j, then Theorem A.1 can be written as
Di(k) ≤ λk2
√ cdx
d̃i . (13)
We then analyze the factors affecting the Smoothing Weight on a specific node vi. From Eq. 13 we know that the nodes with smaller degrees may have larger Di(k). Combined with definition 4.1, we infer that larger Di(k) makes the maxkDi(k) more dominant after the softmax operation, causing that weighted average results depend more on itself and its near neighbors. Inversely, for the nodes with smaller degrees, its result of weighted average depends more equally on all itself, its near neighbors, and its distant neighbors.
At the same time, the Smoothing Weight of the node in a sparser graph (λ2 is positively relative with the sparsity of a graph) decays slower as k increases. Thus, the weighted average result depends more equally on itself, its near neighbors and its distant neighbors. While for the nodes in a denser graph, the weighted average result depends more on its near neighbors and itself.
Prevent over-smoothing. Based on Theorem A.1 we then consider the smoothing distance for weighted averaged embedding features to the station state:
|| K∑ k=0 ωi(k) ∗ [ÂkX]i − [Â∞X]i||2 = || K∑ k=0 ωi(k) ∗ ([ÂkX]i − [Â∞X]i)||2
≤ K∑ k=0 ωi(k) ∗ ||[ÂkX]i − [Â∞X]i||2
= K∑ k=0 ωi(k) ∗Di(k).
(14)
When ωi(k) = 1/(K + 1), which means the average option is non-weighted, we have
lim K→∞ K∑ k=0 ωi(k) ∗Di(k) = lim K→∞ 1 K + 1 K∑ k=0 Di(k)
≤ lim K→∞
1
K + 1 K∑ k=0 λk2
√ cdx
d̃i
= lim K→∞
1
K + 1 1− λK+12 1− λ2
√ cdx
d̃i
= 0,
(15)
causing over-smoothing. When ωi(k) = Di(k)/( ∑K l=0Di(l)), let X k i = [Â
kX]i − [Â∞X]i, we have:
||X̂i − [Â∞X]i||2 = lim K→∞ ωi(0)||X0i + K∑ k=1 ωi(k) ωi(0) Xki ||2 (16)
In the real-world graphs, nodes have different initial features, thus there is little chance that the combination of a node’s neighboring features both lies in the opposite direction and is of the same norm of the node’s initial feature. Suppose that there exits a constant > 0 satisfying min ( Di(0), ||X0i + ∑K k=1 ωi(k) ωi(0) Xki ||2 ) ≥ for node i, we have:
||X̂i − [Â∞X]i||2 ≥ 2 limK→∞ ∑K k=0Di(k) ≥ 2
1 1−λ2 √ cdx d̃i > 0. (17)
Thus we see that even K goes to infinity, we are able to prevent the node representations from reaching the stationary state (the distance bound depends on the node degree di and the initial feature
Xi). Note that in practice the representation at k-th smoothing step X(k) achieves the stationary state much earlier than the infinity-th smoothing step we use in the theoretical analysis, so we use the softmax normalization in Equation 5 to produce a slightly larger bias towards features with longer distances to the stationary state.
Node-adaptive weighting. The above analysis theoretically proved our weighting scheme is able to prevent over-smoothing as the smoothing step goes to infinity. We now show that another advantage of our method is that it can fully leverage the multiple features over different smooth steps in a node-adaptive manner, which is different from the traditional routine of the smooth-blind and fixed scheme. Suppose that the feature of k step is not over-smoothed yet for a specific node i, i.e., the distance of Di(k) = i > 0, we have:
ωi(k) = i limK→∞ ∑K k=0Di(k) ≥ i 1
1−λ2 √ cdx d̃i > 0 . (18)
We see that as long as the feature is not over-smoothed, our method will assign a non-zero weight to the feature. Further, we see that the bound of weight can be affected by the over-smoothing distance of i and degree di of a specific node, implying an adaptive weighting strategy.

heading: B DETAILS ON THE EXPERIMENTS
text: B.1 DATASETS DESCRIPTION
Cora, Citeseer, PubMed, and Wiki are four popular network datasets. The first three (Kipf & Welling, 2016b) are citation networks where nodes stand for research papers and an edge exists between a node pair if one cites the other. Wiki (Yang et al., 2015) is a webpage network where nodes stand for webpages and an edge exists between a node pair if one links the other. Table 3 presents an overview of these four datasets.
B.2 EXPERIMENTAL ENVIRONMENT
The experiments are conducted on a machine with Intel(R) Xeon(R) Gold 5120 CPU @ 2.20GHz, and a single NVIDIA TITAN RTX GPU with 24GB GPU memory. The operating system of the machine is Ubuntu 16.04. For software versions, we use Python 3.6, Pytorch 1.7.1, and CUDA 10.1.

heading: C ADDITIONAL EMPIRICAL RESULTS
text: C.1 PERFORMANCE ON THE NODE CLASSIFICATION TASK
In this subsection, we further assess the quality of the node embeddings generated by NAFS by the evaluation on the node classification task. We follow the linear evaluation protocol, which applies a linear classifier (i.e., Logistic Regression) to the node embeddings to generate final prediction results. We choose GCN (Kipf & Welling, 2016a), JK-Net (Xu et al., 2018), C&S (Huang et al., 2020), SGC (Wu et al., 2019), GAT (Veličković et al., 2017), PPRGo (Bojchevski et al., 2020), APPNP (Klicpera et al., 2018), and DAGNN (Liu et al., 2020) as comparison baselines on the node classification task. Note that C&S in our evaluation adopts MLP as the base model. The evaluation results on three popular datasets, Cora, Citeseer, and PubMed, are provided in Table 4.
The table shows that our method achieves comparable predictive accuracy as one of the state-of-the-art methods DAGNN. However, since the node embeddings are generated before training, our method
avoids performing recursive feature smoothing at each training epoch and storing the entire adjacency matrix on GPU. In this way, our method is more scalable and efficient to apply on large graphs than most state-of-the-art methods like DAGNN.
C.2 SCALABILITY COMPARISON ON SYNTHETIC GRAPHS
A challenging task in today’s graph learning is how to get decent node embedding on large graphs. To test the scalability of our proposed NAFS, we use the Erdős-Rényi graph generator in the Python package NetworkX (Hagberg et al., 2008) to generate artificial graphs of different sizes. The node sizes of the generated artificial graphs vary from 5,000 to 100,000, and the probability of an edge exists between two nodes is set to 0.0001. The experiment settings are altered compared to Sec. 6.4: both GAE and AGE are trained on an NVIDIA TITAN RTX which has 24 GB of memory, since only using CPU to train being unacceptable on large graphs in real-world applications.
The overall experiment results are shown in Figure 6. Figure 6(a) shows a more detailed comparison on relatively small artificial graphs whose node sizes range from 5,000 to 20,000. Besides, Figure 6(b) shows that GAE encounters the out-of-memory problem on the artificial graph composed of 30,000 nodes, while AGE encounters the out-of-memory problem at graph size 80,000. Our proposed NAFS is successfully carried out on all the artificial graphs. On the artificial graph consists of 100,000 nodes, NAFS accomplishes the node embedding generation in 66.1 seconds, which is less than the running time of AGE on the artificial graph consists of 10,000 nodes, 80.9 seconds. The graph size limit of our proposed NAFS is bounded by the CPU memory size. As long as one can successfully execute the multiplication of the sparse adjacency matrix and the feature matrix, NAFS can then be implemented.
C.3 SCALABILITY COMPARISON ON REAL-WORLD GRAPHS.
To verify the scalability advantage of NAFS, we also add the evaluation on real-world datasets. Concretely, we sample subgraphs of different scales from the full ogbn-products graphs via uniformly random node selection and then report the sampled graph size (i.e., number of nodes) and the corresponding runtime (seconds) in Table 5. The experimental results show that NAFS can support
larger graphs (i.e., larger than 30,000 nodes). Besides, it is significantly faster than the compared baselines, especially for large graph datasets.
C.4 VISUALIZATION OF EMBEDDING
To better understand why NAFS is effective, we visualize the node embedding generated by NAFS using T-SNE (Van der Maaten & Hinton, 2008) on the Cora dataset. Moreover, we also visualize the node embedding of ÂkX for reference. All the visualization results are shown in Figure 7. The first row is the results of our proposed NAFS, and the second row is the results of ÂkX.
Figure 7(c) and 7(k) illustrate that at hop 10, the node embedding produced by NAFS and ÂkX are both distinguishable. But as the value of maximal smoothing step becomes larger, the node embedding of ÂkX falls into total disorder, like the situation showed by Figure. 7(n), 7(o) and 7(p). At the same time, NAFS is able to maintain the distinguishable results even when the value of the maximal number of smoothing step increases to 200.

heading: D MORE DETAILS OF NAFS
text: D.1 PSEUDO CODE OF NAFS
Alg. 1 shows the whole pipeline of our proposed NAFS. We first initialize X(0) as the original feature matrix X. Given the normalization parameter rt, we obtained the corresponding normalized adjacency matrix Ârt = D̃
rt−1ÃD̃−rt , which acts as knowledge extractors (line 4). After that, for each node, we use E.q. 4 to calculate the Over-smoothing Distance with all the k ranging from 0 to K (line 6, 7). Then we calculate its Aggregation Weights through Eq. 5 (line 8, 9). After obtaining Aggregation Weights with all k and i, we construct the Aggregation Weight matrix for each k with
Algorithm 1: NAFS pipeline. Input: Smoothing step K, feature matrix X, adjacency matrix A, and {r1, r2, ..., rT }. Output: Graph embedding matrix Z.
1 Initialize the feature matrix X(0) = X; 2 Operation 1: Feature Smoothing 3 for 1 ≤ t ≤ T do 4 Update the normalized adjacency matrix Ârt = D̃
rt−1ÃD̃−rt ; 5 for 1 ≤ i ≤ n do 6 for 0 ≤ k ≤ K do 7 Calculate Di(k) and wi(k) with Eq. 4 and 5, respectively; 8 for 0 ≤ k ≤ K do 9 Construct W(k) with Eq. 6;
10 Smooth the node features X with X̂(t) = K∑ k=0 W(k)ÂkrtX; 11 Operation 2: Feature Ensemble 12 Compute the final embedding Z with Z← ⊕i∈{1,2,...,T}X̂(i).
Eq. 6 (line 10, 11). Next, we compute the NAFS output X̂(t) with X̂(t) = K∑ k=0 W(k)ÂkrtX (line 12). Finally, we compute the final embedding result of all t through Z← ⊕i∈{1,2,...,T}X̂(i) (line 14).
D.2 MOTIVATION OF FEATURE ENSEMBLE
Adopting different smoothing operators in the feature smoothing operation (the normalized adjacency matrix Â in Eq. 1) is equivalent to smoothing features in different manners. Other than the normalized adjacency matrix Â, there are many alternatives that have been proposed recently. For example, GraphSAGE (Hamilton et al., 2017) designs three smoothing operators (i.e., Mean, LSTM and Pooling) to flexibly capture the information of neighboring nodes. SIGN (Frasca et al., 2020) enriches the smoothing operators with Personalized-PageRank-based (Klicpera et al., 2019) and trianglebased (Monti et al., 2018) adjacency matrices. However, these methods are not designed for graph representation learning and different smoothing operators may results in diverse smoothed features. For better node representation, feature ensemble is used to combine the smoothed feature under different smoothing operators.
D.3 RELATIONS WITH OTHER SCALABLE GNN ARCHITECTURES
The scalable GNNs can be roughly classified into two categories: (a) "first propagate then predict"; (b) sampling + ordinary GNN. The first category includes famous GNN models like SGC (Wu et al., 2019) and SIGN (Frasca et al., 2020). They disentangles the coupled propagation and transformation operations in traditional GCN layers (Kipf & Welling, 2016a), and execute all the propagation operations as preprocessing. The most representative GNN model of the second category is GraphSAGE (Hamilton et al., 2017), and many other works have been proposed following it, like FastGCN (Chen et al., 2018) and GraphSAINT (Zeng et al., 2019). The main contribution of these works are effective sampling strategies that can preserve the most valuable information from the original graph. The sampling strategies always act as a plug-and-play module, and can be combined with ordinary GNNs, which allows the latter to perform on large-scale graphs. GNNs belong to the "first propagate then predict" category usually enjoy higher efficiency since they avoid recursively performing propagation in each training epoch. Our proposed NAFS belongs to the "first propagate then predict" category.
D.4 ADVANTAGE IN DISTRIBUTED SETTINGS
NAFS can also be adapted to the distributed environment. The feature smoothing process in NAFS are sparse matrix dense matrix multiplications, which have mature implementations in distributed environments. Besides, this process only need to be pre-computed at once. In contrast, during
each training iteration of GAE and its variants, each node needs to repeatedly pull the intermediate representations of other nodes, which leads to high communication cost.

Review summaries:
The paper presents NAFS (Node-Adaptive Feature Smoothing), a method that constructs node representations by relying on smoothing only, i.e. without parameter learning. To do this, the authors first provide a formulation for the smoothing operator after infinite steps, i.e. when the stationary state is reached. They then define over-smoothing distance as a way to assess how much a node is close to the stationary state after k smoothing steps. Finally, they use over-smoothing distances to calculate a different smoothing weight for each node.  Experiments show that representations obtained by smoothing with these weights, together with feature ensembles obtained by applying different convolution coefficients, provide performances in clustering and link prediction tasks which are comparable, if not better, than many other state-of-the-art approaches.

The authors of this paper took a novel perspective to present the node-adaptive feature smoothing (NAFS) algorithm, which generates node embeddings without explicit training/parameter learning. The method first performs feature smoothing, then combines the smoothed features using adaptive weights which are node-specific. They further enhanced this method by ensembling the smoothed features extracted with different hyper-parameters. The authors have conducted many experiments to validate the model performance, and demonstrate the model's efficiency empirically. 

The paper deals with (unsupervised) learning with graphs, specifically node-level tasks. Inspired by spectral GNNs, specifically Graph Convolutional Networks (GCN), the authors propose a simple neighborhood smoothing technique to capture the graph structure around each node in the given graph. Contrary to GNN, the proposed algorithm is parameter-free and hence scales better than their end-to-end trained GNN counterpart.                      Motivated by the problem of over-smoothing of GCN (Li et al., 2018), they propose the so-called "over-smoothing distance" measuring how close a node's feature is to be over-smoothed, which is simply the row-wise distance between $A^k X$ and $A^{\infty}X$ with regard to the Euclidian distance, see Def. 3.1. Based on this distance they define the smoothing weight matrix which is then used to weight the neighboring node features during neighborhood aggregation. The resulting features $\hat{X}$ are computed as $\hat{X} =  \sum^{K}_{k=0} W(k)\hat{A}X$, where $W(k)$ is the $k$th smoothing weight matrix.  The proposed architecture is evaluated on standard, old, small-scale (unsupervsied) link and node classification (Cora, Citeseer, PubMed) tasks showing somewhat classification performance while showing a considerable speedup in computation time.