p-LAPLACIAN BASED GRAPH NEURAL NETWORKS
Year: 2021

Abstract:
Graph neural networks (GNNs) have demonstrated superior performance for semisupervised node classification on graphs, as a result of their ability to exploit node features and topological information simultaneously. However, most GNNs implicitly assume that the labels of nodes and their neighbors in a graph are the same or consistent, which does not hold in heterophilic graphs, where the labels of linked nodes are likely to differ. Moreover, when the topology is non-informative for label prediction, ordinary GNNs may work significantly worse than simply applying multi-layer perceptrons (MLPs) on each node. To tackle the above problem, we propose a new p-Laplacian based GNN model, termed as GNN, whose message passing mechanism is derived from a discrete regularization framework and could be theoretically explained as an approximation of a polynomial graph filter defined on the spectral domain of p-Laplacians. The spectral analysis shows that the new message passing mechanism works simultaneously as low-pass and high-pass filters, thus making GNNs are effective on both homophilic and heterophilic graphs. Empirical studies on real-world and synthetic datasets validate our findings and demonstrate that GNNs significantly outperform several stateof-the-art GNN architectures on heterophilic benchmarks while achieving competitive performance on homophilic benchmarks. Moreover, GNNs can adaptively learn aggregation weights and are robust to noisy edges.

heading: 1 INTRODUCTION
text: In this paper, we explore the usage of graph neural networks (GNNs) for semi-supervised node classification on graphs, especially when the graphs admit strong heterophily or noisy edges. Semisupervised learning problems on graphs are ubiquitous in a lot of real-world scenarios, such as user classification in social media (Kipf & Welling, 2017), protein classification in biology (Velickovic et al., 2018), molecular property prediction in chemistry (Duvenaud et al., 2015), and many others (Marcheggiani & Titov, 2017; Satorras & Estrach, 2018). Recently, GNNs are becoming the de facto choice for processing graph structured data. They can exploit the node features and the graph topology by propagating and transforming the features over the topology in each layer and thereby learn refined node representations. A series of GNN architectures have been proposed, including graph convolutional networks (Bruna et al., 2014; Henaff et al., 2015; Defferrard et al., 2016; Kipf & Welling, 2017; Wu et al., 2019), graph attention networks (Velickovic et al., 2018; Thekumparampil et al., 2018), and other representatives (Hamilton et al., 2017; Xu et al., 2018; Pei et al., 2020).
Most of the existing GNN architectures work under the homophily assumption, i.e. the labels of nodes and their neighbors in a graph are the same or consistent, which is also commonly used in graph clustering (Bach & Jordan, 2004; von Luxburg, 2007; Liu & Han, 2013) and semi-supervised learning on graphs (Belkin et al., 2004; Hein, 2006; Nadler et al., 2009). However, recent studies (Zhu et al., 2020; 2021; Chien et al., 2021) show that in contrast to their success on homophilic graphs, most GNNs fail to work well on heterophilic graphs, in which linked nodes are more likely to have distinct labels. Moreover, GNNs could even fail on graphs where their topology is not helpful for label prediction. In these cases, propagating and transforming node features over the graph topology could lead to worse performance than simply applying multi-layer perceptrons (MLPs) on each of the nodes independently. Several recent works were proposed to deal with the heterophily issues of GNNs. Zhu et al. (2020) finds that heuristically combining ego-, neighbor, and higher-order embeddings improves GNN performance on heterophilic graphs. Zhu et al. (2021) uses a compatibility matrix to model the graph homophily or heterphily level. Chien et al. (2021) incorporates
the generalized PageRank algorithm with graph convolutions so as to jointly optimize node feature and topological information extraction for both homophilic and heterophilic graphs. However, the problem of GNNs on graphs with non-informative topologies (or noisy edges) remains open.
Unlike previous works, we tackle the above issues of GNNs by proposing the discrete p-Laplacian based message passing scheme, termed as p-Laplacian message passing. It is derived from a discrete regularization framework and is theoretically verified as an approximation of a polynomial graph filter defined on the spectral domain of the p-Laplacian. Spectral analysis of p-Laplacian message passing shows that it works simultaneously as low-pass and high-pass filters 1 and thus is applicable to both homophilic and heterophilic graphs. Moreover, when p ̸= 2, our theoretical results indicate that it can adaptively learn aggregation weights in terms of the variation of node embeddings on edges (measured by the graph gradient (Amghibech, 2003; Zhou & Schölkopf, 2005; Luo et al., 2010)), and work as low-pass or both low-pass and high-pass filters on a node according to the local variation of node embeddings around the node (measured by the norm of graph gradients).
Based on p-Laplacian message passing, we propose a new GNN architecture, called pGNN, to enable GNNs to work with heterophilic graphs and graphs with non-informative topologies. Several existing GNN architectures, including SGC (Wu et al., 2019), APPNP (Klicpera et al., 2019) and GPRGNN (Chien et al., 2021), can be shown to be analogical to pGNN with p = 2. Our empirical studies on real-world benchmark datasets (homophilic and heterophilic datasets) and synthetic datasets (cSBM (Deshpande et al., 2018)) demonstrate that pGNNs obtain the best performance on heterophilic graphs and competitive performance on homophilic graphs against state-of-the-art GNNs. Moreover, experimental results on graphs with different levels of noisy edges show that pGNNs work much more robustly than GNN baselines and even as well as MLPs on graphs with completely random edges. Additional experiments (reported in Appendix F.5) illustrate that intergrating pGNNs with existing GNN architectures (i.e. GCN (Kipf & Welling, 2017), JKNet (Xu et al., 2018)) can significantly improve their performance on heterophilic graphs. In conclusion, our contributions can be summarized as below:
(1) New methodologies. We propose p-Laplacian message passing and pGNN to adapt GNNs to heterophilic graphs and graphs where the topology is non-informative for label prediction. (2) Superior performance. We empirically demonstrate that pGNNs is superior on heterophilic graphs and competitive on homophilic graphs against state-of-the-art GNNs. Moreover, pGNNs work robustly on graphs with noisy edges or non-informative topologies. (3) Theoretical justification. We theoretically demonstrate that p-Laplacian message passing works as both low-pass and high-pass filters and the message passing iteration is guarantee to converge with proper settings. (4) New paradigm of designing GNN architectures. We bridge the gap between discrete regularization framework and GNNs, which could further inspire researchers to develop new graph convolutions or message passing schemes using other regularization techniques with explicit assumptions on graphs. Due to space limit, we defer the discussions on related work and future work and all proofs to the Appendix.

