arXiv CSopen access

The WidthWall: A Strict Expressivity Hierarchy for Hypergraph Neural Networks

arXiv:2605.13690v1 [cs.LG] 13 May 2026

The Width Wall: A Strict Expressivity Hierarchy for Hypergraph Neural Networks

Fengqing Jiang♣ Yuetai Li♣ Yichen Feng♣ Kaiyuan Zheng♣ Luyao Niu♣ Bhaskar Ramasubramanian ♢ Basel Alomair♠♣♡ Linda Bushnell♣ Radha Poovendran♣ ♣ University of Washington ♢ Western Washington University ♠ King Abdulaziz City for Science and Technology ♡ HUMAIN [email protected]

Abstract Hypergraphs provide a natural framework to model higher-order interactions in scientific, social, and biological systems. Hypergraph neural networks (HGNNs) aim to learn from such data, yet it remains unclear which higher-order structures these models can represent. We show that hypergraph expressivity is governed by which small patterns an architecture can detect and count. We formalize this via homomorphism densities, which measure how often a structural motif appears in a hypergraph. Combining classical homomorphism-count completeness with invariant approximation, we show that homomorphism densities generate all continuous hypergraph invariants and organize them into a strict hierarchy indexed by hypertree width. This yields a Width Wall: a fundamental architectural limit beyond which no hidden dimension, training procedure or fixed-depth HGNN can represent invariants requiring wider patterns. Our framework provides a unified characterization of 15 HGNN architectures, precisely identifies information lost by clique expansion, and motivates density-aware models that extend expressivity beyond bounded-width message passing. We experimentally validate this finding on an A PPLICATION N ODE C LASSIFICATION S UITE of real-world hypergraphs, where the Width Wall predicts when graph-reduction baselines fail and when density features help.

1

Introduction

Many real-world systems- such as co-authorship networks, protein complexes, legislative coalitionsare organized by simultaneous interactions among multiple entities. Hypergraphs provide a natural language for modeling these higher-order relations: a single hyperedge can bind an entire group, rather than forcing interaction into pairwise links [6, 7]. This modeling advantage has driven a rapid expansion of hypergraph neural networks (HGNNs) [14, 23, 48, 50]. However, HGNNs are not guaranteed to learn higher-order structures and relations effectively from hypergraphs. The central bottleneck lies in how the architecture summarizes each hyperedge. Some methods collapse group interactions into pairwise messages, while others preserve the incidence structure, but still rely on local, symmetric aggregation over node-hypergraph edge neighborhoods. For example, HGNN and HyperGCN transform each hyperedge into a graph-based propagation operator, so a k-way relation is ultimately processed through pairwise messages [23, 50]. In contrast, native models such as HNHN, UniGNN, AllSet, and ED-HNN retain the incidence structure, yet stil aggregate only local symmetric summaries [14, 20, 28, 48]. Thus, even architectures that preserve Preprint. Under review.

higher-order structure may fail to capture genuinely higher-order dependencies. This raises the question: when can these models actually learn genuine higher-order structure from the input? We address the problem through the lens of information loss. Consider legislative co-sponsorship: a single hyperedge may encode a multi-senator coalition whose political meaning depends on the entire group. If an HGNN first reduces that coalition to pairwise edges, a bipartisan coalition of 8 senators can become indistinguishable from two partisan blocs of 4 after clique expansion, even though their downstream implications may differ sharply. We empirically validate this failure mode on Senate co-sponsorship data [24] using: (1) HyperGCN, which projects hyperedges to pairwise edges and reaches 55% node-classification accuracy, and (2) a model that aggregates over full hyperedge memberships and reaches 93%. This gap validates architectural information loss: once the relevant higher-order pattern is discarded, no increase in hidden dimension or training can recover it. In this paper, we characterize information loss via a notion named hypergraph expressivity. Our key insight is that hypergraph expressivity is governed by which small higher-order patterns an architecture can detect and count. We formalize these counts using homomorphism densities, which measure how often a pattern appears in the input hypergraph. Using classical homomorphism-count completeness, we show that these densities separate non-isomorphic hypergraphs and generate all continuous hypergraph invariants- that is, all continuous permutation-invariant functions on the compact hypergraph-tensor space. Consequently, any continuous structural property of a fixed-size hypergraph can be approximated from pattern densities. This perspective yields a diagnostic basis that makes explicit which information is lost when a model replaces a hypergraph. The density basis induces a graded expressivity scale for HGNNs. Let Hr denote the class of continuous invariants generated by homomorphism densities of patterns with generalized hypertree width at most r, and let Iall denote all continuous hypergraph invariants. We show that these classes form a strict hierarchy, H1 ⊊ H2 ⊊ · · · ⊊ Iall . An HGNN whose aggregation scheme reaches level r can compute invariants generated by patterns of width ≤ r, but cannot represent invariants requiring larger width. We term this obstruction the Width Wall: an architectural ceiling determined by pattern densities accessible to the model, rather than by hidden dimension or optimization. This correspondence gives a principled way to classify existing HGNN architectures within the hierarchy.

2

Homomorphism Features and the Expressivity Framework

Notation and Background. We write [n] = {1, . . . , n}. A k-uniform hypergraph is H = ([n], E),
where E ⊆ [n] k is a set of unordered k-node hyperedges. The adjacency tensor of H is the order-k array A with Ai1 ,...,ik = 1 if (i1 , . . . , ik ) ∈ E, and 0 otherwise. Because hyperedges are unordered, A is symmetric, i.e., invariant under permutations of its k indices. Let Sk,n denote the space of order-k, dimension-n symmetric tensors. Thus, A ∈ Sk,n for any k-uniform hypergraph. The case k = 2 recovers ordinary graphs, where S2,n is the space of symmetric n × n adjacency matrices. k,n Weighted hypergraphs allow entries in [0, 1], yielding a compact domain [0, 1]S . For non-uniform hypergraphs with maximum edge size K, we represent adjacency as A = (A(2) , . . . , A(K) ), where A(j) ∈ Sj,n records hyperedges of size j. Finally, Sym(n) denotes the group of vertex permutations, acting on A by simultaneously relabeling all tensor indices. A hypergraph invariant is a function unchanged under this action. 2.1

Tensor Domain and Hypergraph Invariants

Overview. In graph theory, homomorphism counts provide a canonical coordinate system for finite graphs. They characterize graphs up to isomorphism, underlie dense graph limits, and connect to logical and WL-style expressivity [19, 37, 38]. We extend this viewpoint to hypergraphs and use it to bridge combinatorial patterns and neural architectures. In this section, we characterize the expressivity of an HGNN through the class of hypergraph homomorphism densities it can access, which measure how often a structural pattern appears in the input hypergraph [22, 36, 37, 53]. We then treat these densities as invariant coordinates: each pattern defines one label-independent coordinate of the hypergraph, and different architecture classes expose different subsets of these coordinates. Further algebraic structure appears in Appendix A.1, with an independent invariant-theoretic derivation in Appendix B; concrete examples of the resulting hypergraph invariants are in Appendix I. 2

Our discussion proceeds in three steps. First, we represent hypergraphs as symmetric tensors, making invariance to vertex ordering within each hyperedge. Second, we introduce pattern densities and their aligned-template counterparts. Finally, we show that pattern densities generate all continuous invariant functions. Under this view, expressivity reduces to asking which coordinates an architecture can actually compute. If an architecture cannot access patterns generating a target invariant, increasing hidden dimension alone will not recover that invariant. Tensor Domain and Invariance. To analyze what a model can learn from hypergraph structure, we separate structural information from arbitrary vertex labels. We use the adjacency tensor representation as the common input space for hypergraph invariants, which later serves as features for analyzing HGNN expressivity. A hypergraph invariant must be unchanged under vertex relabeling. For a permutation Π ∈ Sym(n) of the vertex set [n], the relabeling action is (Π · A)i1 ,...,ik = AΠ−1 (i1 ),...,Π−1 (ik ) . Together, tensor symmetry and vertex-relabeling action define the domain of continuous hypergraph k,n invariants. Unweighted hypergraphs correspond to {0, 1} vertices of the compact cube [0, 1]S , while weighted hypergraphs fill the cube. Accordingly, a continuous hypergraph invariant is a continuous function on this tensor domain that is constant on Sym(n)-orbits. For non-uniform hypergraphs with maximum edge size K, the same construction applies component-wise to the tuple (A(2) , . . . , A(K) ), where A(j) ∈ Sj,n records hyperedges of size j. We now formalize hypergraph invariants as label-independent structural functions. Definition 2.1 (Hypergraph invariant). A function f : Sk,n → R is a hypergraph invariant if f (Π · A) = f (A) for all A ∈ Sk,n and all Π ∈ Sym(n). 2.2

Homomorphism Densities and Pattern Alignment

We use two complementary families of invariants. The first is the family of homomorphism densities, which count how often a finite pattern appears in the target hypergraph and provide the countbased algebra underlying our completeness and density results. The second is the family of pattern alignment scores, which compare a target hypergraph to a template under optimal vertex relabeling and provide optimization-based template invariants. Let Fk denote the set of all finite k-uniform hypergraph patterns. Let F be such a pattern and H be a target k-uniform hypergraph. A map ϕ : V (F ) → V (H) places vertices of the pattern into the target, where V (·) denotes the vertex set of a hypergraph. Definition 2.2 (Hypergraph homomorphism [36, 37]). Let F = ([p], EF ) and H = ([n], EH ) be k-uniform hypergraphs. A homomorphism ϕ : V (F ) → V (H) is a map such that hyperedge (ϕ(i1 ), . . . , ϕ(ik )) ∈ EH whenever hyperedge (i1 , . . . , ik ) ∈ EF . The homomorphism count is hom(F, H) = {ϕ : ϕ is a homomorphism from [p] to [n]} . The homomorphism density is the fraction of all maps from [p] to [n] that are homomorphisms: t(F, H) =

hom(F, H) . np k,n

For weighted target hypergraphs with adjacency tensor A ∈ [0, 1]S and a pattern F = ([p], EF ), the homomorphism density becomes a multilinear polynomial [22, 38, 53]: X Y 1 t(F, A) = p Aϕ(e) , (1) n ϕ:[p]→[n] e∈EF

where Aϕ(e) = Aϕ(i1 ),...,ϕ(ik ) for e = (i1 , . . . , ik ). Thus, t(F, A) is a polynomial of degree |EF | in the tensor entries. These features are continuous and permutation-invariant, and their pattern encodes the higher-order interaction being measured. For example, a single edge measures edge density, a star measures local overlap, and more entangled patterns capture higher-order coordination. Thus, the family {t(F, ·) : F ∈ Fk } should be viewed not merely as a collection of statistics, but as a coordinate map that records all pattern measurements of the hypergraph. Homomorphism densities form an algebra under pointwise addition and multiplication.In particular, multiplying two densities is equivalent to counting the disjoint union of the two query patterns 3

t(F1 , A) · t(F2 , A) = t(F1 ⊔ F2 , A) [37, Ch. 5]. We write Ak for the unital algebra generated by {t(F, ·) : F ∈ Fk }, with constants given by the empty pattern t(∅, A) = 1. Further details on algebraic structure appear in Appendix A.1. For the architectural constructions later in the paper, we also introduce pattern alignment as a templatebased analogue of pattern counting. Rather than averaging over all vertex maps, this feature compares a template tensor with the target under the best vertex relabeling. It is therefore invariant to vertex labels and provides an optimization-based invariant. Definition 2.3 (Pattern alignment score). For any template tensor P ∈ Sk,n , the pattern alignment score is ΘP (A) = max ⟨P, Π · A⟩. (2) Π∈Sym(n)

The score ΘP is invariant under vertex relabeling of A by construction. Homomorphism densities provide the count-based algebra used in the completeness and density proofs, while alignment scores provide the optimization-based template features. 2.3

Completeness and Density

The preceding constructions provide concrete invariant coordinates. We now ask whether they capture all structural information available to a continuous invariant learner. Finite densities model structural information. We prove that finite collections of pattern densities can approximate arbitrary continuous invariant targets on the fixed-n tensor domain. k,n

Theorem 2.4 (Pattern-density universality with fixed n). Fix n ≥ 1. Let K = [0, 1]S be the space of weighted k-uniform hypergraphs on [n], and let C(K/Sym(n)) denote the Banach space of continuous Sym(n)-invariant functions on K with the supremum norm. Then for every continuous hypergraph invariant f and every ϵ > 0, there exist patterns F1 , . . . , Fm and a polynomial q : Rm → R such that
sup f (A) − q t(F1 , A), . . . , t(Fm , A) < ϵ. A∈K

Theorem 2.4 shows that any fixed-size structural property that is continuous and invariant to vertex relabeling (e.g., a spectral, cut, or topological statistic) can be approximated by a polynomial in finitely many pattern densities. In this sense, pattern densities form a complete coordinate system for continuous invariants on the fixed tensor domain: once all pattern coordinates are available, no additional continuous label-independent information is required. The fixed n setting matches the neural-network regime, which operates on a tensor of fixed size. Extensions to non-uniform hypergraphs are given in Appendix C.3; extensions to hypergraphon limits are in [22, 53].

3

Hypergraph Invariants Induce a Strict Hierarchy

We have shown that each pattern density t(F, A) defines an invariant coordinate and that these coordinates are complete for continuous invariants on a fixed tensor domain. We now investigate which of these coordinates are accessible to different HGNN architectures. 3.1

Architecture Scope

Before comparing architectures, we fix the scope of our claims. Throughout, we consider label-free local hypergraph neural architectures satisfying the following conditions. Assumption 3.1 (Architecture scope). The architecture satisfies: (C1) permutation equivariance, with no node identifiers or order-dependent operations; (C2) locality, with aggregation only from incident hyperedge neighborhoods; (C3) symmetric aggregation, such as sum, mean, max, DeepSets, or attention over a multiset; (C4) no symmetry-breaking features, such as Laplacian eigenvectors, random IDs, or global positional encodings; and (C5) deterministic structure-dependent preprocessing. Key observation. HGNN architectures can be compared by the pattern densities they can represent. Clique-expansion methods, such as HGNN [23] and HyperGCN [50], factor through the pairwise projection ϕ : Sk,n → Sn , retaining only graph-level structure. Native hyperedge message-passing 4

methods, such as HNHN [20], UniGNN [28], AllSet [14], and ED-HNN [48], preserve incidence but access only patterns that can be assembled through local aggregation. In contrast, the full patterndensity algebra from Theorem 2.4 imposes no such restriction. The hierarchy below formalizes these successive enlargements as nested invariant classes. For an HGNN architecture class F, I(F) is the set of continuous hypergraph invariants that functions in F can approximate. Thus I(F) records the invariant information representable by the architecture. We first formalize the loss induced by the clique-expansion methods. Proposition 3.2 (Clique-expansion limitation). Let ϕ : Sk,n → Sn be the clique-expansion map. Suppose a clique-expansion HGNN has the form fθ (A) = gθ (ϕ(A)), where gθ : Sn → Rd is a graph neural network applied to the clique-expanded graph. If two hypergraphs have the same clique expansion, ϕ(A1 ) = ϕ(A2 ), then fθ (A1 ) = fθ (A2 ) for every parameter choice θ. Proposition 3.2 shows clique-expansion methods discard all information not determined by the pairwise shadow. However, native incidence-based HGNNs remain limited to patterns accessible through local aggregation. 3.2

Width, Invariance Classes, and the Width Wall

We investigate how hyperedge structure must be coordinated jointly in order to count a pattern. We measure this coordination complexity using generalized hypertree width [26]. Definition 3.3 (Generalized hypertree width [26]). A tree decomposition of a hypergraph F is a tree T whose nodes t are assigned vertex sets Bt ⊆ V (F ), called bags, covering all vertices and hyperedges, such that for each v ∈ V (F ) the bags containing v form a connected subtree. The generalized hypertree width Wgh (F ) is the minimum over all tree decompositions of F of the maximum number of hyperedges of F needed to cover any single bag. Intuitively, Wgh (F ) measures how many hyperedges must be tracked jointly when counting F : tree-like patterns have Wgh = 1, while dense or highly connected patterns require larger width. Generalized hypertree width induces a graded structure on the density coordinates: each pattern density t(F, ·) has a width, and increasing the allowed width exposes strictly richer invariant information. The architectural question then becomes which portion of this graded density algebra a given model can access. Clique-expansion models access only graph-level coordinates of the pairwise shadow, native message-passing models access a bounded-width portion of the hypergraph coordinates, and the full density algebra imposes no width restriction. We formalize these three regimes as invariant function classes. Definition 3.4 (Invariant classes via pattern complexity). For a width budget R < ∞, define the following classes of continuous hypergraph invariants: 1. ICE : invariants computable from pattern densities of graph patterns on the clique expansion (i.e., functions factoring through ϕ : Sk,n → Sn ). (R)

2. INH : invariants generated by ICE together with hypergraph pattern densities t(F, ·) for patterns satisfying Wgh (F ) ≤ R. 3. Iall : all continuous hypergraph invariants on Sk,n (generated by the full pattern density algebra Ak , by Theorem 2.4). For a positive integer r, the width-r level consists of pattern densities t(F, ·) with Wgh (F ) ≤ r. Equivalently, level r contains invariants generated by hypergraph patterns whose counting can be organized using at most r covering hyperedges in any bag of a generalized hypertree decomposition. Definition 3.5 (Width Wall). For an architecture class whose invariant information is contained in Hr , the Width Wall at level r is the obstruction formed by invariants in Iall \ Hr . Equivalently, it is the boundary between the pattern densities the architecture can access, those with Wgh (F ) ≤ r, and the wider pattern densities needed to separate some hypergraphs. Thus, the Width Wall is a structural barrier. Once an architecture exposes only width-r coordinates, no increase in hidden dimension, readout complexity, or training can recover invariants with Wgh > r. 5

