# Equivariant Hypergraph Diffusion Neural Operators

Peihao Wang1∗, Shenghao Yang 2, Yunyu Liu 3, Zhangyang Wang 1, Pan Li3
1University of Texas at Austin, 2University of Waterloo, 3Purdue University
<EMAIL_ADDRESS><EMAIL_ADDRESS>
<EMAIL_ADDRESS>Correspondence to: Peihao Wang and Pan Li.

###### Abstract

Hypergraph neural networks (HNNs) using neural networks to encode hypergraphs
provide a promising way to model higher-order relations in data and further
solve relevant prediction tasks built upon such higher-order relations.
However, higher-order relations in practice contain complex patterns and are
often highly irregular. So, it is often challenging to design an HNN that
suffices to express those relations while keeping computational efficiency.
Inspired by hypergraph diffusion algorithms, this work proposes a new HNN
architecture named ED-HNN, which provably approximates any continuous
equivariant hypergraph diffusion operators that can model a wide range of
higher-order relations. ED-HNN can be implemented efficiently by combining
star expansions of hypergraphs with standard message passing neural networks.
ED-HNN further shows great superiority in processing heterophilic hypergraphs
and constructing deep models. We evaluate ED-HNN for node classification on
nine real-world hypergraph datasets. ED-HNN 111Our code is available at:
https://github.com/Graph-COM/ED-HNN. uniformly outperforms the best baselines
over these nine datasets and achieves more than 2%$\uparrow$ in prediction
accuracy over four datasets therein.

##  1 Introduction

Machine learning on graphs has recently attracted great attention in the
community due to the ubiquitous graph-structured data and the associated
inference and prediction problems [1, 2, 3, 4, 5, 6]. Current works primarily
focus on graphs which can model only pairwise relations in data. Emerging
research has shown that higher-order relations that involve more than two
entities often reveal more significant information in many applications [7, 8,
9, 10, 11]. For example, higher-order network motifs build the fundamental
blocks of many real-world networks [12, 13, 14, 15, 16]. Session-based (multi-
step) behaviors often indicate the preferences of web users in more precise
ways [17, 18, 19, 20]. To capture these higher-order relations, hypergraphs
provide a dedicated mathematical abstraction [21]. However, learning
algorithms on hypergraphs are still far underdeveloped as opposed to those on
graphs.

Recently, inspired by the success of graph neural networks (GNNs), researchers
have started investigating hypergraph neural network models (HNNs) [22, 23,
24, 25, 26, 27, 28, 29, 30]. Compared with GNNs, designing HNNs is more
challenging. First, as aforementioned, higher-order relations modeled by
hyperedges could contain very complex information. Second, hyperedges in real-
world hypergraphs are often of large and irregular sizes. Therefore, how to
effectively represent higher-order relations while efficiently processing
those irregular hyperedges is the key question to answer when one is to design
HNNs. Unfortunately, just as other types of neural networks (NNs), HNNs look
like a black box. Given an HNN, one can hardly tell the type of relations that
it can represent and thus often fails to make a suitable tradeoff between the
expressiveness and the computational efficiency of the model.

In this work, inspired by the recently developed hypergraph diffusion
algorithms [31, 32, 33, 34, 35, 36, 37, 38], we design a novel HNN
architecture that holds provable expressiveness to approximate a large class
of hypergraph diffusion while keeping computational efficiency. Hypergraph
diffusion is significant due to its transparency and has been widely applied
to clustering [31, 16, 39], semi-supervised learning [32, 37, 38, 40], ranking
aggregation [16, 41], network analysis [33, 34, 35, 36] and signal processing
[42, 8, 43, 44]. However, traditional hypergraph diffusion needs to first
handcraft potential functions to model higher-order relations and then use
their gradients or some variants of their gradients as the diffusion operators
to characterize the exchange of diffused quantities on the nodes within one
hyperedge. The design of those potential functions often requires significant
insights into the applications, which may not be available in practice.

We observe that the most commonly-used hyperedge potential functions are
permutation invariant, which covers the applications where none of the nodes
in a higher-order relation are treated as inherently special. For such
potential functions, we further show that their induced diffusion operators
must be permutation equivariant. Inspired by this observation, we propose a
NN-parameterized permutation-equivariant architecture that is expressive to
provably represent any continuous hyperedge diffusion operators, whose NN
parameters can be learned in a data-driven way. We also introduce an efficient
implementation of the architecture that can be easily implemented over current
GNN platforms [45, 46]: We just need to combine a bipartite representation (or
star expansion [47, 48], equivalently) of hypergraphs and the standard message
passing neural network (MPNN) [49]. By repeating this architecture by layers
with shared parameters, we finally obtain our model named Equivariant
Diffusion-based HNN (ED-HNN). Fig. 1 shows an illustration of hypergraph
diffusion and the key architecture in ED-HNN.

To the best of our knowledge, we are the first one to establish the connection
between the general class of hypergraph diffusion algorithms and the design of
HNNs. Previous HNNs were either less expressive to represent equivariant
diffusion operators [22, 23, 24, 25, 26, 27] or needed to learn the
representations by adding significantly many auxiliary nodes [28, 29, 30]. We
provide detailed discussion of them in Sec. 3.4. We also show that due to the
equivariant architecture, ED-HNN is by design good at predicting node labels
over heterophilic hypergraphs where hyperedges mix nodes from different
classes. Moreover, ED-HNN can go very deep without much performance decay.

As a by-product, our proof of expressiveness does not use equivariant
polynomials as a bridge, which were used previously in set-function
approximation [50, 51, 52]. Our statement allows precise representations of
continuous equivariant set functions (when each set element has 1-dim feature
channel) by compositing a continuous function and the sum of another
continuous function on each set entry, while previous works have only achieved
an approximation result. This result may be of independent interest for the
community.

We evaluate ED-HNN by performing node classsification over 9 real-world
hypergraph datasets. ED-HNN uniformly outperforms all baseline methods in
averaged performance and achieves significant improvement (more than
2$\%\uparrow$) over 4 datasets therein. ED-HNN also shows super robustness
when going deep. We further carefully design synthetic experiments to
demonstrate the expressiveness of ED-HNN and the superiority of ED-HNN in
processing heterphilic hypergraphs.

Figure 1: Hypergraph diffusion often uses permutation-invariant hyperedge
potentials to model higher-order relations. The gradients or their variants of
those potentials determine the diffusion process and are termed diffusion
operators. Our ED-HNN can universally represent such operators by feeding node
representations into the message from hyperedges to nodes.

##  2 Preliminaries: Hypergraphs and Hypergraph Diffusion

In this section, we formulate the hypergraph diffusion problem, along the way,
introduce the notations.

######  Definition 1 (Hypergraph).

Let $\mathcal{G}=(\mathcal{V},\mathcal{E},\bm{X})$ be an attributed hypergraph
where $\mathcal{V}$ is the node set of size $N$ and $\mathcal{E}$ is the
hyperedge set of size $M$. Each hyperedge $e=\\{v_{1}^{(e)},...,v_{\lvert
e\rvert}^{(e)}\\}$ is a subset of $\mathcal{V}$. Unlike graphs, a hyperedge
may contain more than two nodes.
$\bm{X}=[...,\bm{x}_{v},...]^{T}\in\mathbb{R}^{N}$ denotes the node attributes
and $\bm{x}_{v}$ denotes the attribute of node $v$. Define
$d_{v}=|\\{e\in\mathcal{E}:v\in e\\}|$ as the degree of node $v$. Let
$\bm{D},\,\bm{D}_{e}$ denote the diagonal degree matrix for $v\in\mathcal{V}$
and the sub-matrix for $v\in e$.

Here, we use 1-dim attributes for convenient discussion while our experiments
often have multi-dim attributes. Learning algorithms will combine attributes
and hypergraph structures into (latent) features defined as follows. These
features can be further used to make prediction for downstream tasks.

######  Definition 2 (Latent features).

Let $\bm{h}_{v}\in\mathbb{R}$ denote the (latent) features of node
$v\in\mathcal{V}$. The feature vector
$\bm{H}=\begin{bmatrix}...,\bm{h}_{v},...\end{bmatrix}_{v\in\mathcal{V}}^{\top}\in\mathbb{R}^{N}$
includes node features as entries. Further, collect the features into a edge
feature vector $\bm{H}_{e}=\begin{bmatrix}\bm{h}_{e,1}&\cdots&\bm{h}_{e,\lvert
e\rvert}\end{bmatrix}^{\top}\in\mathbb{R}^{\lvert e\rvert}$, where
$\bm{h}_{e,i}$ corresponds to the feature $\bm{h}_{v_{i}^{(e)}}$ for some node
$v_{i}^{(e)}\in e$. Note that for any $v\in e$, there is one corresponding
index $i\in\\{1,...,\lvert e\rvert\\}$. Later, we may use subscripts $(e,i)$
and $(e,v)$ interchangeably if they cause no confusion.

A widely-used heuristic to generate features $\bm{H}$ is via hypergraph
diffusion algorithms.

######  Definition 3 (Hypergraph Diffusion).

Define node potential functions
$f(\cdot;\bm{x}_{v}):\mathbb{R}\rightarrow\mathbb{R}$ for $v\in\mathcal{V}$
and hyperedge potential functions
$g_{e}(\cdot):\mathbb{R}^{|e|}\rightarrow\mathbb{R}$ for each
$e\in\mathcal{E}$. The hypergrah diffusion combines the node attributes and
the hypergraph structure and asks to solve

$\displaystyle\min_{\bm{H}}\sum_{v\in\mathcal{V}}f(\bm{h}_{v};\bm{x}_{v})+\sum_{e\in\mathcal{E}}g_{e}(\bm{H}_{e}).$
(1)

In practice, $g_{e}$ is often shared across hyperedges of the same size.
Later, we ignore the subscript $e$.

The two potential functions are often designed via heuristics in traditional
hypergraph diffusion literatures. Node potentials often correspond to some
negative-log kernels of the latent features and the initial attributes. For
example, $f(\bm{h}_{v};\bm{x}_{v})$ could be $(\bm{h}_{v}-\bm{x}_{v})_{2}^{2}$
when to compute hypergraph PageRank diffusion [33, 36] or
$-\bm{x}_{v}\bm{h}_{v}$ when to compute hypergraph min-cut [53, 54, 55].
Hyperedge potentials are more significant and complex, as they need to model
those higher-order relations between more than two objects. Here list a few
examples.

######  Example 1 (Hyperedge potentials).

Some practical $g(\bm{H}_{e})$ may be chosen as follows.

  * •

Clique Expansion (CE, hyperedges reduced to cliques) plus pairwise potentials
[31]: $\sum_{u,v\in e}(\bm{h}_{v}-\bm{h}_{u})_{2}^{2}$ or with properly degree
normalization $\sum_{u,v\in
e}(\frac{\bm{h}_{v}}{\sqrt{d_{v}}}-\frac{\bm{h}_{u}}{\sqrt{d_{u}}})_{2}^{2}\,(\triangleq
g(\bm{D}_{e}^{-1/2}\bm{H}_{e}))$.

  * •

Divergence to the mean [38, 56]: $\sum_{v\in
e}(\bm{h}_{v}-\||e|^{-1}\bm{H}_{e}\|_{p})_{2}^{2}$, where $\|\cdot\|_{p}$
computes the $\ell_{p}$-norm.

  * •

Total Variation (TV) [32, 37, 36]: $\max_{u,v\in
e}|\bm{h}_{v}-\bm{h}_{u}|^{p}$, $p\in\\{1,2\\}$.

  * •

Lovász Extension [57] for cardinality-based set functions (LEC) [54, 39, 33,
34]: $\langle\bm{y},\varsigma(\bm{H}_{e})\rangle^{p}$, $p\in\\{1,2\\}$ where
$\bm{y}=\begin{bmatrix}...,\bm{y}_{j},...\end{bmatrix}^{\top}\in\mathbb{R}^{|e|}$
is a constant vector and $\varsigma(\bm{H}_{e})$ sorts the values of
$\bm{H}_{e}$ in a decreasing order.

One may reproduce TV by using LEC and setting $\bm{y}_{1}=-\bm{y}_{|e|}=1$. To
reveal more properties of these hyperedge potentials $g(\cdot)$ later, we give
the following definitions.

######  Definition 4 (Permutation Invariance & Equivariance).

Function $\psi:\mathbb{R}^{K}\rightarrow\mathbb{R}$ is permutation-invariant
if for any $K$-dim permutation matrix $\bm{P}\in\Pi[K]$,
$\psi(\bm{P}\bm{Z})=\psi(\bm{Z})$ for all $\bm{Z}\in\mathbb{R}^{K}$. Function
$\psi:\mathbb{R}^{K}\rightarrow\mathbb{R}^{K}$ is permutation-equivariant if
for any $K$-dim permutation matrix $\bm{P}\in\Pi[K]$,
$\psi(\bm{P}\bm{Z})=\bm{P}\psi(\bm{Z})$ for all $\bm{Z}\in\mathbb{R}^{K}$.

We may easily verify permutation invariance of the hyperedge potentials in
Example 1. The physical meaning behind this observation is that the prediction
goal of an application is often independent of node identities in a hyperedge
so practical $g$’s often keep invariant w.r.t. the node ordering [58].

##  3 Neural Equivariant Diffusion Operators and ED-HNN

Previous design of hyperedge potential functions is tricky. Early works
adopted clique or star expansion by reducing hyperedges into (weighted) edges
[31, 47, 16] and further used traditional graph methods. Later, researchers
proved that those hyperedge reduction techniques cannot well represent higher-
order relations [16, 59]. Therefore, Lovász extensions of set-based cut-cost
functions on hyperedges have been proposed recently and used as the potential
functions [32, 39, 33, 36, 34, 35, 60]. However, designing those set-based cut
cost is practically hard and needs a lot of trials and errors. Other types of
handcrafted hyperedge potentials to model information propagation can also be
found in [61, 62].

Our idea uses data-driven approaches to model such potentials, which naturally
brings us to HNNs. On one hand, we expect to leverage the extreme expressive
power of NNs to learn the desired hypergraph diffusion automatically from the
data. On the other hand, we are interested in having novel hypergraph NN (HNN)
architectures inspired by traditional hypergraph diffusion solvers.