heading: 2 PRELIMINARIES AND BACKGROUND
text: Notation. Let G = (V, E ,W) be an undirected graph, where V = {1, 2, . . . , N} is the set of nodes, E ⊆ V × V is the set of edges, W ∈ RN×N is the adjacency matrix and Wi,j = Wj,i, Wi,j > 0 for [i, j] ∈ E , Wi,j = 0, otherwise. Ni = {j}[i,j]∈E denotes the set of neighbors of node i, D ∈ RN×N = diag(D1,1, . . . DN,N ) denotes the diagonal degree matrix with Di,i = ∑N j=1 Wi,j , for i = 1, . . . , N . f : V → R and g : E → R are functions defined on the vertices and edges of G, respectively. FV denotes the Hilbert space of functions endowed with the inner product ⟨f, f̃⟩FV :=∑
i∈V f(i)f̃(i). Similarly define FE . We also denote by [K] = {1, 2, . . . ,K},∀K ∈ N and we use ∥x∥ = ∥x∥2 = ( ∑d i=1 x 2 i )
1/2,∀x ∈ Rd to denote the Frobenius norm of a vector. Problem Formulation. Given a graph G = (V, E ,W), each node i ∈ V has a feature vector Xi,: which is the i-th row of X and a subset of nodes in G have labels from a label set L = {1, . . . , L}. The goal of semi-supervised node classification on G is to learn a mapping M : V → L and predict the labels of unlabeled nodes.
1Note that if the low frequencies and high frequencies dominate the middle frequencies (the frequencies that are around the cutoff frequency), we say that the filter works both as low-pass and high-pass filters.
Homophily and Heterophily. The homophily or heterophily of a graph is used to describe the relation of labels between linked nodes in the graphs. The level of homophily of a graph can be measured by H(G) = Ei∈V [∣∣{j}j∈Ni,yi=yj ∣∣ /|Ni|] (Pei et al., 2020; Chien et al., 2021), where∣∣{j}j∈Ni,yi=yj ∣∣ denotes the number of neighbors of i ∈ V that share the same label as i and H(G) → 1 corresponds to strong homophily while H(G) → 0 indicates strong heterophily. We say that a graph is a homophilic (heterophilic) graph if it has strong homophily (heterophily).
Graph Gradient. The graph gradient of an edge [i, j], i, j ∈ V is defined to be a measurement of the variation of a function f 2 : V → R on the edge [i, j]. Definition 1 (Graph Gradient). Given a graph G = (V, E) and a function f : V → R, the graph gradient is an operator ∇ : FV → FE defined by
(∇f)([i, j]) := √ Wi,j Dj,j f(j)− √ Wi,j Di,i f(i), for all [i, j] ∈ E . (1)
For [i, j] /∈ E , (∇f)([i, j]) := 0. The graph gradient of a function f at vertex i is defined to be ∇f(i) := ((∇f)([i, 1]), . . . , (∇f)([i,N ])) and its Frobenius norm is given by ∥∇f(i)∥2 := ( ∑N
j=1(∇f)2([i, j]))1/2, which measures the variation of f around node i. We measure the variation of f over the whole graph G by Sp(f) where it is defined to be
Sp(f) := 1
2 N∑ i=1 N∑ j=1 ∥(∇f)([i, j])∥p = 1 2 N∑ i=1 N∑ j=1 ∥∥∥∥∥ √ Wi,j Dj,j f(j)− √ Wi,j Di,i f(i) ∥∥∥∥∥ p , for p ≥ 1, (2)
Note that the definition of Sp here is different with the p-Dirichlet form in Zhou & Schölkopf (2005). Graph Divergence. The graph divergence is defined to be the adjoint of the graph gradient: Definition 2 (Graph Divergence). Given a graph G = (V, E), and functions f : V → R, g : E → R, the graph divergence is an operator div : FE → FV which satisfies ⟨∇f, g⟩ = ⟨f,−divg⟩. (3) The graph divergence can be computed by
(divg)(i) = N∑ j=1 √ Wi,j Di,i (g([i, j])− g([j, i])) . (4)
Fig. 4 in Appendix E.1 gives a tiny example of illustration of graph gradient and graph divergence.
Graph p-Laplacian Operator. By the definitions of graph gradient and graph divergence, we reach the definition of graph p-Laplacian operator as below. Definition 3 (Graph p-Laplacian3). Given a graph G = (V, E) and a function f : V → R, the graph p-Laplacian is an operator ∆p : FV → FV defined by
∆pf := − 1
2 div(∥∇f∥p−2∇f), for p ≥ 1. (5)
where ∥ · ∥p−2 is element-wise, i.e. ∥∇f(i)∥p−2 = (∥(∇f)([i, 1])∥p−2, . . . , ∥(∇f)([i,N ])∥p−2). Substituting Eq. (1) and Eq. (4) into Eq. (5), we obtain
(∆pf)(i) = N∑ j=1 √ Wi,j Di,i ∥(∇f)([j, i])∥p−2 (√ Wi,j Di,i f(i)− √ Wi,j Dj,j f(j) ) (6)
The graph p-Laplacian is semi-definite: ⟨f,∆pf⟩ = Sp(f) ≥ 0 and we have ∂Sp(f) ∂f ∣∣∣∣ i = p(∆pf)(i). (7) When p = 2, ∆2 is refered as the ordinary Laplacian operator and ∆2 = I − D−1/2WD−1/2 and when p = 1, ∆1 is refered as the Curvature operator and ∆1f := − 12div(∥∇f∥
−1∇f). Note that Laplacian ∆2 is a linear operator, while in general for p ̸= 2, p-Laplacian is nonlinear since ∆p(af) ̸= a∆p(f) for a ∈ R.
2f can be a vector function: f : V → Rc for some c ∈ N and here we use f : V → R for better illustration. 3Note that the definition adopted is slightly different with the one used in Zhou & Schölkopf (2005) where ∥ · ∥p−2 is not element-wise and the one used in some literature such as Amghibech (2003); Bühler & Hein (2009), where (∆pf)(i) = ∑N j=1 Wi,j Di,i |f(i)− f(j)|p−2 (f(i)− f(j)) for p > 1 and p = 1 is not allowed.
3 p-LAPLACIAN BASED GRAPH NEURAL NETWORKS
In this section, we derive the p-Laplacian message passing scheme from a p-Laplacian regularization framework and present pGNN, a new GNN architecture developed upon the new message passing scheme. We theoretically characterize how p-Laplacian message passing adaptively learns aggregation weights and profits pGNN for being effective on both homophilic and heterophilic graphs.
3.1 p-LAPLACIAN REGULARIZATION FRAMEWORK
Given an undirected graph G = (V, E) and a signal function with c (c ∈ N) channels f : V → Rc, let X = (X⊤1,:, . . . ,X ⊤ N,:)
⊤ ∈ RN×c with Xi,: ∈ R1×c, i ∈ [N ] denoting the node features of G and F = (F⊤1,:, . . . ,F ⊤ N,:)
⊤ ∈ RN×c be a matrix whose ith row vector Fi,: ∈ R1×c, i ∈ [N ] represents the function value of f at the i-th vertex in G. We present a p-Laplacian regularization problem whose cost function is defined to be
F∗ = argmin F L(F) := argmin F Sp(F) + µ N∑ i=1 ∥Fi,: −Xi,:∥2, (8)
where µ ∈ (0,∞). The first term of the right-hand side in Eq. (8) is a measurement of variation of the signal over the graph based on p-Laplacian. As we will discuss later, different choices of p result in different smoothness constraint on the signals. The second term is the constraint that the optimal signals F∗ should not change too much from the input signal X, and µ provides a trade-off between these two constraints.
Regularization with p = 2. When p = 2, the solution of Eq. (8) satisfies ∆2F∗ + µ(F∗ −X) = 0 and we can obtain the closed form (Zhou et al., 2003; Zhou & Schölkopf, 2005)
F∗ = µ(∆2 + µIN ) −1X. (9)
Then, we could use the following iteration algorithm to get an approximation of Eq. (9):
F(k+1) = αD−1/2WD−1/2F(k) + βX, (10) where k represents the iteration index, α = 11+µ and β = µ 1+µ = 1− α. The iteration converges to a closed-form solution as k goes to infinity (Zhou et al., 2003; Zhou & Schölkopf, 2005). We could relate the the result here with the personalized PageRank (PPR) (Page et al., 1999; Klicpera et al., 2019) algorithm (proof defered to Appendix D.1): Theorem 1 (Relation to personalized PageRank (Klicpera et al., 2019)). µ(∆2 + µIN )−1 in the closed form solution of Eq. (9) is equivalent to the personalized PageRank matrix.
Regularization with p > 1. For p > 1, the solution of Eq. (8) satisfies p∆pF∗ +2µ(F∗ −X) = 0. By Eq. (6) we have that, for all i ∈ [N ],
N∑ j=1 Wi,j√ Di,i
∥(∇f∗)([j, i])∥p−2 (
1√ Di,i F∗i,: − 1√ Dj,j F∗j,:
) + 2µ
p
( F∗i,: −Xi,: ) = 0.
Based on which we can construct a similar iterative algorithm to obtain a solution (Zhou & Schölkopf, 2005):
F (k+1) i,: = α (k) i,i N∑ j=1 M (k) i,j√ Di,iDj,j F (k) j,: + β (k) i,i Xi,:, for all i ∈ [N ], (11)
with M(k) ∈ RN×N , α(k) = diag(α(k)1,1 , . . . , α (k) N,N ), β (k) = diag(β (k) 1,1 , . . . , β (k) N,N ) updated by
M (k) i,j = Wi,j ∥∥∥∥∥ √ Wi,j Di,i F (k) i,: − √ Wi,j Dj,j F (k) j,: ∥∥∥∥∥ p−2 , for all i, j ∈ [N ], (12)
α (k) i,i = 1 / N∑ j=1 M (k) i,j Di,i + 2µ p  , β(k)i,i = 2µp αi,i, for all i ∈ [N ], (13) Note that in Eq. (12), when
∥∥∥√Wi,jDi,i F(k)i,: −√Wi,jDj,j F(k)j,: ∥∥∥ = 0, we set M (k)i,j = 0. It is easy to see that Eq. (10) is the special cases of Eq. (14) with p = 2.
Remark 1 (Discussion on p = 1). For p = 1, when f is a real-valued function (c = 1), ∆1f is a step function, which could make the stationary condition of the objective function Eq. (8) become problematic. Additionally,∆1f is not continuous at ∥(∇f)([i, j])∥ = 0. Therefore, p = 1 is not allowed when f is a real value function. On the other hand, note that there is a Frobenius norm in ∆pf . When f is a vector-valued function (c > 1), the step function in ∆1f only exists on the axes. The stationary condition will be fine if the node embeddings F are not a matrix of vectors that has only one non-zero element, which is true for many graphs. p = 1 may work for these graphs. Overall, we suggest to use p > 1 in practice but p = 1 may work for graphs with multiple channel signals as well. We conduct experiments for p > 1 (e.g., p = 1.5, 2, 2.5) and p = 1 in Sec. 5.
3.2 p-LAPLACIAN MESSAGE PASSING AND pGNN ARCHITECTURE
p-Laplacian Message Passing. Rewrite Eq. (11) in a matrix form we obtain
F(k+1) = α(k)D−1/2M(k)D−1/2F(k) + β(k)X. (14)
Eq. (14) provides a new message passing mechanism, named p-Laplacian message passing.
Remark 2. αD−1/2MD−1/2 in Eq. (14) can be regarded as the learned aggregation weights at each step for message passing. It suggests that p-Laplacian message passing could adaptively tune the aggregation weights during the course of learning, which will be demonstrated theoretically and empirically in the sequel of this paper. βX in Eq. (14) can be regarded as a residual unit, which helps the model escape from the oversmoothing issue (Chien et al., 2021).
We present the following theorem to show the shrinking property of p-Laplacian message passing.
Theorem 2 (Shrinking Property of p-Laplacian Message Passing). Given a graph G = (V, E ,W) with node features X, β(k),F(k),M(k),α(k) are updated accordingly to Equations (11) to (13) for k = 0, 1, . . . ,K and F(0) = X. Then there exist some positive real value µ > 0 which depends on X,G, p and p > 1 such that
Lp(F(k+1)) ≤ Lp(F(k)).
Proof see Appendix D.2. Thm. 2 shows that with some proper positive real value µ and p > 1, the loss of the objective function Eq. (8) is guaranteed to decline after taking one step p-Laplacian message passing. Thm. 2 also demonstrates that the iteration Equations (11) to (13) is guaranteed to converge for p > 1 with some proper µ which is chosen depends on the input graph and p.
pGNN Architecture. We design the architecture of pGNNs using p-Laplacian message passing. Given node features X ∈ RN×c, the number of node labels L, the number of hidden units h, the maximum number of iterations K, and M, α, and β updated by Equations (12) and (13) respectively, we give the pGNN architecture as following:
F(0) = ReLU(XΘ(1)), (15)
F(k+1) = α(k)D−1/2M(k)D−1/2F(k) + β(k)F(0), k = 0, 1, . . . ,K − 1, (16) Z = softmax(F(K)Θ(2)), (17)
where Z ∈ RN×L is the output propbability matrix with Zi,j is the estimated probability that the label at node i ∈ [N ] is j ∈ [L] given the features X and the graph G, Θ(1) ∈ Rc×h and Θ(2) ∈ Rh×L are the first- and the second-layer parameters of the neural network, respectively. Remark 3 (Connection to existing GNN variants). The message passing scheme of pGNNs is different from that of several GNN variants (say, GCN, GAT, and GraphSage), which repeatedly stack message passing layers. In contrast, it is similar with SGC (Wu et al., 2019), APPNP (Klicpera et al., 2019), and GPRGNN (Chien et al., 2021). SGC is an approximation to the closed-form in Eq. (9) (Fu et al., 2020). By Thm. 1, it is easy to see that APPNP, which uses PPR to propagate the node embeddings, is analogical to pGNN with p = 2, termed as 2.0GNN. APPNP and 2.0GNN work analogically and effectively on homophilic graphs. 2.0GNN can also work effectively on heterophilic graphs by letting Θ(2) be negative. However, APPNP fails on heterophilic graphs as its PPR weights
are fixed (Chien et al., 2021). Unlike APPNP, GPRGNN, which adaptively learn the generalized PageRank (GPR) weights, works similarly to 2.0GNN on both homophilic and heterophilic graphs. However, GPRGNN needs more supervised information in order to learn optimal GPR weights. On the contrary, pGNNs need less supervised information to obtain similar results because Θ(2) acts like a hyperplane for classification. pGNNs could work better under weak supervised information. Our analysis is also verified by the experimental results in Sec. 5.
We also provide an upper-bounding risk of pGNNs by Thm. 4 in Appendix C.1 to study the effect of the hyperparameter µ on the performance of pGNNs. Thm. 4 shows that the risk of pGNNs is upperbounded by the sum of three terms: the risk of label prediction using only the original node features X, the norm of p-Laplacian diffusion on X, and the magnitude of the noise in X. µ controls the trade-off between these three terms. The smaller µ, the more weights on the p-Laplacian diffusion term and the noise term and the less weights on the the other term and vice versa.
4 SPECTRAL VIEWS OF p-LAPLACIAN MESSAGE PASSING
In this section, we theoretically demonstrate that p-Laplacian message passing is an approximation of a polynomial graph filter defined on the spectral domain of p-Laplacian. We show by spectral analysis that p-Laplacian message passing works simultaneously as low-pass and high-pass filters.
p-Eigenvalues and p-Eigenvectors of the Graph p-Laplacian. We first introduce the definitions of p-eigenvalues and p-eigenvectors of p-Laplacian. Let ϕp : R → R defined as ϕp(u) = ∥u∥p−2u, for u ∈ R, u ̸= 0. Note that ϕ2(u) = u. For notational simplicity, we denote by ϕp(u) = (ϕp(u1), . . . , ϕp(uN ))
⊤ for u ∈ RN and Φp(U) = (ϕp(U:,1), . . . , ϕp(U:,N )) for U ∈ RN×N and U:,i ∈ RN is the ith column vector of U. Definition 4 (p-Eigenvector and p-Eigenvalue). A vector u ∈ RN is a p-eigenvector of ∆p if it satisfies the equation (∆pu)i = λϕp(ui), for all i ∈ [N ], where λ ∈ R is a real value referred as a p-eigenvalue of ∆p associated with the p-eigenvector u. Definition 5 (p-Orthogonal (Luo, Huang, Ding, and Nie, 2010)). Given two vectors u,v ∈ RN with u,v ̸= 0, we call that u and v is p-orthogonal if ϕp(u)⊤ϕp(v) = ∑N i=1 ϕp(ui)ϕp(vi) = 0.
Luo et al. (2010) demonstrated that the p-eigenvectors of ∆p are p-orthogonal to each other (see Thm. 5 in Appendix C.2 for details). Therefore, the space spanned by the multiple p-eigenvectors of ∆p is p-orthogonal. Additionally, we demonstrate that the p-eigen-decomposition of ∆p is given by: ∆p = Φp(U)ΛΦp(U)⊤ (see Thm. 6 in Appendix C.3 for details), where U is a matrix of p-eigenvectors of ∆p and Λ is a diagonal matrix in which the diagonal is the p-eigenvalues of ∆p.
Graph Convolutions based on p-Laplacian. Based on Thm. 5, the graph Fourier Transform f̂ of any function f on the vertices of G can be defined as the expansion of f in terms of Φ(U) where U is the matrix of p-eigenvectors of ∆p: f̂ = Φp(U)⊤f . Similarly, the inverse graph Fourier transform is then given by: f = Φp(U)f̂ . Therefore, a signal X ∈ RN×c being filtered by a spectral filter gθ can be expressed formally as: gθ ⋆ X = Φp(U)ĝθ(Λ)Φp(U)⊤X, where Λ denotes a diagonal matrix in which the diagonal corresponds to the p-eigenvalues {λl}l=0,...,N−1 of ∆p and ĝθ(Λ) denotes a diagonal matrix in which the diagonal corresponds to spectral filter coefficients. Let ĝθ be a polynomial filter defined as ĝθ = ∑K−1 k=0 θkλ k l , where the parameter θ = [θ0, . . . , θK−1]
⊤ ∈ RK is a vector of polynomial coefficients. By the p-eigen-decomposition of p-Laplacian, we have
gθ ⋆X ≈ K−1∑ k=0 θkΦp(U)Λ kΦp(U) ⊤X = K−1∑ k=0 θk∆ k pX. (18)
Theorem 3. The K-step p-Laplacian message passing is a K-order polynomial approximation to the graph filter given by Eq. (18).
Proof see Appendix D.3. Thm. 3 indicates that p-Laplacian message passing mechanism is implicitly a polynomial spectral filter defined on the spectral domain of p-Laplacian.
Spectral Analysis of p-Laplacian Message Passing. Here, we analyze the spectral propecties of p-Laplacian message passing. We can approximately view p-Laplacian message pasing as a filter of a linear combination of K spectral filters g(Λ)(0), g(Λ)(1), . . . , g(Λ)(K−1) with each spectral filter
defined to be g(Λ)(k) := (αD−1/2MD−1/2)k where Mi,j = Wi,j∥ √
Wi,j Di,i
Fi,: − √
Wi,j Dj,j Fj,:∥p−2
for i, j ∈ [N ] and F is the matrix of node embeddings. We can study the properties of p-Laplacian message passing by studying the spectral properties of αD−1/2MD−1/2 as given below.
Proposition 1. Given a connected graph G = (V, E ,W) with node embeddings F and the pLaplacian ∆p with its p-eigenvectors {u(l)}l=0,1,...,N−1 and the p-eigenvalues {λl}l=0,1,...,N−1. Let gp(λi−1) := αi,i ∑ j D −1/2 i,i Mi,jD −1/2 j,j for i ∈ [N ] be the filters defined on the spectral domain of ∆p, where Mi,j = Wi,j∥∇f([i, j])∥p−2, (∇f)([i, j]) is the graph gradient of the edge between node i and j and ∥∇f(i)∥ is the norm of graph gradient at i. Ni denotes the number of edges connected to i, Nmin = min{Nj}j∈[N ], and k = argmaxj({|u (l) j |/ √ Dl,l}j∈[N ];l=0,...,N−1), then
1. When p = 2, gp(λi−1) works as both low-pass and high-pass filters.
2. When p > 2, if ∥∇f(i)∥ ≤ 2(p−1)/(p−2), gp(λi−1) works as both low-pass and high-pass filters on node i and gp(λi−1) works as low-pass filters on i when ∥∇f(i)∥ ≥ 2(p−1)/(p−2). 3. When 1 ≤ p < 2, if 0 ≤ ∥∇f(i)∥ ≤ 2(2 √ Nk)
1/(p−2), gp(λi−1) works as low-pass filters on node i and gp(λi−1) works as both low-pass and high-pass filters on i when
∥∇f(i)∥ ≥ 2 ( 2 √ Nk )1/(p−2) . Specifically, when p = 1, Nk can be replaced by Nmin.
Proof see Appendix D.7. Proposition 1 shows that when p ̸= 2, p-Laplacian message passing adaptively works as low-pass or both low-pass and high-pass filters on node i in terms of the degree of local node embedding variation around i, i.e. the norm of the graph gradient ∥∇f(i)∥ at node i. When p = 2, p-Laplacian message passing works as both low-pass and high-pass filters on node i regardless of the value of ∥∇f(i)∥. When p > 2, p-Laplacian message passing works as lowpass filters on node i for large ∥∇f(i)∥ and works as both low-pass and high-pass filters for small ∥∇f(i)∥. Therefore, pGNNs with p > 2 can work very effectively on graphs with strong homophily. When 1 ≤ p < 2, p-Laplacian message passing works as low-pass filters for small ∥∇f(i)∥ and works as both low-pass and high-pass filters for large ∥∇f(i)∥. Thus, pGNNs with 1 ≤ p < 2 can work effectively on graphs with low homophily, i.e. heterophilic graphs. The results here confirms our analysis of the aggregation weights of p-Laplacian message passing presented in Thm. 2.