3.3

Expressivity Hierarchy and Architectural Consequences

Invariant classes translate architectural mechanisms into invariant function classes. ICE contains (R) exactly the invariants that depend only on the pairwise shadow. INH is a native-incidence reference class: it retains the CE information and adds hypergraph motifs whose joint coordination complexity is at most R. The parameter R should be read as an architectural width budget, determined by the message-passing mechanism and depth through Proposition 3.8. Finally, Iall is the ideal reference class from Theorem 2.4, with no pattern-complexity restriction imposed. Theorem 3.6 (Strict expressivity hierarchy). Fix a native-incidence reference class with width budget (k) R ≥ 3. For k-uniform hypergraphs with k ≥ 3 on n ≥ max{13, NCFI (R)} vertices, we have (R) (k) ICE ⊊ INH ⊊ Iall , where NCFI (R) denotes the number of vertices in the smallest k-uniform Cai–Furer–Immerman (CFI) lift whose base graph has girth greater than 2R + 1. The proof and explicit separation witnesses are given in Appendices A and F. Theorem 3.6 says that each relaxation of the information bottleneck strictly adds invariant information. Moving from clique expansion to native incidence features recovers higher-order structure lost by the pairwise projection. Moving from any fixed width budget R to the full pattern-density algebra recovers invariants whose witnessing patterns require larger structural complexity. Theorem 3.7 (Infinite expressivity hierarchy). For each r ≥ 1, define Hr as the class of continuous invariants approximable using pattern densities t(F, ·) where F has generalized hypertree width (k) Wgh (F ) ≤ r. Then for every r ≥ 1, k ≥ 3, and n ≥ NCFI (r), the level-r inclusion is strict, i.e., Hr ⊊ Hr+1 . Each strict inclusion is witnessed by explicit families of non-isomorphic hypergraph pairs that are Hr -indistinguishable but Hr+1 -distinguishable. Theorem 3.7 shows that the Width Wall is strict at every level r: the CFI construction provides explicit invariants that cannot be approximated by any architecture restricted to patterns with Wgh (F ) ≤ r, regardless of hidden dimension or training procedure (Corollary A.4 in Appendix A.4). Accessing these invariants requires increasing the available hypertree width. Appendix D summarizes how modern HGNN architectures map into this hierarchy. Appendix K studies related structural frontiers for Sinkhorn relaxations and spectral-width thresholds, with the Steiner-pair verification in Appendix L. Proposition 3.8 (Architecture–pattern width correspondence). Suppose Assumption 3.1 holds. A message-passing HGNN with L layers and aggregate-combine operations over hyperedge neighborhoods is bounded by homomorphism information from patterns of generalized hypertree width at most ck L, where ck ≥ 1 is a constant depending on the arity k and the aggregation scheme (for standard k-ary message passing, ck = 1 [44]). In particular, it cannot distinguish hypergraphs that agree on all such counts. Conversely, sufficiently expressive aggregate-combine maps can realize the corresponding bounded-width homomorphism-count features. An L-layer HGNN uses at most ck L + 1 variables to describe its computation tree, and the resulting logic C ck L+1 captures exactly the homomorphism counts of patterns with Wgh (F ) ≤ ck L. Crucially, L is a fixed architectural parameter that does not grow with n, so the set of accessible patterns has bounded hypertree width regardless of the input size. In the terminology of Definition 3.5, such a model faces the Width Wall at level ck L. We use the known WL–homomorphism correspondence: two k-uniform hypergraphs are r-GWL equivalent if and only if they agree on all homomorphism counts hom(F, ·) from patterns F with Wgh (F ) ≤ r [43, 44]. Quantitative Learnability. The density theorem shows that pattern densities can approximate any continuous invariant, but it does not specify how many samples, features, or parameters are needed. Appendix C collects three quantitative complements that clarify this point. (a) Concentration. For any fixed pattern F , its density t(F, H) can be estimated by sampling random vertex maps. Additive accuracy ϵ with failure probability δ requires O(log(1/δ)/ϵ2 ) samples, independent of the number of vertices n (Proposition C.5). (b) Dimension counting. The Molien series counts the invariant polynomial degrees of freedom. When we restrict to patterns with Wgh (F ) ≤ r, the available polynomial invariant space is strictly smaller; this gives an algebraic measure of the approximation cost imposed by a bounded-width architecture (Proposition B.6). (c) Generalization. For the template-alignment architecture introduced later in Section 4, the Rademacher bound scales p D as O(BP BA BW D log m/N ). Here N is the number of training examples, D and m are the depth 6

and width of the readout MLP, and BP , BA , and BW bound the template, input, and readout weights respectively. The number of learned templates J enters only logarithmically (Proposition C.3).

4

InvNet: An Invariant-Theoretic Architecture

We propose I NV N ET, an HGNN designed to extend beyond local message passing. I NV N ET augments a standard local HGNN with a global branch that learns permutation-invariant pattern-alignment features. The model thus combines two complementary sources of information: local incidence-based messages and global template scores that compare the input hypergraph to higher-order structural patterns. We emphasize that I NV N ET serves primarily as a proof-of-concept architecture illustrating how invariant features can extend expressivity beyond bounded-width message passing. Concretely, I NV N ET computes 
fθ (A) = MLPθ ReadOut h(L) {z } | local branch

  Θ̂P1 (A) − α1 , . . . , Θ̂PJ (A) − αJ , | {z }

(3)

invariant branch

(L)

where h are node embeddings obtained from L rounds of hypergraph message passing. The templates P1 , . . . , PJ are learned global patterns, and each feature Θ̂Pj (A) approximates the optimal alignment between Pj and the input hypergraph over all vertex relabelings. Since exact optimization over Sym(n) is a combinatorial assignment problem, we relax permutation matrices to the Birkhoff polytope Bn and compute approximate solutions via Sinkhorn iterations. See Appendix G for details. This construction illustrates how the invariant hierarchy can guide model design. The local branch retains expressive power of native message-passing HGNNs, while the invariant branch introduces global template features that are not contrained by the bounded-width computation tree of a fixeddepth HGNN. In terms of Definition 3.5, the invariant branch is designed to introduce coordinates beyond the local branch’s Width Wall rather than merely increasing capacity within the same bounded-width class. Theorem 4.1 (InvNet expressivity). I NV N ET strictly extends bounded-depth message passing INH ⊊ ÎInvNet ⊆ Iall while remaining within continuous invariants.

5

Experiments

Our experiments evaluate whether the invariant hierarchy predicts empirical failure modes of HGNN architectures on real-world hypergraph learning tasks. Specifically, we test three claims suggested by the theory: (i) clique-expansion (CE) models can lose higher-order information; (ii) native messagepassing/set-function HGNNs can hit the Width Wall (Definition 3.5) when the useful invariant lies outside the bounded-width summaries exposed by local aggregation; and (iii) density-aware features improve performance when the backbone leaves structural headroom. These experiments directly probe whether architectural performance aligns with invariant classes predicted by the hierarchy. 5.1

Experimental Setup

Models. We compare three architecture classes. The CE tier contains HGNN [23] and HyperGCN [50], which reduce hyperedges to pairwise graph structure. The native tier includes AllDeepSets, AllSetTransformer [14], HNHN [20], and UniGNN [28], which operate directly on incidence structure via bounded-depth message passing or set aggregation. The density-aware tier contains our D ENS N ET-D (PDN features fused with an AllDeepSets backbone; Appendices G–H). D ENS N ET-D is a scalable density-aware HGNN derived from the invariant-feature principle behind I NV N ET. D ENS N ET-D estimates per-node pattern densities by Monte Carlo sampling over each node’s star neighborhood and fuses them with an AllDeepSets backbone through a learned gate. This design preserves the intended role of invariant structural coordinates while avoiding the quadratic Sinkhorn bottleneck on large hypergraphs. For ANCS node classification we additionally report MLP, ED-HNN [48], and SheafHyperGNN [21] where applicable. Task and Datasets. The A PPLICATION N ODE C LASSIFICATION S UITE (ANCS) evaluates transductive node classification on empirical hypergraphs. Nodes carry benchmark-provided features, 7

hyperedges encode observed higher-order relations such as co-citation groups, legislative coalitions, or gene–disease associations, and the model predicts node labels under the standard random-split protocol described in Appendix E. ANCS includes 7 empirical hypergraph node-classification benchmarks, grouped by mean hyperedge size (additional detail in Appendices E.1 and M): • Cora (pairwise, k̄ = 3.03): a co-citation hypergraph with papers as nodes and hyperedges induced by common citing papers; node labels are paper topics and features are bag-of-words vectors [14, 50]. • Citeseer (pairwise, k̄ = 3.20): a co-citation hypergraph constructed analogously to Cora, with document-topic labels and sparse text features [14, 50]. • PubMed (mixed, k̄ = 4.35): a larger co-citation hypergraph over biomedical papers, using the benchmark-provided paper features and disease-category labels [14, 50]. • Cora-CA (mixed, k̄ = 4.28): a co-authorship version of Cora in which each hyperedge groups the co-authors of one publication, with the same paper-topic prediction target [14]. • House (higher-order, k̄ = 34.89): a legislative committee hypergraph from the U.S. House, where hyperedges are committee memberships and task is binary party classification [14, 24]. • Senate Bills (higher-order, k̄ = 7.96): a dense legislative co-sponsorship hypergraph in which each bill forms a hyperedge over sponsoring senators, again evaluated with binary party labels [14, 24]. • Gene-Disease (higher-order, k̄ = 14.0): a biomedical association hypergraph whose hyperedges connect genes to disease contexts, with a 21-class node-classification target under the benchmark protocol [14]. Metrics. All tables report test classification accuracy in percent. For ANCS benchmarks, we report mean accuracy and standard deviations over 5–10 seeds. Bold indicates the best result in each column and underlining denotes the second best. Parameter setup. All ANCS methods use Adam with learning rate 10−3 , weight decay 0, gradient clipping 1.0, label smoothing 0.1, cosine annealing, and 200 epochs (400 for Senate Bills). Baselines use hidden dimension 64 and dropout 0.2. D ENS N ET-D uses hidden dimension 128, dropout 0.5, 200 Monte Carlo density samples, and a profiler-set gate: b = −5.0 for pairwise datasets, b = −2.0 for mixed datasets, and b = 0.0 for higher-order datasets. The baseline settings follow the AllSet benchmark protocol [14]; D ENS N ET-D’s capacity and sampling parameters are chosen by validation grid search, while its gate and freezing schedule are determined by the DATASET P ROFILER from dataset statistics. Full training, hyperparameter-selection, and compute-budget details are in Appendices N and O. 5.2

Experimental Results

Finding 1: architecture class matters most when hyperedges are genuinely higher order. The ANCS results show that performance is primarily determined by architecture class. We observe that models that retain incidence-level or density-aware summaries perform best when the dataset contains label-relevant higher-order structure. This trend appears across most datasets. CE models consistently underperform the best native or density-aware method on Cora, Citeseer, PubMed, CoraCA, House, Senate, and Gene-Disease. The largest gap occurs on Senate Bills, where D ENS N ET-D and AllDeepSets both reach 92.7% while HyperGCN reaches only 55.4%. Here, reducing each bill coalition to pairwise edges appears to discard signals that native set-function models can preserve. Importantly, performance does not follow a simple monotonic trend in mean hyperedge size. HGNN and HyperGCN do not degrade uniformly from pairwise to mixed to higher-order regimes. Instead, the relevant comparison is within each dataset, where architecture choices determine which structural summaries are available. Appendix P analyzes the main anomaly and ceiling cases, including Senate Bills and Gene-Disease. Thus Table 1 supports the overall architectural trend while keeping the claim model-specific: pairwise reductions are often fragile when large empirical hyperedges carry label signal, while density-aware features help primarily when the backbone leaves structural headroom. Finding 2: density features help only when the backbone leaves headroom. We isolate the effect of density features through a strict ablation in Appendix E.4 comparing D ENS N ET-D to 8

Table 1: Node classification accuracy (%), grouped by architecture tier. Best bold, second underlined. MLP is a feature-only reference. ‡ Published results use architecture-specific recipes and exclude Senate/House; all other methods are re-run with the shared ANCS recipe. Baselines use hidden = 64; D ENS N ET-D uses hidden = 128. Pairwise (k̄<4)

Mixed (4≤k̄<7)

Higher-order (k̄≥7)

Method

Cora Citeseer PubMed Cora-CA House Senate Gene-Dis.

MLP

74.1

72.3

87.1

74.6

70.4

61.2

36.4

CE Tier (clique expansion) HGNN 69.9 68.1 HyperGCN 71.0 68.1

71.2 78.3

76.8 74.3

62.8 63.9

87.6 55.4

85.7 81.3

Native Tier (set-function / message passing) HNHN 72.8 68.3 70.0 UniGNN 77.3 72.4 87.6 AllDeepSets 76.2 70.9 87.2 ED-HNN‡ 80.3 73.7 89.0 SheafHyperGNN‡ 81.3 74.7 87.7

76.6 82.9 80.6 84.0 85.5

68.5 70.9 1 67.8 72.5 73.8

79.6 78.8 92.7 — —

87.8 85.4 86.1 — —

Density-Aware Tier (ours) D ENS N ET-D 77.5

81.2

70.3

92.7

86.0

72.6

87.7

an AllDeepSets backbone under identical training settings. Under this controlled comparison, density features improve performance on 6/7 datasets. The largest gains occur on the two highest-k̄ headroom datasets: Gene-Disease rises from 81.7 ± 1.1% to 86.0% (+4.3 points), and House rises from 68.1 ± 2.4% to 70.3 ± 2.0% (+2.2 points). The only exception is a ceiling case: Senate-Bills changes by only −0.3 points because AllDeepSets already reaches 93.0 ± 3.8%. Thus the densityaware tier is not simply adding capacity; its benefit appears when local aggregation leaves structural signal unused, and it disappears when the backbone already saturates the task. Additional experimental results. Appendix E gives the full experimental protocol, dataset statistics, strict density ablation, and node-classification table with error bars. Appendix P collects negative and ceiling cases, including the Senate-Bills ceiling result, the Gene-Disease strict-ablation gain, and normalization contrast between House and Senate Bills. Appendix Q tests whether deeper native models close the Native-Hard gap, showing that depth alone can worsen performance through oversmoothing. Appendix R reports learned gate trajectories, showing when density features are suppressed, active, or indifferent across regimes. Our experimental findings align with the invariant hierarchy: CE models are restricted to pairwise coordinates, while native and density-aware models access progressively richer pattern classes.

6

Conclusion

We studied when hypergraph neural networks can learn genuine higher-order structure, rather than artifacts of the architectural summaries they compute. Our main conclusion is that hypergraph expressivity is governed by pattern densities: these densities form a universal basis for continuous hypergraph invariants, and generalized hypertree width organizes them into a strict hierarchy. This yields the Width Wall, an architectural ceiling determined by which higher-order patterns a model can count, not by hidden dimension or optimization alone. The hierarchy gives a unified map of existing HGNN classes, pinpoints the information lost by clique expansion, and explains when empirical accuracy gaps should appear. Our experiments support this account: constructed witnesses expose the predicted separations, while real node-classification benchmarks show that density-aware features help most when higher-order structure remains beyond the backbone’s native reach. Extending the theory to n → ∞ and quantifying the statistical and computational cost of crossing the Width Wall remain important directions. 9

Limitations, Ethical Statement, Broader Impact. This work provides a structural account of HGNN expressivity through homomorphism densities and generalized hypertree width, but several limitations remain. Our theoretical results primarily concern continuous permutation-invariant functions on fixed-size hypergraphs, leaving a full treatment of asymptotic regimes, statistical sample complexity, and noisy or evolving hypergraphs for future work. Although density-aware features can help models cross the Width Wall, they may also introduce additional computational cost and require careful selection of informative patterns. Improved hypergraph learning can benefit domains where higher-order relations are essential, including biology, scientific collaboration, and social systems analysis. At the same time, these models may be applied to sensitive relational data, where higher expressivity could amplify privacy risks, encode spurious group-level associations, or support misleading inferences about individuals from group memberships. Responsible deployment therefore requires attention to data provenance, privacy protection, fairness evaluation, and domain-specific validation beyond predictive accuracy.

LLM Usage We used large language models (LLMs) to support the preparation of this manuscript, including assistance with writing, editing, and improving clarity of presentation. LLMs were also used to support experimental workflows, such as drafting code, checking implementation details, and organizing analysis. All intellectual contributions, methodological decisions, experimental design, results interpretation, and final manuscript content were reviewed and verified by the authors, who take full responsibility for the accuracy and integrity of the work.

10