To achieve the goals, next, we show by the gradient descent algorithm (GD) or
alternating direction method of multipliers (ADMM) [63], solving objective Eq.
1 amounts to iteratively applying some hyperedge diffusion operators.
Parameterizing such operators using NNs by each step can unfold hypergraph
diffusion into an HNN. We make the key observation that the diffusion
operators are inherently permutation-equivariant. To universally represent
them, we propose a standard form to design equivariant NN-based operators and
use it to build ED-HNN via an efficient implementation.

###  3.1 Emerging Equivariance in Hypergraph Diffusion

We start with discussing the traditional solvers for Eq. 1. If $f$ and $g$ are
both differentiable, one straightforward optimization approach is to adopt
gradient descent. The node-wise update of each iteration can be formulated as
below:

$\displaystyle\bm{h}_{v}^{(t+1)}\leftarrow\bm{h}_{v}^{(t)}-\eta(\nabla
f(\bm{h}_{v}^{(t)};\bm{x}_{v})+\sum_{e:v\in e}[\nabla
g(\bm{H}_{e}^{(t)})]_{v}),\quad\text{for}\;v\in\mathcal{V},$ (2)

where $[\nabla g(\bm{H}_{e})]_{v}$ denotes the gradient w.r.t. $\bm{h}_{v}$
for $v\in e$. We use the supscript $t$ to denote the number of the current
iteration, $\bm{h}_{v}^{(0)}=\bm{x}_{v}$ is the initial features, and $\eta$
is known as the step size.

For general $f$ and $g$, we may adopt ADMM [63]: For each $e\in\mathcal{E}$,
we introduce an auxiliary variable
$\bm{Q}_{e}=\begin{bmatrix}\bm{q}_{e,1}&\cdots&\bm{q}_{e,\lvert
e\rvert}\end{bmatrix}^{\top}\in\mathbb{R}^{\lvert e\rvert}$. We initialize
$\bm{h}_{v}^{(0)}=\bm{x}_{v}$ and $\bm{Q}_{e}^{(0)}=\bm{H}_{e}^{(0)}$. And
then, iterate

$\displaystyle\bm{Q}_{e}^{(t+1)}$
$\displaystyle\leftarrow\operatorname{\textbf{prox}}_{\eta
g}(2\bm{H}_{e}^{(t)}-\bm{Q}_{e}^{(t)})-\bm{H}_{e}^{(t)}+\bm{Q}_{e}^{(t)},\quad\text{for}\;e\in\mathcal{E},$
(3) $\displaystyle\bm{h}_{v}^{(t+1)}$
$\displaystyle\leftarrow\operatorname{\textbf{prox}}_{\eta
f(\cdot;\bm{x}_{v})/d_{v}}(\sum_{e:v\in
e}\bm{q}_{e,v}^{(t+1)}/d_{v}),\quad\text{for}\;v\in\mathcal{V},$ (4)

where
$\operatorname{\textbf{prox}}_{\psi}(\bm{h})\triangleq\arg\min_{\bm{z}}\psi(\bm{z})+\frac{1}{2}\|\bm{z}-\bm{h}\|_{2}^{2}$
is the proximal operator. The detailed derivation can be found in Appendix A.
Note that previous hypergraph diffusion algorithms based on primal-dual hybrid
gradient descent [32] or alternative projection [55, 33] for some specific $f$
and $g$ are close related to the above general form given by ADMM. The
operator $\operatorname{\textbf{prox}}_{\psi}(\cdot)$ has nice properties as
reviewed in Proposition 1, which enables the possibility of NN-based
approximation even when $f$ and $g$ are not differentiable. One noted non-
differentiable example is the LEC case when
$g(\bm{H}_{e})=\langle\bm{y},\varsigma(\bm{H}_{e})\rangle^{p}$.

######  Proposition 1 ([64, 65]).

If $\psi(\cdot):\mathbb{R}^{K}\rightarrow\mathbb{R}$ is a lower semi-
continuous convex function, then $\operatorname{\textbf{prox}}_{\psi}(\cdot)$
is 1-Lipschitz continuous.

The node-side operations, gradient $\nabla f(\cdot;\bm{x}_{v})$ and proximal
gradient $\operatorname{\textbf{prox}}_{\eta f(\cdot;\bm{x}_{v})}(\cdot)$ are
relatively easy to model, while the operations on hyperedges are more
complicated. We name gradient $\nabla
g(\cdot):\mathbb{R}^{|e|}\rightarrow\mathbb{R}^{|e|}$ and proximal gradient
$\operatorname{\textbf{prox}}_{\eta
g}(\cdot):\mathbb{R}^{|e|}\rightarrow\mathbb{R}^{|e|}$ as hyperedge diffusion
operators, since they summarize the collection of node features inside a
hyperedge and dispatch the aggregated information to interior nodes
individually. Next, we reveal one crucial property of those hyperedge
diffusion operators by the following propostion (see the proof in Appendix
B.1):

######  Proposition 2.

Given any permutation-invariant hyperedge potential function $g(\cdot)$,
hyperedge diffusion operators $\operatorname{\textbf{prox}}_{\eta g}(\cdot)$
and $\nabla g(\cdot)$ are permutation equivariant.

It states that an permutation invariant hyperedge potential leads to an
operator that should process different nodes in a permutation-equivariant way.

###  3.2 Building Equivariant Hyperedge Diffusion Operators

To model permutation-equivariant diffusion operators, we may resort to
previous studies on NN architectures to represent equivariant multi-set
functions [50, 51, 52]. However, there are still some gaps when we apply their
approaches in the hypergraph setting. In practice, these works focus on the
case with separate multi-sets, such as prediction over separate sets of images
or 3-D points, while hypergraphs contain a large number of hyperedges that are
overlapped. Moveover, as hyperedges could be large and irregular in practice,
how to build scalable equivariant NN architectures on hyperedges is
challenging. Note that although graphs also have coupled edges, edges in
graphs are of regular size two. Building equivariant edge operation in graphs
is often much less crucial or challenging than that in hypergraphs.

Our design is built upon the following Theorem 1. We leave the proof in
Appendix B.2.

######  Theorem 1.

$\psi(\cdot):[0,1]^{K}\rightarrow\mathbb{R}^{K}$ is a continuous permutation-
equivariant function, if and only if it can be represented as
$[\psi(\bm{Z})]_{i}=\rho(\bm{z}_{i},\sum_{j=1}^{K}\phi(\bm{z}_{j}))$,
$i\in[K]$ for any
$\bm{Z}=\begin{bmatrix}...,\bm{z}_{i},...\end{bmatrix}^{\top}\in[0,1]^{K}$,
where $\rho:\mathbb{R}^{K^{\prime}}\rightarrow\mathbb{R}$,
$\phi:\mathbb{R}\rightarrow\mathbb{R}^{K^{\prime}-1}$ are two continuous
functions, and $K^{\prime}$ could be as small as $K$.

Theorem 1 indicates that any continuous permutation-equivariant function where
each entry of the input $\bm{Z}$ has 1-dim feature channel can be precisely
represented as a composition of a continuous function $\rho$ and the sum of
another continuous function $\phi$ on each input entry. This result
generalizes the representation of permutation-invariant functions by
$\rho(\sum_{i=1}^{K}\phi(\bm{z}_{i}))$ in [50] to the equivariant case. An
architecture with a similar spirit was proposed in [52] for equivariant set
function representations. However, in theory, their proof used the
universality of equivariant polynomials as a bridge so their construction
$\hat{\psi}(\bm{Z})$ only allows approximation, i.e., small
$\|\psi(\bm{Z})-\hat{\psi}(\bm{Z})\|$ instead of precise representation in
Theorem 1. Note that, in practice, $\rho$ and $\phi$ often need to be further
approximated by NNs. So, the precise representaiton in Theorem 1 has more
values of its theoretical implication than its practical usage given the
results in [52].

The above theoretical observation inspires the design of equivariant hyperedge
diffusion operators. Specifically, for an operator
$\hat{\psi}(\cdot):\mathbb{R}^{|e|}\rightarrow\mathbb{R}^{|e|}$ that may
denote either gradient $\nabla g(\cdot)$ or proximal gradient
$\operatorname{\textbf{prox}}_{\eta g}(\cdot)$ for each hyperedge $e$, we may
parameterize it as:

$\displaystyle[\hat{\psi}(\bm{H}_{e})]_{v}=\hat{\rho}(\bm{h}_{v},\sum_{u\in
e}\hat{\phi}(\bm{h}_{u})),\,\text{for $v\in e$,}\,\text{where
$\hat{\rho},\,\hat{\phi}$ are multi-layer perceptions (MLPs). }$ (5)

Intuitively, the inner sum collects the $\hat{\phi}$-encoding node features
within a hyperedge and then $\hat{\rho}$ combines the collection with the
features from each node further to perform separate operation.

The implementation of the above $\hat{\psi}$ is not trivial. A naive
implementation is to generate an auxillary node to representation each
$(v,e)$-pair for $v\in\mathcal{V}$ and $e\in\mathcal{E}$ and learn its
representation as adopted in [29, 28, 30]. However, this may substantially
increase the model complexity.

Our implementation is built upon the widely-used bipartite representations
[24, 25, 26] (or star expansion [48, 47], equivalently) of hypergraphs paired
with the standard message passing NN (MPNN) [49] that can be efficiently
implemented via commonly-used GNN platforms [45, 46]. Specifically, we build a
bipartite graph $\mathcal{\bar{G}}=(\mathcal{\bar{V}},\mathcal{\bar{E}})$. The
node set $\mathcal{\bar{V}}$ contains two parts
$\mathcal{V}\cup\mathcal{V_{\mathcal{E}}}$ where $\mathcal{V}$ is the original
node set while $\mathcal{V_{\mathcal{E}}}$ contains nodes that correspond to
original hyperedges $e\in\mathcal{E}$. Then, add an edge between
$v\in\mathcal{V}$ and $e\in\mathcal{E}$ if $v\in e$. The model ED-HNN is
implemented by following Algorithm 1.

Algorithm 1: ED-HNN
---
Initialization: $\bm{H}^{(0)}=\bm{X}$ and three MLPs
$\hat{\phi},\,\hat{\rho},\,\hat{\varphi}$ (shared across $L$ layers).
For $t=0,1,2,...,L-1$, do:
1\. Designing the messages from $\mathcal{V}$ to $\mathcal{E}$:
$\bm{m}_{u\rightarrow e}^{(t)}=\hat{\phi}(\bm{h}_{u}^{(t)})$, for all
$u\in\mathcal{V}$.
2\. Sum $\mathcal{V}\rightarrow\mathcal{E}$ messages over hyperedges
$\bm{m}_{e}^{(t)}=\sum_{u\in e}\bm{m}_{u\rightarrow e}^{(t)}$, for all
$e\in\mathcal{E}$.
3\. Broadcast $\bm{m}_{e}^{(t)}$ and design the messages from $\mathcal{E}$ to
$\mathcal{V}$: $\bm{m}_{e\rightarrow
v}^{(t)}=\hat{\rho}(\bm{h}_{v}^{(t)},\bm{m}_{e}^{(t)})$, for all $v\in e$.
4\. Update $\bm{h}_{v}^{(t+1)}=\hat{\varphi}(\bm{h}_{v}^{(t)},\sum_{e:u\in
e}\bm{m}_{e\rightarrow u}^{(t)},\bm{x}_{v},d_{v})$, for all $v\in\mathcal{V}$.

The equivariant diffusion operator $\hat{\psi}$ can be constructed via steps
1-3. The last step is to update the node features to accomplish the first two
terms in Eq. 2 or the ADMM update Eq. 4. The initial attributes $\bm{x}_{v}$
and node degrees are included to match the diffusion algorithm by design. As
the diffusion operators are shared across iterations, ED-HNN shares parameters
across layers. Now, we summarize the joint contribution of our theory and
efficient implementation as follows.

######  Proposition 3.

MPNNs on bipartite representations (or star expansions) of hypergraphs are
expressive enough to learn any continuous diffusion operators induced by
invariant hyperedge potentials.

ED-HNN can well approximate the GD algorithm (Eq. 2), while ED-HNN does not
perfectly match ADMM unless we assume $\bm{Q}_{e}^{(t)}=\bm{H}_{e}^{(t)}$ in
Eq. 3 for a practical consideration. This assumption may reduce the
performance of the hypergraph diffusion algorithm while our model ED-HNN has
already achieved good enough empirical performance. Tracking
$\bm{Q}_{e}^{(t)}$ means recording the messages from $\mathcal{E}$ to
$\mathcal{V}$ for every iteration, which is not supported by the current GNN
platforms [45, 46] and may consume more memory. We leave the study of the
algorithms that can track the update of $\bm{Q}_{e}^{(t)}$ as a future study.
A co-design of the algorithm and the system may be needed to guarantee the
scalability of the algorithm.

Extension. Our ED-HNN can be naturally extended to build equivariant node
operators because of the duality between hyperedges and nodes, though this is
not necessary in the traditional hypergraph diffusion problem Eq. 1.
Specifically, we call this model ED-HNNII, which simple revises the step 1 in
ED-HNN as $\bm{m}_{u\rightarrow
e}^{(t)}=\hat{\phi}(\bm{h}_{u}^{(t)},\bm{m}_{e}^{(t-1)})$ for any $e$ such
that $u\in e$. Due to the page limit, we put some experiment results on ED-
HNNII in Appendix C.1.

###  3.3 Discussion: Advantages of ED-HNN in Heterophilic Settings and Deep
Models

Here, we discuss several advantages of ED-HNN due to the design of equivariant
diffusion operators.

Heterophily describes the network phenomenon where nodes with the same labels
and attributes are less likely to connect to each other directly [67].
Predicting node labels in heterophilic networks is known to be more
challenging than that in homophilic networks, and thus has recently become an
important research direction in the community [68, 69, 70, 71]. Heterophily
has been proved as a more common phenomenon in hypergraphs than in graphs
since it is hard to expect all nodes in a giant hyperedge to share a common
label [16, 72]. Moreover, predicting node labels in heterophilic hypergraphs
is more challenging than that in graphs as a hyperedge may consist of the
nodes from multiple categories.