heading: 5 EMPIRICAL STUDIES
text: In this section, we empirically study the effectiveness of pGNNs for semi-supervised node classification using and real-world benchmark and synthetic datasets with heterophily and strong homophily. The experimental results are also used to validate our theoretical findings presented previously.
Datasets and Experimental Setup. We use seven homophilic benchmark datasets: citation graphs Cora, CiteSeer, PubMed (Sen et al., 2008), Amazon co-purchase graphs Computers, Photo, coauthor graphs CS, Physics (Shchur et al., 2018), and six heterophilic benchmark datasets: Wikipedia graphs Chameleon, Squirrel (Rozemberczki et al., 2021), the Actor co-occurrence graph, webpage graphs Wisconsin, Texas, Cornell (Pei et al., 2020). The node classification tasks are conducted in the transductive setting. Following Chien et al. (2021), we use the sparse splitting (2.5%/2.5%/95%) and the dense splitting (60%/20%/20%) to randomly split the homophilic and heterophilic graphs into training/validation/testing sets, respectively. Dataset statistics and their levels of homophily are presented in Appendix E.
Baselines. We compare pGNN with seven models, including MLP, GCN (Kipf & Welling, 2017), SGC (Wu et al., 2019), GAT (Velickovic et al., 2018), JKNet (Xu et al., 2018), APPNP (Klicpera et al., 2019), GPRGNN (Chien et al., 2021). We use the Pytorch Geometric library (Fey & Lenssen, 2019) to implement all baselines except GPRGNN. For GPRGNN, we use the code released by the authors4. The details of hyperparameter settings are deferred to Appendix E.3.
4https://github.com/jianhao2016/GPRGNN
Superior Performance on Real-World Heterophilic Datasets. The results on homophilic benchmark datasets are deferred to Appendix F.1, which show that pGNNs obtains competitive performance against state-of-the-art GNNs on homophilic datasets. Table 1 summarizes the results on heterophilic benchmark datasets. Table 1 shows that pGNNs significantly dominate the baselines and 1.0GNN obtains the best performance on all heterophilic graphs except the Texas dataset. For Texas, 2.0GNN is the best. We also observe that MLP works very well and significantly outperforms most GNN baselines, which indicates that the graph topology is not informative for label prediction on these heterophilic graphs. Therefore, propagating and transforming node features over the graph topology could lead to worse performance than MLP. Unlike ordinary GNNs, pGNNs can adaptively learn aggregation weights and ignore edges that are not informative for label prediction and thus could work better. It confirms our theoretical findings presented in previous sections. Note that GAT can also learn aggregation weights, i.e. the attention weights. However, the aggregation weights learned by GAT are significantly distinct from that of pGNNs, as we will show following.
Computers
Interpretability of the Learned Aggregation Weights of pGNNs. We showcase the interpretability of the learned aggregation weights αi,iD −1/2 i,i Mi,jD −1/2 j,j of
pGNNs by studying its entropy distribution, along with the attention weights of GAT on real-world datasets. Denote {Ai,j}j∈Ni as the aggregation weights of node i and its neighbors. For GAT, {Ai,j}j∈Ni are referred as the attention weights (in the first layer) and for pGNNs are αi,iD −1/2 i,i Mi,jD −1/2 j,j . For any node i, {Ai,j}j∈Ni forms a discrete probability distribution over all its neighbors with the entropy given by H({Ai,j}j∈Ni) = − ∑ j∈Ni Ai,j log(Ai,j). Low entropy means high degree of concentration and vice versa. An entropy of zero means all aggregation weights or attentions are on one source node. The uniform distribution has the highest entropy of log(Di,i). Fig. 1 reports the results on Computers, Wisconsin and we defer more results on other datasets to Appendix F.2 due to space limit. Fig. 1 shows that the aggregation weight entropy distributions of GAT and pGNNs on Computers (homophily) are both similar to the uniform case. It indicates the original graph topology of Computers is very helpful for label prediction and therefore GNNs could work very well on Computers. However, for Wisconsin (heterophily), the entropy distribution of pGNNs is significantly different from that of GAT and the uniform case. Most entropy of pGNNs is around zero, which means that most aggregation weights are on one source node. It states that the original graph topology of Wisconsin is not helpful for label prediction, which explains why MLP works well on Wisconsin.
On the contrary, the entropy distribution of GAT is similar to the uniform case and therefore GAT works similarly to GCN and is significantly worse than pGNNs on Wisconsin. Similar results can be observed on the experiments on more datasets in Appendix F.2.
Results on cSBM Datasets. We exam the performance of pGNNs on heterophilic graphs whose topology is informative for label prediction using synthetic graphs generated by cSBM (Deshpande et al., 2018) with ϕ ∈ {−1,−0.75, . . . , 1}. We use the same settings of cSBM used in Chien et al. (2021). Due to the space limit, we refer the readers to Chien et al. (2021) for more details of cSBM dataset. Fig. 2 reports the results on cSBM using sparse splitting (for results on cSBM with dense splitting see Appendix F.3). Fig. 2 shows that when ϕ ≤ −0.5 (heterophilic graphs), 2.0GNN obtains the best performance and pGNNs and GPRGNN significantly dominate the others. It validates the effectiveness of pGNNs on het-
erophilic graphs. Moreover, 2.0GNN works better than GPRGNN and it again confirms that 2.0GNN is more superior under weak supervision (2.5% training rate), as stated in Remark 3. Note that 1.0GNN and 1.5GNN are not better than 2.0GNN, the reason could be the iteration algorithms Eq. (11) with p = 1, 1.5 are not as stable as the one with p = 2. When the graph topology is almost non-informative for label prediction (ϕ = −0.25, 0), The performance of pGNNs is close to MLP and they outperform the other baselines. Again, it validates that pGNNs can erase noninformative edges and work as well as MLP and confirms the statements in Thm. 4. When the graph is homophilic (ϕ ≥ 0.25), 1.5GNN is the best on weak homophilic graphs (ϕ = 0.25, 0.5) and pGNNs work competitively with all GNN baselines on strong homophilic graphs (ϕ ≥ 0.75). Results on Datasets with Noisy Edges. We conduct experiments to evaluate the performance of pGNNs on graphs with noisy edges by randomly adding edges to the graphs and randomly remove the same number of original edges. We define the random edge rate as r := #random edges#all edges . The experiments are conducted on 4 homophilic datasets (Computers, Photo, CS, Physics) and 2 heterophilic datasets (Wisconsin, Texas) with r = 0.25, 0.5, 1. Fig. 3 reports the results on Computers, Wisconsin and we defer more results to Appendix F.4. Fig. 3 shows that pGNNs significantly outperform all baselines. Specifically, 1.5GNN obtains the best performance on Computers, and 1.5GNN and 2.0GNN even work as well as MLP on Computers with completely random edges (r = 1). For Wisconsin, 1.0GNN is the best, and 1.0GNN and 1.5GNN significantly dominate the others. We also observed that APPNP and GPRGNN, whose architectures are analogical to 2.0GNN, also work better than other GNNs. Nevertheless, they are significantly outperformed by pGNNs overall. Similar results can be observed in the experiments conducted on more datasets as presented in Appendix F.4.
CONCLUSION. We have addressed the problem of generalizing GNNs to heterophilic graphs and graphs with noisy edges. To this end, we derived a novel p-Laplacian message passing scheme from a discrete regularization framework and proposed a new pGNN architecture. We theoretically demonstrate our method works as low-pass and high-pass filters and thereby applicable to both homophilic and heterophilic graphs. We empirically validate our theoretical results and show the advantages of our methods on heterophilic graphs and graphs with non-informative topologies.

heading: REPRODUCIBILITY STATEMENT
text: In order to ensure reproducibility, we have made the efforts in the following respects: (1) Provide a sampled code as the supplementary material; (2) Provide self-contained proofs of the main claims in Appendices C and D; (3) Provide more details on experimental configurations in Appendix E and experimental results in Appendix F. All the datasets are publicly available as described in the main text. We will fully release our training and evaluation code, as well as our train/validation/test splits.

heading: CONTENTS

heading: A Related Work 16

heading: B Discussions and Future Work 17

heading: C Additional Theorems 18
text: C.1 Theorem 4 (Upper-Bounding Risk of pGNN) . . . . . . . . . . . . . . . . . . . . 18
C.2 Theorem 5 (p-Orthogonal Theorem (Luo et al., 2010)) . . . . . . . . . . . . . . . 18
C.3 Theorem 6 (p-Eigen-Decomposition of ∆p) . . . . . . . . . . . . . . . . . . . . . 18
C.4 Theorem 7 (Bounds of p-Eigenvalues) . . . . . . . . . . . . . . . . . . . . . . . . 19

heading: D Proof of Theorems 19
text: D.1 Proof of Theorem 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
D.2 Proof of Theorem 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
D.3 Proof of Theorem 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
D.4 Proof of Theorem 4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
D.5 Proof of Theorem 6 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
D.6 Proof of Theorem 7 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
D.7 Proof of Proposition 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27

heading: E Dataset Statistics and Hyperparameters 29
text: E.1 Illustration of Graph Gradient and Graph Divergence . . . . . . . . . . . . . . . . 29
E.2 Dataset Statistics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
E.3 Hyperparameter Settings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