References [1] P.-A. Absil, Robert Mahony, and Rodolphe Sepulchre. Optimization Algorithms on Matrix Manifolds. Princeton University Press, 2008. [2] Ryan P. Adams and Richard S. Zemel. Ranking via sinkhorn propagation. In arXiv preprint arXiv:1106.1925, 2011. [3] Devanshu Arya, Deepak Kumar Gupta, Stevan Rudinac, and Marcel Worring. HyperSAGE: Generalizing inductive representation learning on hypergraphs. In arXiv preprint arXiv:2010.04558, 2020. [4] Chaewoon Bae, Doyun Choi, Jaehyun Lee, and Jaemin Yoo. Parameter-free hypergraph neural network for few-shot node classification. Advances in Neural Information Processing Systems (NeurIPS), 2025. [5] Peter L. Bartlett, Dylan J. Foster, and Matus J. Telgarsky. Spectrally-normalized margin bounds for neural networks. In Advances in Neural Information Processing Systems (NeurIPS), 2017. [6] Federico Battiston, Giulia Cencetti, Iacopo Iacopini, Vito Latora, Maxime Lucas, Alice Patania, Jean-Gabriel Young, and Giovanni Petri. Networks beyond pairwise interactions: Structure and dynamics. Physics Reports, 874:1–92, 2020. [7] Austin R. Benson, David F. Gleich, and Jure Leskovec. Higher-order organization of complex networks. Science, 353(6295):163–166, 2016. [8] Cristian Bodnar, Fabrizio Frasca, Yu Guang Wang, Nina Otter, Guido Montufar, Pietro Liò, and Michael Bronstein. Weisfeiler and Lehman go topological: Message passing simplicial networks. In International Conference on Machine Learning (ICML), 2021. [9] Jan Böker. Color refinement, homomorphisms, and hypergraphs. In Graph-Theoretic Concepts in Computer Science (WG), volume 11789 of LNCS, pages 338–350, 2019. doi: 10.1007/ 978-3-030-30786-8_26. [10] Michael M. Bronstein, Joan Bruna, Taco Cohen, and Petar Veličković. Geometric deep learning: Grids, groups, graphs, geodesics, and gauges. arXiv preprint arXiv:2104.13478, 2021. [11] Jin-Yi Cai, Martin Fürer, and Neil Immerman. An optimal lower bound on the number of variables for graph identification. Combinatorica, 12(4):389–410, 1992. [12] Davi Castro-Silva, Fernando Mário de Oliveira Filho, Lucas Slot, and Frank Vallentin. A recursive theta body for hypergraphs. Combinatorica, 43:909–938, 2023. [13] T.-H. Hubert Chan, Anand Louis, Zhihao Gavin Tang, and Chenzi Zhang. Spectral properties of hypergraph laplacian and approximation algorithms. Journal of the ACM, 65(3), 2018. [14] Eli Chien, Chao Pan, Jianhao Peng, and Olgica Milenkovic. You are AllSet: A multiset function framework for hypergraph neural networks. In International Conference on Learning Representations (ICLR), 2022. [15] Seongjin Choi, Gahee Kim, and Se-Young Yun. Weisfeiler and Lehman go categorical. arXiv preprint arXiv:2602.06787, 2026. [16] Charles J. Colbourn and Jeffrey H. Dinitz. Handbook of Combinatorial Designs. CRC Press, 2nd edition, 2007. [17] Chandler Davis. All convex invariant functions of Hermitian matrices. Archiv der Mathematik, 8:276–278, 1957. [18] Etienne de Klerk and Renata Sotirov. Exploiting group symmetry in semidefinite programming relaxations of the quadratic assignment problem. Mathematical Programming, 122:225–246, 2010. [19] Holger Dell, Martin Grohe, and Gaurav Rattan. Lovász meets Weisfeiler and Leman. In International Colloquium on Automata, Languages, and Programming (ICALP), 2018. 11

[20] Yihe Dong, Will Sawin, and Yoshua Bengio. HNHN: Hypergraph networks with hyperedge neurons. In ICML Workshop on Graph Representation Learning and Beyond (GRL+), 2020. arXiv:2006.12278. [21] Iulia Duta, Giulia Cassarà, Fabrizio Silvestri, and Pietro Liò. Sheaf hypergraph networks. In Advances in Neural Information Processing Systems (NeurIPS), 2023. [22] Gábor Elek and Balázs Szegedy. A measure-theoretic approach to the theory of dense hypergraphs. Advances in Mathematics, 231(3–4):1731–1772, 2012. [23] Yifan Feng, Haoxuan You, Zizhao Zhang, Rongrong Ji, and Yue Gao. Hypergraph neural networks. In AAAI Conference on Artificial Intelligence, 2019. [24] James H. Fowler. Legislative cosponsorship networks in the US House and Senate. Social Networks, 28(4):454–465, 2006. [25] Noah Golowich, Alexander Rakhlin, and Ohad Shamir. Size-independent sample complexity of neural networks. In Proceedings of the 31st Conference on Learning Theory (COLT), 2018. [26] Georg Gottlob, Nicola Leone, and Francesco Scarcello. Hypertree decompositions and tractable queries. Journal of Computer and System Sciences, 64(3):579–627, 2002. [27] Mustafa Hajij, Ghada Zamzmi, Theodore Papamarkou, Nina Miolane, Aldo Guzmán-Sáenz, Karthikeyan Natesan Ramamurthy, et al. Topological deep learning: Going beyond graph data. arXiv preprint arXiv:2206.00606, 2022. [28] Jing Huang and Jie Yang. UniGNN: a unified framework for graph and hypergraph neural networks. In International Joint Conference on Artificial Intelligence (IJCAI), 2021. [29] Nicolas Keriven and Gabriel Peyré. Universal invariant and equivariant graph neural networks. In Advances in Neural Information Processing Systems (NeurIPS), 2019. [30] Jinwoo Kim, Saeyoon Oh, Sungjun Cho, and Seunghoon Hong. Equivariant hypergraph neural networks. In European Conference on Computer Vision (ECCV), 2022. [31] Diederik P. Kingma and Jimmy Ba. Adam: A method for stochastic optimization. In International Conference on Learning Representations (ICLR), 2015. [32] Tamara G. Kolda and Brett W. Bader. Tensor decompositions and applications. SIAM Review, 51(3):455–500, 2009. [33] Adrian S. Lewis. The convex analysis of unitarily invariant matrix functions. Journal of Convex Analysis, 2:173–183, 1995. [34] Pan Li and Olgica Milenkovic. Submodular hypergraphs: p-Laplacians, Cheeger inequalities and spectral clustering. In International Conference on Machine Learning (ICML), 2018. [35] Xiaoyu Li, Guangyu Tang, and Jiaojiao Jiang. Implicit hypergraph neural networks: A stable framework for higher-order relational learning with provable guarantees. arXiv preprint arXiv:2508.09427, 2025. [36] László Lovász. Operations with structures. Acta Mathematica Hungarica, 18(3–4):321–328, 1967. [37] László Lovász. Large Networks and Graph Limits, volume 60 of Colloquium Publications. American Mathematical Society, 2012. [38] László Lovász and Balázs Szegedy. Limits of dense graph sequences. Journal of Combinatorial Theory, Series B, 96(6):933–957, 2006. [39] Haggai Maron, Heli Ben-Hamu, Hadar Serviansky, and Yaron Lipman. Provably powerful graph networks. In Advances in Neural Information Processing Systems (NeurIPS), 2019. [40] Haggai Maron, Heli Ben-Hamu, Nadav Shamir, and Yaron Lipman. Invariant and equivariant graph networks. In International Conference on Learning Representations (ICLR), 2019. 12

[41] Christopher Morris, Martin Ritzert, Matthias Fey, William L. Hamilton, Jan Eric Lenssen, Gaurav Rattan, and Martin Grohe. Weisfeiler and leman go neural: Higher-order graph neural networks. In AAAI Conference on Artificial Intelligence, 2019. [42] Christopher Morris, Yaron Lipman, Haggai Maron, Bastian Rieck, Nils M. Kriege, Martin Grohe, Matthias Fey, and Karsten Borgwardt. Weisfeiler and Leman go machine learning: The story so far. Journal of Machine Learning Research, 24(333):1–59, 2023. [43] Benjamin Scheidt. On homomorphism indistinguishability and hypertree depth. In International Colloquium on Automata, Languages, and Programming (ICALP), 2024. [44] Benjamin Scheidt and Nicole Schweikardt. Counting homomorphisms from hypergraphs of bounded generalised hypertree width: A logical characterisation. In International Symposium on Mathematical Foundations of Computer Science (MFCS), 2023. [45] Richard Sinkhorn. A relationship between arbitrary positive matrices and doubly stochastic matrices. Annals of Mathematical Statistics, 35:876–879, 1964. [46] Bernd Sturmfels. Algorithms in Invariant Theory. Springer, 2nd edition, 2008. [47] Bohan Tang, Zexi Liu, Keyue Jiang, Siheng Chen, and Xiaowen Dong. Training-free message passing for learning on hypergraphs. In International Conference on Learning Representations (ICLR), 2025. [48] Peihao Wang, Shenghao Yang, Yunyu Liu, Zhangyang Wang, and Pan Li. Equivariant hypergraph diffusion neural operators. In International Conference on Learning Representations (ICLR), 2023. [49] Keyulu Xu, Weihua Hu, Jure Leskovec, and Stefanie Jegelka. How powerful are graph neural networks? In International Conference on Learning Representations (ICLR), 2019. [50] Naganand Yadati, Madhav Nimishakavi, Prateek Yadav, Vikram Niber, Anand Louis, and Partha Talukdar. HyperGCN: A new method for training graph convolutional networks on hypergraphs. In Advances in Neural Information Processing Systems (NeurIPS), 2019. [51] Manzil Zaheer, Satwik Kottur, Siamak Ravanbakhsh, Barnabás Póczos, Ruslan Salakhutdinov, and Alexander J. Smola. Deep sets. In Advances in Neural Information Processing Systems (NeurIPS), 2017. [52] Detian Zhang, Chengqiang Zhang, Yanghui Rao, Qing Li, and Chunjiang Zhu. Improved expressivity of hypergraph neural networks through high-dimensional generalized WeisfeilerLeman algorithms. In International Conference on Machine Learning (ICML), volume 267 of PMLR, pages 76880–76908, 2025. [53] Yufei Zhao. Hypergraph limits: A regularity approach. Random Structures & Algorithms, 47 (2):205–226, 2015.

13

A

Proofs

This appendix collects the proofs of all results stated in the main text. The algebra-structure statement and completeness theorem (§2.3) are proved in §A.A.1–A.2; the strict hierarchy (Theorem 3.6) and its infinite refinement are in §A.4; InvNet’s expressivity claim is in §A.5. A second, algebraically independent proof of the density theorem appears in Appendix B. A.1

Background on Algebra Structure of Pattern Densities

Theorem A.1 (Algebra structure [cf. 37, Ch. 5]). Let Fk denote the set of all finite k-uniform hypergraphs. The collection of functions {t(F, ·) : F ∈ Fk } generates a unital algebra Ak of k,n continuous Sym(n)-invariant functions on [0, 1]S under pointwise multiplication: t(F1 , A) · t(F2 , A) = t(F1 ⊔ F2 , A),

(4)

where F1 ⊔ F2 denotes the disjoint union. The algebra contains constants: t(∅, A) = 1 for the empty hypergraph. A.2

Background on Homomorphism Completeness

Theorem A.2 (Homomorphism completeness [cf. 37, Ch. 5]). Let H1 , H2 be k-uniform hypergraphs on [n]. Then H1 ∼ = H2 if and only if hom(F, H1 ) = hom(F, H2 ) for every k-uniform hypergraph F. Corollary A.3 (Pattern alignment completeness). Let H1 , H2 be non-isomorphic k-uniform hypergraphs with adjacency tensors A1 , A2 ∈ Sk,n . Then there exists P ∈ Sk,n such that ΘP (A1 ) ̸= ΘP (A2 ). Proof. By Theorem A.2, some pattern F satisfies hom(F, H1 ) ̸= hom(F, H2 ); by Möbius inversion, some F ′ has inj(F ′ , H1 ) ̸= inj(F ′ , H2 ). Taking P = AF ′ , the score ΘP differs between A1 and A2 . A.3

Proof of Proposition 3.2

Proof. A clique-expansion HGNN first maps A 7→ ϕ(A) (the clique-expansion adjacency matrix) and then applies a GNN to ϕ(A). Since the GNN’s output depends only on ϕ(A), the composite function fθ (A) = gθ (ϕ(A)) factors through the clique expansion. If ϕ(A1 ) = ϕ(A2 ), then fθ (A1 ) = fθ (A2 ) regardless of θ. A.4

Proof of Theorem 3.6 (R)

Proof. The proof has three parts: (i) ICE ⊆ INH , (ii) the first inclusion is strict, and (iii) the second inclusion is strict. (R)

Part (i): ICE ⊆ INH . This inclusion holds by definition of the native-incidence reference class: (R) INH contains all CE invariants and augments them with bounded-width hypergraph motif counts. (R)

(R)

Part (ii): ICE ⊊ INH . We exhibit an invariant in INH \ ICE . Consider the Steiner pair (H1 , H2 ) from Appendix F.1. Both have identical clique expansions ϕ(A1 ) = ϕ(A2 ) = J − I (the adjacency matrix of K13 ), so any function in ICE assigns the same value to both. However, the two systems have different Pasch counts (13 versus 8; Appendix F.1). Let PPasch be the 4-edge Pasch configuration. (R) Then the polynomial invariant t(PPasch , ·) separates H1 and H2 , and belongs to INH for every R ≥ 3 because the Pasch pattern has generalized hypertree width Wgh = 3 for the representation used here. This proves the strict inclusion for the native function class; it does not assert that every shallow native message-passing implementation computes the Pasch count. (R)

Part (iii): INH ⊊ Iall . Native HGNNs using message passing are bounded by k-GWL for some finite k depending on the depth and architecture [52]. By the homomorphism-logic correspondence [44], k-GWL-equivalent hypergraphs agree on all pattern densities t(F, ·) with Wgh (F ) ≤ k. The Cai–Furer–Immerman (CFI) construction (Appendix F.2) produces pairs (H1CFI , H2CFI ) that are 14

non-isomorphic but indistinguishable by k-GWL for any fixed k (by choosing base graphs with girth greater than 2k + 1). By the completeness theorem (Theorem A.2), some pattern density separates (R) these pairs, yielding an invariant in Iall \ INH . Corollary A.4 (Architectural impossibility). For any r ≥ 1, there exist continuous hypergraph invariants f ∈ C(Sk,n /Sym(n)) such that no architecture limited to patterns of Wgh ≤ r can approximate f to arbitrary precision, regardless of depth, width, or training procedure. Concretely, at each level r, the CFI construction provides explicit functions that witness the gap Hr ⊊ Hr+1 . Remark A.5 (Quantitative vertex bounds). The CFI separation at level r applies whenever n ≥ (k) NCFI (r). This threshold is finite because the construction only requires a base graph of girth > 2r + 1; for the 3-uniform lifts used here, this gives n ≥ Ω(r2 ) vertices [11]. Concretely: the CE-vs. Native separation above requires n ≥ 13 (the Steiner pair on 13 vertices, the smallest v with non-isomorphic STS(v) pairs; larger n follow by adding matched isolated vertices); the CFI refinement gives bounds such as r = 2 requiring n ≥ 15 and r = 3 requiring n ≥ 30 under the stated lift. For practical hypergraph datasets with n ≥ 100, all levels r ≤ 7 of the hierarchy are active. A.5

Proof of Theorem 4.1

We first establish that the Sinkhorn relaxation preserves separation on the Steiner pair. Lemma A.6 (Sinkhorn separation on STS pairs). Let A1 , A2 be the adjacency tensors of two nonisomorphic STS(v) systems. Then the Sinkhorn relaxation satisfies Θ̂SK A1 (A1 ) = v(v − 1)/3 · k! and SK SK SK Θ̂A1 (A2 ) < Θ̂A1 (A1 ). In particular, for v = 13: Θ̂A1 (A1 ) = 156 and Θ̂SK A1 (A2 ) < 156. Proof. The identity permutation achievesQ⟨A1P , A1 ⟩ = 26 · 3! = 156. No D ∈ B13 can exceed 156: each summand satisfies (D⊗3 · A)ijk ≤ ℓ ( jℓ Diℓ jℓ ) = 1 for nonnegative doubly stochastic D, so the sum of 156 such terms is at most 156. Hence Θ̂SK A1 (A1 ) = 156. For separation, suppose for contradiction that some D∗ ∈ B13 achieves ⟨A1 , (D∗ )⊗3 · A2 , =⟩156. Since the 156 summands are each at most 1 and sum to 156, every summand equals 1: for every triple ∗ ∗ ∗ {i, j, k} of STS1 and all a, b, c with Dia , Djb , Dkc > 0, we must have {a, b, c} ∈ STS2 . If some ∗ ∗ ∗ row of D has two positive entries Dia1 , Dia2 > 0 (a1 ̸= a2 ), choose any triple {i, j, k} ∈ STS1 and any b, c in the support of rows j, k; then both {a1 , b, c} and {a2 , b, c} must be triples in STS2 , but each pair {b, c} lies in exactly one triple, so a1 = a2 —contradiction. Hence D∗ is a permutation matrix mapping STS2 to STS1 , contradicting non-isomorphism. In summary, the argument shows that any maximizer in Bv must have singleton row support, and is therefore a vertex of the Birkhoff polytope (a permutation matrix). Proof of Theorem 4.1. We prove each claim in order. Claim 1: ÎInvNet ⊇ INH . I NV N ET contains a full message-passing local branch (Appendix G.2) that subsumes any native HGNN architecture. The MLP output head can ignore the pattern-alignment branch features (by setting their corresponding weights to zero), so I NV N ET can implement any function computable by a native HGNN. Claim 2: ÎInvNet ⊋ INH . Consider the Steiner pair from Appendix F.1. The pattern alignment score ΘA1 (·) separates this pair, and by Lemma A.6, the Sinkhorn relaxation preserves this separation for any finite τ > 0. Meanwhile, any message-passing HGNN bounded by 1-GWL assigns identical colorings to both Steiner triple systems (since every node has degree 6 and every pair co-occurs exactly once). Hence ΘA1 ∈ ÎInvNet \ INH . Claim 3: ÎInvNet ⊆ Iall . The pattern-alignment branch computes the vector zinv = [Θ̂SK P1 (A) − SK SK ⊗k α1 , . . . , Θ̂PJ (A) − αJ ]. Each Θ̂Pj (A) = maxD∈Bn ⟨Pj , D · A⟩ is a continuous function of A (the maximum of a family of continuous functions parameterized over the compact set Bn ). ⊗k For invariance: Θ̂SK · (Π · A)⟩ = maxD ⟨Pj , (DΠ)⊗k · A⟩. Since Pj (Π · A) = maxD∈Bn ⟨Pj , D D 7→ DΠ is a bijection on Bn (the product of a doubly stochastic matrix and a permutation matrix is doubly stochastic), this equals Θ̂SK Pj (A). Thus each component of zinv is a continuous invariant, and 15