Learnable equivariant diffusion operators are expected to be superior in
predicting heterophilic node labels. For example, if a hyperedge $e$ of size 3
is known to cover two nodes $v,u$ from class $\mathcal{C}_{1}$ while one node
$w$ from class $\mathcal{C}_{0}$. We may use the LEC potential
$\langle\bm{y},\varsigma(\bm{H}_{e})\rangle^{2}$ by setting
$\bm{y}=[1,-1,0]^{\top}$. Suppose the three nodes’ attibutes are
$\bm{H}_{e}=[\bm{h}_{v},\bm{h}_{u},\bm{h}_{w}]^{\top}=[0.7,0.5,0.3]^{\top}$,
where $\bm{h}_{v}$, $\bm{h}_{u}$ are close as they are from the same class.
One may check the hyperedge diffusion operator gives $\nabla
g(\bm{H}_{e})=[0.4,-0.4,0]^{\top}$. One-step gradient descent
$\bm{h}_{v}-\eta[\nabla g(\bm{H}_{e})]_{v}$ with a proper step size
($\eta<0.5$) drags $\bm{h}_{v}$ and $\bm{h}_{u}$ closer while keeping
$\bm{h}_{w}$ unchanged. However, invariant diffusion that allocates every node
with the same change does not help with the classification. Moreover, consider
another case when the ground truth is different, say $w,u$ from the same class
and $v$ from the other, using the above $\bm{y}$ may increase the error.
However, learnable operators can solve the problem, e.g., by learning another
$\bm{y}=[0,1,-1]^{\top}$.

Moreover, equivariant hyperedge diffusion operators are also good at building
deep models. GNNs tend to degenerate the performance when going deep, which is
often attributed to their oversmoothing [73, 74] and overfitting problems
[75]. Equivariant operators allocating different messages to different nodes
helps with overcoming the oversmoothing issue. On the other hand, diffusion
operators by sharing parameters across different layers may reduce the risk of
overfitting.

###  3.4 Previous Hypergraph Neural Networks for Representing Hypergraph
Diffusion

ED-HNN is the first HNN inspired by the general class of hypergraph diffusion
and can provably achieve universal approximation of hyperedge diffusion
operators. So, how about previous HNNs representing hypergraph diffusion? We
temporarily ignore the big difference in whether parameters are shared across
layers and give some analysis as follows. HGNN [22] runs graph convolution
[76] on clique expansions of hypergraphs, which directly assumes the hyperedge
potentials follow CE plus pairwise potentials and cannot learn other
operations on hyperedges. HyperGCN [23] essentially leverages the total
variation potential by adding mediator nodes to each hyperedge, which adopts a
transformation technique of the total variantion potential in [77, 78]. A few
works view each hyperedge as a multi-set of nodes and each node as a multi-set
of hyperedges [24, 25, 26, 27, 28, 29, 30, 79]. Among them, HNHN [24], UniGNNs
[25], EHGNN [79] and AllSet [26] build two invariant set-pooling functions on
both the hyperedge and node sides, which cannot represent equivariant
functions. One subtle point is that Uni{GIN,GraphSAGE,GCNII} [25] adopt a jump
link in each layer to pass node features from the former layer directly to the
next and may expect to use a complicated node-side invariant function
$\hat{\psi}^{\prime}(\bm{h}_{v}^{(t)},\sum_{e:v\in e}\bm{m}_{e}^{(t)})$ to
approximate our $\hat{\psi}(\bm{h}_{v}^{(t)},\sum_{e:v\in
e}\bm{m}_{e\rightarrow v}^{(t)})=\hat{\psi}(\bm{h}_{v}^{(t)},\sum_{e:v\in
e}\hat{\rho}(\bm{h}_{v}^{(t)},\bm{m}_{e}^{(t)}))$ in Step 4 of Algorithm 1
directly. This may be doable in theory. However, the jump-link solution may
expect to have a higher dimension of $\bm{m}_{e}^{(t)}$ (dim=$|\cup_{e:v\in
e}e|$) so that the sum pooling $\sum_{e:v\in e}\bm{m}_{e}^{(t)}$ does not lose
anything from the neighboring nodes of $v$ before interacting with
$\bm{h}_{v}^{(t)}$, while in ED-HNN, $\bm{m}_{e}^{(t)}$ needs dim=$|e|$
according to Theorem 1. Our empirical experiments also verify more
expressiveness of ED-HNN. HCHA [27] and UniGAT [30] compute attention weights
as a result of combining node features and hyperedge messages, which looks
like our $\rho$ operation. However, a scalar attention weight is too limited
to represent the potentially complicated $\rho$. HyperSAGE [28], MPNN-R [29]
and LEGCN [30] may learn hyperedge equivariant operators, while they need many
auxillary nodes that represent node-hyperedge pairs. Of course, none of these
works have mentioned any theoretical arguments on universal representation of
hypergraph diffusion as ours.

##  4 Further Discussion on Related Works

We have discussed the most relevant works on HNN architectures and hypergraph
diffusion algorithms in Sec. 3. Here, we review some relevant works on
optimization-inspired NNs. Optimization-inspired NNs often get praised for
their interpretability [80, 81] and certified convergence under certain
conditions of the learned operators [82, 83, 84]. Early works that adopt
optimization-inspired NNs focus on the applications such as compressive
sensing [85, 86, 87, 88], speech processing [89, 90], image denoising [91, 92,
93, 94], partitioning [95, 96], deblurring [97, 98] and so on. Recently,
optimization-inspired NNs have also been applied to graphs [99, 100, 101,
102]. They unroll either learnable [99, 101], fixed $\ell_{1}$-norm [100] or
fixed $\ell_{2}$-norm pair-wise potentials [102] to formulate GNNs. Some works
also view the graph diffusion as an ODE [103, 104], which essentially
corresponds to our gradient descent Eq. 2. As there are only pairwise
relations in graphs, the GNN design is often much simpler than the HNN
counterpart. To construct an expressive operator with for pairwise relations
is trivial. Implicit NNs [105, 106] are to directly parameterize the optimum
of an optimization problem. Implicit NNs have also been applied to graphs
[107, 108, 109, 110] but not to hypergraphs. So, a future study of this work
is to develop implicit models that can directly parameterize the solution to
Eq. 1.

##  5 Experiments

Table 1: Dataset statistics. CE homophily is the homophily score [68] based on
CE of hypergraphs.

| Cora | Citeseer | Pubmed | Cora-CA | DBLP-CA | Congress | Senate | Walmart | House
---|---|---|---|---|---|---|---|---|---
# nodes | 2708 | 3312 | 19717 | 2708 | 41302 | 1718 | 282 | 88860 | 1290
# edges | 1579 | 1079 | 7963 | 1072 | 22363 | 83105 | 315 | 69906 | 340
# classes | 7 | 6 | 3 | 7 | 6 | 2 | 2 | 11 | 2
avg. $|e|$ | 1.748 | 1.541 | 1.963 | 1.929 | 2.213 | 8.654 | 9.531 | 3.461 | 8.056
CE Homophily | 0.897 | 0.893 | 0.952 | 0.803 | 0.869 | 0.555 | 0.498 | 0.530 | 0.509

###  5.1 Results on Benchmarking Datasets

Experiment Setting. In this subsection, we evaluate ED-HNN on nine real-world
benchmarking hypergraphs. We focus on the semi-supervised node classification
task. The nine datasets include co-citation networks (Cora, Citeseer, Pubmed),
co-authorship networks (Cora-CA, DBLP-CA) [23], Walmart [111], House [112],
Congress and Senate [113, 114]. More details of these datasets can be found in
Appendix C.2. Since the last four hypergraphs do not contain node attributes,
we follow the method of [26] to generate node features from label-dependent
Gaussian distribution [115].

As we show in Table 1, these datasets already cover sufficiently diverse
hypergraphs in terms of scale, structure, and homo-/heterophily. We compare
our method with top-performing models on these benchmarks, including HGNN
[22], HCHA [27], HNHN [24], HyperGCN [23], UniGCNII [25], AllDeepSets [26],
and AllSetTransformer [26]. All the hyperparameters for baselines follow from
[26] and we fix the learning rate, weight decay and other training recipes
same with the baselines. Other model specific hyperparameters are obtained via
grid search (see Appendix C.3). We randomly split the data into
training/validation/test samples using 50%/25%/25% splitting percentage by
following [26]. We choose prediction accuracy as the evaluation metric. We run
each model for ten times with different training/validation splits to obtain
the standard deviation.

Performance Analysis. Table 2 shows the results. Our ED-HNN uniformly
outperforms all the compared models on all the datasets. We observe that the
top-performing baseline models are AllSetTransformer, AllDeepSets and
UniGCNII. As having been analyzed in Sec. 3.4, they model invariant set
functions on both node and hyperedge sides. UniGCNII also adds initial and
jump links, which accidentally resonates with our design principle (the step 4
in Algorithm 1). However, their performance has large variation across
different datasets. For example, UniGCNII attains promising performance on
citation networks, however, has subpar results on Walmart dataset. In
contrast, our model achieves stably superior results, surpassing AllSet models
by 12.9% on Senate and UniGCNII by 12.5% on Walmart. We owe our empirical
significance to the theoretical design of exact equivariant function
representation. Regarding the computational efficiency, we report the wall-
clock times for training (for 100 epochs) and testing of all models over the
largest hypergraph Walmart, where all models use the same hyperparameters that
achieve the reported performance on the left. Our model achieves efficiency
comparable to AllSet models [26] and is much faster than UniGCNII. So, the
implementation of equivariant computation in ED-HNN is still efficient.

Table 2: Prediction Accuracy (%). Bold font† highlights when ED-HNN
significantly (difference in means > 0.5 $\times$ std) outperforms all
baselines. The best baselines are underlined. The training and testing times
are test on Walmart by using the same server with one GPU NVIDIA RTX A6000.

| Cora | Citeseer | Pubmed | Cora-CA | DBLP-CA |  Training Time ($10^{-1}$ s)
---|---|---|---|---|---|---
HGNN [25] | 79.39 $\pm$ 1.36 |  72.45 $\pm$ 1.16 |  86.44 $\pm$ 0.44 |  82.64 $\pm$ 1.65 |  91.03 $\pm$ 0.20 |  0.24 $\pm$ 0.51
HCHA [27] |  79.14 $\pm$ 1.02 |  72.42 $\pm$ 1.42 |  86.41 $\pm$ 0.36 |  82.55 $\pm$ 0.97 |  90.92 $\pm$ 0.22 | 0.24 ± 0.01
HNHN [24] |  76.36 $\pm$ 1.92 |  72.64 $\pm$ 1.57 |  86.90 $\pm$ 0.30 |  77.19 $\pm$ 1.49 |  86.78 $\pm$ 0.29 |  0.30 $\pm$ 0.56
HyperGCN [23] |  78.45 $\pm$ 1.26 |  71.28 $\pm$ 0.82 |  82.84 $\pm$ 8.67 |  79.48 $\pm$ 2.08 |  89.38 $\pm$ 0.25 |  0.42 $\pm$ 1.51
UniGCNII [25] |  78.81 $\pm$ 1.05 |  73.05 $\pm$ 2.21 |  88.25 $\pm$ 0.40 | 83.60 $\pm$ 1.14 | 91.69 $\pm$ 0.19 |  4.36 $\pm$ 1.18
AllDeepSets [26] |  76.88 $\pm$ 1.80 |  70.83 $\pm$ 1.63 | 88.75 $\pm$ 0.33 |  81.97 $\pm$ 1.50 |  91.27 $\pm$ 0.27 |  1.23 $\pm$ 1.09
AllSetTransformer [26] |  78.58 $\pm$ 1.47 | 73.08 $\pm$ 1.20 |  88.72 $\pm$ 0.37 |  83.63 $\pm$ 1.47 |  91.53 $\pm$ 0.23 |  1.64 $\pm$ 1.63
ED-HNN (ours) | 80.31 $\pm$ 1.35† | 73.70 $\pm$ 1.38† | 89.03 $\pm$ 0.53† |  83.97 $\pm$ 1.55 | 91.90 $\pm$ 0.19† |  1.71 $\pm$ 1.13
| Congress | Senate | Walmart | House | Avg. Rank |  Inference Time ($10^{-2}$ s)
HGNN [25] |  91.26 $\pm$ 1.15 |  48.59 $\pm$ 4.52 |  62.00 $\pm$ 0.24 |  61.39 $\pm$ 2.96 | 4.89 |  1.01 $\pm$ 0.04
HCHA [27] |  90.43 $\pm$ 1.20 |  48.62 $\pm$ 4.41 |  62.35 $\pm$ 0.26 |  61.36 $\pm$ 2.53 | 5.44 |  1.54 $\pm$ 0.18
HNHN [24] |  53.35 $\pm$ 1.45 |  50.93 $\pm$ 6.33 |  47.18 $\pm$ 0.35 |  67.80 $\pm$ 2.59 | 7.11 |  6.11 $\pm$ 0.05
HyperGCN [23] |  55.12 $\pm$ 1.96 |  42.45 $\pm$ 3.67 |  44.74 $\pm$ 2.81 |  48.32 $\pm$ 2.93 | 7.33 |  0.87 $\pm$ 0.06
UniGCNII [25] | 94.81 $\pm$ 0.81 |  49.30 $\pm$ 4.25 |  54.45 $\pm$ 0.37 |  67.25 $\pm$ 2.57 | 3.67 |  21.22 $\pm$ 0.13
AllDeepSets [26] |  91.80 $\pm$ 1.53 |  48.17 $\pm$ 5.67 |  64.55 $\pm$ 0.33 |  67.82 $\pm$ 2.40 | 4.89 |  5.35 $\pm$ 0.33
AllSetTransformer [26] |  92.16 $\pm$ 1.05 | 51.83 $\pm$ 5.22 | 65.46 $\pm$ 0.25 | 69.33 $\pm$ 2.20 | 2.67 |  6.06 $\pm$ 0.67
ED-HNN (ours) |  95.00 $\pm$ 0.99 | 64.79 $\pm$ 5.14† | 66.91 $\pm$ 0.41† | 72.45 $\pm$ 2.28† | 1.00 |  5.87 $\pm$ 0.36

###  5.2 Results on Synthetic Heterophilic Hypergraph Dataset

Table 3: Prediction Accuracy (%) over Synthetic Hypergraphs with Controlled
Heterophily $\alpha$.