heading: F Additional Experiments 31
text: F.1 Experimental Results on Homophilic Benchmark Datasets . . . . . . . . . . . . . 31
F.2 Experimental Results of Aggregation Weight Entropy Distribution . . . . . . . . . 31
F.3 Experimental Results on cSBM . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
F.4 Experimental Results on Graphs with Noisy Edges . . . . . . . . . . . . . . . . . 34
F.5 Experimental Results of Intergrating pGNNs with GCN and JKNet . . . . . . . . . 36
F.6 Experimental Results of pGNNs on PPI Dataset for Inductive Learning . . . . . . . 37
F.7 Experimental Results of pGNNs on OGBN arXiv Dataset . . . . . . . . . . . . . . 37
F.8 Running Time of pGNNs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
F.9 Experimental Results on Benchmark Datasets for 64 Hidden Units . . . . . . . . . 38
F.10 Training Curves for p = 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
F.11 Visualization Results of Node Embeddings . . . . . . . . . . . . . . . . . . . . . 40

heading: A RELATED WORK
text: Graph Neural Networks. Graph neural networks (GNNs) are a variant of neural networks for graph-structured data, which can propagate and transform the node features over the graph topology and exploit the information in the graphs. Graph convolutional networks (GCNs) are one type of GNNs whose graph convolution mechanisms or the message passing schemes were mainly inspired by the field of graph signal processing. Bruna et al. (2014) defined a nonparametric graph filter using the Fourier coefficients. Defferrard et al. (2016) introduced Chebyshev polynomial to avoid computational expensive eigen-decomposition of Laplacian and obtain localized spectral filters. GCN (Kipf & Welling, 2017) used the first-order approximation and reparameterized trick to simplify the spectral filters and obtain the layer-wise graph convolution. SGC (Wu et al., 2019) further simplify GCN by removing non-linear transition functions between each layer. Chen et al. (2018) propose importance sampling to design an efficient variant of GCN. Xu et al. (2018) explored a jumping knowledge architecture that flexibly leverages different neighborhood ranges for each node to enable better structure-aware representation. Atwood & Towsley (2016); Liao et al. (2019); Abu-El-Haija et al. (2019) exploited multi-scale information by diffusing multi-hop neighbor information over the graph topology. Wang & Leskovec (2020) used label propagation to improve GCNs. Klicpera et al. (2019) incorporated personalized PageRank with GCNs. Liu et al. (2021) introduced a l1 norm-based graph smoothing term to enhance the local smoothnesss adaptivity of GNNs. Hamilton et al. (2017); Zeng et al. (2020) proposed sampling and aggregation frameworks to extent GCNs to inductive learning settings. Another variant of GNNs is graph attention networks (Velickovic et al., 2018; Thekumparampil et al., 2018; Abu-El-Haija et al., 2018), which use attention mechanisms to adaptively learn aggregation weights based on the nodes features. There are many other works on GNNs (Pei et al., 2020) (Ying et al., 2018; Xinyi & Chen, 2019; Velickovic et al., 2019; Zeng et al., 2020), we refer to Zhou et al. (2020); Battaglia et al. (2018); Wu et al. (2021) for a comprehensive review. Most GNN models implicitly assume that the labels of nodes and their neighbors should be the same or consistent, while it does not hold for heterophilic graphs. Zhu et al. (2020) investigated the issues of GNNs on heterophilic graphs and proposed to separately learn the embeddings of ego-node and its neighborhood. Zhu et al. (2021) proposed a framework to model the heterophily or homophily levels of graphs. Chien et al. (2021) incorporated generalized PageRank with graph convolution to adapt GNNs to heterophilic graphs.
There are also some works on the interpretability of GNNs proposed recently. Li et al. (2018); Ying et al. (2019); Fu et al. (2020) showed that spectral graph convolutions work as conducting Laplacian smoothing on the graph signals and Wu et al. (2019); NT & Maehara (2019) demonstrated that GCN, SGC work as low-pass filters. Gama et al. (2020) studied the stability properties of GNNs. Xu et al. (2019); Oono & Suzuki (2020); Loukas (2020) studied the expressiveness of GNNs. Verma & Zhang (2019); Garg et al. (2020) work on the generalization and representation power of GNNs.
Graph based Semi-supervised Learning. Graph-based semi-supervised learning works under the assumption that the labels of a node and its neighbors shall be the same or consistent. Many methods have been proposed in the last decade, such as Smola & Kondor (2003); Zhou et al. (2003); Belkin et al. (2004) use Laplacian regularization techniques to force the labels of linked nodes to be the same or consistent. Zhou & Schölkopf (2005) introduce discrete regularization techniques to impose different regularizations on the node features based on p-Laplacian. Lable propagation (Zhu et al., 2003) recursively propagates the labels of labeled nodes over the graph topology and use the convergence results to make predictions. To mention but a few, we refer to Zhou & Schölkopf (2005); van Engelen & Hoos (2020) for a more comprehensive review.

heading: B DISCUSSIONS AND FUTURE WORK
text: In this section, we discuss the future work of pGNNs. Our theoretical results and experimental results could lead to several potential extensions of pGNNs.
New Paradigm of Designing GNN Architectures. We bridge the gap between discrete regularization framework, graph-based semi-supervised learning, and GNNs, which provides a new paradigm of designing new GNN architectures. Following the new paradigm, researchers could introduce more regularization techniques, e.g., Laplacian regularization (Smola & Kondor, 2003; Belkin et al., 2004), manifold regularization (Sindhwani et al., 2005; Belkin et al., 2005; Niyogi, 2013), highorder regularization (Zhou & Belkin, 2011), Bayesian regularization (Liu et al., 2014), entropy regularization (Grandvalet & Bengio, 2004), and consider more explicit assumptions on graphs, e.g. the homophily assumption, the low-density region assumption (i.e. the decision boundary is likely to lie in a low data density region), manifold assumption (i.e. the high dimensional data lies on a low-dimensional manifold), to develop new graph convolutions or message passing schemes for graphs with specific properties and generalize GNNs to a much broader range of graphs. Moreover, the paradigm also enables us to explicitly study the behaviors of the designed graph convolutions or message passing schemes from the theory of regularization (Belkin & Niyogi, 2008; Niyogi, 2013; Slepcev & Thorpe, 2017).
Applications of pGNNs to learn on graphs with noisy topologies. The empirical results (as shown in Fig. 3 and Tables 6 and 7) on graphs with noisy edges show that pGNNs are very robust to noisy edges, which suggests the applications of p-Laplacian message passing and pGNNs on the graph learning scenarios where the graph topology could potentially be seriously intervened.
Integrating with existing GNN architectures. As shown in Table 9, the experimental results on heterophilic benchmark datasets illustrate that integrating GCN, JKNet with pGNNs can significantly improve their performance on heterophilic graphs. It shows that pGNN could be used as a plug-and-play component to be integrated into existing GNN architectures and improve their performance on real-world applications.
Inductive learning for pGNNs. pGNNs are shown to be very effective for inductive learning on PPI datasets as reported in Table 10. pGNNs even outperforms GAT on PPI, while using much fewer parameters than GAT. It suggests the promising extensions of pGNNs to inductive learning on graphs.

heading: C ADDITIONAL THEOREMS
text: C.1 THEOREM 4 (UPPER-BOUNDING RISK OF pGNN)
Theorem 4 (Upper-bounding risks of pGNNs). Given a graph G = (V, E ,W) with N nodes, let X ∈ RN×c be the node features and y ∈ RN be the node labels and M(k),α(k),βk,Fk are updated accordingly by Equations (12) to (14) for k = 0, 1, . . . ,K − 1 and F(0) = X, K ∈ N. Assume that G is d-regular and the ground-truth node features X∗ = X + ϵ, where ϵ ∈ RN×c represents the noise in the node features and there exists a L-Lipschitz function σ : RN×c → RN such that σ(X∗) = y. let ỹ(k+1) = α(k)D−1/2M(k)D−1/2σ(F(k)) + β(k)σ ( F(0) ) , we have
1
N N∑ i=1 |yi − ỹi| ≤ 1 N N∑ i=1 β (K−1) i,i |yi − σ(Xi,:)|
+ L
N N∑ i=1 α (K−1) i,i ∥∥∥∥∥∥∆(K−1)p F(K−1)i,: + K−2∑ k=0 K−2∏ l=k  N∑ j=1 M (l) i,j d ∆(k)p Xi,: ∥∥∥∥∥∥
+ L
N N∑ i=1 (1− β(K−1)i,i ) ∥ϵi,:∥ .
Proof see Appendix D.4. Thm. 4 shows that the risk of pGNNs is upper-bounded by the sum of three terms: The first term of the r.h.s in the above inequation represents the risk of label prediction using only the original node features X, the second term is the norm of p-Laplacian diffusion on the node features X, and the third term is the magnitude of the noise in the node features. αi,i and βi,i control the trade-off between these three terms and they are related to the hyperparameter µ in Eq. (10). The smaller µ, the smaller βi,i and larger αi,i, thus the more important of the p-Laplacian diffusion term but also the more effect from the noise. Therefore, for graphs whose topological information is not helpful for label prediction, we could impose more weights on the first term by using a large µ so that pGNNs work more like MLPs which simply learn on node features. While for graphs whose topological information is helpful for label prediction, we could impose more weights on the second term by using a small µ so that pGNNs can benefit from p-Laplacian smoothing on node features.
In practice, to choose a proper value of µ one may first simply apply MLPs on the node features to have a glance at the helpfulness of the node features. If MLPs work very well, there is not much space for the graph’s topological information to further improve the prediction performance and we may choose a large µ. Otherwise, there could be a large chance for the graph’s topological information to further improve the performance and we should choose a small µ.
C.2 THEOREM 5 (p-ORTHOGONAL THEOREM (LUO ET AL., 2010))
Theorem 5 (p-Orthogonal Theorem (Luo et al., 2010)). If u(l) and u(r) are two eigenvectors of p-Laplacian ∆p associated with two different non-zero eigenvalues λl and λr, W is symmetric and p ≥ 1, then u(l) and u(r) are p-orthogonal up to the second order Taylor expansion.
Thm. 5 implies that ϕp(u)(l)⊤ϕp(u(r)) ≈ 0, for all l, r = 0, . . . , N − 1 and λl ̸= λr. Therefore, the space spanned by the multiple eigenvectors of the graph p-Laplacian is p-orthogonal.
C.3 THEOREM 6 (p-EIGEN-DECOMPOSITION OF ∆p)
Theorem 6 (p-Eigen-Decomposition of ∆p). Given the p-eigenvalues {λl ∈ R}l=0,1,...,N−1, and the p-eigenvectors {u(l) ∈ RN}l=0,1,...,N−1 of p-Laplacian ∆p and ∥u(l)∥p = ( ∑N i=1(u (l) i ) p)1/p = 1, let U be a matrix of p-eigenvectors with U = (u(0),u(1), . . . ,u(N−1)) and Λ be a diagonal matrix with Λ = diag(λ0, λ1, . . . , λN−1), then the p-eigen-decomposition of p-Laplacian ∆p is given by
∆p = Φp(U)ΛΦp(U) ⊤.
When p = 2, it reduces to the standard eigen-decomposition of the Laplacian matrix.
Proof see Appendix D.5.
C.4 THEOREM 7 (BOUNDS OF p-EIGENVALUES)
Theorem 7 (Bounds of p-Eigenvalues). Given a graph G = (V, E ,W), if G is connected and λ is a p-eigenvalue associated with the p-eigenvector u of ∆p, let Ni denotes the number of edges connected to node i, Nmin = min{Ni}i=1,2,...,N , and k = argmax({|ui|/ √ Di,i}i=1,2,...,N ), then
1. for p ≥ 2, 0 ≤ λ ≤ 2p−1; 2. for 1 < p < 2, 0 ≤ λ ≤ 2p−1 √ Nk; 3. for p = 1, 0 ≤ λ ≤ √ Nmin.
Proof see Appendix D.6.

heading: D PROOF OF THEOREMS

heading: D.1 PROOF OF THEOREM 1
text: Proof. Let i be the one-hot indicator vector whose i-th element is one and the other elements are zero. Then, we can obtain the personalized PageRank on node i, denoted as πPPR(i), by using the recurrent equation (Klicpera et al., 2019):
π (k+1) PPR (i) = αD −1/2WD−1/2π (k) PPR(i) + βi,
where k is the iteration step, 0 < α < 1 and β = (1− α) represents the restart probability. Without loss of generality, suppose π(0)PPR(i) = i. Then we have,
π (k) PPR(i) = αD −1/2WD−1/2π (k−1) PPR (i) + βi = αD−1/2WD−1/2 ( αD−1/2WD−1/2π
(k−2) PPR (i) + βi
) + βi
= ( αD−1/2WD−1/2 )2 π
(t−2) PPR (i) + βαD −1/2WD−1/2i+ βi
= ( αD−1/2WD−1/2 )k π
(0) PPR(i) + β k−1∑ t=0 ( αD−1/2WD−1/2 )t i
= ( αD−1/2WD−1/2 )k i+ β k−1∑ t=0 ( αD−1/2WD−1/2 )t i
Since 0 < α < 1 and the eigenvalues of D−1/2WD−1/2 in [−1, 1], we have
lim k→∞
( αD−1/2WD−1/2 )k = 0,
and we also have
lim k→∞ k−1∑ t=0 ( αD−1/2WD−1/2 )t = ( IN − αD−1/2WD−1/2 )−1 .
Therefore,
πPPR(i) = lim k→∞
π (k) PPR(i) = β ( IN − αD−1/2WD−1/2 )−1 i
= β (α∆2 + (1− α)IN )−1 i = µ(∆2 + µIN ) −1i,
where we let α = 11+µ and β = µ 1+µ , µ > 0. Then the fully personalized PageRank matrix can be obtained by substituting i with IN :
ΠPPR = µ(∆2 + µIN ) −1.