compositions of continuous invariants with continuous functions (here MLPout ) remain continuous invariants, so ÎInvNet ⊆ Iall .

B

Orbit Separation via Symmetric Tensor Invariant Theory

This appendix provides a second, completely independent proof that hypergraph invariant functions are generated by pattern-counting features, using classical algebraic invariant theory. This algebraic route—Hilbert’s finiteness theorem and the Reynolds operator—complements the combinatorial (Lovász–Stone–Weierstrass) approach of Section 2 and provides additional structural insight, particularly for the k ≥ 3 case where matrix invariant theory no longer suffices. B.1

The Invariant Ring of Symmetric Tensors

Let Sk,n denote the space of real symmetric order-k tensors on Rn , viewed as a finite-dimensional real vector space. The symmetric group Sym(n) acts linearly on Sk,n by (Π · A)i1 ,...,ik = AΠ−1 (i1 ),...,Π−1 (ik ) . We work over R throughout; since Sym(n) is a finite group, all invarianttheoretic results (finite generation, orbit separation) hold over R without requiring algebraic closure [46]. We consider the ring of Sym(n)-invariant polynomial functions: R[Sk,n ]Sym(n) = {p ∈ R[Sk,n ] : p(Π · A) = p(A) for all Π ∈ Sym(n), A ∈ Sk,n }. Theorem B.1 (Finite generation — Hilbert). The invariant ring R[Sk,n ]Sym(n) is finitely generated as an R-algebra. That is, there exist finitely many invariant polynomials p1 , . . . , pN such that every invariant polynomial is a polynomial combination of p1 , . . . , pN . Proof. Sym(n) is a finite group acting linearly on the finite-dimensional space Sk,n . By Hilbert’s finiteness theorem for reductive groups [46], the invariant ring is finitely generated. More directly, since Sym(n) is finite, Noether’s degree bound guarantees that generators can be chosen with degree ≤ |Sym(n)| = n!. B.2

The Reynolds Operator and Explicit Generators

The Reynolds operator R : R[Sk,n ] → R[Sk,n ]Sym(n) is defined by averaging over the group: X 1 R(p)(A) = p(Π · A). n!

(5)

Π∈Sym(n)

This is a surjective projection: R(p) = p if and only if p is already invariant, and R(p) is invariant for any p. Tensor contraction traces as generators. For any multigraph Γ on vertex set [r] with each vertex labeled by a k-tensor index, define the tensor contraction trace: X Y TΓ (A) = Aiv,1 ,...,iv,k , (6) i1 ,...,ir ∈[n] v∈V (Γ) indices contracted per Γ

where the contraction pattern is specified by Γ: indices sharing an edge in Γ are identified. Each TΓ is manifestly invariant under Sym(n) (it is a sum over all index assignments, which is permutationinvariant). Proposition B.2 (Generators are pattern counts). Every tensor contraction trace TΓ is a linear combination of homomorphism counts hom(F, H) for appropriate k-uniform patterns F , and vice versa. P Q Proof. The contraction trace TΓ (A) = i1 ,...,ir v Aiv,1 ,...,iv,k sums over all maps ϕ : V (Γ) → [n] (with ϕ possibly non-injective) weighted by the product of tensor entries. The dual hypergraph FΓ is defined as follows: the vertex set of FΓ is the set of distinct summation indices in TΓ , and each vertex v ∈ V (Γ) (representing one tensor factor Aiv,1 ,...,iv,k ) yields a k-hyperedge {iv,1 , . . . , iv,k } in FΓ ; two tensor factors sharing a contracted index become hyperedges sharing a vertex. With this identification, TΓ (A) = hom(FΓ , H) exactly. Conversely, every hom(F, H) can be written as a tensor contraction trace by choosing Γ as the dual of F (vertices of Γ correspond to hyperedges of F , with edges in Γ encoding shared vertices). 16

B.3

Algebraic Completeness: Orbits Separated by Polynomial Invariants

Theorem B.3 (Orbit separation). Let A1 , A2 ∈ Sk,n lie in distinct Sym(n)-orbits (i.e., H1 ∼ ̸ H2 ). = Then there exists an invariant polynomial p ∈ R[Sk,n ]Sym(n) with p(A1 ) ̸= p(A2 ). Proof. Since Sym(n) is a finite group, every orbit O(A) = {Π · A : Π ∈ Sym(n)} is a finite set, hence Zariski-closed (and Euclidean-closed; for finite sets in Rd , Zariski and Euclidean closure coincide). Two distinct orbits O(A1 ) and O(A2 ) are disjoint Zariski-closed sets. Since R[Sk,n ] separates points (as Sk,n ∼ = Rd ), there exists a polynomial q with q(A1 ) ̸= q(A′2 ) for all A′2 ∈ O(A2 ). Applying the Reynolds operator: R(q)(A1 ) = q(A1 ) (since q is already correct at A1 after averaging), and R(q)(A2 ) is the average of q over O(A2 ). More precisely, since O(A1 ) and O(A2 ) are finite disjoint sets in Rd , there exists a polynomial q taking the value 1 on O(A1 ) and 0 on O(A2 ) (by Lagrange interpolation in high enough degree). Then p = R(q) satisfies p(A1 ) = 1 and p(A2 ) = 0, and p is invariant. B.4

From Algebraic Completeness to the Density Theorem

Combining the results above: Corollary B.4 (Algebraic density). The invariant ring R[Sk,n ]Sym(n) is dense in C(K/Sym(n)) for any compact K ⊆ Sk,n . Moreover, this ring is generated by tensor contraction traces (equivalently, by homomorphism counts). Proof. By Theorem B.3, invariant polynomials separate orbits. By Theorem B.1, the invariant ring is finitely generated. Polynomial invariants form a separating subalgebra of C(K/Sym(n)): they separate orbits (Theorem B.3), contain constants (the polynomial 1), and are closed under multiplication and addition. The passage from polynomial to continuous density then follows directly from the Stone–Weierstrass theorem: any separating subalgebra of C(X) (for compact Hausdorff X) containing constants is dense. By Proposition B.2, the generators are homomorphism counts. This provides a second, algebraic proof of Theorem 2.4, arrived at through a completely different route: Hilbert’s finiteness theorem and the Reynolds operator rather than Lovász’s homomorphism completeness. The two proofs reinforce each other: the combinatorial proof (main text) gives the explicit connection to homomorphism counting and hypertree width, while the algebraic proof gives finite generation and the connection to the invariant ring structure. B.5

Computational Hardness of Pattern Alignment

Theorem B.5 (Computational hardness gap). Computing the pattern alignment score ΘP (A) = maxΠ∈Sym(n) ⟨P, Π · A⟩ is equivalent to the multidimensional assignment problem (MAP) for k ≥ 3, which is NP-hard and APX-hard (no polynomial-time approximation scheme exists unless P = NP). For k = 2, the problem reduces to the quadratic assignment problem (QAP): ΘP (A) =

max ⟨P, ΠAΠ⊤ ⟩, Π∈Sym(n)

which is also NP-hard but admits polynomial-time spectral relaxations with the closed-form solution ΛP (A) = λ(P )⊤ λ(A) [17, 33]. For k ≥ 3, no such closed-form spectral relaxation exists: the natural relaxation over O(n) requires solving a non-convex optimization on the Stiefel manifold. This hardness gap is what motivates the Sinkhorn and Stiefel relaxation layers in I NV N ET (§G.1). B.6

The Hilbert Series and Approximation Dimension

The Molien series counts the dimension of invariant polynomial spaces by degree: X X 1 1 Sym(n)
d M (t) = dim R[Sk,n ]d t = , n! det(I − t ρ(Π)) d≥0

Π∈Sym(n)

17

where ρ(Π) denotes the matrix of the linear action of Π on Sk,n . For k = 2 (symmetric matrices), the invariant polynomials of degree d are spanned by the power-sum symmetric functions of the eigenvalues, giving the familiar Molien series with generating function related to the partition function. For k ≥ 3, the Molien series is computable but lacks a simple closed form, reflecting the richer structure of tensor invariants. The Molien series provides a quantitative bridge between the qualitative universality of Theorem 2.4 and concrete approximation complexity. Proposition B.6 (Approximation space dimension). Let I≤d denote the space of Sym(n)-invariant polynomials of degree ≤ d in the entries of A ∈ Sk,n . Pd (a) dim(I≤d ) = i=0 [ti ] M (t), where [ti ]M (t) is the i-th Molien coefficient. (b) Every basis element of Id corresponds to a homomorphism count hom(F, ·) for a pattern F with |V (F )| ≤ d (by Proposition B.2). (c) Forp k = 2 (graphs): √ dim(Id ) = p(d), the number of integer partitions of d, giving dim(I≤d ) ∼ exp(π 2d/3)/(4d 3) by the Hardy–Ramanujan formula. (d) For the ghtw-restricted subspace I≤d,w (patterns with ghtw(F ) ≤ w): dim(I≤d,w ) ≤ {F : |V (F )| ≤ d, ghtw(F ) ≤ w} /∼iso , which is strictly smaller than dim(I≤d ) for d sufficiently large. This connects the Molien series to approximation power: to ϵ-approximate a target invariant f by a degree-d polynomial in pattern densities, one needs at most dim(I≤d ) basis functions. The ghtw restriction (part d) quantifies the cost of architectural limitations: bounded-depth message-passing models operate in the strictly smaller space I≤d,w , and the dimension gap grows with d.

C

Quantitative Analysis

The universality theorem (Theorem 2.4) is qualitative: it guarantees approximation but provides no convergence rates. This appendix develops three quantitative complements: a generalization bound for I NV N ET via Rademacher complexity (Appendix C.1), a concentration inequality for estimating pattern densities from finite samples (Appendix C.2), and the Molien series connection to approximation dimension. C.1

Generalization Bound for InvNet

We analyze the sample complexity of learning with I NV N ET by bounding the Rademacher complexity of its hypothesis class. Definition C.1 (InvNet hypothesis class). Let FJ,D denote the class of functions f : Sk,n → R computed by an I NV N ET with: • J Sinkhorn template layers with templates satisfying ∥Pj ∥F ≤ BP and temperature τ > 0; • a D-layer MLP readout with spectral norm ≤ BW per layer and width ≤ w; applied to inputs satisfying ∥A∥F ≤ BA . The key technical ingredient is the Lipschitz continuity of the Sinkhorn relaxation: Lemma C.2 (Lipschitz continuity of Sinkhorn scores). The entropy-regularized pattern alignment score sτj (A) = maxD∈Bn [⟨Pj , D⊗k · A⟩ + τ H(D)] satisfies: (i) Boundedness: |sτj (A)| ≤ BP BA + τ log n. (ii) Lipschitz continuity: |sτj (A) − sτj (A′ )| ≤ BP ∥A − A′ ∥F . Proof. (i) By Cauchy–Schwarz, |⟨Pj , D⊗k · A⟩| ≤ ∥Pj ∥F ∥D⊗k · A∥F ≤ BP BA , where the second inequality uses that D⊗k is a contraction in Frobenius norm for doubly stochastic D (since 18

∥D∥op ≤ 1 and the tensor product preserves spectral norm bounds). The entropy term satisfies 0 ≤ H(D) ≤ log n. (ii) For any fixed D ∈ Bn , the map A 7→ ⟨Pj , D⊗k · A⟩ is linear with operator norm ≤ BP . The pointwise supremum of BP -Lipschitz functions (parameterized by D ∈ Bn ) is itself BP -Lipschitz: |sτj (A) − sτj (A′ )| ≤ sup |⟨Pj , D⊗k · (A − A′ )⟩| ≤ BP ∥A − A′ ∥F . D∈Bn

Proposition C.3 (Generalization bound for InvNet). Let ℓ : R × R → [0, C] be a ρ-Lipschitz loss function. For N i.i.d. training samples {(Ai , yi )}N i=1 , with probability at least 1 − δ, every f ∈ FJ,D satisfies: r p D 2ρ BP BA · BW 2D log(2w) log(2/δ) √ R(f ) ≤ R̂(f ) + , + 3C 2N N P where R(f ) = E[ℓ(f (A), y)] is the population risk and R̂(f ) = N1 i ℓ(f (Ai ), yi ) is the empirical risk. Proof sketch. By Lemma C.2, the Sinkhorn feature map A 7→ (sτ1 (A), . . . , sτJ (A)) is BP -Lipschitz with outputs bounded by BP BA + τ log n. By the spectral-norm Rademacher bound for deep networks [5, 25], p the MLP readout √applied to bounded inputs has empirical Rademacher complexity D R̂N ≤ Bin · BW 2D log(2w)/ N , where Bin ≤ BP BA (absorbing the local branch bound and logarithmic entropy terms into constants). Composing with the ρ-Lipschitz loss and applying the standard Rademacher-to-generalization conversion [5] yields the stated bound. Remark C.4 (Architectural implications). The bound makes the dependence on I NV N ET’s design parameters explicit: • Number of templates J: enters only logarithmically (through the input dimension to the MLP). Adding templates is “cheap” in sample complexity. • Template norm BP : enters linearly. Large templates require proportionally more data. D • MLP depth D: enters exponentially through BW . Shallow readouts generalize better, consistent with InvNet’s expressivity coming from the Sinkhorn branch rather than the MLP depth.

• Temperature τ : does not appear in the leading term. The Sinkhorn relaxation incurs no additional sample complexity cost relative to the hard alignment (though it affects approximation quality via the smoothing gap). C.2

Concentration of Pattern Density Estimators

For a fixed pattern F , the homomorphism density can be estimated from random vertex-maps. Proposition C.5 (Pattern density concentration). Let F be a k-uniform pattern with |V (F )| = v, and let H be a k-uniform hypergraph on n nodes with adjacency tensor A. The normalized density t(F, H) = hom(F, H)/nv can be estimated by drawing N independent uniformly random maps ϕi : V (F ) → [n] and computing: t̂(F, H) =

N 1 X Y Aϕi (e) . N i=1 e∈E(F )

Then E[t̂] = t(F, H) and Pr[|t̂ − t(F, H)| ≥ ϵ] ≤ 2 exp(−2N ϵ2 ). Q k Proof. Each sample Xi = e∈E(F ) Aϕi (e) ∈ [0, 1] for unweighted hypergraphs (A ∈ {0, 1}n ). P Q By linearity, E[Xi ] = n−v ϕ:V (F )→[n] e∈E(F ) Aϕ(e) = t(F, H). The result follows from Hoeffding’s inequality applied to the i.i.d. bounded random variables X1 , . . . , XN . To estimate any single pattern density to additive accuracy ϵ with confidence 1 − δ requires N = O(log(1/δ)/ϵ2 ) random samples—independent of both n and |V (F )|. However, the variance grows 19

with sparsity: for rare patterns (t(F, H) ≪ 1), multiplicative estimation requires N = Ω(1/t(F, H)) samples, which can be exponential in v for sparse hypergraphs. Combined with Proposition B.6, this yields a concrete (if loose) sample complexity pipeline: (i) choose target accuracy ϵ; (ii) determine the polynomial degree d(ϵ) via Jackson-type bounds (open for k ≥ 3; for graphs, d = O(1/ϵ2 ) suffices for cut-metric approximation [37]); (iii) enumerate the dim(I≤d ) basis patterns; (iv) estimate each density to accuracy ϵ/ dim(I≤d ) using O(dim(I≤d )2 /ϵ2 ) samples. C.3

Non-Uniform Density

Corollary C.6 (Non-uniform density). For non-uniform hypergraphs with maximum arity K and LK j,n direct-sum adjacency A = (A(2) , . . . , A(K) ) ∈ , the algebra generated by pattern j=2 S densities {t(F, ·) : F is j-uniform, 2 ≤ j ≤ K} together with cross-order products t(F1 , ·) · t(F2 , ·) (where F1 , F2 have different arities) is dense in the continuous Sym(n)-invariant functions on L Sj,n . j [0, 1] LK j,n Proof. The product space j=2 [0, 1]S is compact Hausdorff under the product topology, and Sym(n) acts diagonally. For each arity j, the j-uniform pattern densities separate j-uniform components (Theorem A.2). Cross-order products separate configurations that agree per-order but differ in joint structure. The combined algebra separates all orbits, contains constants, and is closed under multiplication; Stone–Weierstrass applies.

D

Architecture Classification of Modern HGNNs

Four-way correspondence diagram. The following diagram summarizes the structural correspondence governing the entire hierarchy: Hypertree width r

Pattern densities ghtw(F ) ≤ r

[44]

Thm. 2.4

Counting logic C r+1

HGNN class Hr