| Homophily | Heterophily
---|---|---
$\alpha=1$ | $\alpha=2$ | $\alpha=3$ | $\alpha=4$ | $\alpha=6$ | $\alpha=7$
HGNN [25] |  99.86 $\pm$ 0.05 |  90.42 $\pm$ 1.14 |  66.60 $\pm$ 1.60 |  57.90 $\pm$ 1.43 |  50.77 $\pm$ 1.94 |  48.68 $\pm$ 1.21
HGNN + JumpLink [25] |  99.29 $\pm$ 0.19 |  92.61 $\pm$ 0.97 |  79.59 $\pm$ 2.32 |  64.39 $\pm$ 1.55 |  55.39 $\pm$ 2.01 |  51.14 $\pm$ 2.05
AllDeepSets [26] |  99.86 $\pm$ 0.07 |  99.13 $\pm$ 0.31 |  92.39 $\pm$ 1.68 |  80.12 $\pm$ 2.31 |  56.42 $\pm$ 0.85 |  52.14 $\pm$ 2.88
AllDeepSets + JumpLink [25] |  99.95 $\pm$ 0.07 |  99.08 $\pm$ 0.24 |  93.52 $\pm$ 1.31 |  83.45 $\pm$ 2.64 |  59.68 $\pm$ 1.27 |  52.21 $\pm$ 2.36
AllSetTransformer [26] |  98.60 $\pm$ 0.48 |  95.89 $\pm$ 0.90 |  88.38 $\pm$ 1.81 |  83.34 $\pm$ 1.33 |  54.13 $\pm$ 4.83 |  50.15 $\pm$ 1.73
ED-HNN (ours) |  99.95 $\pm$ 0.08 | 99.68 $\pm$ 0.13† | 95.37 $\pm$ 1.23† | 90.35 $\pm$ 0.83† | 65.58 $\pm$ 0.76† | 54.31 $\pm$ 2.40†

#### Experiment Setting.

As discussed, ED-HNN is expected to perform well on heterophilic hypergraphs.
We evaluate this point by using synthetic datasets with controlled
heterophily. We generate data by using contextual hypergraph stochastic block
model [115, 116, 117]. Specifically, we draw two classes of $2,500$ nodes each
and then randomly sample 1,000 hyperedges. Each hyperedge consists of 15
nodes, among which $\alpha_{i}$ many are sampled from class $i$. We use
$\alpha=\min\\{\alpha_{1},\alpha_{2}\\}$ to denote the heterophily level.
Afterwards, we generate label-dependent Gaussian node features. We test both
homophilic ($\alpha=1,2$ or CE homophily$\geq 0.7$) and heterophilic
($\alpha=4\sim 7$ or CE homophily $\leq 0.7$) cases. We compare ED-HNN with
HGNN [22], AllSet models [26] and their variants with jump links. We follow
previous 50%/25%/25% data splitting methods and repeated 10 times the
experiment.

Results. Table 3 shows the results. On homophilic datasets, all the models can
achieve good results, while ED-HNN keeps slightly better than others. Once
$\alpha$ surpasses 3, i.e., entering the heterophilic regime, the superiority
of ED-HNN is more obvious. The jump-link trick indeed also helps, while
building equivariance as ED-HNN does directly provides more significant
improvement.

###  5.3 Benefits in Deepening Hypergraph Neural Networks

We also demonstrate that by using diffusion models and parameter tying, ED-HNN
can benefit from deeper architectures, while other HNNs cannot. Fig. 2
illustrates the performance of different models versus the number of network
layers. We compare with HGNN [22], AllSet models [26], and UniGCNII [25].
UniGCNII inherits from [118] which is known to be effective to counteract
oversmoothness. The results reveal that AllSet models suffer from going deep.
HGNN working through more lightweight mechanism has better tolerance to depth.
However, none of them can benefit from deepening. On the contrary, ED-HNN
successfully leverages deeper architecture to achieve higher accuracy. For
example, adding more layers boosts ED-HNN by ~1% in accuracy on Pubmed and
House, while elevating ED-HNN from 58.01% to 64.79% on Senate dataset.

Figure 2: Comparing Deep Achitectures across Different Models: Test Acc. v.s.
# Layers Figure 3: Comparing the Powers to Represent known Diffusion: MAE v.s.
Latent Dimensions

###  5.4 Expressiveness Justification on the Synthetic Diffusion Dataset

We are to evaluate the ability of ED-HNN to express given hypergraph
diffusion. We generate semi-synthetic diffusion data using the Senate
hypergraph [112] and synthetic node features. The data consists of 1,000 pairs
$(\bm{H}^{(0)},\bm{H}^{(1)})$. The initial node features $\bm{H}^{(0)}$ are
sampled from 1-dim Gaussian distributions. To obtain $\bm{H}^{(1)}$ we apply
the gradient step in Eq. 2. For non-differential cases, we adopt subgradients
for convenient computation. We fix node potentials as
$f(\bm{h}_{v};\bm{x}_{v})=(\bm{h}_{v}-\bm{x}_{v})^{2}$ and consider 3
different edge potentials in Example 1 with varying complexities: a) CE, b) TV
($p$=2) and c) LEC. The goal is to let one-layer models
$\mathcal{V}\rightarrow\mathcal{E}\rightarrow\mathcal{V}$ to recover
$\bm{H}^{(1)}$. We compare ED-HNN with our implemented baseline (Invariant)
with parameterized invariant set functions on both the node and hyperedge
sides, and AllSet models [26] that also adopt invariant set functions. We keep
the scale of all models almost the same to give fair comparison, of which more
details are given in Appendix C.5. The results are reported in Fig. 3. The CE
case gives invariant diffusion so all models can learn it well. The TV case
mostly shows the benefit of the equivariant architecture, where the error
almost does not increase even when the dimension decreases to 32. The LEC case
is challenging for all models, though ED-HNN is still the best. The reason, we
think, is that learning the sorting operation in LEC via the sum pooling in
Eq. 5 is empirically challenging albeit theoretically doable. A similar
phenomenon has been observed in previous literatures [119, 66].

##  6 Conclusion

This work introduces a new hypergraph neural network ED-HNN that can model
hypergraph diffusion process. We show that any hypergraph diffusion with
permutation-invariant potential functions can be represented by iterating
equivariant diffusion operators. ED-HNN provides an efficient way to model
such operators based on the commonly-used GNN platforms. ED-HNN shows
superiority in processing heterophilic hypergraphs and constructing deep
models. For future works, ED-HNN can be applied to other tasks such as
hypergraph regression tasks. Or implicit HNNs can be developed to directly
optimize the optima of hypergraph diffusion problems.