heading: D.2 PROOF OF THEOREM 2
text: Proof. By the definition of Lp(f) in Eq. (8), we have for some positive real value µ, µ > 0
Lp(F) = 1
2 N∑ i=1 N∑ j=1 ∥∥∥∥∥ √ Wi,j Di,i Fi,: − √ Wi,j Dj,j Fj,: ∥∥∥∥∥ p + µ N∑ i=1 ∥Fi,: −Xi,:∥2.
and by Eq. (12),
M (k) i,j := Wi,j ∥∥∥∥∥ √ Wi,j Di,i F (k) i,: − √ Wi,j Dj,j F (k) j,: ∥∥∥∥∥ p−2
Then, we have
∂Lp(F(k)) ∂F
(k) i,:
= p N∑ j=1 √ Wi,j Di,i ∥∥∥∥∥ √ Wi,j Di,i F (k) i,: − √ Wi,j Dj,j F (k) j,: ∥∥∥∥∥ p−2(√ Wi,j Di,i F (k) i,: − √ Wi,j Dj,j F (k) j,: ) + 2µ(F (k) i,: −Xi,:)
= p  N∑ j=1 M (k) i,j Di,i F (k) i,: − N∑ j=1 M (k) i,j√ Di,iDj,j F (k) j,: + 2µ(F(k)i,: −Xi,:) = p
 N∑ j=1 M (k) i,j Di,i + 2µ p F(k)i,: −  N∑ j=1 M (k) i,j√ Di,iDj,j F (k) j,: + 2µ p Xi,:  = p
α (k) i,i
F(k)i,: − α(k)i,i N∑
j=1
M (k) i,j√
Di,iDj,j F
(k) j,: + β (k) i,i Xi,:  = p
α (k) i,i
( F
(k) i,: − F (k+1) i,:
) ,
which indicates that
F (k) i,: − F (k+1) i,: =
α (k) i,i
p · ∂Lp(F (k))
∂F (k) i,:
.
For all i, j ∈ [N ],v ∈ R1×c, denote by
∂Lp(F(k)i,: ) := ∂Lp(F(k)) ∂F
(k) i,:
,
M ′(k) i,j := Wi,j ∥∥∥∥∥ √ Wi,j Di,i (F (k) i,: + v)− √ Wi,j Dj,j F (k) j,: ∥∥∥∥∥ p−2 ,
α ′(k) i,i := 1 / N∑ j=1 M ′(k) i,j Di,i + 2µ p  , β ′(k) i,i := 2µ
p α ′(k) i,i
F ′(k+1) i,: := α ′(k) i,i N∑ j=1 M ′(k) i,j√ Di,iDj,j F (k) j,: + β ′ i,iXi,:.
Then
∥∥∥∂Lp(F(k)i,: + v)− ∂Lp(F(k)i,: )∥∥∥ =
∥∥∥∥∥ pα′(k)i,i ( F (k) i,: + v − F ′(k+1) i,: ) − p α (k) i,i ( F (k) i,: − F (k+1) i,: )∥∥∥∥∥ ≤ p
α ′(k) i,i
∥v∥+ ∥∥∥∥∥ pα′(k)i,i ( F (k) i,: − F ′(k) i,: ) − p α (k) i,i ( F (k) i,: − F (k+1) i,: )∥∥∥∥∥ = p
α ′(k) i,i
∥v∥+ ∥∥∥∥∥ ( p
α ′(k) i,i − p α (k) i,i
) F
(k) i,: −
p
α ′(k) i,i
F ′(k+1) i,: +
p
α (k) i,i
F (k+1) i,: ∥∥∥∥∥ = p
α ′(k) i,i
∥v∥+ p ∥∥∥∥∥∥  N∑
j=1
M ′(k) i,j
Di,i − N∑ j=1 M (k) i,j Di,i F(k)i,: − N∑ j=1 M ′(k) i,j√ Di,iDj,j F (k) j,: − 2µ p Xi,: + N∑ j=1 M (k) i,j√ Di,iDj,j F (k) j,: + 2µ p Xi,: ∥∥∥∥∥∥ = p
α ′(k) i,i
∥v∥+ p ∥∥∥∥∥∥  N∑
j=1
M ′(k) i,j
Di,i − N∑ j=1 M (k) i,j Di,i F(k)i,: − N∑ j=1 M ′(k) i,j√ Di,iDj,j F (k) j,: + N∑ j=1 M (k) i,j√ Di,iDj,j F (k) j,: ∥∥∥∥∥∥ =
p N∑ j=1 M (k) i,j Di,i + 2µ  ∥v∥+ p N∑ j=1 M ′(k) i,j −M (k) i,j Di,i ∥v∥+ p ∥∥∥∥∥∥ N∑ j=1 M ′(k) i,j −M (k) i,j Di,i F (k) i,: − N∑ j=1 M ′(k) i,j −M (k) i,j√ Di,iDj,j F (k) j,: ∥∥∥∥∥∥ = p
 N∑ j=1 M (k) i,j Di,i + 2µ p + o (p,v,X,G)  ∥v∥.
Therefore, there exists some real positive value µ ∈ o (p,v,X,G) > 0 such that
∥∥∥∂Lp(F(k)i,: + v)− ∂Lp(F(k)i,: )∥∥∥ ≤ p( N∑ j=1 M (k) i,j Di,i + 2µ p ) ∥v∥ = p α (k) i,i ∥v∥. (19)
Let γ = (γ1, . . . , γN )⊤ ∈ RN and η ∈ RN×c. By Taylor’s theorem, we have:
Lp(F(k)i,: + γiηi,:) = Lp(F(k)i,: ) + γi ∫ 1 0 ⟨∂Lp(F(k)i,: + ϵγiηi,:),ηi,:⟩dϵ
= Lp(F(k)i,: ) + γi⟨ηi,:, ∂Lp(F (k) i,: )⟩+ γi ∫ 1 0 ⟨∂Lp(F(k)i,: + ϵγiηi,:)− ∂Lp(F (k) i,: ),ηi,:⟩dϵ
≤ Lp(F(k)i,: ) + γi⟨ηi,:, ∂Lp(F (k) i,: )⟩+ γi ∫ 1 0 ∥∂Lp(F(k)i,: + ϵγiηi,:)− ∂Lp(F (k) i,: )∥∥ηi,:∥dϵ
≤ Lp(F(k)i,: ) + γi⟨ηi,:, ∂Lp(F (k) i,: )⟩+
p
2α (k) i,i
γ2i ∥ηi,:∥2
Let η = −∇Lp(F(k)) and choose some positive real value µ which depends on X,G, p and p > 1, i.e. µ ∈ o (p,X,G). By Eq. (19), we have for all i ∈ [N ],
Lp(F(k)i,: − γ∂Lp(F (k) i,: )) ≤ Lp(F (k) i,: )− ⟨γi∂Lp(F (k) i,: ), ∂Lp(F (k) i,: )⟩+
p
2α (k) i,i
γ2i ∥∂Lp(F (k) i,: )∥ 2
= Lp(F(k)i,: )− p
2α (k) i,i
( 2α
(k) i,i γi p − γ2i
) ∥∂Lp(Fi,:)∥2
= Lp(F(k)i,: )− p
2α (k) i,i
 ( α (k) i,i )2 p2 − ( γi − α (k) i,i p )2 ∥∂Lp(F(k)i,: )∥2. Then for all i ∈ [N ], when 0 ≤ γi ≤ 2α (k) i,i p , we have Sp(F (k) i,: − γi∂Sp(F (k) i,: )) ≤ Sp(F (k) i,: ) and
γi = α
(k) i,i
p minimizes Sp(F (k) i,: − γi∂Sp(F (k) i,: )). Therefore,
Lp(F(k+1)) = Lp(F(k) − 1
p ·α(k)∇Lp(F(k))) ≤ Lp(F(k)).

heading: D.3 PROOF OF THEOREM 3
text: Proof. Without loss of generality, suppose F(0) = X. Denote M̃(k) = D−1/2M(k)D−1/2, by Eq. (14), we have for K ≥ 2,
F(K) = α(K−1)D−1/2M(K−1)D−1/2F(K−1) + β(K−1)X
= α(K−1)M̃K−1F(K−1) + β(K−1)X = α(K−1)M̃K−1 ( α(K−2)M̃K−2F(K−2) + β(K−2)X ) + β(K−1)X
= α(K−1)α(K−2)M̃K−1M̃K−2F(K−2) +α(K−1)M̃K−1β(K−2)X+ β(K−1)X
= ( K−1∏ k=0 α(k) )( K−1∏ k=0 M̃(k) ) F(0) + K−1∑ k=1 ( K−1∏ l=K−k α(l)M̃(l) ) β(K−1−k)X+ β(K−1)X
= ( K−1∏ k=0 α(k) )( K−1∏ k=0 M̃(k) ) X+ K−1∑ k=1 ( K−1∏ l=K−k α(l)M̃(l) ) β(K−1−k)X+ β(K−1)X.
(20)
Recall Equations (12) and (13), we have
M̃ (k) i,j = Wi,j√ Di,iDj,j ∥∥∥∥∥ √ Wi,j Di,i F (k) i,: − √ Wi,j Dj,j F (k) j,: ∥∥∥∥∥ p−2 , for all i, j = 1, 2, . . . , N,
and α (k) i,i =
1∑N j=1 M (k) i,j Di,i + 2µp , for all i = 1, 2, . . . , N.
Note that the eigenvalues of M̃ are not infinity and 0 < αi,i < 1 for all i = 1, . . . , N . Then we have
lim K→∞ K−1∏ k=0 α(k) = 0,
and
lim K→∞ ( K−1∏ k=0 α(k) )( K−1∏ k=0 M̃(k) ) = 0.
Therefore,
lim K→∞ F(K) = lim K→∞ ( K−1∑ k=1 ( K−1∏ l=K−k α(l)M̃(l) ) β(K−1−k)X+ β(K−1)X ) . (21)
By Equations (6) and (12), we have
∆pf(i) = N∑ j=1 Wi,j Di,i ∥(∇f)([j, i])∥p−2f(i)− N∑ j=1 Wi,j√ Di,iDj,j ∥(∇f)([j, i])∥p−2f(j)
= N∑ j=1 Mi,j Di,i f(i)− N∑ j=1 Mi,j√ Di,iDj,j f(j). (22)
By Eq. (13), we have N∑ j=1 M (k) i,j Di,i = 1 α (k) i,i − 2µ p . (23)
Equations (22) and (23) show that
∆(k)p =
(( α(k) )−1 − 2µ
p IN
) − M̃(k), (24)
which indicates α(k)M̃(k) = IN − 2µ
p α(k) −α(k)∆(k)p . (25)
Eq. (25) shows that α(k)M̃(k) is linear w.r.t ∆p and therefore can be expressed by a linear combination in terms of ∆p: α(k)M̃(k) = θ′(k)∆p, (26) where θ′ = diag(θ′0, θ ′ 1, . . . , θ ′ N−1) are the parameters. Therefore, we have
lim K→∞ F(K) = lim K→∞ ( K−1∑ k=1 ( K−1∏ l=K−k α(l)M̃(l) ) β(K−1−k)X+ β(K−1)X )
= lim K→∞ ( K−1∑ k=1 ( K−1∏ l=K−k θ′(l)∆p ) β(K−1−k)X+ β(K−1)X )
= lim K→∞ ( K−1∑ k=1 β(K−1−k) ( K−1∏ l=K−k θ′(l) ) ∆kpX+ β (K−1)X )
= lim K→∞ K−1∑ k=0 θ′′(k)∆kpX,
where θ′′(k) = diag(θ′′(k)1 , θ ′′(k) 2 , . . . , θ ′′(k) N ) defined as θ ′′(0) = β(K−1) and
θ ′′(k) i = β (K−1−k) i,i K−1∏ l=K−k θ ′(l) i , for k = 1, 2, . . . ,K − 1.
Let θ = (θ0, θ1, . . . , θK−1) defined as θk = ∑N i=1 θ ′′(k) i for all k = 0, 1, . . . ,K − 1, then
lim K→∞ F(K) = lim K→∞ ( K−1∑ k=1 θk∆ k pX+ θ0X )
= lim K→∞ K−1∑ k=0 θk∆ k pX.
Therefore complete the proof.