Equivariant higher-order networks. Equivariant higher-order networks [39, 29] process k-tensor features directly at O(nk ) cost, achieving maximal expressivity in principle. These are included in our hierarchy at the Iall level; their impractical scaling motivates the bounded-width alternatives studied in the main text. Example D.1 (Worked instance: L = 2, k = 3). A 2-layer native HGNN on 3-uniform hypergraphs uses c3 · 2 + 1 = 3 variables (with c3 = 1 for standard message passing). We count all patterns F with ghtw(F ) ≤ 2. Concretely: (i) triangle patterns (three vertices, one hyperedge; ghtw = 1); (ii) diamond patterns (four vertices, two hyperedges sharing a (k−1)-face; ghtw = 2, since two hyperedges can be covered by bags of size 4 in a star decomposition). However, the complete
(3) (3) 3-uniform hypergraph K5 (five vertices, 53 = 10 hyperedges) has ghtw(K5 ) = 3, requiring 4 variables and hence L ≥ 3 layers. Thus H2 detects pairwise overlap of hyperedges (diamond patterns) but cannot resolve the global arrangement of five mutually overlapping hyperedges. Remark D.2 (Two faces of the same structure). Pattern densities and the WL hierarchy are two views of the same object: WL captures structural indistinguishability (a decision problem), while homomorphism densities capture structural quantification (a function approximation problem). The density theorem shows that the latter is complete—the algebra of pattern densities generates all continuous invariants—whereas WL with any fixed number of variables is not. D.1

Extended Classification Table

The table below summarizes HGNN architectures. 20

Table 2: Extended classification of modern HGNN architectures. ✓ = maps cleanly under (C1)–(C5), △ = partial, — = out of scope. Model

Operator

Class / Clean

Key consideration

Clique-expansion / graph reduction −1/2

−1/2

Dv HW De−1 H ⊤ Dv ; factors through clique expansion by Prop. 3.2 Mediator selection via Fiedler vector is a global spectral heuristic; learned component is CE, but preprocessing may inject information beyond ICE

HGNN [23]

Star expansion

ICE / ✓

HyperGCN [50]

Mediator GCN

ICE / △

Native incidence message passing HNHN [20] Nonlinear n↔e

INH / ✓

UniGNN [28]

Unified MP

INH / ✓

AllDeepSets [14]

DeepSets/edge

INH / ✓

AllSetTransformer [14]

Attention/edge

INH / ✓

HyperSAGE [3]

SAGE-style MP

INH / ✓

Diffusion / operator learning Equiv. diffusion ED-HNN [48]

INH / ✓

Per-step operator is local (hypergraph Laplacian); continuous-time formulation ≈ deep native MP

Directly implements the ghtw hierarchy at level k; validates our framework Operates on sheaf data (restriction maps); richer than bare hypergraph ⇒ out of scope Higher-order topological structure beyond hypergraphs

Degree-weighted bipartite MP; satisfies (C1)– (C5) GCN/GAT/GIN variants on incidence structure Symmetric multiset aggregation within each hyperedge Attention is an equivariant multiset function; does not increase accessible pattern width Inductive neighborhood sampling on incidence structure; satisfies (C1)–(C5)

Higher-order / topological k-GWL HNNs [52]

k-dim WL

Hk / ✓

SheafHNN [21]

Sheaf Laplacian

—/—

TopoNN [27]

Cell complex

—/—

Training-free MP Closed-form linear

INH / △

IHGNN [35]

Implicit equilib.

INH / ✓

InvNet (ours)

Pattern align.

ÎInvNet / ✓

Preprocessing-heavy TF-HNN [47] ZEN [4]

D.2

INH / △

Expressivity lives in multi-hop preprocessing; effective depth ̸= learned depth Parameter-free; closed-form weight matrix from linearized incidence MP; expressivity bounded by propagation depth Fixed-point of incidence MP; infinite effective depth mitigates oversmoothing; satisfies (C1)– (C5) Sinkhorn relaxation accesses patterns beyond bounded ghtw (Thm. 4.1)

Conditions for the Architecture–Width Correspondence

The correspondence between message-passing depth and pattern complexity holds under the following conditions: (C1) Permutation equivariance. The model’s node-level computations commute with permutations of the node set. No node identifiers, positional indices, or order-dependent operations are used. 21

(C2) Locality. Each layer aggregates information only from the incident neighborhood (nodes sharing a hyperedge). Global all-pairs attention or full-graph pooling within a single layer can effectively increase the accessible pattern width beyond ck L. (C3) Symmetric aggregation. Neighborhood aggregation functions are permutation-invariant over the neighbor multiset. This includes sum, mean, max, DeepSets, and attention-weighted aggregation (since attention weights are computed from the multiset, not from any ordering). (C4) No symmetry-breaking features. The model does not use Laplacian eigenvectors, random node features, or other global positional encodings that can distinguish structurally identical nodes (see Remark D.3 below). (C5) Deterministic structure-dependent preprocessing. If the architecture includes preprocessing steps (e.g., HyperGCN’s mediator selection), these must be deterministic functions of the hypergraph structure. Stochastic or learned preprocessing requires separate analysis. All “✓” entries satisfy conditions (C1)–(C5). The “△” entries (HyperGCN, TF-HNN) satisfy (C1)– (C3) but require care with (C4)–(C5): HyperGCN’s mediator selection is a spectral heuristic that may inject global information beyond ICE , and TF-HNN’s expressivity depends on the preprocessing depth rather than the learned model depth. Remark D.3 (Positional encodings and the hierarchy). Positional encodings (PEs) such as random node features [42], Laplacian eigenvectors, or random walk profiles augment the input with globallycomputed node identifiers that provide information beyond the hypergraph structure. With sufficiently informative PEs, architectures can in principle distinguish arbitrary non-isomorphic hypergraphs— analogous to the graph setting, where random features elevate 1-WL to maximal distinguishing power without increasing the WL dimension [42]. Our hierarchy isolates a complementary question: what can the architecture compute from structure alone? This is the relevant measure in three common scenarios: (i) when canonical node identifiers are unavailable or unreliable (e.g., molecular or biological hypergraphs without fixed atom orderings), (ii) when one seeks guarantees about structural pattern detection that hold regardless of how PEs are constructed, and (iii) when designing PEs themselves—understanding which structural invariants the base architecture already captures reveals which PEs provide genuinely new information rather than redundant signal. The pattern-complexity grading of the hierarchy is orthogonal to the information injected by PEs: PEs break the symmetry assumption (C4) rather than enabling the model to compute homomorphism counts of higher-width patterns.

E

Additional Experimental Details

This appendix provides the dataset statistics, training details, and ablation studies supporting §5 of the main text. All non-‡ results in Table 1 are reproduced under the identical recipe described here; the density ablation (Table 4) isolates the contribution of the density branch. E.1

Dataset Statistics

Table 3: Dataset statistics for the seven node classification benchmarks in Table 1. Expanded statistics including edge-size percentiles and regime classification are in Table 9. Dataset

Nodes

Hyperedges

Features

Classes

Avg. HE size

Cora Citeseer Pubmed Cora-CA House Senate Bills Gene-Disease

2,708 3,312 19,717 2,708 1,290 294 5,012

1,579 1,079 7,963 1,072 340 29,157 2,009

1,433 3,703 500 1,433 100 100 100

7 6 3 7 2 2 21

3.0 3.2 4.3 4.3 34.9 8.0 14.0

For House, Senate Bills, and Gene-Disease, which lack published node features, we follow the AllSet protocol [14]: one-hot label vectors concatenated with N (0, 1) noise to produce 100-dimensional features. 22

E.2

Separation Hypergraph Construction Details

Steiner pair verification. We verify computationally that the two STS(13) systems in Appendix F.1 have identical clique expansions: for both H1 and H2 , the co-occurrence matrix Cij = |{e ∈ E : i ∈ e, j ∈ e}| equals 1 for all i ̸= j, confirming ϕ(A1 ) = ϕ(A2 ) = J − I (= K13 ). Non-isomorphism is verified via Pasch count: H1 contains 13 Pasch configurations while H2 contains 8. Since Pasch count is an isomorphism invariant, this provides a definitive proof of non-isomorphism without requiring automorphism group computation. CFI pair construction. For the CFI construction, we use as base graphs the incidence graphs of projective planes PG(2, q) for q ∈ {3, 4, 5}, which are (q + 1)-regular bipartite graphs with girth 6. Each gadget replaces an edge with 4 auxiliary hyperedges on 2 fresh nodes. The resulting hypergraphs have n = |VG | + 2|EG | nodes and 4|EG | hyperedges. E.3

Training Details

Optimization. All models are trained with the Adam optimizer [31]. For IWS experiments (Exp. 1, 2): learning rate 10−3 , weight decay 5 × 10−4 (applied only to weight matrices, excluding biases and normalization parameters), batch size 32, trained for 300 epochs with early stopping (patience 50). For node classification (Exp. 3): learning rate 10−3 , weight decay 0, cosine-annealing schedule, label smoothing 0.1, 200 epochs (400 for Senate-Bills) with patience 50 (150 for Senate-Bills). Baselines use hidden dimension 64, dropout 0.2 (AllSet protocol); D ENS N ET-D uses hidden 128, dropout 0.5 (see Appendix N for full details). Weight decay is applied via parameter group splitting: normalization layer parameters (γ, β) and all bias terms are always excluded from weight decay regularization. I NV N ET-specific hyperparameters. Number of templates J ∈ {4, 8, 16}; Sinkhorn iterations Titer ∈ {5, 10, 20}; temperature τ annealed from 1.0 to 0.01 over training. Templates Pj are initialized with orthogonal initialization and identity bias. Compute. All experiments are run on NVIDIA H100 (80GB) GPUs. IWS experiments complete in < 1 hour each. Node classification uses 8 GPUs in parallel with per-job isolation (no DDP); full sweep completes in < 12 hours. E.4

Density Ablation

Table 4: Strict density ablation: AllDeepSets with identical training recipe (hidden = 128, dropout = 0.5, label smoothing = 0.1, cosine schedule) vs. D ENS N ET-D. The only difference is the density branch. Density features add +0.6 to +4.3 points on 6/7 datasets; the ceiling dataset shows ≤ 0.3 change. Dataset

k̄

AllDeepSets

DensNet-D

Density lift

Senate-Bills Cora Citeseer PubMed Cora-CA House Gene-Disease

8 3 3 4.3 4.3 35 14

93.0±3.8 76.1±1.7 70.8±1.7 87.1±0.6 80.4±1.2 68.1±2.4 81.7±1.1

92.7±4.1 77.5±1.2 72.6±0.5 87.7 81.2 70.3±2.0 86.0

−0.3 (ceiling) +1.4 +1.8 +0.6 +0.8 +2.2 +4.3

E.5

Node Classification Results with Error Bars

E.6

Additional Results

We provide ablation studies on: (i) the number of templates J, (ii) the Sinkhorn temperature τ , (iii) the number of Sinkhorn iterations Titer , and (iv) the relative contribution of the local vs. invariant branches. Full results with confidence intervals are provided in the supplementary material. 23

Density lift (%)

Gene-Dis. 5 Headroom 4 House Ceiling 3 Citeseer Cora 2 Cora-CA 1 PubMed 0 −0.5 3 4 8 Senate14 35 Mean hyperedge size k̄

Figure 1: Density lift vs. mean hyperedge size k̄ (strict ablation: identical recipe, density branch on/off). Filled: headroom datasets; open: ceiling. The density branch provides the largest lift on Gene-Disease (+4.3%) and House (+2.2%), the two highest-k̄ datasets, and is correctly suppressed on the Senate ceiling case. Table 5: Node classification accuracy (%) with standard deviations over 5–10 seeds (full version of Table 1). ‡ Published results with published error bars. Pairwise (k̄<4) Method MLP HGNN HyperGCN HNHN UniGNN AllDeepSets

Cora

Mixed (4≤k̄<7)

Higher-order (k̄≥7)

Citeseer PubMed Cora-CA

74.1±1.0 72.3±0.3 87.1±0.3 69.9±6.7 68.1±0.7 71.2±9.7 71.0±5.9 68.1±1.7 78.3±5.4 72.8±5.1 68.3±0.7 70.0±11.5 77.3±0.7 72.4±0.7 87.6±0.5 76.2±1.5 70.9±1.1 87.2±0.4

74.6±1.0 76.8±1.0 74.3±0.6 76.6±0.8 82.9±0.8 80.6±1.1

House

Senate Gene-Dis.

70.4 61.2±6.0 62.8±2.3 87.6±2.4 63.9±1.5 55.4±8.6 68.5 79.6±5.0 70.9±2.2 78.8±4.6 67.8±1.9 92.7±2.8

ED-HNN‡ 80.3±1.2 73.7±1.5 89.0±0.6 84.0±1.0 72.5±1.8 SheafHyperGNN‡ 81.3±1.1 74.7±1.3 87.7±0.5 85.5±0.9 73.8±1.6 D ENS N ET-D

— —

36.4±0.9 85.7±0.9 81.3±0.9 87.8±1.0 85.4±1.3 86.1±0.8 — —

77.5±1.2 72.6±0.5 87.7±0.5 81.2±1.0 70.3±2.3 92.7±4.1 86.0±1.0

Cells without ±: House MLP and HNHN are published values without reported variance. All other cells report 10-seed means ± standard deviation.

Pattern density separation across benchmarks. To understand which invariants separate each benchmark pair, we estimate the homomorphism density t(F, H) for 12 patterns F (ranging from single edges to CFI gadgets) on both hypergraphs in each IWS pair, using 10,000 Monte Carlo samples per density (Table 6). Table 6: Pattern density gaps |t(F, H1 ) − t(F, H2 )| for each benchmark pair. Only non-zero entries shown; all other pattern–benchmark combinations have zero gap. The separating patterns align with the hierarchy: CE-Hard requires Pasch-family patterns (ghtw = 3), Native-Hard is separated by even single-edge densities (ghtw = 1), and Multi-Order is separated by no pattern density—only the Sinkhorn alignment branch distinguishes this pair. Benchmark

Separating Pattern

Gap

Pattern ghtw −5

CE-Hard (STS pair)

Pasch / AntiPasch CFI gadgets (T0 , T1 )

9.7 × 10 9.7 × 10−5

3 3

Native-Hard (CFI pair)

SingleEdge DisjointPair

1.0 × 10−3 2.0 × 10−4

1 1

Multi-Order (Fano pair)

(none)

0

—

The results confirm the hierarchy prediction. CE-Hard requires patterns of ghtw ≥ 3 (Pasch configurations) for separation—simple edge or pair statistics are blind, explaining why cliqueexpansion methods fail. Native-Hard is separated by elementary statistics (ghtw = 1), consistent with the CFI construction targeting higher-order symmetry rather than density differences. MultiOrder shows identical densities for all 12 patterns, confirming that pattern density estimation alone cannot distinguish this pair; the Sinkhorn alignment branch (which tests structural alignment beyond counting) is essential. This three-way complementarity validates the combined InvNet architecture: different branches activate for different separation mechanisms. 24

F

Separation Constructions

This appendix provides the explicit constructions of hypergraph pairs that witness the strict inclusions in Theorem 3.6. Each pair serves a dual purpose: a theoretical proof artifact (witnessing the separation in the proof of Theorem 3.6, Appendix A.4) and an I NVARIANT W ITNESS S UITE expressivity benchmark. F.1

The Steiner Pair (CE vs. Native)

We construct two non-isomorphic 3-uniform hypergraphs on n = 13 nodes with identical clique expansions using Steiner triple systems. A Steiner triple system STS(v) is a collection of 3-element
subsets of [v] such that every 2-element
subset is contained in exactly one triple. An STS(v) has v2 / 32 = v(v − 1)/6 triples, and its clique expansion is Kv . The number of non-isomorphic STS(v) systems is 1 for v ∈ {7, 9} and 2 for v = 13 [16]. Thus v = 13 is the smallest order admitting a non-isomorphic STS pair. Explicit construction. Let H1 be the cyclic STS(13), constructed from base blocks {0, 1, 4} and {0, 2, 7} under the action of Z13 :   H1 = {i, i+1, i+4} : i ∈ Z13 ∪ {i, i+2, i+7} : i ∈ Z13 giving 26 triples. Let H2 be obtained by a Pasch trade on H1 : replacing the four triples {3, 4, 7}, {3, 5, 10}, {4, 10, 12}, {5, 7, 12} with {3, 4, 10}, {3, 5, 7}, {4, 7, 12}, {5, 10, 12}. The trade preserves all pairwise coverages (the same 12 pairs on 6 points are covered). Properties. 1. Identical clique expansions: Both cover every pair exactly once, so ϕ(A1 ) = ϕ(A2 ) = J − I (= K13 ). 2. Non-isomorphic: H1 has 13 Pasch configurations and automorphism group Z13 ⋊ Z3 (order 39); H2 has 8 Pasch configurations and a smaller automorphism group. Since the Pasch count is an isomorphism invariant, the two systems are provably non-isomorphic. 3. Separated by elementary invariant: ΘA1 (A1 ) = 156 > ΘA1 (A2 ), since no permutation maps H2 to H1 . F.2

The CFI Pair (Native vs. All)

We adapt the Cai–Fürer–Immerman construction [11] to 3-uniform hypergraphs. Let G = (VG , EG ) be a connected 3-regular base graph with girth > 2k + 1. For each edge e = {u, v}, introduce auxiliary nodes a1e , a2e and hyperedges: Type 0:

{u, a1e , a2e },

{v, a1e , a2e }

Type 1:

{u, v, a1e },

{u, v, a2e }

Define H1CFI using Type 0 for all edges and H2CFI using Type 1 for one edge. Proposition F.1. For base graphs of sufficient girth, H1CFI and H2CFI are non-isomorphic, indistinguishable by k-GWL (for k bounded by the girth), and separated by ΘP for appropriate P . F.3

Multi-Order Separation

For non-uniform hypergraphs in S2,n ⊕ S3,n , on n = 6: H1 : Edges: {1, 2}, {3, 4}, {5, 6}; H2 : Edges: {1, 2}, {3, 4}, {5, 6};