## References

  * [1] X. Zhu, _Semi-supervised learning with graphs_. Carnegie Mellon University, 2005.
  * [2] X. Zhu, Z. Ghahramani, and J. D. Lafferty, “Semi-supervised learning using gaussian fields and harmonic functions,” in _International Conference on Machine Learning (ICML)_ , 2003.
  * [3] D. Zhou, O. Bousquet, T. N. Lal, J. Weston, and B. Schölkopf, “Learning with local and global consistency,” in _Advances in Neural Information Precessing Systems (NeurIPS)_ , 2004.
  * [4] W. L. Hamilton, “Graph representation learning,” 2020.
  * [5] M. Nickel, K. Murphy, V. Tresp, and E. Gabrilovich, “A review of relational machine learning for knowledge graphs,” _Proceedings of the IEEE_ , 2015\.
  * [6] I. Chami, S. Abu-El-Haija, B. Perozzi, C. Ré, and K. Murphy, “Machine learning on graphs: A model and comprehensive taxonomy,” _arXiv preprint arXiv:2005.03675_ , 2020.
  * [7] A. R. Benson, D. F. Gleich, and D. J. Higham, “Higher-order network analysis takes off, fueled by classical ideas and new data,” _arXiv preprint arXiv:2103.05031_ , 2021.
  * [8] M. T. Schaub, Y. Zhu, J.-B. Seby, T. M. Roddenberry, and S. Segarra, “Signal processing on higher-order networks: Livin’on the edge… and beyond,” _Signal Processing_ , 2021.
  * [9] F. Battiston, G. Cencetti, I. Iacopini, V. Latora, M. Lucas, A. Patania, J.-G. Young, and G. Petri, “Networks beyond pairwise interactions: structure and dynamics,” _Physics Reports_ , 2020.
  * [10] R. Lambiotte, M. Rosvall, and I. Scholtes, “From networks to optimal higher-order models of complex systems,” _Nature physics_ , 2019.
  * [11] Y. Liu, J. Ma, and P. Li, “Neural predicting higher-order patterns in temporal networks,” in _Proceedings of the ACM Web Conference 2022_ , 2022, pp. 1340–1351.
  * [12] S. Mangan and U. Alon, “Structure and function of the feed-forward loop network motif,” _National Academy of Sciences_ , 2003.
  * [13] A. R. Benson, D. F. Gleich, and J. Leskovec, “Higher-order organization of complex networks,” _Science_ , 2016.
  * [14] C. E. Tsourakakis, J. Pachocki, and M. Mitzenmacher, “Scalable motif-aware graph clustering,” in _Proceedings of the Web Conference_ , 2017, pp. 1451–1460.
  * [15] P. Li, H. Dau, G. Puleo, and O. Milenkovic, “Motif clustering and overlapping clustering for social network analysis,” in _INFOCOM 2017-IEEE Conference on Computer Communications, IEEE_. IEEE, 2017, pp. 1–9.
  * [16] P. Li and O. Milenkovic, “Inhomogeneous hypergraph clustering with applications,” in _Advances in Neural Information Processing Systems_ , 2017, pp. 2305–2315.
  * [17] X. Xia, H. Yin, J. Yu, Q. Wang, L. Cui, and X. Zhang, “Self-supervised hypergraph convolutional networks for session-based recommendation,” in _Proceedings of the AAAI Conference on Artificial Intelligence_ , vol. 35, no. 5, 2021, pp. 4503–4511.
  * [18] J. Wang, K. Ding, Z. Zhu, and J. Caverlee, “Session-based recommendation with hypergraph attention networks,” in _Proceedings of the 2021 SIAM International Conference on Data Mining (SDM)_. SIAM, 2021, pp. 82–90.
  * [19] J. Wang, K. Ding, L. Hong, H. Liu, and J. Caverlee, “Next-item recommendation with sequential hypergraphs,” in _Proceedings of the 43rd international ACM SIGIR conference on research and development in information retrieval_ , 2020, pp. 1101–1110.
  * [20] N. Wang, S. Wang, Y. Wang, Q. Z. Sheng, and M. A. Orgun, “Exploiting intra-and inter-session dependencies for session-based recommendations,” _World Wide Web_ , vol. 25, no. 1, pp. 425–443, 2022.
  * [21] C. Berge, _Hypergraphs: combinatorics of finite sets_. Elsevier, 1984, vol. 45.
  * [22] Y. Feng, H. You, Z. Zhang, R. Ji, and Y. Gao, “Hypergraph neural networks,” in _Proceedings of the AAAI Conference on Artificial Intelligence_ , vol. 33, no. 01, 2019, pp. 3558–3565.
  * [23] N. Yadati, M. Nimishakavi, P. Yadav, V. Nitin, A. Louis, and P. Talukdar, “Hypergcn: A new method for training graph convolutional networks on hypergraphs,” _Advances in neural information processing systems_ , vol. 32, 2019.
  * [24] Y. Dong, W. Sawin, and Y. Bengio, “Hnhn: Hypergraph networks with hyperedge neurons,” _arXiv preprint arXiv:2006.12278_ , 2020.
  * [25] J. Huang and J. Yang, “Unignn: a unified framework for graph and hypergraph neural networks,” _arXiv preprint arXiv:2105.00956_ , 2021.
  * [26] E. Chien, C. Pan, J. Peng, and O. Milenkovic, “You are allset: A multiset function framework for hypergraph neural networks,” in _International Conference on Learning Representations_ , 2022.
  * [27] S. Bai, F. Zhang, and P. H. Torr, “Hypergraph convolution and hypergraph attention,” _Pattern Recognition_ , vol. 110, p. 107637, 2021.
  * [28] D. Arya, D. K. Gupta, S. Rudinac, and M. Worring, “Hypersage: Generalizing inductive representation learning on hypergraphs,” _arXiv preprint arXiv:2010.04558_ , 2020.
  * [29] N. Yadati, “Neural message passing for multi-relational ordered and recursive hypergraphs,” _Advances in Neural Information Processing Systems_ , vol. 33, pp. 3275–3289, 2020.
  * [30] C. Yang, R. Wang, S. Yao, and T. Abdelzaher, “Hypergraph learning with line expansion,” _arXiv preprint arXiv:2005.04843_ , 2020.
  * [31] D. Zhou, J. Huang, and B. Schölkopf, “Learning with hypergraphs: Clustering, classification, and embedding,” in _Advances in Neural Information Processing Systems_ , 2007, pp. 1601–1608.
  * [32] M. Hein, S. Setzer, L. Jost, and S. S. Rangapuram, “The total variation on hypergraphs-learning on hypergraphs revisited,” in _Advances in Neural Information Processing Systems_ , 2013, pp. 2427–2435.
  * [33] P. Li, N. He, and O. Milenkovic, “Quadratic decomposable submodular function minimization: Theory and practice.” _Journal of Machine Learning Research (JMLR)_ , 2020.
  * [34] M. Liu, N. Veldt, H. Song, P. Li, and D. F. Gleich, “Strongly local hypergraph diffusions for clustering and semi-supervised learning,” in _Proceedings of the Web Conference 2021_ , 2021, pp. 2092–2103.
  * [35] K. Fountoulakis, P. Li, and S. Yang, “Local hyper-flow diffusion,” _Advances in Neural Information Processing Systems_ , vol. 34, 2021.
  * [36] Y. Takai, A. Miyauchi, M. Ikeda, and Y. Yoshida, “Hypergraph clustering based on pagerank,” in _Proceedings of the 26th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining_, 2020, pp. 1970–1978.
  * [37] C. Zhang, S. Hu, Z. G. Tang, and T. H. Chan, “Re-revisiting learning on hypergraphs: confidence interval and subgradient method,” in _International Conference on Machine Learning_ , 2017, pp. 4026–4034.
  * [38] F. Tudisco, A. R. Benson, and K. Prokopchik, “Nonlinear higher-order label spreading,” in _Proceedings of the Web Conference_ , 2021.
  * [39] P. Li and O. Milenkovic, “Submodular hypergraphs: p-laplacians, cheeger inequalities and spectral clustering,” in _International Conference on Machine Learning (ICML)_ , 2018.
  * [40] P. Li, N. He, and O. Milenkovic, “Quadratic decomposable submodular function minimization,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2018.
  * [41] U. Chitra and B. Raphael, “Random walks on hypergraphs with edge-dependent vertex weights,” in _International Conference on Machine Learning (ICML)_ , 2019.
  * [42] S. Zhang, Z. Ding, and S. Cui, “Introducing hypergraph signal processing: Theoretical foundation and practical applications,” _IEEE Internet of Things Journal_ , 2019.
  * [43] N. Komodakis and N. Paragios, “Beyond pairwise energies: Efficient optimization for higher-order mrfs,” in _IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , 2009.
  * [44] T. Werner, “High-arity interactions, polyhedral relaxations, and cutting plane algorithm for soft constraint optimisation (map-mrf),” in _IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , 2008.
  * [45] M. Fey and J. E. Lenssen, “Fast graph representation learning with PyTorch Geometric,” in _International Conference on Learning Representations (ICLR) Workshop on Representation Learning on Graphs and Manifolds_ , 2019.
  * [46] M. Wang, D. Zheng, Z. Ye, Q. Gan, M. Li, X. Song, J. Zhou, C. Ma, L. Yu, Y. Gai, T. Xiao, T. He, G. Karypis, J. Li, and Z. Zhang, “Deep graph library: A graph-centric, highly-performant package for graph neural networks,” _arXiv preprint arXiv:1909.01315_ , 2019.
  * [47] S. Agarwal, K. Branson, and S. Belongie, “Higher order learning with graphs,” in _International Conference on Machine Learning (ICML)_ , 2006.
  * [48] J. Y. Zien, M. D. Schlag, and P. K. Chan, “Multilevel spectral hypergraph partitioning with arbitrary vertex sizes,” _IEEE Transactions on computer-aided design of integrated circuits and systems_ , 1999.
  * [49] J. Gilmer, S. S. Schoenholz, P. F. Riley, O. Vinyals, and G. E. Dahl, “Neural message passing for quantum chemistry,” in _International Conference on Machine Learning (ICML)_ , 2017.
  * [50] M. Zaheer, S. Kottur, S. Ravanbakhsh, B. Poczos, R. R. Salakhutdinov, and A. J. Smola, “Deep sets,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2017.
  * [51] H. Maron, H. Ben-Hamu, N. Shamir, and Y. Lipman, “Invariant and equivariant graph networks,” in _International Conference on Learning Representations (ICLR)_ , 2018.
  * [52] N. Segol and Y. Lipman, “On universal equivariant set networks,” in _International Conference on Learning Representations (ICLR)_ , 2020.
  * [53] P. Stobbe and A. Krause, “Efficient minimization of decomposable submodular functions,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2010.
  * [54] S. Jegelka, F. Bach, and S. Sra, “Reflection methods for user-friendly submodular optimization,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2013.
  * [55] P. Li and O. Milenkovic, “Revisiting decomposable submodular function minimization with incidence relations,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2018.
  * [56] F. Tudisco, K. Prokopchik, and A. R. Benson, “A nonlinear diffusion method for semi-supervised learning on hypergraphs,” _arXiv preprint arXiv:2103.14867_ , 2021.
  * [57] L. Lovász, “Submodular functions and convexity,” in _Mathematical Programming The State of the Art_ , 1983.
  * [58] N. Veldt, A. R. Benson, and J. Kleinberg, “Hypergraph cuts with general splitting functions,” _SIAM Review_ , 2021.
  * [59] I. E. Chien, H. Zhou, and P. Li, “$hs^{2}$: Active learning over hypergraphs with pointwise and pairwise queries,” in _International Conference on Artificial Intelligence and Statistics (AISTATS)_ , 2019.
  * [60] Y. Yoshida, “Cheeger inequalities for submodular transformations,” in _Proceedings of the Thirtieth Annual ACM-SIAM Symposium on Discrete Algorithms_. SIAM, 2019, pp. 2582–2601.
  * [61] L. Neuhäuser, R. Lambiotte, and M. T. Schaub, “Consensus dynamics and opinion formation on hypergraphs,” in _Higher-Order Systems_ , 2022.
  * [62] L. Neuhäuser, M. T. Schaub, A. Mellor, and R. Lambiotte, “Opinion dynamics with multi-body interactions,” in _International Conference on Network Games, Control and Optimization_ , 2021.
  * [63] S. Boyd, N. Parikh, E. Chu, B. Peleato, J. Eckstein _et al._ , “Distributed optimization and statistical learning via the alternating direction method of multipliers,” _Foundations and Trends® in Machine learning_ , 2011.
  * [64] N. Parikh and S. Boyd, “Proximal algorithms,” _Foundations and Trends in optimization_ , 2014.
  * [65] N. G. Polson, J. G. Scott, and B. T. Willard, “Proximal algorithms in statistics and machine learning,” _Statistical Science_ , 2015.
  * [66] E. Wagstaff, F. Fuchs, M. Engelcke, I. Posner, and M. A. Osborne, “On the limitations of representing functions on sets,” in _International Conference on Machine Learning (ICML)_ , 2019.
  * [67] E. M. Rogers, _Diffusion of innovations_ , 2010.
  * [68] H. Pei, B. Wei, K. C.-C. Chang, Y. Lei, and B. Yang, “Geom-gcn: Geometric graph convolutional networks,” in _International Conference on Learning Representations (ICLR)_ , 2020.
  * [69] J. Zhu, Y. Yan, L. Zhao, M. Heimann, L. Akoglu, and D. Koutra, “Beyond homophily in graph neural networks: Current limitations and effective designs,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2020.
  * [70] E. Chien, J. Peng, P. Li, and O. Milenkovic, “Adaptive universal generalized pagerank graph neural network,” in _International Conference on Learning Representations (ICLR)_ , 2021.
  * [71] D. Lim, F. Hohne, X. Li, S. L. Huang, V. Gupta, O. Bhalerao, and S. N. Lim, “Large scale learning on non-homophilous graphs: New benchmarks and strong simple methods,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2021.
  * [72] N. Veldt, A. R. Benson, and J. Kleinberg, “Higher-order homophily is combinatorially impossible,” _SIAM Review (SIREV)_ , 2022.
  * [73] Q. Li, Z. Han, and X.-M. Wu, “Deeper insights into graph convolutional networks for semi-supervised learning,” in _AAAI conference on artificial intelligence (AAAI)_ , 2018.
  * [74] K. Oono and T. Suzuki, “Graph neural networks exponentially lose expressive power for node classification,” in _International Conference on Learning Representations (ICLR)_ , 2019.
  * [75] W. Cong, M. Ramezani, and M. Mahdavi, “On provable benefits of depth in training graph convolutional networks,” _Advances in Neural Information Processing Systems (NeurIPS)_ , 2021.
  * [76] T. N. Kipf and M. Welling, “Semi-supervised classification with graph convolutional networks,” in _International Conference on Learning Representations (ICLR)_ , 2017.
  * [77] T.-H. H. Chan, A. Louis, Z. G. Tang, and C. Zhang, “Spectral properties of hypergraph laplacian and approximation algorithms,” _Journal of the ACM (JACM)_ , 2018.
  * [78] T.-H. H. Chan and Z. Liang, “Generalizing the hypergraph laplacian via a diffusion process with mediators,” _Theoretical Computer Science_ , 2020\.
  * [79] J. Jo, J. Baek, S. Lee, D. Kim, M. Kang, and S. J. Hwang, “Edge representation learning with hypergraphs,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2021.
  * [80] V. Monga, Y. Li, and Y. C. Eldar, “Algorithm unrolling: Interpretable, efficient deep learning for signal and image processing,” _IEEE Signal Processing Magazine_ , 2021.
  * [81] T. Chen, X. Chen, W. Chen, H. Heaton, J. Liu, Z. Wang, and W. Yin, “Learning to optimize: A primer and a benchmark,” _arXiv preprint arXiv:2103.12828_ , 2021.
  * [82] E. Ryu, J. Liu, S. Wang, X. Chen, Z. Wang, and W. Yin, “Plug-and-play methods provably converge with properly trained denoisers,” in _International Conference on Machine Learning (ICML)_ , 2019.
  * [83] A. M. Teodoro, J. M. Bioucas-Dias, and M. A. T. Figueiredo, “Scene-adapted plug-and-play algorithm with convergence guarantees,” in _International Workshop on Machine Learning for Signal Processing (MLSP)_ , 2017.
  * [84] S. H. Chan, X. Wang, and O. A. Elgendy, “Plug-and-play admm for image restoration: Fixed-point convergence and applications,” _IEEE Transactions on Computational Imaging_ , vol. 3, no. 1, pp. 84–98, 2016.
  * [85] K. Gregor and Y. LeCun, “Learning fast approximations of sparse coding,” in _International Conference on Machine Learning (ICML)_ , 2010.
  * [86] B. Xin, Y. Wang, W. Gao, D. Wipf, and B. Wang, “Maximal sparsity with deep networks?” _Advances in Neural Information Processing Systems_ , vol. 29, 2016.
  * [87] J. Liu and X. Chen, “Alista: Analytic weights are as good as learned weights in lista,” in _International Conference on Learning Representations (ICLR)_ , 2019.
  * [88] X. Chen, J. Liu, Z. Wang, and W. Yin, “Theoretical linear convergence of unfolded ista and its practical weights and thresholds,” 2018.
  * [89] J. R. Hershey, J. L. Roux, and F. Weninger, “Deep unfolding: Model-based inspiration of novel deep architectures,” _arXiv preprint arXiv:1409.2574_ , 2014.
  * [90] Z.-Q. Wang, J. L. Roux, D. Wang, and J. R. Hershey, “End-to-end speech separation with unfolded iterative phase reconstruction,” in _Conference of the International Speech Communication Association (INTERSPEECH)_ , 2018.
  * [91] K. Zhang, W. Zuo, S. Gu, and L. Zhang, “Learning deep cnn denoiser prior for image restoration,” in _Conference on Computer Vision and Pattern Recognition (CVPR)_ , 2017.
  * [92] T. Meinhardt, M. Moller, C. Hazirbas, and D. Cremers, “Learning proximal operators: Using denoising networks for regularizing inverse imaging problems,” in _International Conference on Computer Vision (ICCV)_ , 2017\.
  * [93] J. H. R. Chang, C.-L. Li, B. Poczos, B. V. Kumar, and A. C. Sankaranarayanan, “One network to solve them all–solving linear inverse problems using deep projection models,” in _International Conference on Computer Vision (ICCV)_ , 2017.
  * [94] Y. Chen and T. Pock, “Trainable nonlinear reaction diffusion: A flexible framework for fast and effective image restoration,” _IEEE transactions on pattern analysis and machine intelligence_ , vol. 39, no. 6, pp. 1256–1272, 2016.
  * [95] S. Zheng, S. Jayasumana, B. Romera-Paredes, V. Vineet, Z. Su, D. Du, C. Huang, and P. H. Torr, “Conditional random fields as recurrent neural networks,” in _International Conference on Computer Vision (ICCV)_ , 2015.
  * [96] Z. Liu, X. Li, P. Luo, C. C. Loy, and X. Tang, “Deep learning markov random field for semantic segmentation,” _IEEE transactions on pattern analysis and machine intelligence_ , vol. 40, no. 8, pp. 1814–1828, 2017.
  * [97] C. J. Schuler, M. Hirsch, S. Harmeling, and B. Schölkopf, “Learning to deblur,” _IEEE transactions on pattern analysis and machine intelligence_ , vol. 38, no. 7, pp. 1439–1451, 2015.
  * [98] Y. Li, M. Tofighi, J. Geng, V. Monga, and Y. C. Eldar, “Efficient and interpretable deep blind image deblurring via algorithm unrolling,” _IEEE Transactions on Computational Imaging_ , vol. 6, pp. 666–681, 2020\.
  * [99] Y. Yang, T. Liu, Y. Wang, J. Zhou, Q. Gan, Z. Wei, Z. Zhang, Z. Huang, and D. Wipf, “Graph neural networks inspired by classical iterative algorithms,” in _International Conference on Machine Learning (ICML)_ , 2021\.
  * [100] X. Liu, W. Jin, Y. Ma, Y. Li, H. Liu, Y. Wang, M. Yan, and J. Tang, “Elastic graph neural networks,” in _International Conference on Machine Learning (ICML)_ , 2021.
  * [101] S. Chen, Y. C. Eldar, and L. Zhao, “Graph unrolling networks: Interpretable neural networks for graph signal denoising,” _IEEE Transactions on Signal Processing_ , vol. 69, pp. 3699–3713, 2021.
  * [102] J. Klicpera, A. Bojchevski, and S. Günnemann, “Predict then propagate: Graph neural networks meet personalized pagerank,” in _International Conference on Learning Representations_ , 2019.
  * [103] B. Chamberlain, J. Rowbottom, M. I. Gorinova, M. Bronstein, S. Webb, and E. Rossi, “Grand: Graph neural diffusion,” in _International Conference on Machine Learning_. PMLR, 2021, pp. 1407–1418.
  * [104] M. Thorpe, T. M. Nguyen, H. Xia, T. Strohmer, A. Bertozzi, S. Osher, and B. Wang, “Grand++: Graph neural diffusion with a source term,” in _International Conference on Learning Representations_ , 2022.
  * [105] S. Bai, J. Z. Kolter, and V. Koltun, “Deep equilibrium models,” _Advances in Neural Information Processing Systems (NeurIPS)_ , 2019.
  * [106] L. El Ghaoui, F. Gu, B. Travacca, A. Askari, and A. Tsai, “Implicit deep learning,” _SIAM Journal on Mathematics of Data Science_ , vol. 3, no. 3, pp. 930–958, 2021.
  * [107] H. Dai, Z. Kozareva, B. Dai, A. Smola, and L. Song, “Learning steady-states of iterative algorithms over graphs,” in _International Conference on Machine Learning (ICML)_ , 2018.
  * [108] J. Liu, K. Kawaguchi, B. Hooi, Y. Wang, and X. Xiao, “Eignn: Efficient infinite-depth graph neural networks,” _Advances in Neural Information Processing Systems (NeurIPS)_ , 2021.
  * [109] F. Gu, H. Chang, W. Zhu, S. Sojoudi, and L. El Ghaoui, “Implicit graph neural networks,” _Advances in Neural Information Processing Systems (NeurIPS)_ , 2020.
  * [110] Y. Yang, Y. Wang, Z. Huang, and D. Wipf, “Implicit vs unfolded graph neural networks,” _arXiv preprint arXiv:2111.06592_ , 2021.
  * [111] I. Amburg, N. Veldt, and A. R. Benson, “Clustering in graphs and hypergraphs with categorical edge labels,” in _Proceedings of The Web Conference_ , 2020\.
  * [112] P. S. Chodrow, N. Veldt, and A. R. Benson, “Hypergraph clustering: from blockmodels to modularity,” _Science Advances_ , 2021.
  * [113] J. H. Fowler, “Legislative cosponsorship networks in the us house and senate,” _Social networks_ , vol. 28, no. 4, pp. 454–465, 2006.
  * [114] ——, “Connecting the congress: A study of cosponsorship networks,” _Political Analysis_ , vol. 14, no. 4, pp. 456–487, 2006.
  * [115] Y. Deshpande, S. Sen, A. Montanari, and E. Mossel, “Contextual stochastic block models,” in _Advances in Neural Information Processing Systems (NeurIPS)_ , 2018.
  * [116] D. Ghoshdastidar and A. Dukkipati, “Consistency of spectral partitioning of uniform hypergraphs under planted partition model,” _Advances in Neural Information Processing Systems (NeurIPS)_ , 2014.
  * [117] I. C. C.-Y. Lin and I.-H. Wang, “Community detection in hypergraphs: Optimal statistical limit and efficient algorithms,” in _International Conference on Artificial Intelligence and Statistics (AISTATS)_ , 2018.
  * [118] T. Chen, S. Kornblith, M. Norouzi, and G. Hinton, “A simple framework for contrastive learning of visual representations,” in _International Conference on Machine Learning (ICML)_ , 2020.
  * [119] R. L. Murphy, B. Srinivasan, V. Rao, and B. Ribeiro, “Janossy pooling: Learning deep permutation-invariant functions for variable-size inputs,” in _International Conference on Learning Representations (ICLR)_ , 2018.