heading: D.4 PROOF OF THEOREM 4
text: Proof. The first-order Taylor expansion with Peano’s form of remainder for σ at X∗i,: is given by:
σ(F (K−1) j,: ) = σ(X ∗ i,:) +
∂σ(X∗i,:)
∂X
( F
(K−1) j,: −X ∗ i,: )⊤ + o(∥F(K−1)j,: −X ∗ i,:∥).
Note that in general the output non-linear layer σ(·) is simple. Here we assume that it can be well approximated by the first-order Taylor expansion and we can ignore the Peano’s form of remainder. For all i = 1, . . . , N , Di,i = Dj,j = d, we have αi,i ∑N j=1 D −1/2 i,i M (K−1) i,j D −1/2 j,j + β (K−1) i,i = 1. Then∣∣∣yi − ỹ(K)i ∣∣∣ =
∣∣∣∣∣∣yi − α(K−1)i,i N∑ j=1 M (K−1) i,j√ Di,iDj,j σ ( FK−1j,: ) − β(K−1)i,i σ (Xi,:) ∣∣∣∣∣∣ =
∣∣∣∣∣∣yi − β(K−1)i,i σ (Xi,:)− α(K−1)i,i N∑ j=1 M (K−1) i,j√ Di,iDj,j ( σ ( X∗i,: ) + ∂σ ( X∗i,: ) ∂X ( F (K−1) j,: −X ∗ i,: )⊤)∣∣∣∣∣∣ =
∣∣∣∣∣∣yi − α(K−1)i,i N∑ j=1 M (K−1) i,j√ Di,iDj,j yi − β(K−1)i,i σ(Xi,:)− α (K−1) i,i N∑ j=1 M (K−1) i,j√ Di,iDj,j ( ∂σ(X∗i,:) ∂X ( F (K−1) j,: −X ∗ i,: )⊤)∣∣∣∣∣∣ =
∣∣∣∣∣∣β(K−1)i,i (yi − σ(Xi,:))− α(K−1)i,i N∑ j=1 M (K−1) i,j√ Di,iDj,j ( ∂σ(X∗i,:) ∂X ( F (K−1) j,: −Xi,: − ϵi,: )⊤)∣∣∣∣∣∣ ≤ β(K−1)i,i |yi − σ(Xi,:)|+ α (K−1) i,i ∣∣∣∣∣∣ N∑ j=1 M (K−1) i,j√ Di,iDj,j ∂σ(X∗i,:) ∂X ( F (K−1) j,: −Xi,: )⊤∣∣∣∣∣∣+ (1− β(K−1)i,i ) ∣∣∣∣∂σ(X∗i,:)∂X ϵ⊤i,:
∣∣∣∣ ≤ β(K−1)i,i |yi − σ(Xi,:)|+ α (K−1) i,i ∥∥∥∥∂σ(X∗i,:)∂X ∥∥∥∥ ∥∥∥∥∥∥ N∑ j=1 M (K−1) i,j√ Di,iDj,j ( F (K−1) j,: −Xi,: )∥∥∥∥∥∥+ (1− β(K−1)i,i ) ∥∥∥∥∂σ(X∗i,:)∂X
∥∥∥∥ ∥ϵi,:∥ ≤ β(K−1)i,i |yi − σ(Xi,:)|+ α (K−1) i,i L ∥∥∥∥∥∥ N∑ j=1 M (K−1) i,j d ( F (K−1) j,: −Xi,:
)∥∥∥∥∥∥+ (1− β(K−1)i,i )L ∥ϵi,:∥ = β
(K−1) i,i |yi − σ(Xi,:)|+ α (K−1) i,i L ∥∥∥∥∥∥ N∑ j=1 M (K−1) i,j d ( F (K−1) j,: − F (K−1) i,: + F (K−1) i,: −Xi,: )∥∥∥∥∥∥ + (1− β(K−1)i,i )L ∥ϵi,:∥
= β (K−1) i,i |yi − σ(Xi,:)|+ α (K−1) i,i L ∥∥∥∥∥∥∆pF(K−1)i,: + N∑ j=1 M (K−1) i,j d · N∑ j=1 M (K−2) i,j d ( F (K−2) j,: −Xi,: )∥∥∥∥∥∥ + (1− β(K−1)i,i )L ∥ϵi,:∥
= β (K−1) i,i |yi − σ(Xi,:)|+ α (K−1) i,i L ∥∥∥∥∥∥∆(K−1)p F(K−1)i,: + K−2∑ k=0 K−2∏ l=k  N∑ j=1 M (l) i,j d ∆(k)p Xi,: ∥∥∥∥∥∥+ (1− β(K−1)i,i )L ∥ϵi,:∥
= β (K−1) i,i |yi − σ(Xi,:)|+ α (K−1) i,i L ∥∥∥∥∥∥∆(K−1)p F(K−1)i,: + K−2∑ k=0 K−2∏ l=k  N∑ j=1 M (l) i,j d ∆(k)p Xi,: ∥∥∥∥∥∥+ (1− β(K−1)i,i )L ∥ϵi,:∥
Therefore,
1
N N∑ i=1 |yi − ỹi| ≤ 1 N N∑ i=1 β (K−1) i,i |yi − σ(Xi,:)|
+ L
N N∑ i=1 α (K−1) i,i ∥∥∥∥∥∥∆(K−1)p F(K−1)i,: + K−2∑ k=0 K−2∏ l=k  N∑ j=1 M (l) i,j d ∆(k)p Xi,: ∥∥∥∥∥∥
+ L
N N∑ i=1 (1− β(K−1)i,i ) ∥ϵi,:∥

heading: D.5 PROOF OF THEOREM 6
text: Proof. Note that
ϕp(u) ⊤u = N∑ i=1 ϕp(ui)ui = N∑ i=1 ∥ui∥p−2u2i = N∑ i=1 ∥ui∥p = N∑ i=1 |ui|p = ∥u∥pp = 1,
then we have ∆pU = Φp(U)Λ = Φp(U)ΛΦ(U) ⊤U.
Therefore, ∆p = Φp(U)ΛΦp(U)⊤.
When p = 2, by Φ2(U) = U, we get ∆2 = Φ2(U)ΛΦ2(U)⊤ = UΛU⊤.

heading: D.6 PROOF OF THEOREM 7
text: Proof. By the definition of graph p-Laplacian, we have for all i = 1, 2, . . . , N ,
(∆pu)i = N∑ j=1 √ Wi,j Di,i ∥∥∥∥∥ √ Wi,j Di,i ui − √ Wi,j Dj,j uj ∥∥∥∥∥ p−2(√ Wi,j Di,i ui − √ Wi,j Dj,j uj ) = λϕp(ui).
Then, for all i = 1, 2, . . . , N ,
λ = 1
ϕp(ui) N∑ j=1 √ Wi,j Di,i ∥∥∥∥∥ √ Wi,j Di,i ui − √ Wi,j Dj,j uj ∥∥∥∥∥ p−2(√ Wi,j Di,i ui − √ Wi,j Dj,j uj )
= 1
∥ui∥p−2ui N∑ j=1 √ Wi,j Di,i ∥∥∥∥∥ √ Wi,j Di,i ui − √ Wi,j Dj,j uj ∥∥∥∥∥ p−2(√ Wi,j Di,i ui − √ Wi,j Dj,j uj )
= N∑ j=1 √ Wi,j Di,i ∥∥∥√Wi,jDi,i ui −√Wi,jDj,j uj∥∥∥p−2 ∥ui∥p−2 (√ Wi,j Di,i − √ Wi,j Dj,j uj ui )
= N∑ j=1 ( Wi,j Di,i − Wi,j√ Di,iDj,j uj ui ) ∥∥∥√Wi,jDi,i ui −√Wi,jDj,j uj∥∥∥ ∥ui∥ p−2
= N∑ j=1 ( Wi,j Di,i − Wi,j√ Di,iDj,j uj ui )∥∥∥∥∥∥ √ Wi,j Di,i ui − √ Wi,j Dj,j uj ui ∥∥∥∥∥∥ p−2
= N∑ j=1 ( Wi,j Di,i − Wi,j√ Di,iDj,j uj ui )∥∥∥∥∥ √ Wi,j Di,i − √ Wi,j Dj,j uj ui ∥∥∥∥∥ p−2
Let l = argmax{∥ui∥}i=1,2,...,N , the above equation holds for all i = 1, 2, . . . , N , then
λ = N∑ j=1 ( Wl,j Dl,l − Wl,j√ Dl,lDj,j uj ul )∥∥∥∥∥ √ Wl,j Dl,l − √ Wl,j Dj,j uj ul ∥∥∥∥∥ p−2
≥ N∑ j=1 ( Wl,j Dl,l − Wl,j√ Dl,lDj,j ∣∣∣∣ujul ∣∣∣∣ )∥∥∥∥∥ √ Wl,j Dl,l − √ Wl,j Dj,j uj ul ∥∥∥∥∥ p−2 ≥ N∑ j=1 ( Wl,j Dl,l − Wl,j√ Dl,lDj,j )∥∥∥∥∥ √ Wl,j Dl,l − √ Wl,j Dj,j uj ul ∥∥∥∥∥ p−2
≥ 0.
When p = 1,
λ = N∑ j=1 ( Wi,j Di,i − Wi,j√ Di,iDj,j uj ui )∥∥∥∥∥ √ Wi,j Di,i − √ Wi,j Dj,j uj ui ∥∥∥∥∥ −1
≤ N∑ j=1 √ Wi,j Di,i ≤ √√√√Ni N∑ j=1 Wi,j Di,i = √ Ni,
where the last inequality holds by using the Cauchy-Schwarz inequality. The above inequality holds for all i = 1, 2, . . . , N , therefore,
λ ≤ N∑ j=1 √ Wi,j Di,i ≤ √ Nmin.
When p > 1, we have for i = 1, 2, . . . , N ,
λ = N∑ j=1 ( Wi,j Di,i − Wi,j√ Di,iDj,j uj ui )∥∥∥∥∥ √ Wi,j Di,i − √ Wi,j Dj,j uj ui ∥∥∥∥∥ p−2
≤ N∑ j=1 √ Wi,j Di,i ∥∥∥∥∥ √ Wi,j Di,i − √ Wi,j Dj,j uj ui ∥∥∥∥∥ p−1 ≤ N∑ j=1 √ Wi,j Di,i ∥∥∥∥∥ √ Wi,j Di,i + √ Wi,j Dj,j ∣∣∣∣ujui ∣∣∣∣ ∥∥∥∥∥ p−1
= N∑ j=1 (√ Wi,j Di,i )p ∥∥∥∥∥1 + √ Di,i Dj,j ∣∣∣∣ujui ∣∣∣∣ ∥∥∥∥∥ p−1 .
Without loss of generality, let k = argmax({|ui|/ √ Di,i}i=1,2,...,N ). Because the above inequality holds for all i = 1, 2, . . . , N , then we have
λ ≤ N∑ j=1 (√ Wk,j Dk,k )p( 1 + √ Dk,k Dj,j ∣∣∣∣ujuk ∣∣∣∣ )p−1
≤ 2p−1 N∑ j=1 (√ Wk,j Dk,k )p .
For p ≥ 2,
λ ≤ 2p−1 N∑ j=1 (√ Wk,j Dk,k )p ≤ 2p−1 N∑ j=1 (√ Wk,j Dk,k )2 = 2p−1.
For 1 < p < 2,
λ ≤ 2p−1 N∑ j=1 (√ Wk,j Dk,k )p ≤ 2p−1 N∑ j=1 √ Wk,j Dk,k ≤ 2p−1 √√√√Nk N∑ j=1 Wk,j Dk,k = 2p−1 √ Nk.