Triples: {1, 3, 5}, {2, 4, 6} Triples: {1, 4, 5}, {2, 3, 6}

Both have identical per-order statistics but different cross-order coupling. A cross-order template P = (P (2) , P (3) ) separates them; architectures processing each order independently cannot. 25

G

InvNet Architecture Details

This appendix expands on the I NV N ET architecture introduced in Section 4, providing the full algorithm (Algorithm 1), relaxation layer details, and computational complexity analysis. The ⊗k Sinkhorn relaxation of the pattern alignment score is Θ̂SK · A⟩, where Bn P (A) = maxD∈Bn ⟨P, D is the Birkhoff polytope of doubly stochastic matrices. The invariant branch cost is O(J · (Titer · n2 + |E| · k 2 ))—comparable to a single spectral computation for sparse hypergraphs. Algorithm 1, relaxation layers, and template parameterization are detailed below. G.1

Relaxed Elementary Invariant Layers

Computing ΘP (A) = maxΠ∈Sym(n) ⟨P, Π · A⟩ exactly is NP-hard (multidimensional assignment). We use two tractable relaxations. Sinkhorn relaxation.

Relax over the Birkhoff polytope Bn : ⊗k Θ̂SK · A⟩, P (A) = max ⟨P, D

(7)

D∈Bn

solved via the Sinkhorn operator [45, 2] with temperature τ . ⊗k Spectral relaxation. Relax over O(n): Θ̂spec · A⟩. For k = 2: closedP (A) = maxV ∈O(n) ⟨P, V ⊤ form λ(P ) λ(A) [17]. For k ≥ 3: Riemannian gradient ascent on the Stiefel manifold [1]. spec The relaxations satisfy ΘP ≤ Θ̂SK since Sym(n) ⊂ Bn and Sym(n) ⊂ O(n). P ≤ Θ̂P

G.2

Full Architecture