##  Appendix A Derivation of Eq. 3 and 4

We derive iterative update Eq. 3 and 4 as follows. Recall that our problem is
defined as:

$\displaystyle\min_{\bm{H}}\sum_{v\in\mathcal{V}}f(\bm{h}_{v};\bm{x}_{v})+\sum_{e\in\mathcal{E}}g_{e}(\bm{H}_{e})$

where $\bm{H}\in\mathbb{R}^{N}$ is the node feature matrix.
$\bm{H}_{e}=\begin{bmatrix}\bm{h}_{e,1}&\cdots&\bm{h}_{e,\lvert
e\rvert}\end{bmatrix}^{\top}\in\mathbb{R}^{\lvert e\rvert}$ collects the
associated node features inside edge $e$, where $\bm{h}_{e,i}$ corresponds to
the node features of the $i$-th node in edge $e$. We use $\bm{x}_{v}$ to
denote the (initial) attribute of node $v$. Here, neither
$f(\cdot;\bm{x}_{v})$ nor $g_{e}(\cdot)$ is necessarily continuous. To solve
such problem, we borrow the idea from ADMM [63]. We introduce an auxiliary
variable $\bm{R}_{e}=\bm{H}_{e}$ for every $e\in\mathcal{E}$ and reformulate
the problem as a constrained optimization problem:

$\displaystyle\min_{\bm{H}}\sum_{v\in\mathcal{V}}f(\bm{h}_{v};\bm{x}_{v})+\sum_{e\in\mathcal{E}}g_{e}(\bm{R}_{e})$
$\displaystyle\text{subject to }\bm{R}_{e}-\bm{H}_{e}=\bm{0},\,\forall
e\in\mathcal{E}$

By the Augmented Lagrangian Method (ALM), we assign a Lagrangian multiplier
$\bm{S}_{e}$ (scaled by $1/\lambda$) for each edge, then the objective
function becomes:

$\displaystyle\max_{\\{\bm{S}_{e}\\}_{e\in\mathcal{E}}}\min_{\bm{H},\\{\bm{R}_{e}\\}_{e\in\mathcal{E}}}\sum_{v\in\mathcal{V}}f(\bm{h}_{v};\bm{x}_{v})+\sum_{e\in\mathcal{E}}g(\bm{R}_{e})+\frac{\lambda}{2}\sum_{e\in\mathcal{E}}\lVert\bm{R}_{e}-\bm{H}_{e}+\bm{S}_{e}\|_{F}^{2}-\frac{\lambda}{2}\sum_{e\in\mathcal{E}}\lVert\bm{S}_{e}\rVert_{F}^{2}.$
(6)

We can iterate the following primal-dual steps to optimize Eq. 6. The primal
step can be computed in a block-wise sense:

$\displaystyle\bm{R}_{e}^{(t+1)}$
$\displaystyle\leftarrow\operatorname*{arg\,min}_{\bm{R}_{e}}\left\\{\frac{\lambda}{2}\left\lVert\bm{R}_{e}-\bm{H}_{e}^{(t)}+\bm{S}_{e}^{(t)}\right\rVert_{F}^{2}+g(\bm{R}_{e})\right\\}$
$\displaystyle=\operatorname{\textbf{prox}}_{g/\lambda}\left(\bm{H}_{e}^{(t)}-\bm{S}_{e}^{(t)}\right),\forall
e\in\mathcal{E},$ (7) $\displaystyle\bm{h}_{v}^{(t+1)}$
$\displaystyle\leftarrow\operatorname*{arg\,min}_{\bm{h}_{v}}\left\\{\frac{\lambda}{2}\sum_{e:v\in
e}\left\lVert\bm{h}_{v}-\bm{S}_{e,v}^{(t)}-\bm{R}_{e,v}^{(t+1)}\right\rVert_{F}^{2}+f(\bm{h}_{v};\bm{x}_{v})\right\\}$
$\displaystyle=\operatorname{\textbf{prox}}_{f(\cdot;\bm{x}_{v})/\lambda
d_{v}}\left(\frac{\sum_{e:v\in
e}(\bm{S}_{e,v}^{(t)}+\bm{R}_{e,v}^{(t+1)})}{d_{v}}\right),\forall
v\in\mathcal{V}.$ (8)

The dual step can be computed as:

$\displaystyle\bm{S}_{e}^{(t+1)}\leftarrow\bm{S}_{e}^{(t)}+\bm{R}_{e}^{(t+1)}-\bm{H}_{e}^{(t+1)},\forall
e\in\mathcal{E}.$ (9)

Denote $\bm{Q}_{e}^{(t+1)}:=\bm{S}_{e}^{(t)}+\bm{R}_{e}^{(t+1)}$ and
$\eta:=1/\lambda$, then the iterative updates become:

$\displaystyle\bm{Q}_{e}^{(t+1)}$
$\displaystyle=\operatorname{\textbf{prox}}_{\eta
g}(2\bm{H}_{e}^{(t)}-\bm{Q}_{e}^{(t)})+\bm{Q}_{e}^{(t)}-\bm{H}_{e}^{(t)},\
\mbox{for}\ e\in\mathcal{E},$ (10) $\displaystyle\bm{h}_{v}^{(t+1)}$
$\displaystyle=\operatorname{\textbf{prox}}_{\eta
f(\cdot;\bm{x}_{v})/d_{v}}\left(\frac{\sum_{e:e\in
v}\bm{Q}_{e,v}^{(t+1)}}{d_{v}}\right),\ \mbox{for}\ v\in\mathcal{V}.$ (11)

##  Appendix B Deferred Proofs

###  B.1 Proof of Proposition 2

###### Proof.

Define $\pi:[K]\rightarrow[K]$ be an index mapping associated with the
permutation matrix $\bm{P}\in\Pi(K)$ such that
$\bm{P}\bm{Z}=\begin{bmatrix}\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(K)}\end{bmatrix}^{\top}$.
To prove that $\nabla g(\cdot)$ is permutation equivariant, we show by using
the definition of partial derivatives. For any
$\bm{Z}=\begin{bmatrix}\bm{z}_{1}&\cdots&\bm{z}_{K}\end{bmatrix}^{\top}$, and
permutation $\pi$, we have:

$\displaystyle[\nabla g(\bm{P}\bm{Z})]_{i}$
$\displaystyle=\lim_{\delta\rightarrow
0}\frac{g(\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(i)}+\delta,\cdots,\bm{z}_{\pi(K)})-g(\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(i)},\cdots,\bm{z}_{\pi(K)})}{\delta}$
$\displaystyle=\lim_{\delta\rightarrow
0}\frac{g(\bm{z}_{1},\cdots,\bm{z}_{\pi(i)}+\delta,\cdots,\bm{z}_{K})-g(\bm{z}_{1},\cdots,\bm{z}_{\pi(i)},\cdots,\bm{z}_{K})}{\delta}$
$\displaystyle=[\nabla g(\bm{Z})]_{\pi(i)},$ (12)

where the second equality Eq. B.1 is due to the permutation invariance of
$g(\cdot)$. To prove Proposition 2 for the proximal gradient, we first define:
$\bm{H}^{*}=\operatorname{\textbf{prox}}_{g}(\bm{Z})=\operatorname*{arg\,min}_{\bm{H}}g(\bm{H})+\frac{1}{2}\lVert\bm{H}-\bm{Z}\rVert_{F}^{2}$
for some $\bm{Z}$. For arbitrary permutation matrix $\bm{P}\in\Pi(K)$, we have

$\displaystyle\operatorname{\textbf{prox}}_{g}(\bm{P}\bm{Z})$
$\displaystyle=\operatorname*{arg\,min}_{\bm{H}}g(\bm{H})+\frac{1}{2}\lVert\bm{H}-\bm{P}\bm{Z}\rVert_{F}^{2}=\operatorname*{arg\,min}_{\bm{H}}g(\bm{H})+\frac{1}{2}\lVert\bm{P}(\bm{P}^{\top}\bm{H}-\bm{Z})\rVert_{F}^{2}$
$\displaystyle=\operatorname*{arg\,min}_{\bm{H}}g(\bm{H})+\frac{1}{2}\lVert\bm{P}^{\top}\bm{H}-\bm{Z}\rVert_{F}^{2}=\operatorname*{arg\,min}_{\bm{H}}g(\bm{P}^{\top}\bm{H})+\frac{1}{2}\lVert\bm{P}^{\top}\bm{H}-\bm{Z}\rVert_{F}^{2}$
(13) $\displaystyle=\bm{P}\bm{H}^{*},$

where Eq. 13 is due to the permutation invariance of $g(\cdot)$. ∎

###  B.2 Proof of Theorem 1

###### Proof.

To prove Theorem 1, we first summarize one of our key results in the following
Lemma 2.

######  Lemma 2.