heading: D.7 PROOF OF PROPOSITION 1
text: Proof. We proof Proposition 1 based on the bounds of p-eigenvalues as demonstrated in Thm. 7.
By Eq. (6) and Eq. (12), we have
∆pf(i) = N∑ j=1 Wi,j Di,i ∥(∇f)([j, i])∥p−2f(i)− N∑ j=1 Wi,j√ Di,iDj,j ∥(∇f)([j, i])∥p−2f(j)
= N∑ j=1 Mi,j Di,i f(i)− N∑ j=1 Mi,j√ Di,iDj,j f(j). (27)
By Eq. (13), we have N∑ j=1 M (k) i,j Di,i = 1 α (k) i,i − 2µ p . (28)
Equations (27) and (28) show that
∆(k)p =
(( α(k) )−1 − 2µ
p IN
) −D−1/2M(k)D−1/2, (29)
which indicates
α(k)D−1/2M(k)D−1/2 = IN − 2µ
p α(k) −α(k)∆(k)p . (30)
For i = 1, 2, . . . , N , let α̃ := (α̃1, . . . , α̃N ), α̃i := 1/ ∑N
j=1 Mi,j Di,i , then
αi,i N∑ j=1 Mi,j√ Di,iDj,j = (1− 2µ p αi,i)− αi,iλi
=
∑N j=1
Mi,j Di,i∑N
j=1 Mi,j Di,i + 2µp 1− 1∑N j=1 Mi,j Di,i λi  = 1
1 + 2µα̃ip (1− α̃iλi) , (31)
Recall the Eq. (12) that
Mi,j = Wi,j ∥∥∥∥∥ √ Wi,j Di,i Fi,: − √ Wi,j Dj,j Fj,: ∥∥∥∥∥ p−2 = Wi,j∥(∇f)([i, j])∥p−2,
1. When p = 2, for all i = 1, . . . , N , α̃i = 1 and 0 ≤ λi−1 ≤ 2, g2(λi−1) works as both low-pass and high-pass filters.
2. When p > 2, by Thm. 7 we have for all i = 1, . . . , N , 0 ≤ λi−1 ≤ 2p−1. If 0 ≤ α̃i ≤ 21−p, then 0 ≤ 1 − α̃iλi ≤ 1, which indicates that gp(λi−1) works as a low-pass filter; If
α̃i > 2 1−p, then gp(λi−1) works as both low-pass and high-pass filters. Since
N∑ j=1 Mi,j Di,i = N∑ j=1 Wi,j∥(∇f)([i, j])∥p−2 Di,i
≤ √√√√ N∑ j=1 ( Wi,j Di,i )2 N∑ j=1 ∥(∇f)([i, j])∥2(p−2)
≤ √√√√ N∑ j=1 ∥(∇f)([i, j])∥2(p−2)
≤ ∥∇f(i)∥p−2, which indicates that α̃i ≥ ∥∇f(i)∥2−p. 0 ≤ α̃i ≤ 21−p directly implies that 0 ≤ ∥∇f(i)∥2−p ≤ 21−p, i.e. ∥∇f(i)∥ ≥ 2(p−1)/(p−2) and when ∥∇f(i)∥2−p ≥ 21−p, i.e. ∥∇f(i)∥ ≤ 2(p−1)/(p−2), α̃i ≥ 21−p always holds. Therefore, if ∥∇f(i)∥ ≤ 2(p−1)/(p−2), gp(λi−1) works as both low-pass and high-pass filters on node i; If gp(λi−1) works as a low-pass filter, ∥∇f(i)∥ ≥ 2(p−1)/(p−2). 3. When 1 ≤ p < 2, by Thm. 7 we have for all i = 1, . . . , N , 0 ≤ λi−1 ≤ 2p−1 √ Nk. If
0 ≤ α̃i ≤ 21−p/ √ Nk, 0 ≤ 1− α̃iλi ≤ 1, which indicates that gp(λi−1) work as low-pass filters; If α̃i ≥ 21−p/ √ Nk, gp(λi−1) work as both low-pass and high-pass filters. By
N∑ j=1 Wi,j Di,i
1
∥(∇f)([i, j])∥p−2
/ 1∑N
j=1 Wi,j Di,i ∥(∇f)([i, j])∥p−2
= N∑ j=1 Wi,j Di,i
1
∥(∇f)([i, j])∥p−2 · N∑ j=1 Wi,j Di,i ∥(∇f)([i, j])∥p−2
≥ N∑ j=1 Wi,j Di,i
( 1
∥(∇f)([i, j])∥p−2 · ∥(∇f)([i, j])∥p−2 ) = 1,
we have
α̃i = 1∑N
i=1 Mi,j Di,i
= 1∑N
j=1 Wi,j Di,i ∥(∇f)([i, j])∥p−2
≤ N∑ j=1 Wi,j Di,i
1
∥(∇f)([i, j])∥p−2
= N∑ j=1 Wi,j Di,i ∥(∇f)([i, j])∥2−p
≤ ∥∇f(i)∥2−p. α̃i ≥ 21−p/ √ Nk directly implies that ∥∇f(i)∥2−p ≥ 21−p/ √ Nk, i.e. ∥∇f(i)∥ ≥ 2(2 √ Nk) 1/(p−2) and when 0 ≤ ∥∇f(i)∥2−p ≤ 21−p/ √ Nk, i.e. 0 ≤ ∥∇f(i)∥ ≤ 2(2 √ Nk) 1/(p−2), 0 ≤ α̃i ≤ 21−p/ √ Nk always holds. Therefore, if 0 ≤ ∥∇f(i)∥ ≤ 2(2 √ Nk) 1/(p−2), gp(λi−1) work as low-pass filters; If gp(λi−1) work as both low-pass
and high-pass filters, ∥∇f(i)∥ ≥ 2 ( 2 √ Nk )1/(p−2) .
Specifically, when p = 1, by Thm. 7 we have for all i = 1, . . . , N , 0 ≤ λi−1 ≤ 2p−1 √ Nmin. Following the same derivation above we attain if 0 ≤ ∥∇f(i)∥ ≤ 2(2 √ Nmin) 1/(p−2), gp(λi−1) work as low-pass filters; If gp(λi−1) work as both low-pass
and high-pass filters, ∥∇f(i)∥ ≥ 2 ( 2 √ Nmin )1/(p−2) .

heading: E DATASET STATISTICS AND HYPERPARAMETERS
text: E.1 ILLUSTRATION OF GRAPH GRADIENT AND GRAPH DIVERGENCE

heading: E.2 DATASET STATISTICS
text: Table 2 summarizes the dataset statistics and the levels of homophily H(G) of all benchmark datasets. Note that the homophily scores here is different with the scores reported by Chien et al. (2021). There is a bug in their code when computing the homophily scores (doing division with torch integers) which caused their homophily scores to be smaller.

heading: E.3 HYPERPARAMETER SETTINGS
text: We set the number of layers as 2, the maximum number of epochs as 1000, the number for early stopping as 200, the weight decay as 0 or 0.0005 for all models. The other hyperparameters for each model are listed as below:
• 1.0GNN, 1.5GNN, 2.0GNN, 2.5GNN: – Number of hidden units: 16 – Learning rate: {0.001, 0.01, 0.05} – Dropout rate: {0, 0.5} – µ: {0.01, 0.1, 0.2, 1, 10} – K: 4, 6, 8
• MLP: – Number of hidden units: 16 – Learning rate: {0.001, 0.01} – Dropout rate: {0, 0.5}
• GCN: – Number of hidden units: 16 – Learning rate: {0.001, 0.01} – Dropout rate: {0, 0.5}
• SGC: – Number of hidden units: 16 – Learning rate: {0.2, 0.01} – Dropout rate: {0, 0.5} – K: 2
• GAT: – Number of hidden units: 8 – Number of attention heads: 8 – Learning rate: {0.001, 0.005} – Dropout rate: {0, 0.6}
• JKNet: – Number of hidden units: 16 – Learning rate: {0.001, 0.01} – Dropout rate: {0, 0.5} – K: 10 – α: {0.1, 0.5, 0.7, 1} – The number of GCN based layers: 2 – The layer aggregation: LSTM with 16 channels and 4 layers
• APPNP: – Number of hidden units: 16 – Learning rate: {0.001, 0.01} – Dropout rate: {0, 0.5} – K: 10 – α: {0.1, 0.5, 0.7, 1}
• GPRGNN: – Number of hidden units: 16 – Learning rate: {0.001, 0.01, 0.05} – Dropout rate: {0, 0.5} – K: 10 – α: {0, 0.1, 0.2, 0.5, 0.7, 0.9, 1} – dprate: {0, 0.5, 0.7}

heading: F ADDITIONAL EXPERIMENTS

heading: F.1 EXPERIMENTAL RESULTS ON HOMOPHILIC BENCHMARK DATASETS
text: Competitive Performance on Real-World Homophilic Datasets. Table 3 summarizes the averaged accuracy (the micro-F1 score) and standard deviation of semi-supervised node classification on homophilic benchmark datasets. Table 3 shows that the performance of pGNN is very close to APPNP, JKNet, GCN on Cora, CiteSeer, PubMed datasets and slightly outperforms all baselines on Computers, Photo, CS, Physics datasets. Moreover, we observe that pGNNs outperform GPRGNN on all homophilic datasets, which confirms that pGNNs work better under weak supervised information (2.5% training rate) as discussed in Remark 3. We also see that all GNN models work significantly better than MLP on all homophilic datasets. It illustrates that the graph topological information is helpful for the label prediction tasks. Notably, 1.0GNN is slightly worse than the other pGNNs with larger p, which suggests to use p ≈ 2 for homophilic graphs. Overall, the results of Table 3 indicates that pGNNs obtain competitive performance against all baselines on homophilic datasets.

heading: F.2 EXPERIMENTAL RESULTS OF AGGREGATION WEIGHT ENTROPY DISTRIBUTION
text: Here we present the visualization results of the learned aggregation weight entropy distribution of pGNNs and GAT on all benchmark datasets. Fig. 5 and Fig. 6 show the results obtained on homophilic and heterophilic benchmark datasets, respectively.
We observe from Fig. 5 that the aggregation weight entropy distributions learned by pGNNs and GAT on homophilic benchmark datasets are similar to the uniform cases, which indicates that aggregating and transforming node features over the original graph topology is very helpful for label prediction. It explains why pGNNs and GNN baselines obtained similar performance on homophilic benchmark datasets and all GNN models significantly outperform MLP.
Contradict to the results on homophilic graphs shown in Fig. 5, Fig. 6 shows that the aggregation weight entropy distributions of pGNNs on heterophilic benchmark datasets are very different from that of GAT and the uniform cases. We observe from Fig. 6 that the entropy of most of the aggregation weights learned by pGNNs are around zero, which means that most aggregation weights are on one source node. It indicates that the graph topological information in these heterophilic benchmark graphs is not helpful for label prediction. Therefore, propagating and transforming node features over the graph topology could lead to worse performance than MLPs, which validates the results in Table 3 that the performance of MLP is significantly better most GNN baselines on all heterophilic graphs and closed to pGNNs.
0 2 4 6 8 entropy
0
200
400
600
800
1000
1200
1400
# of
n od
es uniform
0 2 4 6 8 entropy
1.0GNN (Acc: 78.43%)
0 2 4 6 8 entropy
1.5GNN (Acc: 78.58%)
0 2 4 6 8 entropy
2.5GNN (Acc: 78.39%)
0 2 4 6 8 entropy
GAT (Acc: 73.86%)
Cora
1 0 1 2 3 4 5 entropy
0
250
500
750
1000
1250
1500
1750
2000
# of
n od
es
uniform
1 0 1 2 3 4 5 entropy
1.0GNN (Acc: 62.95%)
1 0 1 2 3 4 5 entropy
1.5GNN (Acc: 64.87%)
1 0 1 2 3 4 5 entropy
2.5GNN (Acc: 63.32%)
1 0 1 2 3 4 5 entropy
GAT (Acc: 61.59%)
Citeseer
0 2 4 6 8 entropy
0
2000
4000
6000
8000
10000
12000
14000
# of
n od
es
uniform
0 2 4 6 8 entropy
1.0GNN (Acc: 83.2%)
0 2 4 6 8 entropy
1.5GNN (Acc: 83.73%)
0 2 4 6 8 entropy
2.5GNN (Acc: 84.35%)
0 2 4 6 8 entropy
GAT (Acc: 83.85%)
Pubmed
0 2 4 6 8 10 12 14 entropy
0
1000
2000
3000
4000
5000
6000
7000
8000
# of
n od
es
uniform
0 2 4 6 8 10 12 14 entropy
1.0GNN (Acc: 85.24%)
0 2 4 6 8 10 12 14 entropy
1.5GNN (Acc: 85.3%)
0 2 4 6 8 10 12 14 entropy
2.5GNN (Acc: 83.68%)
0 2 4 6 8 10 12 14 entropy
GAT (Acc: 82.73%)
Computers
Photo