Algorithm 1 I NV N ET: Invariant-Theoretic Hypergraph Neural Network Require: Adjacency tensor A ∈ Sk,n , node features X ∈ Rn×d0 Require: Learnable templates {Pj }Jj=1 ⊂ Sk,n , offsets {αj }Jj=1 Require: Message-passing layers {W (ℓ) }L ℓ=1 , Sinkhorn temperature τ 1: // Local branch: hyperedge message passing 2: h(0) ← X 3: for ℓ = 1, . . . , L do
(ℓ−1) (ℓ) : v ∈ e} 4: me ← ϕe {hv
(ℓ) (ℓ−1) (ℓ) , {me : v ∈ e} 5: hv ← ϕv hv 6: end for 7: zlocal ← ReadOut(h(L) ) 8: 9: // Invariant branch: relaxed elementary invariants 10: for j = 1, . . . , J do 11: Sj ← MLPj (A) 12: Dj ← Sinkhorn(Sj /τ, Titer ) 13: θ̂j ← ⟨Pj , Dj⊗k · A⟩ − αj 14: end for 15: zinv ← [θ̂1 , . . . , θ̂J ] 16: 17: // Combine and predict 18: return MLPout (zlocal ∥zinv )

▷ Hyperedge aggregation ▷ Node update ▷ e.g., mean/sum pooling

▷ Score matrix from tensor ▷ Doubly stochastic relaxation

Computational complexity. Local branch: O(L · |E| · k · d). Invariant branch: O(J · (Titer · n2 + |E| · k 2 ))—O(n2 ) for sparse hypergraphs. Template parameterization. For IWS experiments (fixed n, k = 3), we use explicit learnable
templates over all nk potential hyperedges. For ANCS tasks, we parameterize templates implicitly 26

via a small MLP mapping edge-level embeddings to scalar weights, producing sparse template values over the input edge set. Remark G.1 (Towards universality). As J → ∞ and τ → 0, I NV N ET accesses pattern alignment scores for arbitrarily many templates with vanishing relaxation gap. However, two gaps prevent a formal universality claim: (i) the Sinkhorn relaxation gap is nonzero for fixed τ > 0; (ii) finite template families may not suffice for uniform approximation. Characterizing the approximation gap as a function of J, τ , and the target invariant’s hypertree width is an important open problem.

H

DensNet-D Architecture Details

This appendix provides the full architecture specification of D ENS N ET-D, the density-aware node classifier evaluated in Table 1 of the main text. D ENS N ET-D combines an AllDeepSets backbone [14] with per-node pattern density features via concatenation fusion and a learned per-node gate, implementing the density branch of the hierarchy (Theorem 2.4) in a practical, scalable form. H.1

Architecture

Given a hypergraph H = (V, E) with node features X ∈ Rn×d , D ENS N ET-D computes: 1. AllDeepSets backbone.

A two-layer message-passing network with set-function aggregation:

h(0) = Win xv , P (ℓ−1) m(ℓ) , e = u∈e hu

(8) P (ℓ) (ℓ) h̃v = e∋v me ,

(9)

(ℓ) (ℓ) (ℓ−1) h(ℓ) + MLP(ℓ) v = LN hv v (MLPe (h̃v )) ,

(10)

where each MLP is a two-layer network with ReLU and dropout, and LN denotes LayerNorm. For (ℓ) (ℓ) datasets with large hyperedges (k̄ ≥ 20, e.g., House), we use mean normalization: me ← me /|e| (ℓ) (ℓ) and h̃v ← h̃v / deg(v). 2. Density branch. For each node v, we estimate local pattern densities ρv = [t̂(F1 , star(v)), . . . , t̂(FP , star(v))] via Monte Carlo sampling over the star neighborhood star(v) = {e ∈ E : v ∈ e}, using 200 Monte Carlo samples per node (selected from grid {100, 200, 500}; Table 10). The pattern library {F1 , . . . , FP } is selected adaptively by the DATASET P ROFILER based on hypergraph statistics (edge size distribution, mean degree). Density features are projected: zv = LN(Wρ ρv + bρ ). Gradients are not propagated through the MC sampling; the density features are detached. 3. Per-node gate.

A learned gate controls the density contribution per node: gv = σ(w⊤ h(L) v + bg ),

(11)

(L)

where σ is the sigmoid function, hv is the backbone embedding, and bg is initialized to gate_init— a hyperparameter set by the profiler: • Higher-order (k̄ ≥ 7): bg = 0.0

(density active from epoch 0)

• Mixed (4 ≤ k̄ < 7): bg = −2.0

(density partially suppressed)

• Pairwise (k̄ < 4): bg = −5.0

(density strongly suppressed; σ(−5) ≈ 0.007)

The gate weight w is initialized to zero, so the initial gate value depends only on bg . 4. Concatenation fusion (cf. equation (11)). The gated density embedding is concatenated with the backbone output:
ŷv = MLPout h(L) ∥ gv · zv , (12) v where MLPout : R2dh → RC is a two-layer classifier. 27

5. Staged unfreezing. To prevent early noise from the density branch corrupting the backbone, we freeze the gate and density projection for an initial period: higher-order datasets unfreeze at epoch 0; mixed at epoch 30; pairwise at epoch 50. During the frozen phase, the model trains as a pure AllDeepSets backbone. H.2

Training Recipe

All datasets share the same core hyperparameters: hidden dimension 128, learning rate 10−3 , dropout 0.5, 2 message-passing layers, gradient clipping at 1.0, label smoothing 0.1, cosine learning rate schedule. Epochs and patience are dataset-dependent: 400 epochs with patience 150 for Senate-Bills; 200 epochs with patience 50 for all other datasets (Table 11). The only remaining dataset-dependent settings are gate_init and mean normalization, both determined automatically by the profiler. H.3

Relationship to PDN

PDN (used in IWS experiments) and D ENS N ET-D share the same pattern density estimation mechanism. The key differences are: (i) PDN operates on graph-level classification with global pooling, while D ENS N ET-D performs node classification; (ii) D ENS N ET-D uses the AllDeepSets backbone (which we validated as the strongest set-function baseline) rather than a generic local branch; (iii) D ENS N ET-D uses concatenation fusion (matching the theoretical prescription of equation (11)) rather than additive fusion.

I

Examples of Hypergraph Invariants

This appendix catalogs important hypergraph invariants and classifies each by the pattern complexity required for detection within the hierarchy of §3. The catalog provides a practical reference: given a target invariant, Table 7 identifies the minimum architecture tier needed to compute it. Let A ∈ Sk,n be the adjacency tensor of a k-uniform hypergraph. (k)

Local invariants (detectable by small patterns). Number of hyperedges: t(Kk , H) · nk /k!. Maximum weighted degree: detectable from star-pattern counts. Maximum k-clique density: ΘP (A) for P equal to the complete sub-hypergraph tensor. Spectral invariants (requiring global pattern information). Clique-expansion spectral invariants: any graph invariant of the clique expansion AGH , but this map loses information. Hodge Laplacian invariants: eigenvalues of the Hodge Laplacian Lj capture topological information invisible to clique expansion, requiring patterns of hypertree width > 1. Cut and optimization invariants. Hypergraph MAXCUT: continuous, NP-hard. Cheeger constant [13, 34]: characterizes expansion. Relaxation-based invariants. Tensor theta body [12]: SDP relaxation of independence number. Spectral relaxation ΛP (A) = maxV ∈O(n) ⟨P, V ⊗k · A⟩: for k = 2, equals inner product of sorted eigenvalues; for k ≥ 3, no closed form. Table 7: Hypergraph invariants and the pattern complexity required for detection. Invariant

Pattern width

Order

Tractable

Detectable by

Number of hyperedges Maximum degree Max clique density CE spectrum Hodge spectral gap Hypergraph MAXCUT Cheeger constant ΛP

1 (single edge) 1 (star) ≥2 1 (graph) >1 Unbounded Unbounded Varies

Any k Any k k-uniform Any k Simplicial k-uniform k-uniform k-uniform

Yes Yes No Yes Yes No No k=2: Yes

All HGNNs All HGNNs Pattern alignment CE HGNNs Native HGNNs — — Spectral layers

28

J

Related Work

We situate our contributions within four threads of related work. The key distinction is that prior work addresses either homomorphism counting or neural network expressivity, but not the two together in the hypergraph setting. Homomorphism counting and GNN expressivity. Lovász [36] proved that graph homomorphism counts determine graphs up to isomorphism; Lovász [37] developed this into the theory of graph limits. Dell et al. [19] connected homomorphism counts to WL expressivity. For hypergraphs, Scheidt and Schweikardt [44] established that homomorphism indistinguishability over patterns of bounded generalized hypertree width corresponds to bounded-variable counting logic, and Scheidt [43] refined this to hypertree depth. Our work synthesizes these into a complete HGNN expressivity characterization. GNN expressivity and the WL hierarchy. Xu et al. [49] and Morris et al. [41] established the connection between message-passing GNN expressivity and 1-WL. Higher-order networks match k-WL [39, 40]. For hypergraphs, Böker [9] extended color refinement, and Zhang et al. [52] developed k-dimensional generalized WL, proving that (k+1)-GWL is strictly more expressive than k-GWL (and that k-GWL reduces to k-WL on simple graphs). Their k-HNN architecture validates the hierarchy experimentally on hypergraph-classification benchmarks (IMDB, SteamPlayer, TwitterFriend). Our framework subsumes and explains the k-GWL hierarchy: by the known WL– homomorphism correspondence, each GWL level corresponds to the pattern-density class Hk indexed by generalized hypertree width, and our contribution is showing this hierarchy is induced by a single mechanism—pattern counting—that yields a complete characterization of HGNN expressivity (Theorem 3.6). Choi et al. [15] extended WL to categorical structures. Hypergraph neural networks. The landscape has evolved from clique-expansion methods [23, 50] to native hyperedge processing [20, 28, 14, 3] and equivariant architectures [30, 48]. Zaheer et al. [51] established the Deep Sets foundation. Duta et al. [21] introduced sheaf-theoretic structure. Topological and geometric deep learning. Bodnar et al. [8] extended message passing to simplicial complexes; Hajij et al. [27] proposed a general topological framework. Bronstein et al. [10] proposed a group-theoretic design framework. Our work instantiates this for the symmetric group on hypergraph tensors. Tensor methods. Tensor decomposition [32] does not exploit permutation invariance. Our pattern alignment scores bridge tensor optimization and hypergraph structure.

K

Structural Frontiers: Sinkhorn Integrality and Width Thresholds

The strict hierarchy of Theorem 3.6 raises two natural structural questions that we develop here as research programs. (i) When is the Sinkhorn relaxation exact? That is, when does convex optimization over the Birkhoff polytope recover the combinatorial pattern count? (ii) What is the minimum generalized hypertree width required to approximate the Hodge spectral gap? Both questions identify frontiers of the hierarchy: when convex relaxation preserves combinatorial invariants, and when bounded-width pattern counting suffices to approximate topological spectra. A computational verification of the Sinkhorn separation on the Steiner pair follows in Appendix L. K.1

The Sinkhorn Relaxation Gap

For a k-uniform template P and adjacency tensor A, define the Sinkhorn relaxation gap: ⊗k gapSK (P ; A) := Θ̂SK · A⟩ − P (A) − ΘP (A) = max ⟨P, D D∈Bn

max ⟨P, Π⊗k · A⟩ ≥ 0. Π∈Sym(n)

(13) For the entropy-regularized objective used in the architecture, the temperature-dependent gap is   gapτ (P ; A) := Θ̂τP (A) − ΘP (A), where Θ̂τP (A) = max ⟨P, D⊗k · A⟩ + τ H(D) . D∈Bn

29

Since Sym(n) ⊂ Bn , both gaps are nonnegative. For k = 2 (graphs), gapSK reduces to the well-studied QAP relaxation gap: Bn replaces Sym(n) in maxΠ ⟨P, ΠAΠ⊤ ⟩ [18]. For k ≥ 3, the gap involves the tensor Birkhoff polytope {D⊗k : D ∈ Bn }, a strict subset of the set of doubly stochastic k-tensors—a structure with no known analogue of the Birkhoff–von Neumann theorem. Understanding when these gaps vanish—or can be structurally bounded—is central to quantifying the fidelity of the relaxation. K.2

Integrality by Linearity

Lemma A.6 proves that for Steiner triple system templates, every maximizer over Bn has singleton row support and is therefore a permutation matrix, giving gapSK = 0. The proof exploits a combinatorial forcing mechanism: (1) each alignment term is bounded by 1; (2) optimality forces all terms to equal 1; (3) this forces singleton row support; (4) hence D∗ is a permutation matrix. The key structural property is the uniqueness of pairwise overlaps between triples. We isolate the abstract condition and separate the analytic ingredient (termwise bounds and tightness) from the combinatorial one (shadow uniqueness): Definition K.1 (Unique (k−1)-shadow property). A k-uniform hypergraph P satisfies the unique (k−1)-shadow property if every (k−1)-subset of vertices belongs to at most one hyperedge: |{e ∈ E(P ) : S ⊂ e}| ≤ 1 for all |S| = k − 1. This property generalizes linearity of Steiner systems: an STS is linear precisely because every pair appears in exactly one triple. Partial Steiner systems, Steiner quadruple systems, and more generally all 1-designs satisfy this condition. We define the co-degree δ(P ) := max|S|=k−1 |{e ∈ E(P ) : S ⊂ e}|, so the unique shadow property is precisely δ(P ) = 1. Theorem K.2 (Integrality under compatible unique-shadow structure (sketch)). Let P, A ∈ Sk,n be nonnegative k-uniform adjacency tensors satisfying the following conditions: (C1) Nonnegativity: P ≥ 0, A ≥ 0. (C2) Termwise unit bound with tightness: For every k-tuple i ∈ supp(P ), the contribution (D⊗k · A)i ≤ 1 for all D ∈ Bn , and there exists a permutation Π∗ ∈ Sym(n) achieving equality on every supported tuple (i.e., the hard optimum saturates every term). (C3) Compatible shadows: Both supp(P ) and supp(A) satisfy the unique (k−1)-shadow property. Then every maximizer D∗ of Θ̂SK P (A) over Bn is a permutation matrix, and gapSK (P ; A) = 0. Proof sketch (mirrors Lemma A.6 in k-uniform form). Step 1 (termwise upper bound): By (C1) and doubly-stochasticity, each “edge contribution” in ⟨P, D⊗k · A⟩ satisfies (D⊗k · A)i ≤ 1 for i ∈ supp(P ). This is the “156 summands each ≤ 1” step of Lemma A.6. Step 2 (tightness forces support preservation): By (C2), the hard optimum achieves value |E(P )| · k! (or the analogous saturated value). Since Θ̂SK P (A) ≥ ΘP (A) = |E(P )| · k! and each of the |E(P )| terms is ≤ 1, the maximizer D∗ must achieve equality on every supported term. Concretely: whenever Di∗1 a1 · · · Di∗k ak > 0 and (i1 , . . . , ik ) ∈ supp(P ), we must have (a1 , . . . , ak ) ∈ supp(A). Step 3 (unique shadow forces singleton support): Suppose row i of D∗ has two positive entries ∗ ∗ Dia , Dia > 0 with a1 ̸= a2 . Choose any hyperedge e ∈ E(P ) containing vertex i, and let 1 2 b2 , . . . , bk be vertices in the support of the remaining rows. Step 2 forces both {a1 , b2 , . . . , bk } and {a2 , b2 , . . . , bk } to be hyperedges in supp(A). These two edges share the (k−1)-subset {b2 , . . . , bk }, violating (C3). Hence every row has singleton support, and D∗ is a permutation matrix. Scope of the conditions. Condition (C2) is the most restrictive: it requires a permutation achieving simultaneous saturation of all alignment terms. This holds automatically for self-template alignment (P = A1 , evaluating Θ̂SK A1 (A1 ) where the identity permutation saturates) and for any A isomorphic to P . For general non-isomorphic pairs where the optimal permutation does not align all edges, (C2) may fail and the gap can be nonzero even with δ = 1. This precisely delineates the theorem’s applicability: it guarantees that the Sinkhorn relaxation correctly identifies isomorphic copies, which is the separation-relevant setting (cf. Lemma A.6). 30

When the gap is nonzero. For dense templates with δ(P ) > 1, doubly stochastic matrices can “spread” mass across partial overlaps, yielding gapSK > 0. The forcing argument breaks because two edges of A may share a (k−1)-subset when δ(A) > 1, so the contradiction in Step 3 does not arise. K.3

Gap vs. Width: A Structural Frontier

The integrality conjecture isolates a combinatorial sufficient condition (δ(P ) = 1) for exactness. We conjecture a quantitative relationship between the gap and template structure: Conjecture K.3 (Sinkhorn gap bound). For k-uniform templates P with |V (P )| = v vertices and ∥A∥F ≤ BA :
gapSK (P ; A) ≤ Ck · δ(P ) − 1 · v k−2 · BA , where Ck depends only on k. In particular, δ(P ) = 1 (the unique shadow property) implies gapSK = 0. Evidence. For k = 3, δ = 1 (STS): Theorem K.2 gives gap = 0, consistent with the bound vanishing. For k = 2 (graphs): the QAP relaxation gap scales as O(∆(P ) · n) where ∆(P ) is the maximum degree [18], consistent with δ = ∆, v = n. A broader structural question connects the gap to the hierarchy itself: Question K.4 (Gap vs. hypertree width). Does there exist a function f such that sup∥A∥F ≤1 gapSK (P ; A) ≤ f (ghtw(P ))? In particular: does ghtw(P ) = 1 imply exactness? Do large-width templates admit provably positive tensor-Birkhoff gaps? A positive answer would identify a structural frontier beyond which convex relaxation can “fractionally cheat.” Combined with the generalization bound (Proposition C.3), this suggests an expressivity– generalization sweet spot: templates with low δ(P ) and moderate ghtw achieve small relaxation gap and favorable sample complexity—quantifying I NV N ET’s practical reach into the hierarchy. K.4

Width Threshold for the Hodge Spectral Gap

We now turn to a second structural frontier: the minimal hierarchy level required to approximate the Hodge spectral gap. Setup. For a k-uniform hypergraph H viewed as a (k−1)-dimensional simplicial complex (taking ⊤ the downward closure), the Hodge j-Laplacian is Lj = ∂j+1 ∂j+1 + ∂j⊤ ∂j , where ∂j is the j-th boundary operator. The Hodge spectral gap λ1 (L1 ) (the smallest nonzero eigenvalue of L1 ) controls edge random walk mixing and topological connectivity [8, 27]. For k = 2, L0 is the graph Laplacian and λ1 (L0 ) has pattern width 1. For k ≥ 3, L1 captures genuinely higher-order structure invisible to the clique expansion. Why r∗ > 1. Two hypergraphs with the same clique expansion can have different Hodge spectra: the Steiner pair (§F.1) has identical clique expansions (both are K13 weighted by pair-co-occurrence) but potentially different λ1 (L1 ). Therefore λ1 (L1 ) ∈ / ICE , confirming pattern width > 1 (Table 7). K.5

The Spectral Moments Route

The key technical tool for determining the width threshold is the expansion of spectral moments as pattern counts, combined with polynomial approximation of spectral functions. P Moments as pattern densities. For any symmetric matrix M , tr(M m ) = i λm i . Each entry of L1 encodes adjacency relations between simplices. Each term in Lm 1 corresponds to a length-m closed walk in the incidence/boundary structure, inducing a finite simplicial pattern with at most O(m) simplices and ghtw = O(m). Therefore: X tr(Lm cF · t(F, H), (14) 1 ) = F ∈Pm

where Pm is a finite set of patterns with ghtw(F ) = O(m), and cF are integer combinatorial coefficients. Concretely: m = 1 yields degree statistics (ghtw = 1); m = 2 involves the squared 31

boundary operator ∂2 ∂2⊤ , whose (e, e′ )-entry counts triangles sharing edges e and e′ —determined by the diamond pattern (4 vertices, 2 hyperedges, ghtw = 2). Since ∂1⊤ ∂1 depends only on the clique expansion (ghtw = 1), the full tr(L21 ) requires ghtw ≤ 2. From moments to spectral gap approximation. We formalize the width–accuracy relationship: Proposition K.5 (Moment-based spectral gap approximation (informal)). Let L1 ⪰ 0 be the Hodge 1-Laplacian of a 3-uniform hypergraph H with ∥L1 ∥op ≤ Λ and nullity b1 = dim ker(L1 ). Fix
PL ϵ ∈ (0, 1). There exists a polynomial pL (x) = ℓ=0 cℓ xℓ of degree L = O (Λ/ϵ) log(1/ϵ) such that: • pL (0) = 0 and |pL (x)| ≤ ϵ/10 for x ∈ [0, λ1 − ϵ]; • |pL (x) − 1| ≤ ϵ/10 for x ∈ [λ1 , Λ]. ℓ The quantity d1 tr(pL (L1 )) (where d is the matrix dimension) can be computed from moments P tr(L1 ), ℓ ≤ L, and determines λ1 to within O(ϵ) provided moment estimation error satisfies η ℓ |cℓ | ≤ ϵ/10.

Proof sketch. Construct pL as a Chebyshev polynomial approximation to a smoothed step function on [0, Λ] that transitions from 0 to 1 in an ϵ-neighborhood of λ1 . Standard P Chebyshev approximation theory gives degree L = O((Λ/ϵ) log(1/ϵ)) with coefficient growth |cℓ | quasi-polynomial in L. 1 Then d1 tr(pL (L1 )) ≈ d−b if λ1 > ϵ, and the transition region determines λ1 via binary search. d Width–accuracy curve. Since each moment tr(Lℓ1 ) decomposes into pattern densities of ghtw = O(ℓ) via (14), approximation at hierarchy level Hr can estimate moments up to degree L ≈ c · r. Combining with Theorem K.5:   Λ (up to logarithmic factors). width r =⇒ spectral gap accuracy ϵ ≈ O r This is exactly the Jackson-type width–accuracy relationship flagged as open in §6, instantiated for the Hodge spectral gap. It connects directly to the approximation pipeline of Proposition B.6: the polynomial degree d(ϵ) from the Molien series determines the pattern enumeration, and (14) shows these patterns have bounded hypertree width. K.6

The Minimal Width Conjecture

Conjecture K.6 (Hodge gap width threshold). For 3-uniform hypergraphs: (a) r∗ ≥ 2: width-1 pattern densities cannot detect glued 2-simplices sharing edges. (b) r∗ = 2 suffices for bounded-degree families: the diamond pattern at ghtw = 2 captures the dominant spectral structure. (c) In the worst case, r∗ may scale as O(log(1/ϵ)), since the number of moment terms m(ϵ) needed to resolve the spectral gap can grow logarithmically. More generally, for k-uniform hypergraphs, r∗ ≤ k − 1. Geometric evidence for r∗ = 2. The Hodge 1-Laplacian on a 3-uniform hypergraph is L1 = ∂2 ∂2⊤ + ∂1⊤ ∂1 . The term ∂1⊤ ∂1 depends only on the 1-skeleton (ghtw = 1). The term ∂2 ∂2⊤ counts shared faces between triangles—precisely the diamond patterns with ghtw = 2. Since eigenvalues are continuous in matrix entries (Weyl’s perturbation theorem), λ1 (L1 ) is a continuous function of pattern densities up to ghtw = 2. Lower bound strategy. To prove r∗ ≥ 2, it suffices to construct hypergraph pairs that agree on all width-1 pattern densities but differ in Hodge spectral gap. Such constructions would parallel CFI-type separations (§F.2) but target spectral topology rather than isomorphism power. The Steiner pair is a natural candidate: both STS(13) systems have identical degree sequences and clique expansions (ghtw-1 information); verifying computationally that they differ in λ1 (L1 ) would establish r∗ = 2 for k = 3. 32

The general case. For k-uniform hypergraphs, Lk−2 involves boundary operators mapping (k−1)⊤ chains to (k−2)-chains. The squared boundary operator ∂k−1 ∂k−1 counts shared faces between hyperedges—patterns with ghtw ≤ k − 1, yielding the general bound r∗ ≤ k − 1. K.7

Stress Tests and Scope of the Conjectures

We explicitly examine where the conjectures could fail, distinguishing worst-case from distributional settings. Stress test 1: when does integrality fail? Theorem K.2 requires condition (C2): the hard optimum saturates every alignment term. Consider a partial Steiner system P with δ(P ) = 1 but only 5 triples on 13 vertices. If A has 26 triples and the optimal permutation aligns only 5 of 26 edges of A with P , the remaining 21 terms are not saturated. The optimizer could then “spread” mass to improve unsaturated terms, potentially breaking integrality even with δ = 1. This confirms that (C2) is not merely technical—it delineates when the forcing mechanism operates. The theorem is strongest in the separation-relevant regime (testing A against itself or its isomorphic copies), which is precisely the setting of Lemma A.6. Stress test 2: worst-case vs. distributional necessity of r∗ ≥ 2. On certain random ensembles (e.g., Erdős–Rényi 3-uniform H3 (n, p)), the Hodge spectral gap concentrates and correlates strongly with degree statistics. A width-1 model could achieve small average error simply by learning the typical value λ1 ≈ g(n, p). On random ensembles, all architectures (including CE-based ones) can achieve low average error by learning the concentrated typical value, since λ1 concentrates around a function of the edge density. Therefore the “r∗ ≥ 2 is necessary” claim should be understood as a worst-case statement: there exist specific instance pairs where width-1 provably fails, even though width-1 may suffice in distribution. Stress test 3: could r∗ = 2 fail uniformly? Width-2 patterns capture the diamond structure (two triangles sharing an edge), which determines ∂2 ∂2⊤ entrywise. However, the spectral gap may depend on global arrangements of diamonds—e.g., existence of long homological 1-cycles in the simplicial complex—that require resolving many moments. By Theorem K.5, the number of moments scales as O((Λ/ϵ) log(1/ϵ)), and each moment at order m requires ghtw = O(m). For complexes with near-degenerate spectra (eigenvalues clustered near λ1 ), many moments are needed, pushing the effective width beyond 2. The safest summary: r∗ = 2 is the first plausible threshold, and likely suffices for bounded-degree families, but uniform approximation across all 3-uniform hypergraphs may require r∗ = O(log(1/ϵ))—exactly Conjecture K.6(c). Stress test 4: computational verification via the Steiner pair. The lower-bound strategy (§K.6) hinges on the Steiner pair having different Hodge spectra. This is a concrete computational question: construct the Hodge 1-Laplacian L1 for both STS(13) systems and compare λ1 (L1 ). If the spectra coincide, the Steiner pair is not a valid separator for the Hodge gap and a different construction is needed—potentially hypergraphs with the same width-1 pattern densities but different simplicial topology (e.g., a triangulated cylinder vs. a triangulated Möbius strip with matching local statistics). K.8

Unified Perspective

Both structural programs ask related questions at different levels of the hierarchy: • Sinkhorn integrality: When does convex relaxation over Bn preserve the combinatorial optimum over Sym(n)? The answer depends on the co-degree δ(P ) and condition (C2): combinatorial sparsity forces integrality. • Hodge width threshold: When does bounded-width pattern counting capture spectral topology? The answer depends on the walk structure of Lm 1 : spectral moments decompose into pattern densities of bounded hypertree width. The unifying theme is that combinatorial sparsity and bounded overlap enforce both integrality and spectral detectability; high overlap and high width enable fractional phenomena and topological invariants invisible at low levels of the hierarchy. Both frontiers are controlled by the same generalized hypertree width that stratifies our expressivity hierarchy, suggesting a deeper structural theory 33

connecting relaxation gaps, spectral approximation, and architectural expressivity. Identifying these frontiers precisely—through formal proofs, counterexamples, or computational verification—would transform the qualitative hierarchy of Section 3 into a quantitative structural theory with direct architectural implications.

L

Steiner Pair: Sinkhorn Separation Verification

Appendix K develops the theory of the Sinkhorn relaxation gap; here we verify it computationally. We directly confirm the separation guarantee of Lemma A.6 on the two non-isomorphic STS(13) Steiner triple systems (n=13, |E|=26 triples each). We compute the Sinkhorn alignment score θ̂SK using co-occurrence profile similarity matrices and 200 Sinkhorn iterations, sweeping the temperature τ from 0.01 to 1.0. Table 8: Sinkhorn alignment scores on the STS(13) pair. Scores are edge-normalized (max = |E| = 26, corresponding to 26 × 3! = 156 in the tensor inner product of Lemma A.6). The gap θ̂SK (A1 , A1 ) − θ̂SK (A1 , A2 ) > 0 at every τ , confirming the theoretical separation. τ

Self θ̂SK (A1 , A1 )

Cross θ̂SK (A1 , A2 )

Gap

0.01 0.05 0.10 0.30 0.50 1.00

26.000 26.000 25.958 8.962 1.525 0.247

22.000 22.000 21.964 7.583 1.291 0.209

4.000 4.000 3.993 1.379 0.235 0.038

best: 15.0

mean: 2.36

Random perm. (100K samples):

(exact: D∗ = I)

Table 8 confirms the prediction of Lemma A.6: self-alignment strictly exceeds cross-alignment at every temperature. At low τ (≤ 0.05), Sinkhorn converges to a near-permutation matrix and the self-score reaches the theoretical maximum of 26 (= 156/3!); the cross-score saturates at 22, reflecting the best achievable alignment between the non-isomorphic pair. A reverse check confirms symmetry: θ̂SK (A2 , A2 ) − θ̂SK (A2 , A1 ) = 3.99 at τ = 0.1. The random-permutation baseline (mean cross-score 2.36 over 100K samples, best 15.0) shows that the Sinkhorn alignment score is far more discriminative than random matching, validating the practical utility of the Sinkhorn relaxation for hypergraph isomorphism testing. As τ increases, both scores decay (the doubly-stochastic matrix spreads mass), but the gap remains strictly positive—consistent with the proof that any maximizer in B13 achieving equality must be a permutation matrix (Lemma A.6). This provides direct empirical evidence that I NV N ET’s Sinkhorn branch can distinguish CE-indistinguishable hypergraphs at practical temperature settings (τ ∈ [0.05, 0.3]).

M

Dataset Statistics and Preprocessing

Table 9 summarizes the structural characteristics of all seven ANCS benchmarks used in Table 1. The regime classification follows the DATASET P ROFILER: pairwise (k̄ < 4), mixed (4 ≤ k̄ < 7), and higher-order (k̄ ≥ 7). This classification determines the theory-guided gate initialization for density-aware models (Appendix N). Data provenance. Cora, Citeseer, and PubMed are co-citation hypergraphs from the AllSet benchmark [14], originally constructed by Yadati et al. [50]; each hyperedge groups the papers co-cited by a common paper. Cora-CA is a co-authorship hypergraph from the same benchmark, where each hyperedge groups the co-authors of a single publication. Senate-Bills and House-Committees are legislative co-sponsorship hypergraphs from the Cornell SocioPatterns collection [24]; we follow the binary label protocol of Chien et al. [14] with Gaussian noise features (d=100, noise=1.0), matching the evaluation of AllSet, ED-HNN [48], and SheafHyperGNN [21]. All datasets are publicly available; data-loading scripts and preprocessed splits are included in the supplementary code. 34

Data splits. All experiments use 50/25/25 train/validation/test random splits. We generate 10 deterministic splits using seeds 0, . . . , 9 and report mean ± standard deviation across all 10. Table 9: Dataset statistics. n: nodes, m: hyperedges, d: input features, C: classes, k̄: mean edge ¯ mean node degree. cardinality, k95 : 95th-percentile edge size, kmax : largest edge, deg: Dataset

n

m

d

C

k̄

k95

kmax

¯ deg

Regime

Cora Citeseer

2,708 3,312

1,579 1,079

1,433 3,703

7 6

3.03 3.20

5 7

5 26

1.8 1.0

pairwise pairwise

PubMed Cora-CA Senate House Gene-Disease

19,717 2,708 294 1,290 5,012

7,963 1,072 29,157 340 2,009

500 1,433 100 100 100

3 7 2 2 21

4.35 4.28 7.96 34.89 14.0

12 13 28 66 38

171 43 99 82 130

1.8 1.7 789.6 9.2 5.6

mixed mixed higher-order higher-order higher-order

Two structural extremes are worth noting. Senate-Bills has extreme density: 294 nodes but 29,157 hyperedges, giving a mean degree of ≈790 (each senator appears in ≈790 bills); the heavy tail (k95 = 28, kmax = 99) confirms that co-sponsorship networks are genuinely higher-order, not clique expansions of pairwise data. House has the largest mean edge size (k̄ ≈ 35; committees average ≈35 members) but far fewer edges (m=340), representing the opposite higher-order extreme: sparse-large-edges versus dense-many-edges. These contrasting structures exercise different aspects of the density-estimation pipeline (normalization, MC sample count) and motivate the adaptive configuration in Appendix N.

N

Hyperparameter Configuration and Reproducibility

We describe the full training protocol for all models in Table 1. All models share a common backbone and differ only in their aggregation mechanism and (for density-aware models) in theory-guided hyperparameters derived from the DATASET P ROFILER (Appendix M). Shared training protocol. All models use: Adam optimizer, learning rate 0.001, weight decay 0, gradient clipping at 1.0, label smoothing 0.1, cosine-annealing learning rate schedule, and early stopping with patience 50 (150 for Senate-Bills due to high variance). Training runs for 200 epochs (400 for Senate-Bills). All experiments use 10 random seeds per dataset. Baseline models. MLP, HGNN, HyperGCN, HNHN, UniGNN, AllDeepSets, and AllSetTransformer all use hidden dimension 64, dropout 0.2, and 2 message-passing layers. These settings follow the AllSet benchmark protocol [14]. DensNet-D. D ENS N ET-D uses hidden dimension 128, dropout 0.5, 2 AllDeepSets backbone layers, and a single-layer gate network gv = σ(w⊤ hv + b) that modulates the density branch per node. The gate initialization b is set by the DATASET P ROFILER according to regime: b = 0.0 (higher-order, k̄ ≥ 7), b = −2.0 (mixed, 4 ≤ k̄ < 7), b = −5.0 (pairwise, k̄ < 4). The density branch uses ℓ2 -normalization of density features only when k̄ ≥ 20 (House); on smaller-edge datasets, raw features perform better (Appendix P). Staged unfreezing freezes the gate and density projection for 0 epochs (higher-order), 30 (mixed), or 50 (pairwise), allowing the backbone to converge before the density branch begins learning. Table 10 summarizes the model-specific hyperparameters; Table 11 gives the per-dataset D ENS N ET-D configuration derived automatically by the DATASET P ROFILER. Hyperparameter selection. For baselines, we adopt the AllSet benchmark protocol [14] without modification to ensure fair comparison. For D ENS N ET-D, we performed a grid search over hidden ∈ {64, 128, 256}, dropout ∈ {0.2, 0.5}, and MC samples ∈ {100, 200, 500} on the Senate-Bills validation set. Increasing hidden to 256 added 3.6× parameters for +0.3% accuracy; increasing MC beyond 200 did not change accuracy by more than 0.1%. We verified that these choices transfer to all other datasets without re-tuning. All regime-dependent hyperparameters (gate init, freeze epochs, 35

Table 10: Model-specific hyperparameters. All models share: Adam optimizer, lr=0.001, weight decay=0, gradient clipping=1.0, label smoothing=0.1, cosine-annealing LR schedule (min lr=0). “Gate init” is b in σ(b); “MC” is the number of Monte Carlo samples for density estimation. Hidden

Dropout

Layers

Gate init

MC

Params†

MLP HGNN HyperGCN HNHN UniGNN AllDeepSets AllSetTransformer

64 64 64 64 64 64 64

0.2 0.2 0.2 0.2 0.2 0.2 0.2

2 2 2 2 2 2 2

— — — — — — —

— — — — — — —

96,135 96,135 96,135 96,263 96,135 100,487 108,839

D ENS N ET-D

128

0.5

2

by regime

200

∼180K††

Model

†

Param counts on Cora (d=1,433). †† DensNet-D count is for Senate (d=100); on Cora the input projection alone adds ∼183K parameters.

Table 11: D ENS N ET-D per-dataset configuration, derived automatically by the DATASET P ROFILER from dataset statistics (Table 9). Gate init b controls the initial density-branch contribution σ(b); freeze epochs delay density-branch learning until the backbone converges; ℓ2 -norm applies to density features before fusion. Dataset Cora Citeseer PubMed Cora-CA Senate House Gene-Disease

Regime

Gate init b

σ(b)

Freeze ep.

ℓ2 -norm

Epochs

Patience

pairwise pairwise mixed mixed higher-order higher-order higher-order

−5.0 −5.0 −2.0 −2.0 0.0 0.0 0.0

0.007 0.007 0.119 0.119 0.500 0.500 0.500

50 50 30 30 0 0 0

No No No No No Yes No

200 200 200 200 400 200 200

50 50 50 50 150 50 50

normalization) are derived automatically from dataset statistics by the DATASET P ROFILER, requiring zero manual configuration per dataset. Software and hardware. Baseline sweeps and InvNet tuning run on 8× NVIDIA H100 NVL GPUs (95,830 MiB each) with PyTorch 2.6.0, CUDA 12.4, Python 3.10 (conda hginv environment). D ENS N ET-D density precomputation (a one-time CPU cost) runs on an Apple M4 Max (128 GB unified memory): <1 minute for most datasets, ≈2 hours for PubMed (n=19,717), and ≈7 hours for Gene-Disease (n=5,012, k̄=14) due to larger star neighborhoods; model training runs on the same H100 GPUs as all other models. Configuration files for all experiments and the exact launch commands are included in the supplementary code.

O

Compute Budget

Table 12 reports the computational cost of the main experimental campaigns. The total cost is ≈200 H100 GPU-hours, of which approximately 40% (≈80 GPU-hours) produced the final reported results and 60% (≈120 GPU-hours) was spent on development iterations, debugging, and hyperparameter exploration. The dominant cost is invariant-based models (I NV N ET), which account for all five timeout events (3-hour cap per job) in the 8-GPU sweeps due to the Sinkhorn iteration bottleneck; density-based models (D ENS N ET-D) and all baselines complete in under 5 minutes per seed on every dataset. At approximately 700W per H100, the total energy cost is ≈140 kWh, corresponding to ≈56 kg CO2 e at the US average grid intensity of 0.4 kg/kWh. This is modest relative to the scale of the evaluation (approximately 860 individual training runs across the full sweep alone) and reflects the efficiency of the 8-GPU parallel scheduling: most jobs complete in under 30 minutes, with only InvNet-based models requiring hour-scale compute. 36

Table 12: Compute budget across experimental campaigns. “Wall-clock” is end-to-end time including job scheduling overhead. Peak GPU memory: I NV N ET on Senate-Bills uses ≈3,000 MiB; all other models use < 500 MiB. Campaign

Hardware

GPUs

Wall-clock

GPU-hours

Full sweep (8 models × 14 tasks × 10 seeds) Replication run (independent verification) InvNet tuning (Senate, House, Cora-CA) Development iterations (canary/ramp/full cycles) DensNet-D evaluation (7 datasets × 10 seeds)

H100 NVL H100 NVL H100 NVL H100 NVL M4 Max

8 8 1 8 CPU

∼8h ∼8h ∼12h ∼15h ∼24h

∼40 ∼40 ∼12 ∼100 — ≈ 200

Total H100

P

When Density Does Not Help

Not every design choice improves performance, and not every dataset benefits from density features. We document seven negative results that shaped the final D ENS N ET-D design (summarized in Table 13); each is informative about when and why the Width Wall matters. Table 13: Summary of negative results. ∆: accuracy change relative to the best configuration. Design variant

Dataset

∆ acc.

MEAN normalization MEAN normalization Senate config → Cora-CA Additive fusion Additive fusion Density branch Hierarchy scaling

Senate House Cora-CA Senate House Gene-Dis. IWS

−5.4% +4.5% −0.3% ≈−2% ≈−3% 0.0% flat at 100%

Root cause Sum-pooling preserves magnitude on small edges Scale imbalance on large edges (k̄=35) Width Wall gap too narrow (k̄=4.3) Shared embedding space limits expressivity Same Node features dominate (d=100) Local MP solves all benchmarks

1. MEAN normalization hurts small-edge datasets. We expected normalizing density features by edge count (replacing sum-pooling with mean-pooling over MC samples) to improve performance by removing scale dependence. Instead, it reduces Senate-Bills accuracy by 5.4 percentage points (90.5% → 85.1%) but improves House by 4.5 points. The explanation is structural: Senate has many small-to-medium edges (k̄ ≈ 8) where sum-pooling preserves discriminative magnitude variation, while House has few very large edges (k̄ ≈ 35) where unnormalized sums create scale imbalance across nodes with different neighborhood sizes. The final design normalizes only when k̄ ≥ 20. 2. Cora-CA tuning does not transfer from Senate. Applying the Senate-Bills recipe (hidden=128, J=4, dropout=0.2) to Cora-CA produced 79.7 ± 0.7%, marginally below the default configuration (80.0 ± 1.1%). A grid search over hidden ∈ {128, 256} and J ∈ {4, 8} yielded at best 80.1%. The lesson: Cora-CA’s small edges (k̄ = 4.28) leave little room for invariant-based features; the gap between the Width Wall and the baseline is too narrow for tuning to exploit. 3. Additive fusion underperforms concatenation. An early D ENS N ET-D variant used additive fusion (hv + gv · zdensity ) instead of the concatenation-based design of Equation 5 in the main text. v Additive fusion constrains density features to the same embedding space as the backbone, limiting the expressivity of the combined representation. Switching to concatenation with a learned projection improved Senate-Bills by ≈2 points and House by ≈3 points. 4. AllDeepSets and AllSetTransformer collapse on structured data. In an early baseline sweep (v3 code, different hyperparameters), AllDeepSets achieved only 50.2% on House, near chance for a binary task; the final tuned result in Table 1 is higher (67.8%). AllSetTransformer similarly underperformed on the Multi-Order IWS task: 69.7% (vs. 100% for HNHN/UniGNN) with high variance (±23.7%). These failures highlight that set-function aggregations require careful tuning for node classification on heterogeneous graphs, despite performing well on Senate-Bills (92.7%) where dense co-sponsorship structure favors global aggregation. 37

5. Gene-Disease: density lift revealed by strict ablation (+4.3%). Gene-Disease (n=5,012, m=2,009, k̄=14, 21 classes) was initially classified as “indifferent” based on the original ablation (AllDeepSets at hidden=64 scoring 86.1% vs. D ENS N ET-D at hidden=128 scoring 86.0%). The strict ablation (Table 4) reveals the true picture: AllDeepSets at the identical DensNet-D recipe (hidden=128, dropout=0.5, label smoothing, cosine schedule) drops to 81.7 ± 1.1%—the larger model overfits on this 21-class dataset with label-noise features. D ENS N ET-D recovers to 86.0%, a density lift of +4.3 points. The density branch is providing genuine topological signal that compensates for the backbone’s capacity-induced overfitting on this challenging multi-class task. How the anomaly cases should be read. The main ANCS table should therefore be interpreted within each dataset rather than as a monotone trend in k̄. Senate-Bills is a ceiling case: AllDeepSets already reaches 93.0 ± 3.8% under the strict recipe, so adding density features changes accuracy by only −0.3 points. Gene-Disease is different: the tier-wise ordering is not clean because HGNN remains competitive and HNHN is best, but the strict ablation still reveals that density features recover +4.3 points under the identical AllDeepSets recipe. These cases explain why the main text treats the hierarchy as a predictor of failure modes and headroom, not as a claim that accuracy should increase monotonically with mean hyperedge size. 6. Hierarchy scaling curve is flat. Running PDN with pattern libraries restricted to vmax ∈ {3, 4, 5, 6} on all three IWS benchmarks yields 100% at every level. On CE-Hard this is expected (any native architecture suffices). On Native-Hard, AllDeepSets alone reaches only 92.4%, so the density branch is contributing—but even the simplest patterns (vmax =3) suffice to bridge the gap. The curve is flat not because density is unused, but because the density signal needed for these benchmarks does not require complex patterns. Isolating the dependence on pattern complexity requires benchmarks with higher-width separating invariants; designing such benchmarks is an open problem. 7. Sinkhorn alignment is prohibitive at scale. I NV N ET PDN on Gene-Disease (n=5,012) achieves only 64.7% (single seed), far below every baseline including MLP. The Sinkhorn alignment branch requires O(n2 ) pairwise distance computation and iterative normalization; at this scale, the optimization landscape becomes intractable and the invariant features degenerate. This motivates the D ENS N ET-D design: density estimation via Monte Carlo sampling scales linearly in n, whereas Sinkhorn alignment does not. Takeaways. Five patterns emerge: (i) density features help on 6/7 datasets when a strict ablation (identical recipe, density on/off) is used, with the largest lift on the highest-k̄ datasets (GeneDisease +4.3%, House +2.2%); (ii) normalization choices must match the edge-size distribution, not be applied uniformly; (iii) fusion architecture matters—concatenation with a learned projection outperforms additive fusion because it decouples the backbone and density embedding spaces; (iv) the learned gate is a reliable diagnostic—when it stays at initialization (Gene-Disease) or drifts toward suppression (Cora, Citeseer), density features are uninformative; (v) Sinkhorn-based invariant features do not scale beyond n ≈ 5,000 nodes, favoring Monte Carlo density estimation for large datasets. These findings suggest that future density-aware models should adapt their normalization and fusion strategy to dataset structure automatically, a direction the DATASET P ROFILER partially addresses.

Q

Depth Scaling: Can Deeper Native Models Close the Gap?

A natural question is whether increasing L from 2 to 4 or 8 closes the Native-Hard gap in practice. The results reveal a clear pattern: depth degrades performance rather than improving it. HNHN and UniGNN maintain 100% at L=4 but collapse to ≈45–53% at L=8—below chance level (50% for binary classification), indicating the models learn degenerate representations that anti-correlate with the label. AllDeepSets and AllSetTransformer degrade monotonically from L=2. This is consistent with oversmoothing in hypergraph message passing: at L=8, node representations converge to a shared equilibrium, destroying the discriminative signal. The theoretical benefit of deeper architectures (access to higher-ghtw patterns) is overwhelmed by the practical cost of optimization collapse. 38

Table 14: Depth scaling on Native-Hard: accuracy (%) for native-tier models at L ∈ {2, 4, 8} layers. 10 seeds, 200 pairs, 300 epochs, hidden = 64. All models collapse at L=8; depth worsens performance rather than improving it. Model AllDeepSets HNHN UniGNN AllSetTransformer

L=2

L=4

L=8

84.6 ± 21.6 100.0 ± 0.0 100.0 ± 0.0 88.0 ± 19.9

47.0 ± 4.5 100.0 ± 0.0 100.0 ± 0.0 58.7 ± 18.3

46.8 ± 4.4 45.1 ± 2.3 52.5 ± 12.9 46.8 ± 4.4

The implication is architectural: to access invariants beyond the Width Wall, one cannot simply add layers. Density estimation (PDN, D ENS N ET-D) achieves 100% on Native-Hard at L=2 by accessing higher-order patterns through a fundamentally different mechanism—Monte Carlo pattern counting—that does not suffer from oversmoothing.

R

Gate Dynamics: Theory Predicts Model Behavior

The DATASET P ROFILER classifies each dataset by regime and sets the gate initialization accordingly (Table 11). A natural question is whether the learned gate confirms the profiler’s prediction—does training push the gate toward suppression on pairwise datasets, as the theory suggests? We investigate this with a controlled experiment: on Cora and Citeseer (both pairwise, k̄ < 4), we initialize the density gate at g0 = 0 (σ(g0 ) = 0.5, giving equal weight to local backbone and density branch) and let gradient descent decide the allocation. Figure 2(a) plots the gate trajectory σ(g) over training epochs for all six runs (three splits per dataset). On both datasets, every run drifts monotonically below 0.5 and converges to ≈0.497: the model consistently allocates less weight to the density branch. The drift is small in magnitude (≈0.003) but its direction is unambiguous—across all six runs, the gate never sustains a value above initialization. The Citeseer traces show a brief upward excursion in epochs 10–20 (the backbone has not yet converged, so the density branch temporarily provides a noisy gradient signal) before settling into the same downward trend. Figure 2(b) contrasts the final learned gate across regimes. Three distinct gate behaviors emerge: (i) Active (IWS higher-order witnesses: CE-Hard, Native-Hard, Multi-Order): gate stabilizes at σ(g) ≈ 0.11, maintaining meaningful density-branch contribution sufficient for 100% accuracy. (ii) Suppression (pairwise ANCS datasets: Cora, Citeseer): gate converges to ≈0.497—density features carry no discriminative signal beyond the backbone. (iii) Indifference (Gene-Disease, k̄=14): gate remains frozen at σ(g) = 0.5000 across all 10 seeds. Unlike suppression (where the gate drifts below 0.5), indifference means gradient descent finds no signal through the gate. Nonetheless, the strict ablation (Table 4) reveals a +4.3% density lift—the density features contribute through the concatenation pathway even when the gate does not actively modulate them. The theory-guided initialization (b = −5.0 for pairwise in production, Table 11) amplifies the natural suppression by starting the gate at σ(b) ≈ 0.007, ensuring the density branch does not interfere with backbone learning during early training.

39

Init. σ(0)=0.5 Cora (k̄=3.0) Citeseer (k̄=3.2)

0.5

0.5 Final σ(g)

Learned gate σ(g)

0.502

0.498 0.496

init0.5 σ(0) 0.5 init0.5 σ(0)

0.4 0.3 0.2

0.11 0.11 init0.11 σ(−2)

0.1 0.494

0

10

20

30 40 Epoch

50

60

70

(a) Gate trajectories on pairwise datasets (neutral init). All six runs (3 splits × 2 datasets; light traces show additional seeds) drift below 0.5: gradient descent suppresses the density branch.

0

-H t-H -O ora Cite ne-D CE Na M C Ge (b) Final learned gate σ(g) across regimes (±1σ). Active (left 3): ∼11%. Suppressed (Cora/Cite): ∼50%. Indifferent (Gene-D): frozen at 0.500.

Figure 2: Gate dynamics validate the profiler’s regime prediction. (a) On pairwise datasets (k̄ < 4), starting from a neutral initialization σ(0) = 0.5, the density gate drifts monotonically below 0.5: gradient descent finds no useful density signal and suppresses the branch. (b) Three gate regimes emerge: active (σ(g) ≈ 0.11, IWS higher-order witnesses—density contributes), suppression (σ(g) ≈ 0.497, Cora/Citeseer—density uninformative), and indifference (σ(g) = 0.500, GeneDisease—zero gradient signal, node features dominate). Different initializations reflect the profiler’s theory-guided configuration (dashed lines); the gate remains near or below initialization except when density features carry genuine discriminative signal.

40