$\psi(\cdot):[0,1]^{K}\rightarrow\mathbb{R}^{K}$ is a permutation-equivariant
function if and only if there is a function
$\rho(\cdot):[0,1]^{K}\rightarrow\mathbb{R}$ that is permutation invariant to
the last $K-1$ entries, such that
$[\psi(\bm{Z})]_{i}=\rho(\bm{z}_{i},\underbrace{\bm{z}_{i+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{i-1}}_{K-1})$
for any $i$.

###### Proof.

(Sufficiency) Define $\pi:[K]\rightarrow[K]$ be an index mapping associated
with the permutation matrix $\bm{P}\in\Pi(K)$ such that
$\bm{P}\bm{Z}=\begin{bmatrix}\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(K)}\end{bmatrix}^{\top}$.
Then
$[\psi(\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(K)})]_{i}=\rho(\bm{z}_{\pi(i)},\bm{z}_{\pi(i+1)},\cdots,\bm{z}_{\pi(K)},\cdots,\bm{z}_{\pi(i-1)})$.
Since $\rho(\cdot)$ is invariant to the last $K-1$ entries,
$[\psi(\bm{P}\bm{Z})]_{i}=\rho(\bm{z}_{\pi(i)},\bm{z}_{\pi(i)+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{\pi(i)-1})=[\psi(\bm{Z})]_{\pi(i)}$.

(Necessity) Given a permutation-equivariant function
$\psi:[0,1]^{K}\rightarrow\mathbb{R}^{K}$, we first expand it to the following
form: $[\psi(\bm{Z})]_{i}=\rho_{i}(\bm{z}_{1},\cdots,\bm{z}_{K})$.
Permutation-equivariance means
$\rho_{\pi(i)}(\bm{z}_{1},\cdots,\bm{z}_{K})=\rho_{i}(\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(K)})$.
Suppose given an index $i$, consider any permutation $\pi:[K]\rightarrow[K]$,
where $\pi(i)=i$. Then, we have
$\rho_{i}(\bm{z}_{1},\cdots,\bm{z}_{i},\cdots,\bm{z}_{K})=\rho_{\pi(i)}(\bm{z}_{1},\cdots,\bm{z}_{i},\cdots,\bm{z}_{K})=\rho_{i}(\bm{z}_{\pi(1)},\cdots,\bm{z}_{i},\cdots,\bm{z}_{\pi(K)})$,
which implies $\rho_{i}:\mathbb{R}^{K}\rightarrow\mathbb{R}$ must be invariant
to the $K-1$ elements other than the $i$-th element. Now, consider a
permutation $\pi$ where $\pi(1)=i$. Then
$\rho_{i}(\bm{z}_{1},\bm{z}_{2},\cdots,\bm{z}_{K})=\rho_{\pi(1)}(\bm{z}_{1},\bm{z}_{2},\cdots,\bm{z}_{K})=\rho_{1}(\bm{z}_{\pi(1)},\bm{z}_{\pi(2)},\cdots,\bm{z}_{\pi(K)})=\rho_{1}(\bm{z}_{i},\bm{z}_{i+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{i-1})$,
where the last equality is due to our previous argument. This implies two
results. First, for all $i$,
$\rho_{i}(\bm{z}_{1},\bm{z}_{2},\cdots,\bm{z}_{i},\cdots,\bm{z}_{K}),\forall
i\in[K]$ should be written in terms of
$\rho_{1}(\bm{z}_{i},\bm{z}_{i+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{i-1})$.
Moreover, $\rho_{1}$ is permutation invariant to its last $K-1$ entries.
Therefore, we just need to set $\rho=\rho_{1}$ and broadcast it accordingly to
all entries. We conclude the proof. ∎

To proceed the proof, we bring in the following mathematical tools [50]:

######  Definition 5.

Given a vector
$\bm{z}=\begin{bmatrix}\bm{z}_{1},\cdots,\bm{z}_{K}\end{bmatrix}^{\top}\in\mathbb{R}^{K}$,
we define power mapping $\phi_{M}:\mathbb{R}\rightarrow\mathbb{R}^{M}$ as
$\phi_{M}(z)=\begin{bmatrix}z&z^{2}&\cdots&z^{M}\end{bmatrix}^{\top}$, and
sum-of-power mapping $\Phi_{M}:\mathbb{R}^{K}\rightarrow\mathbb{R}^{M}$ as
$\Phi_{M}(\bm{z})=\sum_{i=1}^{K}\phi_{M}(\bm{z}_{i})$, where $M$ is the
largest degree.

######  Lemma 3.

Let
$\mathcal{X}=\\{\begin{bmatrix}\bm{z}_{1},\cdots,\bm{z}_{K}\end{bmatrix}^{\top}\in[0,1]^{K}\text{
such that }\bm{z}_{1}<\bm{z}_{2}<\cdots<\bm{z}_{K}\\}$. We define mapping
$\tilde{\phi}:\mathbb{R}\rightarrow\mathbb{R}^{K+1}$ as
$\tilde{\phi}(z)=\begin{bmatrix}z^{0}&z^{1}&z^{2}&\cdots&z^{K}\end{bmatrix}^{\top}$,
and mapping $\tilde{\Phi}:\mathbb{R}^{K}\rightarrow\mathbb{R}^{K+1}$ as
$\tilde{\Phi}(\bm{z})=\sum_{i=1}^{K}\tilde{\phi}(\bm{z}_{i})$, where $M$ is
the largest degree. Then $\tilde{\Phi}$ restricted on $\mathcal{X}$, i.e.,
$\tilde{\Phi}:\mathcal{X}\rightarrow\mathbb{R}^{K+1}$, is a homeomorphism.

###### Proof.

Proved in Lemma 6 in [50]. ∎

We note that Definition 5 is slightly different from the mappings defined in
Lemma 3 [50] as it removes the constant (zero-order) term. Combining with
Lemma 3 [50] and results in [66], we have the following result:

######  Lemma 4.

Let
$\mathcal{X}=\\{\begin{bmatrix}\bm{z}_{1},\cdots,\bm{z}_{K}\end{bmatrix}^{\top}\in[0,1]^{K}\text{
such that }\bm{z}_{1}<\bm{z}_{2}<\cdots<\bm{z}_{K}\\}$, then there exists a
homeomorphism $\Phi_{M}:\mathcal{X}\rightarrow\mathbb{R}^{M}$ such that
$\Phi_{M}(\bm{z})=\sum_{i=1}^{K}\phi_{M}(\bm{z}_{i})$ where
$\phi_{M}:\mathbb{R}\rightarrow\mathbb{R}^{M}$ if $M\geq K$.

###### Proof.

For $M=K$, we choose $\phi_{K}$ and $\Phi_{K}$ to be the power mapping and
power sum with largest degree $K$ defined in Definition 5. We note that
$\tilde{\Phi}(\bm{z})=\begin{bmatrix}K&\Phi_{K}(\bm{z})^{\top}\end{bmatrix}^{\top}$.
Since $K$ is a constant, there exists a homeomorphism between the images of
$\tilde{\Phi}(\bm{z})$ and $\Phi_{K}(\bm{z})$. By Lemma 3,
$\tilde{\Phi}:\mathcal{X}\rightarrow\mathbb{R}^{K+1}$ is a homeomorphism,
which implies $\Phi_{K}(\bm{z}):\mathcal{X}\rightarrow\mathbb{R}^{K}$ is also
a homeomorphism.

For $M>K$, we first pad every input $\bm{z}\in\mathcal{X}$ with a constant
$k>1$ to be an $M$-dimension $\hat{\bm{z}}\in\mathbb{R}^{M}$. Note that
padding is homeomorphic since $k$ is a constant. All such $\hat{\bm{z}}$ form
a subset
$\mathcal{X}^{\prime}\subset\\{\begin{bmatrix}\bm{z}_{1},\cdots,\bm{z}_{M}\end{bmatrix}^{\top}\in[0,k]^{M}\text{
such that }\bm{z}_{1}<\bm{z}_{2}<\cdots<\bm{z}_{M}\\}$. We choose
$\hat{\phi}:\mathbb{R}\rightarrow\mathbb{R}^{M}$ to be power mapping and
$\hat{\Phi}:\mathcal{X}^{\prime}\rightarrow\mathbb{R}^{M}$ to be sum-of-power
mapping restricted on $\mathcal{X}^{\prime}$, respectively. Following [66], we
construct $\Phi_{M}(\bm{z})$ as below:

$\displaystyle\hat{\Phi}_{M}(\hat{\bm{z}})$
$\displaystyle=\sum_{i=1}^{K}\hat{\phi}_{M}(\bm{z}_{i})+\sum_{i=K+1}^{M}\hat{\phi}_{M}(k)=\sum_{i=1}^{K}\hat{\phi}_{M}(\bm{z}_{i})+\sum_{i=1}^{M}\hat{\phi}_{M}(k)-\sum_{i=1}^{K}\hat{\phi}_{M}(k)$
(14)
$\displaystyle=\sum_{i=1}^{K}(\hat{\phi}_{M}(\bm{z}_{i})-\hat{\phi}_{M}(k))+\sum_{i=1}^{M}\hat{\phi}_{M}(k)=\sum_{i=1}^{K}(\hat{\phi}_{M}(\bm{z}_{i})-\hat{\phi}_{M}(k))+M\hat{\phi}_{M}(k),$
(15)

which induces
$\sum_{i=1}^{K}(\hat{\phi}_{M}(\bm{z}_{i})-\hat{\phi}_{M}(k))=\hat{\Phi}_{M}(\hat{\bm{z}})-M\hat{\phi}_{M}(k)$.
Let $\phi_{M}(\bm{z})=\hat{\phi}_{M}(\bm{z}_{i})-\hat{\phi}_{M}(k)$, and
$\Phi_{M}(\bm{z})=\sum_{i=1}^{K}\phi_{M}(\bm{z}_{i})$. Since $[0,k]$ is
naturally homeomorphic to $[0,1]$, by our argument for $M=K$,
$\hat{\Phi}_{M}:\mathcal{X}^{\prime}\rightarrow\mathbb{R}^{M}$ is a
homeomorphism. This implies $\Phi_{M}:\mathcal{X}\rightarrow\mathbb{R}^{M}$ is
also a homeomorphism. ∎

It is straightforward to show the sufficiency of Theorem 1 by verifying that
for arbitrary permutation $\pi:[K]\rightarrow[K]$,
$[\psi(\bm{z}_{\pi(1)},\cdots,\bm{z}_{\pi(K)})]_{i}=\rho(\bm{z}_{\pi(i)},\sum_{j=1}^{K}\phi(\bm{z}_{j}))=[\psi(\bm{z}_{1},\cdots,\bm{z}_{K})]_{\pi(i)}$.
With Lemma 2 and 4, we can conclude the necessity of Theorem 1 by the
following construction:

  1. 1.

By Lemma 2, any permutation equivariant function
$\psi(\bm{z}_{1},\cdots,\bm{z}_{K})$ can be written as
$[\psi(\cdot)]_{i}=\tau(\bm{z}_{i},\bm{z}_{i+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{i-1})$
such that $\tau(\cdot)$ is invariant to the last $K-1$ elements.

  2. 2.

By Lemma 4, we know that there exists a homeomorphism mapping $\Phi_{K}$ which
is continuous, invertible, and invertibly continuous. For arbitrary
$\bm{z}\in[0,1]^{K}$, the difference between $\bm{z}$ and
$\Phi_{K}^{-1}\circ\Phi_{K}(\bm{z})$ is up to a permutation.

  3. 3.

Since $\tau(\cdot)$ is permutation invariant to the last $K-1$ elements, we
can construct the function
$\tau(\bm{z}_{i},\bm{z}_{i+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{i-1})=\tau(\bm{z}_{i},\Phi_{K-1}^{-1}\circ\Phi_{K-1}(\bm{z}_{i+1},\cdots,\bm{z}_{K},\cdots,\bm{z}_{i-1}))=\tilde{\rho}(\bm{z}_{i},\sum_{j\neq
i}\phi_{K-1}(\bm{z}_{j}))=\rho(\bm{z}_{i},\sum_{j=1}^{K}\phi_{K-1}(\bm{z}_{j}))$,
where $\tilde{\rho}(x,\bm{y})=\tau(x,\Phi_{K-1}^{-1}(\bm{y}))$ and
$\rho(x,\bm{y})=\tilde{\rho}(x,\bm{y}-\phi_{K-1}(x))$.

Since $\tau$, $\Phi_{K-1}^{-1}$, and $\phi_{K-1}$ are all continuous
functions, their composition $\rho$ is also continuous. ∎

##  Appendix C Additional Experiments and Details

###  C.1 Experiments on ED-HNNII

As we described in Sec. 3.2, an extension to our ED-HNN is to consider
$\mathcal{V}\rightarrow\mathcal{E}$ message passing and
$\mathcal{E}\rightarrow\mathcal{V}$ message passing as two equivariant set
functions. The detailed algorithm is illustrated in Algorithm 2, where the red
box highlights the major difference with ED-HNN.

Algorithm 2: ED-HNNII
---
Initialization: $\bm{H}^{(0)}=\bm{X}$ and three MLPs
$\hat{\phi},\,\hat{\rho},\,\hat{\varphi}$ (shared across $L$ layers).
For $t=0,1,2,...,L-1$, do:
1\. Designing the messages from $\mathcal{V}$ to $\mathcal{E}$:
$\bm{m}_{u\rightarrow e}^{(t)}=\hat{\phi}(\bm{m}_{u\rightarrow
e}^{(t-1)},\bm{h}_{u}^{(t)})$, for all $u\in\mathcal{V}$.
2\. Sum $\mathcal{V}\rightarrow\mathcal{E}$ messages over hyperedges
$\bm{m}_{e}^{(t)}=\sum_{u\in e}\bm{m}_{u\rightarrow e}^{(t)}$, for all
$e\in\mathcal{E}$.
3\. Broadcast $\bm{m}_{e}^{(t)}$ and design the messages from $\mathcal{E}$ to
$\mathcal{V}$: $\bm{m}_{e\rightarrow
v}^{(t)}=\hat{\rho}(\bm{h}_{v}^{(t)},\bm{m}_{e}^{(t)})$, for all $v\in e$.
4\. Update $\bm{h}_{v}^{(t+1)}=\hat{\varphi}(\bm{h}_{v}^{(t)},\sum_{e:u\in
e}\bm{m}_{e\rightarrow u}^{(t)},\bm{x}_{v},d_{v})$, for all $v\in\mathcal{V}$.

In our tested datasets, hyperedges do not have initial attributes. In our
implementation, we assign a common learnable vector for every hyperedge as
their first-layer features $\bm{m}_{u\rightarrow e}^{(0)}$. The performance of
ED-HNNII and the comparison with other baselines are presented in Table 4. Our
finding is that ED-HNNII can outperform ED-HNN on datasets with relatively
larger scale or more heterophily. And the average accuracy improvement by ED-
HNNII is around 0.5%. We argue that ED-HNNII inherently has more complex
computational mechanism, and thus tends to overfit on small datasets (e.g.,
Cora, Citeseer, etc.). Moreover, superior performance on heterophilic datasets
also implies that injecting more equivariance benefits handling heterophilic
data. From Table 4, the measured computational efficiency is also comparable
to ED-HNN and other baselines.

Table 4: Additional experiments and updated leaderboard with ED-HNNII.
Prediction accuracy (%). Bold font highlights when ED-HNNII outperforms the
original ED-HNN. Other details are kept consistent with Table 2.

| Cora | Citeseer | Pubmed | Cora-CA | DBLP-CA |  Training Time ($10^{-1}$ s)
---|---|---|---|---|---|---
HGNN [25] |  79.39 $\pm$ 1.36 |  72.45 $\pm$ 1.16 |  86.44 $\pm$ 0.44 |  82.64 $\pm$ 1.65 |  91.03 $\pm$ 0.20 |  0.24 $\pm$ 0.51
HCHA [27] |  79.14 $\pm$ 1.02 |  72.42 $\pm$ 1.42 |  86.41 $\pm$ 0.36 |  82.55 $\pm$ 0.97 |  90.92 $\pm$ 0.22 | 0.24 ± 0.01
HNHN [24] |  76.36 $\pm$ 1.92 |  72.64 $\pm$ 1.57 |  86.90 $\pm$ 0.30 |  77.19 $\pm$ 1.49 |  86.78 $\pm$ 0.29 |  0.30 $\pm$ 0.56
HyperGCN [23] |  78.45 $\pm$ 1.26 |  71.28 $\pm$ 0.82 |  82.84 $\pm$ 8.67 |  79.48 $\pm$ 2.08 |  89.38 $\pm$ 0.25 |  0.42 $\pm$ 1.51
UniGCNII [25] |  78.81 $\pm$ 1.05 |  73.05 $\pm$ 2.21 |  88.25 $\pm$ 0.40 |  83.60 $\pm$ 1.14 |  91.69 $\pm$ 0.19 |  4.36 $\pm$ 1.18
AllDeepSets [26] |  76.88 $\pm$ 1.80 |  70.83 $\pm$ 1.63 |  88.75 $\pm$ 0.33 |  81.97 $\pm$ 1.50 |  91.27 $\pm$ 0.27 |  1.23 $\pm$ 1.09
AllSetTransformer [26] |  78.58 $\pm$ 1.47 |  73.08 $\pm$ 1.20 |  88.72 $\pm$ 0.37 |  83.63 $\pm$ 1.47 |  91.53 $\pm$ 0.23 |  1.64 $\pm$ 1.63
ED-HNN |  80.31 $\pm$ 1.35 |  73.70 $\pm$ 1.38 |  89.03 $\pm$ 0.53 |  83.97 $\pm$ 1.55 |  91.90 $\pm$ 0.19 |  1.71 $\pm$ 1.13
ED-HNNII |  78.47 $\pm$ 1.62 |  72.65 $\pm$ 1.56 | 89.56 $\pm$ 0.62 |  82.17 $\pm$ 1.68 | 91.93 $\pm$ 0.29 |  1.85 $\pm$ 1.09
| Congress | Senate | Walmart | House | Avg. Rank |  Inference Time ($10^{-2}$ s)
HGNN [25] |  91.26 $\pm$ 1.15 |  48.59 $\pm$ 4.52 |  62.00 $\pm$ 0.24 |  61.39 $\pm$ 2.96 | 5.67 |  1.01 $\pm$ 0.04
HCHA [27] |  90.43 $\pm$ 1.20 |  48.62 $\pm$ 4.41 |  62.35 $\pm$ 0.26 |  61.36 $\pm$ 2.53 | 6.22 |  1.54 $\pm$ 0.18
HNHN [24] |  53.35 $\pm$ 1.45 |  50.93 $\pm$ 6.33 |  47.18 $\pm$ 0.35 |  67.80 $\pm$ 2.59 | 7.11 |  8.11 $\pm$ 0.05
HyperGCN [23] |  55.12 $\pm$ 1.96 |  42.45 $\pm$ 3.67 |  44.74 $\pm$ 2.81 |  48.32 $\pm$ 2.93 | 8.33 |  0.87 $\pm$ 0.06
UniGCNII [25] |  94.81 $\pm$ 0.81 |  49.30 $\pm$ 4.25 |  54.45 $\pm$ 0.37 |  67.25 $\pm$ 2.57 | 4.33 |  21.22 $\pm$ 0.13
AllDeepSets [26] |  91.80 $\pm$ 1.53 |  48.17 $\pm$ 5.67 |  64.55 $\pm$ 0.33 |  67.82 $\pm$ 2.40 | 5.89 |  5.35 $\pm$ 0.33
AllSetTransformer [26] |  92.16 $\pm$ 1.05 |  51.83 $\pm$ 5.22 |  65.46 $\pm$ 0.25 |  69.33 $\pm$ 2.20 | 3.33 |  6.06 $\pm$ 0.67
ED-HNN |  95.00 $\pm$ 0.99 |  64.79 $\pm$ 5.14 |  66.91 $\pm$ 0.41 |  72.45 $\pm$ 2.28 | 1.56 |  5.87 $\pm$ 0.36
ED-HNNII | 95.19 $\pm$ 1.34 |  63.81 $\pm$ 6.17 | 67.24 $\pm$ 0.45 | 73.95 $\pm$ 1.97 | 2.56 |  6.07 $\pm$ 0.40

###  C.2 Details of Benchmarking Datasets

Our benchmark datasets consist of existing seven datasets (Cora, Citeseer,
Pubmed, Cora-CA, DBLP-CA, Walmart, and House) from [26], and two newly
introduced datasets (Congress [114] and Senate [113]). For existing datasets,
we downloaded the processed version by [26]. For Congress [114] and
Senate[113], we used the same setting as [72]. In Congress dataset, nodes are
US Congresspersons and hyperedges are comprised of the sponsor and co-sponsors
of legislative bills put forth in both the House of Representatives and the
Senate. In Senate dataset, nodes are US Congresspersons and hyperedges are
comprised of the sponsor and co-sponsors of bills put forth in the Senate.
Each node in both datasets is labeled with political party affiliation. Both
datasets were from James Fowler’s data [113, 114]. We also list more detailed
statistical information on the tested datasets in Table 5.

Table 5: More dataset statistics. CE homophily is the homophily score [68]
based on CE of hypergraphs.

| Cora | Citeseer | Pubmed | Cora-CA | DBLP-CA | Congress | Senate | Walmart | House
---|---|---|---|---|---|---|---|---|---
# nodes | 2708 | 3312 | 19717 | 2708 | 41302 | 1718 | 282 | 88860 | 1290
# hyperedges | 1579 | 1079 | 7963 | 1072 | 22363 | 83105 | 315 | 69906 | 340
# features | 1433 | 3703 | 500 | 1433 | 1425 | 100 | 100 | 100 | 100
# classes | 7 | 6 | 3 | 7 | 6 | 2 | 2 | 11 | 2
avg. $d_{v}$ | 2.767 | 2.043 | 2.756 | 2.693 | 3.411 | 427.251 | 20.177 | 6.184 | 10.180
avg. $|e|$ | 1.748 | 1.541 | 1.963 | 1.929 | 2.213 | 8.654 | 9.531 | 3.461 | 8.056
CE Homophily | 0.897 | 0.893 | 0.952 | 0.803 | 0.869 | 0.555 | 0.498 | 0.530 | 0.509

###  C.3 Hyperparameters for Benchmarking Datasets

For a fair comparison, we use the same training recipe for all the models. For
baseline models, we precisely follow the hyperparameter settings from [26].
For ED-HNN, we adopt Adam optimizer with fixed learning rate=0.001 and weight
decay=0.0, and train for 500 epochs for all datasets. The standard deviation
is reported by repeating experiments on ten different data splits. We fix the
input dropout rate to be 0.2, and dropout rate to be 0.3. For internal MLPs,
we add a LayerNorm for each layer similar to [26]. Other parameters regarding
model sizes are obtained by grid search, which are enumerated in Table 6. The
search range of layer number is {1, 2, 4, 6, 8} and the hidden dimension is
{96, 128, 256, 512}. We find the model size is proportional to the dataset
scale, and in general heterophilic data need deeper architecture. For ED-
HNNII, due to the inherent model complexity, we need to prune model depth and
width to fit each dataset.

Algorithm 3: Contextual Hypergraph Stochastic Block Model
---
Initialization: Empty hyperedge set $\mathcal{E}=\emptyset$. Draw vertex set
$\mathcal{V}_{1}$ of 2,500 nodes with class 1.
Draw vertex set $\mathcal{V}_{2}$ of 2,500 nodes with class 2.
For $i=0,1,2,...,1,000$, do:
1\. Sample a subset $e_{1}$ with $\alpha_{1}$ nodes from $\mathcal{V}_{1}$.
2\. Sample a subset $e_{2}$ with $\alpha_{2}$ nodes from $\mathcal{V}_{2}$.
3\. Construct the hyperedge $\mathcal{E}\leftarrow\mathcal{E}\cup\\{e_{1}\cup
e_{2}\\}$.
Table 6: Choice of hyperparameters for ED-HNN. # is short for “number of”, hd.
stands for hidden dimension, cls. means the classifier. When number of MLP
layers equals to 0, the MLP boils down to be an identity mapping.

ED-HNN | Cora | Citeseer | Pubmed | Cora-CA | DBLP-CA | Congress | Senate | Walmart | House
---|---|---|---|---|---|---|---|---|---
# iterations | 1 | 1 | 8 | 1 | 1 | 4 | 8 | 6 | 8
# layers of $\hat{\phi}$ | 0 | 0 | 2 | 1 | 1 | 2 | 2 | 2 | 2
# layers of $\hat{\rho}$ | 1 | 1 | 2 | 1 | 1 | 2 | 2 | 2 | 2
# layers of $\hat{\varphi}$ | 1 | 1 | 1 | 1 | 1 | 2 | 2 | 2 | 2
# layer cls. | 1 | 1 | 2 | 2 | 2 | 2 | 2 | 2 | 2
$\operatorname{MLP}$ hd. | 256 | 256 | 512 | 128 | 128 | 156 | 512 | 512 | 256
cls. hd. | 256 | 256 | 256 | 96 | 96 | 128 | 256 | 256 | 128
ED-HNNII | Cora | Citeseer | Pubmed | Cora-CA | DBLP-CA | Congress | Senate | Walmart | House
# iterations | 1 | 1 | 4 | 1 | 1 | 4 | 4 | 4 | 4
# layers of $\hat{\phi}$ | 1 | 1 | 2 | 1 | 1 | 2 | 2 | 2 | 2
# layers of $\hat{\rho}$ | 1 | 1 | 2 | 1 | 1 | 2 | 2 | 2 | 2
# layers of $\hat{\varphi}$ | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1
# layer cls. | 1 | 1 | 2 | 2 | 2 | 2 | 2 | 2 | 2
$\operatorname{MLP}$ hd. | 128 | 128 | 512 | 128 | 128 | 156 | 512 | 512 | 256
cls. hd. | 96 | 96 | 256 | 96 | 96 | 128 | 256 | 256 | 128

###  C.4 Synthetic Heterophilic Datasets

We use the contextual hypergraph stochastic block model [115, 116, 117] to
synthesize data with controlled heterophily. The generated graph contains
5,000 nodes and two classes in total, and 2,500 nodes for each class. We
construct hyperedges by randomly sampling $\alpha_{1}$ nodes from class 1, and
$\alpha_{2}$ nodes from class 2 without replacement. Each hyperedge has a
fixed cardinality $|e|=\alpha_{1}+\alpha_{2}=15$. We draw 1,000 hyperedges in
total. The detailed data synthesis pipeline is summarized in Algorithm 3. We
use $\alpha=\min\\{\alpha_{1},\alpha_{2}\\}$ to characterize the heterophily
level of the hypergraph. For a more intuitive illustration, we list the CE
homophily corresponding to different $\alpha$ in Table 7.

Table 7: Correspondence between Heterophily $\alpha$ and CE Homophily. $\alpha$ | $1$ | $2$ | $3$ | $4$ | $6$ | $7$
---|---|---|---|---|---|---
CE Homophily | 0.875 | 0.765 | 0.672 | 0.596 | 0.495 | 0.474

Experiments on synthetic heterophilic datasets fix the training
hyperparameters and hidden dimension=256 to guarantee a fair parameter budget
(~1M). Baseline HNNs are all of one-layer architecture as they are not
scalable with the depth shown in Sec. 5.3. Since our ED-HNN adopts parameter
sharing scheme, we can easily repeat the diffusion layer twice to achieve
better results without parameter overheads.

###  C.5 Synthetic Diffusion Datasets and Additional Experiments

In order to evaluate the ability of ED-HNN to express given hypergraph
diffusion, we generate semi-synthetic diffusion data using the Senate
hypergraph [112] and synthetic node features. The data consists of 1,000 pairs
$(\bm{H}^{(0)},\bm{H}^{(1)})$. The initial node features $\bm{H}^{(0)}$ are
sampled from 1-dim Gaussian distributions with mean 0 and variance randomly
drawn between 1 and 100. That is, to generate a single instance of
$\bm{H}^{(0)}$, we first pick $\sigma$ uniformly from [1,10], and then sample
the coordinate entries as $\bm{h}^{(0)}_{v}\sim N(0,\sigma^{2})$. Then we
apply the gradient step in Eq. 2 to obtain the corresponding $\bm{H}^{(1)}$ .
For non-differentiable node or edge potentials, we adopt subgradients for
convenient computation. We fix the node potential as
$f(\bm{h}_{v};\bm{x}_{v})=(\bm{h}_{v}-\bm{x}_{v})^{2}$ where
$\bm{x}_{v}\equiv\bm{h}^{(0)}_{v}$. We consider 3 different edge potentials
from Example 1 with varying complexities: a) CE, b) TV ($p=2$) and c) LEC
($p=2$). For LEC, we set $\bm{y}$ as follows: if $|e|$ is even, then
$y_{i}=2/|e|$ if $i\leq|e|/2$ and $-2/|e|$ otherwise; if $|e|$ is odd, then
$y_{i}=2/(|e|-1)$ if $i\leq(|e|-1)/2$, $y_{i}=0$ if $i=(|e|+1)/2$, and
$y_{0}=-2/(|e|-1)$ otherwise. In order to apply the gradient step in Eq. 2 we
need to specify the learning rate $\eta$. We choose $\eta$ in a way such that
Var($\bm{H}^{(1)}$)$/$Var($\bm{H}^{(0)}$) does not vary too much among the
three different edge potentials. Specifically, we set $\eta=0.5$ for CE,
$\eta=0.02$ for TV and $\eta=0.1$ for LEC.

Beyond semi-synthetic diffusion data generated from the gradient step Eq. 2,
we also considered synthetic diffusion data obtained from the proximal
operators Eq. 3 and 4. We generated a random uniform hypergraph with 1,000
nodes and 1,000 hyperedges of constant hyperedge size 20. The diffusion data
on this hypergraph consists of 1,000 pairs $(\bm{H}^{(0)},\bm{H}^{(1)})$. The
initial node features $\bm{H}^{(0)}$ are sampled in the same way as before. We
apply the updates given by Eq. 3 and 4 to obtain $\bm{H}^{(1)}$. We consider
the same node potential and 2 edge potentials TV ($p=2$) and LEC ($p=2$). We
set $\eta=1/2$ for both cases. We show the results in Figure 4. The additional
results resonate with our previous results in Figure 3. Again, our ED-HNN
outperforms other baseline HNNs by a significant margin when hidden dimension
is limited.

Figure 4: Comparing the Powers to Represent known Diffusion (using ADMM with
the proximal operators in Eq. 3 and 4): MAE v.s. Latent Dimensions

###  C.6 More Complexity Analysis

In this section, we provide more efficiency comparison with HyperSAGE [28] and
LEGCN [30]. As we discussed in 3.4, HyperSAGE and LEGCN may be able to learn
hyperedge equivariant operator. However, both HyperSAGE and LEGCN need to
build node-hyperedge pair, which is memory consuming and unfriendly to
efficient message passing implementation. To demonstrate this we compare the
model our model with HyperSAGE and LEGCN in terms of training and inference
efficiency. We reuse the official code provided in OpenReview to benchmark
HyperSAGE. Since LEGCN has not released code, we reimplement it using Pytorch
Geometric Library. We can only test the speed of these models on Cora dataset,
because HyperSAGE and LEGCN cannot scale up as there is no efficient
implementation yet for their computational prohibitive preprocessing
procedures. The results are presented in Table 8. We note that HyperSAGE
cannot easily employ message passing between $\mathcal{V}$ and $\mathcal{E}$,
since neighbor aggregation step in HyperSAGE needs to rule out the central
node. Its official implementation adopts a naive “for”-loop for each forward
pass which is unable to fully utilize the GPU parallelism and significant
deteriorates the speed. LEGCN can be implemented via message passing, however,
the graphs expanded over node-edge pairs are much denser thus cannot scale up
to larger dataset.

Table 8: Performance and Efficiency Comparison with HyperSAGE and LEGCN. The prediction accuracy of HyperSAGE is copied from the original manuscript [28]. |  HyperSAGE [28] |  LEGCN [30] | ED-HNN | ED-HNNII
---|---|---|---|---
Training Time ($10^{-1}$ s) |  43.93 $\pm$ 2.15 |  0.56 $\pm$ 0.71 |  0.15 $\pm$ 0.68 |  0.25 $\pm$ 0.48
Inference Time ($10^{-2}$ s) |  297.57 $\pm$ 30.57 |  0.42 $\pm$ 0.06 |  0.15 $\pm$ 0.03 |  0.20 $\pm$ 0.08
Prediction Accuracy (%) |  69.30 $\pm$ 2.70 |  73.34 $\pm$ 1.06 |  80.31 $\pm$ 1.35 |  78.47 $\pm$ 1.62