heading: Physics
text: 1 0 1 2 3 4 5 entropy
0
250
500
750
1000
1250
1500
1750
2000
# of
n od
es
uniform
1 0 1 2 3 4 5 entropy
1.0GNN (Acc: 48.8%)
1 0 1 2 3 4 5 entropy
1.5GNN (Acc: 48.8%)
1 0 1 2 3 4 5 entropy
2.5GNN (Acc: 49.23%)
1 0 1 2 3 4 5 entropy
GAT (Acc: 44.64%)
Chameleon
1 0 1 2 3 4 5 6 entropy
0
1000
2000
3000
4000
5000
# of
n od
es
uniform
1 0 1 2 3 4 5 6 entropy
1.0GNN (Acc: 32.08%)
1 0 1 2 3 4 5 6 entropy
1.5GNN (Acc: 28.15%)
1 0 1 2 3 4 5 6 entropy
2.5GNN (Acc: 33.33%)
1 0 1 2 3 4 5 6 entropy
GAT (Acc: 28.72%)
Squirrel
1 0 1 2 3 4 5 entropy
0
1000
2000
3000
4000
5000
6000
# of
n od
es
uniform
1 0 1 2 3 4 5 entropy
1.0GNN (Acc: 41.21%)
1 0 1 2 3 4 5 entropy
1.5GNN (Acc: 40.63%)
1 0 1 2 3 4 5 entropy
2.5GNN (Acc: 39.46%)
1 0 1 2 3 4 5 entropy
GAT (Acc: 31.12%)
Actor
1 0 1 2 3 4 5 entropy
0
25
50
75
100
125
150
175
200
# of
n od
es
uniform
1 0 1 2 3 4 5 entropy
1.0GNN (Acc: 98.15%)
1 0 1 2 3 4 5 entropy
1.5GNN (Acc: 96.3%)
1 0 1 2 3 4 5 entropy
2.5GNN (Acc: 84.26%)
1 0 1 2 3 4 5 entropy
GAT (Acc: 62.96%)
Wisconsin
1 0 1 2 3 4 5 entropy
0
20
40
60
80
100
120
140
# of
n od
es
uniform
1 0 1 2 3 4 5 entropy
1.0GNN (Acc: 86.25%)
1 0 1 2 3 4 5 entropy
1.5GNN (Acc: 83.75%)
1 0 1 2 3 4 5 entropy
2.5GNN (Acc: 86.25%)
1 0 1 2 3 4 5 entropy
GAT (Acc: 48.75%)
Texas

heading: Cornell

heading: F.3 EXPERIMENTAL RESULTS ON CSBM
text: In this section we present the experimental results on cSBM using sparse splitting and dense splitting, respectively. We used the same settings in Chien et al. (2021) in which the number of nodes n = 5000, the number of features f = 2000, ϵ = 3.25 for all experiments. Table 4 reports the results on cSBM with sparse splitting setting, which also are presented in Fig. 2 and discussed in Sec. 5. Table 5 reports the results on cSBM with dense splitting settings.
1.0GNN 95.75±1.21 93.06±1.13 77.39±4.21 61.38±0.39 61.80±0.29 65.73±2.11 85.85±3.24 96.80±0.87 97.40±1.10 1.5GNN 95.90±3.01 94.10±4.57 73.08±2.59 61.44±0.30 61.77±0.35 66.01±1.88 90.57±0.71 97.38±0.43 97.76±0.86 2.0GNN 98.37±0.78 96.32±1.50 84.93±0.39 61.13±0.51 61.79±0.34 63.55±1.73 88.55±1.05 97.56±0.16 97.94±0.39 2.5GNN 97.74±0.99 96.78±0.44 83.21±2.12 61.30±0.41 61.74±0.34 62.88±2.31 79.64±2.15 95.71±0.34 97.25±0.58
1.0GNN 98.19±0.28 94.38±0.44 86.40±1.00 80.57±0.43 80.21±0.42 87.32±0.47 92.42±0.62 97.52±0.33 98.37±0.26 1.5GNN 98.88±0.16 95.62±0.21 86.87±1.22 80.70±0.71 80.28±0.31 86.29±0.43 92.40±0.24 97.56±0.25 98.24±0.32 2.0GNN 99.21±0.09 96.91±0.16 92.96±0.31 80.83±0.61 80.04±0.49 84.96±0.60 91.18±0.27 97.41±0.14 98.45±0.14 2.5GNN 99.21±0.14 96.94±0.16 93.28±0.37 80.93±0.44 80.28±0.38 83.83±0.70 86.10±0.39 96.28±0.43 97.76±0.18
Table 5 shows that pGNNs obtain the best performance on weak homophilic graphs (ϕ = 0, 0.25) while competitive performance against GPRGNN on strong heterophilic graphs (ϕ = −0.75,−1) and competitive performance with state-of-the-art GNNs on strong homophilic graphs (ϕ = 0.75, 1). We also observe that GPRGNN is slightly better than pGNNs on weak heterophilic graphs (ϕ = −0.25,−0.5), which suggests that GPRGNN could work very well using strong supervised information (60% training rate and 20% validation rate). However, as shown in Table 4, pGNNs work better than GPRGNN under weak supervised information (2.5% training rate and 2.5%) on all heterophilic graphs. The result is reasonable, as discussed in Remark 3 in Sec. 3.2, GPRGNN can adaptively learn the generalized PageRank (GPR) weights and it works similarly to 2.0GNN on both homophilic and heterophilic graphs. However, it needs more supervised information in order to learn optimal GPR weights. On the contrary, pGNNs need less supervised information to obtain similar results because Θ(2) acts like a hyperplane for classification. Therefore, pGNNs can work better under weak supervised information.

heading: F.4 EXPERIMENTAL RESULTS ON GRAPHS WITH NOISY EDGES
text: Here we present more experimental results on graph with noisy edges. Table 6 reports the results on homophilic graphs (Computers, Photo, CS, Physics) and Table 7 reports the results on heterophilic graphs (Wisconsin Texas). We observe from Tables 6 and 7 that pGNNs dominate all baselines. Moreover, pGNNs even slightly better than MLP when the graph topologies are completely random, i.e. the noisy edge rate r = 1. We also observe that the performance of GCN, SGC, JKNet on homophilic graphs dramatically degrades as the noisy edge rate r increases while they do not change
a lot for the cases on heterophilic graphs. It is reasonable since the original graph topological information is very helpful for label prediction on these homophilic graphs. Adding noisy edges and remove the same number of original edges could significantly degrade the performance of ordinary GNNs. On the other hand, since we find that the original graph topological information in Wisconsin and Texas is not helpful for label prediction. Therefore, adding noisy edges and removing original edges on these heterophilic graphs would not affect too much their performance.
F.5 EXPERIMENTAL RESULTS OF INTERGRATING pGNNS WITH GCN AND JKNET
Here we further conduct experiments to study whether pGNNs can be intergrated into existing GNN architectures and improve their performance on heterophilic graphs. We use two popular GNN architectures: GCN (Kipf & Welling, 2017) and JKNet (Xu et al., 2018).
To incorporate pGNNs with GCN, we use the pGNN layers as the first layer of the combined models, termed as pGNN + GCN, and GCN layer as the second layer. Specifically, we use the aggregation weights αD−1/2MD−1/2 learned by the pGNN in the first layer as the input edge weights of GCN layer in the second layer. To combine pGNN with JKNet, we use the pGNN layer as the GNN layers in the JKNet framework, termed as pGNN + JKNet. Tables 8 and 9 report the experimental results on homophilic and heterophilic benchmark datasets, respectively.
We observe from Table 8 that intergrating pGNNs with GCN and JKNet does not improve their performance on homophilic graphs. The performance of GCN slightly degrade after incorporating pGNNs. The performance of JKNet also slightly degrade on Cora, CiteSeer, and PubMed but is improved on Computers, Photo, CS, Physics. It is reasonable since GCN and JKNet can predict well on these homophilic benchmark datasets based on their original graph topology.
However, for heterophilic benchmark datasets, Table 9 shows that there are significant improvements over GCN, and JKNet after intergrating with pGNNs. Moreover, pGNNs + JKNet obtain advanced performance on all heterophilic benchmark datasets and even better than pGNNs on Squirrel. The results of Table 9 demonstrate that intergrating pGNNs with GCN and JKNet can sigificantly improve their performance on heterophilic graphs.
F.6 EXPERIMENTAL RESULTS OF pGNNS ON PPI DATASET FOR INDUCTIVE LEARNING
Additionally, we conduct comparison experiments of pGNNs against GAT on PPI dataset (Zitnik & Leskovec, 2017) using the inductive learning settings as in Velickovic et al. (2018) (20 graphs for training, 2 graphs for validation, 2 graphs for testing). We use three layers of GAT architecture with 256 hidden units, use 1 attention head for GAT (1 head) and 4 attention heads for GAT (4 heads). We use three pGNN layers and a MLP layer as the first layer for pGNNs, set µ = 0.01, K = 1, and use 256 hidden units for pGNN-256 and 512 hidden units for pGNN-512. The experimental results are reported in Table 10.
From Appendix F.6 we observe that the results of 2.5GNN on PPI slightly better than GAT with 4 attention heads and other pGNNs are very close to it. Moreover, all results of pGNNs significantly outperform GAT with one attention head. The results of pGNNs on PPI is impressive. pGNNs have much less parameters than GAT with 4 attention heads while obtain very completitive performance on PPI. When we use more hidden units, 512 hidden units, pGNNs-512 significantly outperform GAT, while pGNNs-512 still have less parameters. It illustrates the superior potential of applying pGNNs to inducting learning on graphs.
F.7 EXPERIMENTAL RESULTS OF pGNNS ON OGBN ARXIV DATASET
Here we present the experimental of pGNNs on OGBN arXiv dataset (Hu et al., 2020). We use the official data split setting of OGBN arXiv. We use three layers pGNN architecture and 256 hidden units with µ = 0.5, K = 2. We also combine pGNNs with correct and smooth model (C&S) (Huang et al., 2021) and introduce residual units. The results of MLP, GCN, JKNet, DeepGCN (Li et al.,
2019), GCN with residual units, C&S model are extracted from the leaderboard for OGBN arXiv. dataset5. Table 11 summaries the results of pGNNs against the baselines.
We observe from Table 11 that pGNNs outperform MLP, GCN, JKNet, and DeepGCN. The performance of pGNNs can be further improved by combining it with C&S model and residual units and 3GNN + residual (6 layers) + C&S obtains the best performance against the baselines.
F.8 RUNNING TIME OF pGNNS
Tables 12 and 13 report the averaged running time of pGNNs and baselines on homophilic and heterophilic benchmark datasets, respectively.
GPRGNN 86.5 ms / 12.42s 195.8 ms / 12.6s 88.6 ms / 12.59s 93.3 ms / 15.98s 86.7 ms / 12.65s 92 ms / 17.8s 217.1 ms / 26.33s 1.0GNN 96 ms / 20.12s 98.1 ms / 19.81s 100.2 ms / 21.74s 151.4 ms / 64.08s 121.3 ms / 34.07s 109.7 ms / 25.03s 122.9 ms / 49.59s 1.5GNN 98.2 ms / 20.19s 97 ms / 20.26s 100.2 ms / 22.6s 140.3 ms / 64.08s 120 ms / 34.22s 112.3 ms / 25.11s 127.9 ms / 49.54s 2.0GNN 98.1 ms / 20.11s 96.3 ms / 19.97s 99.3 ms / 22.17s 141 ms / 64.04s 129.3 ms / 34.14s 104.7 ms / 24.93s 124.6 ms / 49.35s 2.5GNN 96.6 ms / 20.12s 92.9 ms / 20.16s 103 ms / 22.17s 141.6 ms / 64.01s 128.1 ms / 34.22s 110.8 ms / 25.07s 124 ms / 49.39s

heading: F.9 EXPERIMENTAL RESULTS ON BENCHMARK DATASETS FOR 64 HIDDEN UNITS
text: F.10 TRAINING CURVES FOR p = 1
F.11 VISUALIZATION RESULTS OF NODE EMBEDDINGS

Review summaries:
This paper proposes a p-Laplacian based GNN to handle heterophilic graphs and graphs with non-informative topologies. Both the above cases are assumed in most existing GNN architectures and hence this work breaks away from the norm. This work proposes a discrete p-Laplacian message passing scheme which is derived from a discrete regularization framework. The authors do a spectral analysis of their novel p-Laplacian message passing scheme and show that it works as both a low-pass and high-pass filter, which is then applicable to both homophilic and heterophilic graphs. More specifically, they show that p-GNN with p>2, works well for graphs which exhibit strong homophily, while for p \in [1,2) p-GNN works effectively on heterophilic graphs. The empirical results support the theoretical justifications and outperform the baselines quite significantly especially on heterophilic and non-informative topology bearing graphs.

The paper derives the p-Laplacian message passing formula under the p-Laplacian based regularization framework and further proposes p-GNN architecture. Authors further justify the relations of p-Laplacian message passing with low and high-pass filters and the upper bound of one layer risk of p-GNNs. Experiments show the superiority of p-GNNs on both heterophilic and homophilic settings. However, we still have some concerns for the paper before further evaluation.