利用图on与神经逆算子构建复杂网络的生成模型

# A Generative Model of Complex Networks Using Graphons and Neural Inverse Operators
Source: [https://arxiv.org/html/2610.02439](https://arxiv.org/html/2610.02439)
Wooseong Choi††thanks:Equal contribution\.Affiliation:Department of Biomedical EngineeringAffiliation:University of Southern CaliforniaAffiliation:Los Angeles, CA 90089Email:[wooseong@usc\.edu](mailto:)Italo’Ivo Lima Dias Pinto11footnotemark:1Affiliation:Department of Electrical EngineeringAffiliation:University of Southern CaliforniaAffiliation:Los Angeles, CA 90089Email:[italoivo@gmail\.com](mailto:)Chen SunAffiliation:Department of Biomedical EngineeringAffiliation:University of Southern CaliforniaAffiliation:Los Angeles, CA 90089Email:[csun8109@usc\.edu](mailto:)Gaurav Gupta††thanks:Work done outside of the author’s role at Amazon\.Dong Song††thanks:Co\-corresponding authors\.Affiliation:Department of Biomedical EngineeringAffiliation:University of Southern CaliforniaAffiliation:Los Angeles, CA 90089Email:[dsong@usc\.edu](mailto:)Paul Bogdan33footnotemark:3Affiliation:Department of Electrical EngineeringAffiliation:University of Southern CaliforniaAffiliation:Los Angeles, CA 90089Email:[pbogdan@usc\.edu](mailto:)

###### Abstract

Generative graph models are central to understanding and simulating complex networks\. However, existing approaches have complementary strengths and limitations\. Mechanistic models offer interpretability but rely on instance\-specific estimation methods\. Deep generative models, on the other hand, offer amortized inference at the cost of interpretability and are largely limited to graph sizes seen during training\. Scientific applications motivate a framework that retains the strengths of both paradigms\. We bridge them by formulating both the generative model and parameter recovery in function space\. A multifractal step graphon extends standard step graphons with a recursive construction that compactly parameterizes complex networks\. This formulation admits a neural inverse operator to recover its parameters, enabling inference on unseen graph sizes\. We evaluate our model, trained only on synthetic multifractal step graphon realizations, against both paradigms\. Against a graph foundation model pretrained on empirical networks, our method achieves the best average performance on three of four metrics in a zero\-shot graph\-generation benchmark, indicating that the model transfers to real\-world graphs\. We also apply our method to single\-observation networks, a regime largely inaccessible to deep models that require training corpora, where it performs comparably to an instance\-specific method that optimizes on each graph\. In a multi\-subject EEG case study, the inferred parameters track a reversible change in brain state more sensitively than traditional network statistics\. Together, these results indicate that mechanistic interpretability and amortized inference can be effectively unified in a generative graph model to enhance our understanding of complex networks\.

## 1Introduction

Complex networks model biological\([La Rocca et al\., 2018](https://arxiv.org/html/2610.02439#bib.bib32);[Yin et al\., 2020](https://arxiv.org/html/2610.02439#bib.bib33);[Yang et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib5)\), social\([Lee and Jung, 2006](https://arxiv.org/html/2610.02439#bib.bib48);[Liu et al\., 2017](https://arxiv.org/html/2610.02439#bib.bib47)\), and infrastructure systems\([Long and Chen, 2021](https://arxiv.org/html/2610.02439#bib.bib34);[Batty and Longley, 1994](https://arxiv.org/html/2610.02439#bib.bib35)\)\. Summary statistics of such networks aggregate features across scales and can obscure differences in the underlying structure\([Xue and Bogdan, 2017](https://arxiv.org/html/2610.02439#bib.bib4)\)\. In addition, these metrics can vary spuriously across different partial observations of the same network\([Lee et al\., 2006](https://arxiv.org/html/2610.02439#bib.bib49);[Stumpf et al\., 2005](https://arxiv.org/html/2610.02439#bib.bib50);[Achlioptas et al\., 2005](https://arxiv.org/html/2610.02439#bib.bib51)\)\. This motivates going beyond descriptive statistics and recovering the compact generative mechanisms that reproduce the observed network\. Such mechanisms are particularly valuable for revealing hidden differences underlying complex networks\([Yang and Bogdan, 2020](https://arxiv.org/html/2610.02439#bib.bib52);[Akarca et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib53)\)\. As such, recovering these generative processes is essential for enhancing our understanding of complex systems\.

Generative graph models synthesize realistic networks by capturing meaningful representations of their underlying structure\. Two paradigms have emerged, each with its own approach\. Mechanistic models define explicit wiring rules that are algorithmically implemented to produce graphs with prescribed properties\([Erdős and Rényi, 1960](https://arxiv.org/html/2610.02439#bib.bib1);[Barabási and Albert, 1999](https://arxiv.org/html/2610.02439#bib.bib2)\)\. In contrast, deep generative models learn graph distributions directly from data\([Li et al\., 2018](https://arxiv.org/html/2610.02439#bib.bib29);[Kipf and Welling, 2016](https://arxiv.org/html/2610.02439#bib.bib31)\)\.

Classical mechanistic models capture specific regimes, from homogeneous random graphs\([Erdős and Rényi, 1960](https://arxiv.org/html/2610.02439#bib.bib1)\)to small\-world\([Watts and Strogatz, 1998](https://arxiv.org/html/2610.02439#bib.bib3)\)and preferential\-attachment networks\([Barabási and Albert, 1999](https://arxiv.org/html/2610.02439#bib.bib2)\)\. However, these frameworks are restricted to modeling a prescribed class of networks\. To overcome these constraints, recursive models such as SKG\([Leskovec et al\., 2010](https://arxiv.org/html/2610.02439#bib.bib14)\), multifractal network generators\([Palla et al\., 2010](https://arxiv.org/html/2610.02439#bib.bib13)\), and WMGM\([Yang et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib5)\)established the foundation for using recursive rules to compactly describe networks observed across many domains\.

However, recovering these models from data remains challenging\. SKGs\([Leskovec et al\., 2010](https://arxiv.org/html/2610.02439#bib.bib14)\)rely on maximum likelihood estimation \(MLE\) with permutation sampling\. Multifractal network generators\([Palla et al\., 2010](https://arxiv.org/html/2610.02439#bib.bib13)\)use simulated annealing, while the weighted multifractal generative model \(WMGM\)\([Yang et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib5)\)relies on expectation\-maximization \(EM\)\. MLE and EM are susceptible to local optima\([Wu, 1983](https://arxiv.org/html/2610.02439#bib.bib8);[Blei et al\., 2017](https://arxiv.org/html/2610.02439#bib.bib9)\), whereas simulated annealing suffers from prohibitively slow convergence\([Hajek, 1988](https://arxiv.org/html/2610.02439#bib.bib26)\)\. These methods also search over latent node\-to\-community assignments, a combinatorial step whose cost grows rapidly with graph size\([Leskovec et al\., 2010](https://arxiv.org/html/2610.02439#bib.bib14);[Daudin et al\., 2008](https://arxiv.org/html/2610.02439#bib.bib55);[Brandes et al\., 2008](https://arxiv.org/html/2610.02439#bib.bib54)\)\. Moreover, these methods are instance\-specific, so the aforementioned costs recur with every observation\.

Deep generative graph models\([Simonovsky and Komodakis, 2018](https://arxiv.org/html/2610.02439#bib.bib28);[You et al\., 2018](https://arxiv.org/html/2610.02439#bib.bib23);[Liao et al\., 2019](https://arxiv.org/html/2610.02439#bib.bib24);[Vignac et al\., 2022](https://arxiv.org/html/2610.02439#bib.bib27)\)amortize generation from data, but their learned representations are often less interpretable and typically do not extrapolate to graph sizes outside the training distribution\. Autoregressive models\([Li et al\., 2018](https://arxiv.org/html/2610.02439#bib.bib29);[You et al\., 2018](https://arxiv.org/html/2610.02439#bib.bib23);[Liao et al\., 2019](https://arxiv.org/html/2610.02439#bib.bib24)\)cannot extrapolate beyond the training distribution due to architectural constraints \(e\.g\., GraphRNN’s BFS window and GRAN’s maximum node capacity\)\. Variational autoencoders and diffusion\-based models similarly generate graphs whose sizes are effectively limited to those observed during training\([Simonovsky and Komodakis, 2018](https://arxiv.org/html/2610.02439#bib.bib28);[Niu et al\., 2020](https://arxiv.org/html/2610.02439#bib.bib25);[Vignac et al\., 2022](https://arxiv.org/html/2610.02439#bib.bib27)\)\.

This work\.Taken together, there is a strong need for a framework that retains the strengths of both paradigms, such as interpretability and amortized inference, which are jointly necessary for generative models to serve in scientific applications\. Generalization to unseen graph sizes further broadens its applicability\. We address this need by formulating both the generative model and its inference in function space:

- •Compact Function\-Space Representation of Complex Networks:We introduce the multifractal step graphon: a recursively\-refined, piecewise\-constant kernel parameterized by a small set of base coefficientsθ\(1\)=\(l\(1\),p\(1\)\)\\theta^\{\(1\)\}=\(l^\{\(1\)\},p^\{\(1\)\}\)that directly encode hierarchical community structure\.
- •Neural Inverse Operator Formulation:We cast the parametric graphon recovery as an inverse operator\-learning problem\. This formulation enables amortized inference that generalizes to unseen graph sizes\.

![Refer to caption](https://arxiv.org/html/2610.02439v1/mwnot_fig1_v3.png)Figure 1:Discovery of Hidden Structures in Real\-World Networks\(a\)The base graphon encodes hidden community structures underlying complex networks\.\(b\)Kronecker expansion of base parameters produces the multifractal step graphon\.\(c\)Finite realization of the graphon produces networks that share similar properties with those seen in various empirical domains222Figure[2](https://arxiv.org/html/2610.02439#footnote2)incorporates a brain tractography image by[Yeh et al\. \(2018\)](https://arxiv.org/html/2610.02439#bib.bib46)\[[CC BY\-SA 4\.0](https://creativecommons.org/licenses/by-sa/4.0/)\] and a protein interaction network by[Chen et al\. \(2019\)](https://arxiv.org/html/2610.02439#bib.bib45)\[[CC BY 4\.0](https://creativecommons.org/licenses/by/4.0/)\] via Wikimedia Commons\. The citation network by citationgraph\.org is Public Domain \(CC0\)\. This derived schematic is distributed under CC BY\-SA 4\.0\.\. Recovering the generative rule reveals\(d\)hidden differences between networks with the same statistics and\(e\)the shared latent structure of a network from different partial observations, where subsampling can produce different statistics\.
## 2Related Work

##### Graphons\.

A graphon is a bounded, measurable functionW:\[0,1\]2→\[0,1\]W:\[0,1\]^\{2\}\\to\[0,1\]encoding connection probabilities between latent node positions\([Lovász and Szegedy, 2006](https://arxiv.org/html/2610.02439#bib.bib17);[Borgs et al\., 2008](https://arxiv.org/html/2610.02439#bib.bib18)\)\. Step graphons approximate this continuous space as piecewise\-constant surfaces, admitting an interpretable parametrization in which latent communities and connection rules correspond to domain partitions and uniform probabilities within each block, respectively\([Airoldi et al\., 2013](https://arxiv.org/html/2610.02439#bib.bib20)\)\. Prior work on graphon estimation has been limited to instance\-specific statistical methods\([Airoldi et al\., 2013](https://arxiv.org/html/2610.02439#bib.bib20);[Chan and Airoldi, 2014](https://arxiv.org/html/2610.02439#bib.bib30)\)\. Moreover, standard step graphons implicitly assume non\-hierarchical topologies and cannot represent the recursive, self\-similar structure present in many real\-life networks\([Xue and Bogdan, 2017](https://arxiv.org/html/2610.02439#bib.bib4)\)\.

##### Neural Operators\.

Recently, neural operators have emerged as powerful frameworks for learning the mapping between function spaces, with applications primarily in solving forward partial differential equations\([Lu et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib16);[Kovachki et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib6);[Li et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib7);[Gupta et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib10);[Xiao et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib11)\)\.*Inverse*operator learning, in which the network recovers underlying parameters from functions, has been applied to PDE inverse problems\([Molinaro et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib36)\)\. However, no other work has applied it directly to graphon estimation to the best of our knowledge\.

##### Graph Generation\.

The complementary limitations of mechanistic and deep generative graph models are discussed in Section 1\. Existing graphon\-based generative models offer either interpretability, through a linear factorization over learned graphon factors in a graphon autoencoder \(GNAE\)\([Xu et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib56)\)or flexibility across graph sizes by representing graphons as implicit neural representations at arbitrary resolution\([Xia et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib57);[Azizpour et al\., 2024](https://arxiv.org/html/2610.02439#bib.bib58)\)\. Our approach provides both\. The method’s parameters define an explicit wiring rule with a closed\-form link to network properties \(Appendix[B](https://arxiv.org/html/2610.02439#A2)\), whereas GNAE’s interpretability lies in its learned latent factors\. Its function\-space formulation further admits a neural inverse operator that generalizes across graph sizes\.

## 3Multifractal Step Graphon Formulation

### 3\.1Overview

Here, we extend step graphons\([Airoldi et al\., 2013](https://arxiv.org/html/2610.02439#bib.bib20)\)to better represent hierarchical structure observed in real\-life networks\([Xue and Bogdan, 2017](https://arxiv.org/html/2610.02439#bib.bib4)\)\. We formalize this through a recursively refined piecewise\-constant kernel on the unit square, which gives a continuum model parameterized byθ\(1\)=\(l\(1\),p\(1\)\)\\theta^\{\(1\)\}=\(l^\{\(1\)\},p^\{\(1\)\}\)\. The parametersl\(1\)l^\{\(1\)\}andp\(1\)p^\{\(1\)\}yield heterogeneous node attribution and weighted connectivity, respectively\.MMdefines the number of base communities\. A Kronecker recursion expands these base parameters to arbitrary depthKK, partitioning the latent space intoMKM^\{K\}communities whose sizes and interaction strengths factorize multiplicatively across scales\. We useM=2M=2throughout, which keeps theMKM^\{K\}communities populated in tractable training graphs\. Because two recursion steps together form a validM=4M=4base rule, this choice still covers the subset of larger\-MMmodels whose base rule is Kronecker\-factorizable \(Appendix[A](https://arxiv.org/html/2610.02439#A1)\)\.

We model generation of graphs from the multifractal step graphon as realizations of a jointly exchangeable array\. By the Aldous representation theorem\([Aldous, 1981](https://arxiv.org/html/2610.02439#bib.bib19)\), any such array admits a generative representation as a measurable functionW:\[0,1\]2→\[0,1\]W:\[0,1\]^\{2\}\\to\[0,1\]from which finite graphs are sampled\. This provides the foundation for graphon\-based generative modeling\. Thus, our use of the graphon follows the Aldous tradition—as a generative kernel from which finite networks are sampled—rather than the Lovász–Szegedy tradition of graph limits\.

#### 3\.1\.1Continuous Latent Space

A graphon is a bounded, Lebesgue\-measurable functionW:\[0,1\]2→\[0,1\]W:\[0,1\]^\{2\}\\to\[0,1\]belonging to the Hilbert spaceL2​\(\[0,1\]2\)L^\{2\}\(\[0,1\]^\{2\}\)\([Lovász and Szegedy, 2006](https://arxiv.org/html/2610.02439#bib.bib17);[Borgs et al\., 2008](https://arxiv.org/html/2610.02439#bib.bib18)\)\. While traditionally defined as symmetric \(W⁡\(x,y\)=W⁡\(y,x\)W\(x,y\)=W\(y,x\)\) for undirected graphs, we follow theoretical\([Aldous, 1981](https://arxiv.org/html/2610.02439#bib.bib19)\)and empirical\([Airoldi et al\., 2013](https://arxiv.org/html/2610.02439#bib.bib20)\)precedents and generalize to asymmetric kernels \(W⁡\(x,y\)≠W⁡\(y,x\)W\(x,y\)\\neq W\(y,x\)\) to capture directed topologies\. In this formulation,W⁡\(x,y\)W\(x,y\)represents the probability of a directed edge from a node at latent coordinatexxto a node atyy\.

#### 3\.1\.2Step Graphon Parameterization

A step graphon approximates this continuous space as a piecewise\-constant surface\. We partition the unit interval\[0,1\]\[0,1\]intoMMdisjoint subsets\{I1,…,IM\}\\\{I\_\{1\},\\dots,I\_\{M\}\\\}\. A probability vectorl∈ΔM−1l\\in\\Delta^\{M\-1\}defines the length of each subset such thatλ⁡\(Ii\)=li\\lambda\(I\_\{i\}\)=l\_\{i\}\. The connection probability between subsets is defined by a blockwise affinity matrixp∈\[0,1\]M×Mp\\in\[0,1\]^\{M\\times M\}\. The step graphon is then:

W⁡\(x,y\)=∑i=1M∑j=1Mpi​j​𝕀Ii​\(x\)​𝕀Ij​\(y\)W\(x,y\)=\\sum\_\{i=1\}^\{M\}\\sum\_\{j=1\}^\{M\}p\_\{ij\}\\,\\mathbb\{I\}\_\{I\_\{i\}\}\(x\)\\,\\mathbb\{I\}\_\{I\_\{j\}\}\(y\)\(1\)where𝕀\\mathbb\{I\}is the indicator function\.

#### 3\.1\.3Multifractal Connectivity via Kronecker Recursion

Standard step graphons implicitly assume non\-hierarchical topologies\. However, real\-world networks exhibit varying levels of heterogeneity across different scales\([Xue and Bogdan, 2017](https://arxiv.org/html/2610.02439#bib.bib4)\)\. This heterogeneity is formalized by the Lipschitz–Hölder exponentα\\alpha, which measures how the local connectivity scales with resolution\. A*monofractal*object has the same scaling behavior everywhere \(a singleα\\alpha\), while a*multifractal*object exhibits a distribution of scaling behaviors, characterized by the spectrumf⁡\(α\)f\(\\alpha\)as its width \(Δ​α=αmax−αmin\\Delta\\alpha=\\alpha\_\{\\max\}\-\\alpha\_\{\\min\}\) and shape\.

To capture this self\-similar structure, we recursively expand the base parametersθ\(1\)=\(l\(1\),p\(1\)\)\\theta^\{\(1\)\}=\(l^\{\(1\)\},p^\{\(1\)\}\), wherel\(1\)∈int​\(ΔM−1\)l^\{\(1\)\}\\in\\text\{int\}\(\\Delta^\{M\-1\}\)andp\(1\)∈\(0,1\)M×Mp^\{\(1\)\}\\in\(0,1\)^\{M\\times M\}, to a multiscale depthKKvia the Kronecker product:

l\(K\)=l\(1\)⊗⋯⊗l\(1\)⏟K​times,p\(K\)=p\(1\)⊗⋯⊗p\(1\)⏟K​timesl^\{\(K\)\}=\\underbrace\{l^\{\(1\)\}\\otimes\\cdots\\otimes l^\{\(1\)\}\}\_\{K\\text\{ times\}\},\\quad p^\{\(K\)\}=\\underbrace\{p^\{\(1\)\}\\otimes\\cdots\\otimes p^\{\(1\)\}\}\_\{K\\text\{ times\}\}\(2\)This recursive expansion partitions the space intoMKM^\{K\}contiguous sub\-intervals indexed byq=\(q1,…,qK\)q=\(q\_\{1\},\\ldots,q\_\{K\}\)\. The multifractal step graphon is defined as:

WK​\(x,y\)=∑q∑hpq​h\(K\)​𝕀Iq​\(x\)​𝕀Ih​\(y\)W\_\{K\}\(x,y\)=\\sum\_\{q\}\\sum\_\{h\}p^\{\(K\)\}\_\{qh\}\\,\\mathbb\{I\}\_\{I\_\{q\}\}\(x\)\\,\\mathbb\{I\}\_\{I\_\{h\}\}\(y\)\(3\)where the interaction amplitude factorizes across scales aspq​h\(K\)=∏k=1Kpqk​hk\(1\)p^\{\(K\)\}\_\{qh\}=\\prod\_\{k=1\}^\{K\}p^\{\(1\)\}\_\{q\_\{k\}h\_\{k\}\}\. This multiplicative factorization enables varying levels of heterogeneity across scales:log⁡pq​h\(K\)\\log p^\{\(K\)\}\_\{qh\}is a sum ofKKterms fromlog⁡p\(1\)\\log p^\{\(1\)\}, so different index paths yield different local scaling exponents\. Uniform entries inp\(1\)p^\{\(1\)\}collapse to a single exponent \(Δ​α=0\\Delta\\alpha=0, monofractal\), while heterogeneous entries widen the spectrum\. Thus,θ\(1\)\\theta^\{\(1\)\}compactly controls whereWKW\_\{K\}lies on the monofractal–multifractal continuum\. The Kronecker structure also yields closed\-form expressions for several network statistics, such as the multifractal spectrum, derived in Appendix[B](https://arxiv.org/html/2610.02439#A2)\.

![Refer to caption](https://arxiv.org/html/2610.02439v1/mwnot_fig2.png)Figure 2:Neural inverse operator architecture\.\(a\)The input adjacency matrix is decomposed acrossLLscales by a 2D multiwavelet transform, yielding hierarchical coefficients that are tokenized and processed by a branch network \(Keys and Values\)\.\(b\)Learnable queries forℓ\(1\)\\ell^\{\(1\)\}andp\(1\)p^\{\(1\)\}pass through a trunk network and are combined with the branch via cross\-attention\. Scale\-specific outputs are decoded by MLPs to predict the base parameters\.
#### 3\.1\.4Sampling Networks from the Multifractal Step Graphon

A finite graphG=\(V,E\)G=\(V,E\)is constructed by sampling from the multifractal step graphonWKW\_\{K\}\. By the Aldous representation, each nodeuuis associated with a latent coordinatexu∼Uniform​\(0,1\)x\_\{u\}\\sim\\text\{Uniform\}\(0,1\), which falls in the intervalIquI\_\{q\_\{u\}\}with probabilitylqu\(K\)l^\{\(K\)\}\_\{q\_\{u\}\}\. BecauseWKW\_\{K\}is piecewise\-constant on these intervals, this is equivalent to drawing community assignments directly:

qu∼Multinomial​\(l\(K\)\),q\_\{u\}\\sim\\text\{Multinomial\}\\bigl\(l^\{\(K\)\}\\bigr\),\(4\)which is the form used in implementation\. The continuous interaction intensity for pair\(u,v\)\(u,v\)is then the constant value ofWKW\_\{K\}on the blockIqu×IqvI\_\{q\_\{u\}\}\\times I\_\{q\_\{v\}\}:

Wu​v=WK​\(xu,xv\)=pqu​qv\(K\)\.W\_\{uv\}=W\_\{K\}\(x\_\{u\},x\_\{v\}\)=p^\{\(K\)\}\_\{q\_\{u\}q\_\{v\}\}\.\(5\)To map this continuous amplitude to a discrete edge weightru​v∈ℕ0r\_\{uv\}\\in\\mathbb\{N\}\_\{0\}, we sample from a conditional geometric distribution:

P⁡\(ru​v=r∣Wu​v\)=\(1−Wu​v\)​\(Wu​v\)r\.P\(r\_\{uv\}=r\\mid W\_\{uv\}\)=\(1\-W\_\{uv\}\)\(W\_\{uv\}\)^\{r\}\.\(6\)Higher connection probabilities naturally produce heavier\-tailed edge weights while maintaining support on all non\-negative integers\.

#### 3\.1\.5Identifiability

Inverse problems are typically ill\-posed\([Molinaro et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib36)\)\. We show identifiability \(Proposition 1\) of our formulation, which ensures that the inverse map is well\-defined and distinct parameters produce distinct graphons\.

###### Proposition 1\(Identifiability of the multifractal step graphon\)\.

Assume additionally that the expanded block profiles ofp\(K\)p^\{\(K\)\}are pairwise distinct, so that no two blocks have identical incoming and outgoing interaction profiles\. Under this generic non\-degeneracy condition, letM≥2M\\geq 2andK≥2K\\geq 2\. Consider two base parameterizationsθ1=\(l1\(1\),p1\(1\)\)\\theta\_\{1\}=\(l\_\{1\}^\{\(1\)\},p\_\{1\}^\{\(1\)\}\)andθ2=\(l2\(1\),p2\(1\)\)\\theta\_\{2\}=\(l\_\{2\}^\{\(1\)\},p\_\{2\}^\{\(1\)\}\), whereli∈int⁡\(ΔM−1\)l\_\{i\}\\in\\mathrm\{int\}\(\\Delta^\{M\-1\}\)andpi∈\(0,1\)M×Mp\_\{i\}\\in\(0,1\)^\{M\\times M\}\. LetWK\(1\)W\_\{K\}^\{\(1\)\}andWK\(2\)W\_\{K\}^\{\(2\)\}denote the corresponding multifractal step graphons obtained viaKK\-fold Kronecker expansion\. IfWK\(1\)=WK\(2\)W\_\{K\}^\{\(1\)\}=W\_\{K\}^\{\(2\)\}almost everywhere on\[0,1\]2\[0,1\]^\{2\}, then there exists a permutationσ∈SM\\sigma\\in S\_\{M\}such thatl1\(1\)=σ⁡\(l2\(1\)\)l\_\{1\}^\{\(1\)\}=\\sigma\(l\_\{2\}^\{\(1\)\}\)andp1\(1\)=σ⁡\(p2\(1\)\)p\_\{1\}^\{\(1\)\}=\\sigma\(p\_\{2\}^\{\(1\)\}\), whereσ\\sigmaacts by reindexing:σ​\(l\)i=lσ⁡\(i\)\\sigma\(l\)\_\{i\}=l\_\{\\sigma\(i\)\}andσ​\(p\)i​j=pσ⁡\(i\)​σ​\(j\)\\sigma\(p\)\_\{ij\}=p\_\{\\sigma\(i\)\\sigma\(j\)\}\.

The proof, given in Appendix[C](https://arxiv.org/html/2610.02439#A3), proceeds in three steps: \(i\) the piecewise\-constant graphon determines its block structure up to permutation via its discontinuities\. \(ii\) the Kronecker factorization of the expanded parameters into base parameters is unique, by rank\-1 symmetric tensor decomposition\. \(iii\) the coupling betweenl\(K\)l^\{\(K\)\}andp\(K\)p^\{\(K\)\}through shared block indices forces the permutation to be consistent across both\.

## 4Inference via Neural Operators

Recoveringθ\(1\)=\(l\(1\),p\(1\)\)\\theta^\{\(1\)\}=\(l^\{\(1\)\},p^\{\(1\)\}\)from an observed graph requires learning the inverse operatorℱ:WK↦θ\(1\)\\mathcal\{F\}:W\_\{K\}\\mapsto\\theta^\{\(1\)\}\. The observed adjacency matrix is treated as a finite\-sample realization of a graph generated fromWK∈L2​\(\[0,1\]2\)W\_\{K\}\\in L^\{2\}\(\[0,1\]^\{2\}\)\. BecauseW1W\_\{1\}is fully parameterized byθ\(1\)\\theta^\{\(1\)\}, recovering these parameters fully specifies the multiscale graphonWKW\_\{K\}\.

![Refer to caption](https://arxiv.org/html/2610.02439v1/expressivity_panel.png)Figure 3:Compact Parameterization of Diverse Topologies\.\(a\)Degree distributions span Poisson\-like, log\-normal\-like, power\-law\-like \(with exponential\-tail\), and bimodal regimes that are present in real\-life networks\.\(b\)Sweepingwstructw\_\{\\mathrm\{struct\}\}allows expression of varying levels of assortativity and disassortativity\.\(c\)wstructw\_\{\\mathrm\{struct\}\}also controls the full range of modularity \(left axis\) and spectral gapλ2\\lambda\_\{2\}\(right axis\)\.\(d\)Multifractal spectra forms a continuous envelope from monofractal to multifractal\.Multifractal step graphons are low\-regularity, piecewise\-constant objects, for which multiwavelet bases are known to provide efficient representations\([Alpert, 1993](https://arxiv.org/html/2610.02439#bib.bib21)\)\. We verify empirically that the multiwavelet basis resolves the recursive structure ofWKW\_\{K\}: the proportion of energy in the detail\-coefficients grows exponentially with recursion depth \(Appendix[E](https://arxiv.org/html/2610.02439#A5)\)\. We instantiate the inverse operator as a Multiwavelet Neural Operator Transformer \(MWNOT\) \(Fig\.[2](https://arxiv.org/html/2610.02439#S3.F2)\), combining a 2D multiwavelet decomposition\([Gupta et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib10);[Xiao et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib11)\)with cross\-attention for operator approximation\([Shih et al\., 2025](https://arxiv.org/html/2610.02439#bib.bib12)\)\. We use MWNOT as one architecture compatible with the proposed formulation, not as a claim of architectural optimality\.

### 4\.1Architecture

First, a 2D discrete multiwavelet transform\([Gupta et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib10)\)decomposesAAacrossLLspatial scales\. The multiwavelet basis is constructed from shifted Legendre polynomials of orderkkand orthogonalized via Gram–Schmidt so that detail coefficients at each scale are orthogonal to the coarse approximation\. At each scaleℓ\\ell, the transform yields coarse coefficientss\(ℓ\)s^\{\(\\ell\)\}capturing macroscopic connectivity and detail coefficientsd\(ℓ\)d^\{\(\\ell\)\}capturing high\-frequency variation\.

Second, we adopt a branch–trunk cross\-attention architecture\([Lu et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib16);[Shih et al\., 2025](https://arxiv.org/html/2610.02439#bib.bib12)\)to map these multiscale coefficients toθ\(1\)\\theta^\{\(1\)\}\. The branch network encodes the input function through its wavelet tokens, producing keys and values\. The trunk network processes learnable query embeddings forl\(1\)l^\{\(1\)\}andp\(1\)p^\{\(1\)\}\. Cross\-attention approximates the operator through a linear combination of input\-derived features weighted by target\-dependent queries\([Shih et al\., 2025](https://arxiv.org/html/2610.02439#bib.bib12)\)\.

Third, scale\-specific attention outputs are concatenated acrossℓ\\elland decoded by MLPs to predictl^\(1\)\\hat\{l\}^\{\(1\)\}andp^\(1\)\\hat\{p\}^\{\(1\)\}\. Full architectural details are provided in Appendix[D](https://arxiv.org/html/2610.02439#A4)\.

### 4\.2Training

Training details, including data generation, loss functions, optimization schedule, and hyperparameter settings, are provided in Appendix[F](https://arxiv.org/html/2610.02439#A6)\.

## 5Experiments

The central claim of this work is that the multifractal step graphon, paired with a neural inverse operator trained only on its synthetic realizations, recovers interpretable parameters of real\-world networks through amortized inference that generalizes across graph sizes\. We support this claim through a series of experiments evaluating expressivity \(Sec\.[5\.1](https://arxiv.org/html/2610.02439#S5.SS1)\), parameter recovery on unseen graph sizes \(Sec\.[5\.2](https://arxiv.org/html/2610.02439#S5.SS2)\), and zero\-shot transfer to empirical data of various domains and sizes \(Sec\.[5\.3\.1](https://arxiv.org/html/2610.02439#S5.SS3.SSS1)\)\. Finally, we apply our framework to single\-observation empirical networks \(Sec\.[5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2)\), where the inferred parameters reconstruct each network comparably to per\-instance fitting and track a reversible change in brain state in a multi\-subject EEG case study\.

### 5\.1Compact Parameterization of Diverse Topologies

We empirically validate that the multifractal step graphon family covers various topologies observed in real\-world networks\. We sweep two axes of𝐩\(1\)\\mathbf\{p\}^\{\(1\)\}: a heterogeneity parameterwhet∈\[0,1\]w\_\{\\mathrm\{het\}\}\\in\[0,1\]controlling block\-density variation, and a structure parameterwstruct∈\[−1,1\]w\_\{\\mathrm\{struct\}\}\\in\[\-1,1\]interpolating from disassortative to assortative community organization\. Sweep details are in Appendix[G](https://arxiv.org/html/2610.02439#A7)\.

Figure[3](https://arxiv.org/html/2610.02439#S4.F3)confirms that varying𝜽\(1\)\\boldsymbol\{\\theta\}^\{\(1\)\}alone produces networks with several qualitatively distinct regimes of \(a\) degree distributions\([Broido and Clauset, 2019](https://arxiv.org/html/2610.02439#bib.bib37)\), \(b\) assortativity, \(c\) modularity with inverse relationship to spectral gap, and \(d\) multifractal spectrum width, spanning from monofractal to multifractal\. Together, these results suggest that the model provides a broad synthetic pretraining distribution that may enable transfer to real\-world networks in various domains\.

Figure 4:Model Accuracy on Synthetic Graphs\.\(a\)The model generalizes to unseen graph sizes and\(b\)shows stable performance across varying levels of hierarchical structure\.
### 5\.2Zero\-Shot Transfer to Unseen Size

We next verify that the neural inverse operator recovers𝜽\(1\)\\boldsymbol\{\\theta\}^\{\(1\)\}accurately across the conditions encountered at deployment: varying graph sizesNNand recursive depthsKK\. A single MWNOT is trained onN∈\{128,512,1024\}N\\in\\\{128,512,1024\\\}andK∈\{2,4,6\}K\\in\\\{2,4,6\\\}, then evaluated on 100 held\-out test samples per setting where performance is averaged acrossK∈\{2,…,6\}K\\in\\\{2,\\ldots,6\\\}for eachNN\. We test our model’s ability to generalize to unseen graph sizes\. MAE for𝐩\(1\)\\mathbf\{p\}^\{\(1\)\}and𝐥\(1\)\\mathbf\{l\}^\{\(1\)\}remains stable as size is changed fromN=128N=128toN=8192N=8192\(Fig\.[4](https://arxiv.org/html/2610.02439#S5.F4)a\)\. Consequently, our model allows the user to train on small graphs for inference on larger ones, which may enhance the broad applicability of the model\. Recovery remains stable across hierarchical depths as well\.

Table 1:Zero\-shot Performance ComparisonDEG, CC, Spec, Orb: MMD of degree, clustering coefficient, Laplacian eigenvalues and orbit counts\. Ours and Ours\-Mono are means over 10 seeds\. Full table with SDs in Table[S1](https://arxiv.org/html/2610.02439#A8.T1)\)\. On ECO, the orbit kernel saturates, giving the same value for any set of distinct dense samples\. This cell is not used to mark the best method\.

### 5\.3Zero\-Shot Transfer to Empirical Data

#### 5\.3\.1Graph Generation Benchmark

The preceding sections show that the model spans diverse topologies \(Sec\.[5\.1](https://arxiv.org/html/2610.02439#S5.SS1)\) and that its parameters are recovered accurately in graph sizes not seen during training \(Sec\.[5\.2](https://arxiv.org/html/2610.02439#S5.SS2)\)\. We now test whether these combine to transfer to real\-world graphs, using a model trained only on synthetic graphon realizations\. We use the zero\-shot graph\-generation benchmark and baselines from[Wang et al\. \(2025a\)](https://arxiv.org/html/2610.02439#bib.bib38)\. LGGM\-X is pretrained on over5,0005\{,\}000empirical graphs, with the target domainXXheld out\. We use the same code, test sets, and metrics, and report baseline results from[Wang et al\. \(2025a\)](https://arxiv.org/html/2610.02439#bib.bib38)\. The complete protocol is in Appendix[H](https://arxiv.org/html/2610.02439#A8)\. As an ablation,Ours\-Monois trained on graphons withp\(1\)≡0\.5p^\{\(1\)\}\\equiv 0\.5, removing the topological diversity of the prior \(Sec\.[5\.1](https://arxiv.org/html/2610.02439#S5.SS1)\)\.

Our model achieves the best average performance on three of four metrics \(“ALL” in Table[1](https://arxiv.org/html/2610.02439#S5.T1)\) and wins in all metrics in 4/12 domains\. These results are notable given the contrasting training paradigms\. LGGM\-X is trained on thousands of real\-world networks across many domains, whereas our model is trainedexclusivelyon synthetic realizations of the multifractal step graphon\. In addition, our method uses∼\\sim22×\\timesfewer parameters \(∼\\sim205K vs\. 4\.61M\)\. LGGM\-X achieves stronger performance on some domains, such as infrastructure networks, where our method may struggle as it was trained on synthetic graphs without geometric priors\. Averaged across domains, the ablation increases MMD on all four metrics, indicating that the topological diversity of the prior contributes to the transfer\. Taken together, the results suggest that the generative rule and the inference of its parameters transfer to real\-world networks in various domains\.

#### 5\.3\.2Applications to Single\-Observation Empirical Data

Recent advances in connectomics have produced datasets with unprecedented scale, such as the complete connectome ofDrosophila\([Dorkenwald et al\., 2024](https://arxiv.org/html/2610.02439#bib.bib39)\)\. However, these data often exist as a single observation\. Whole\-graph generation on these datasets has been inaccessible to existing deep generative models since they require training corpora, and foundation graph models are still in development\([Wang et al\., 2025b](https://arxiv.org/html/2610.02439#bib.bib44)\)\. In order to bridge this gap, we assess our model’s ability to infer meaningful structures in such datasets, with sizes ranging fromN=163N=163toN=6301N=6301\.

Evaluation for single\-graph datasets\.Standard graph generation evaluation uses MMD between distributions of graph statistics computed across many graphs in test and generated sets\([You et al\., 2018](https://arxiv.org/html/2610.02439#bib.bib23);[Vignac et al\., 2022](https://arxiv.org/html/2610.02439#bib.bib27);[Wang et al\., 2025a](https://arxiv.org/html/2610.02439#bib.bib38)\)\. This metric is ill\-defined atN=1N=1since no test distribution exists with a single empirical graph\. We instead compare intra\-graph distributions between empirical and synthetic networks using 1\-Wasserstein distances on degree distributions and Laplacian spectra, as well as mean squared error on the multifractal spectra \(Appendix[J](https://arxiv.org/html/2610.02439#A10)\)\.

Reference comparison with per\-instance fitting\.We calibrate against WMGM\([Yang et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib5)\), which is to our knowledge the only generative\-rule estimator for self\-similar networks with weighted edges and node heterogeneity\. WMGM optimizes on each network separately via EM\. WMGM thus serves as a reference for what direct fitting achieves, not as a competing baseline\. Our zero\-shot method is better than or tied with WMGM on 10 of 18 comparisons \(2 better, 8 tied, 8 worse\), while running 9–550×\\timesfaster \(Table[2](https://arxiv.org/html/2610.02439#S5.T2)\)\. Differences are largest in spectral error, where WMGM is better on four of six networks\. This indicates that our zero\-shot method reconstructs the multifractal and degree structure of real networks comparably to direct fitting, with a remaining gap in spectral structure\.

Table 2:Reconstruction of Single\-Instance NetworksMFrac: multifractal MSE; DegW: degree Wasserstein; SpecW: spectral Wasserstein\. Error metrics are means over 10 seeds, with the full table with SD values reported in Table[S2](https://arxiv.org/html/2610.02439#A9.T2)\. Times are defined in Appendices[F](https://arxiv.org/html/2610.02439#A6)and[I](https://arxiv.org/html/2610.02439#A9)\. Bold marks the lower value where the two methods’ mean±\\pmSD intervals do not overlap\.

Table 3:Revealing Hidden Differences in EEG Networks Under Sedation\.Values are mean±\\pmSD overn=20n=20sessions per state; the same 20 subjects contribute to every state\.ddis Cohen'sdd\(pooled SD\) for the first state minus the second\. Bold marks values whose two\-sided Mann–WhitneyUUtest survives Holm correction over the four contrasts \(p<0\.05p<0\.05\)\. The correctedppis shown in parentheses only for significant values\.

Inferred parameters track brain state\.To demonstrate the utility of our approach in scientific applications, we next ask whether the inferred parameters reflect reversible changes in network organization\. We analyze alpha\-band EEG connectivity from 20 subjects at baseline, under mild and moderate propofol sedation, and after recovery\([Chennu et al\., 2016](https://arxiv.org/html/2610.02439#bib.bib40)\)\(Appendix[K](https://arxiv.org/html/2610.02439#A11)\)\. We focus on the alpha band \(8–15 Hz\) functional network, where propofol produces its most prominent EEG changes on this frequency band\([Purdon et al\., 2013](https://arxiv.org/html/2610.02439#bib.bib60)\)\.

Fromp^\(1\)\\hat\{p\}^\{\(1\)\}we derive two summaries, the assortativity index and multifractal width \(Δ​α\\Delta\\alpha\) \(Appendix[K](https://arxiv.org/html/2610.02439#A11)\)\. Both metrics separate moderate sedation from every other state \(Cohen’s\|d\|=0\.8\|d\|=0\.8–0\.90\.9, the difference in means divided by the pooled SD, Holm\-correctedp<0\.05p<0\.05\), while values at baseline, mild sedation and recovery are similar \(Table[3](https://arxiv.org/html/2610.02439#S5.T3)\)\. Under moderate sedation the rule becomes less disassortative and its multifractal spectrum narrows \(Δ​α\\Delta\\alpha:3\.4→1\.73\.4\\to 1\.7\), consistent with reported changes in multiscale complexity\([Casali et al\., 2013](https://arxiv.org/html/2610.02439#bib.bib41);[Wilsenach et al\., 2025](https://arxiv.org/html/2610.02439#bib.bib43)\)and integration\([Luppi et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib42)\)under altered consciousness\. Critically, both metrics return to levels similar to baseline on recovery, suggesting that they reflect a reversible, state\-dependent change rather than drift over the recording session\. In contrast, traditional network statistics show smaller effects \(\|d\|=0\.2\|d\|=0\.2–0\.60\.6\), none significant after correction\. The two summaries also decode moderate sedation in unseen subjects \(AUC0\.690\.69–0\.750\.75\), suggesting that they carry state information that is shared between subjects \(Fig\.[S2](https://arxiv.org/html/2610.02439#A11.F2)A\)\. Moreover, the fit does not differ across states, indicating that the difference in inferred parameters during moderate sedation is not confounded by changes in the model fit \(Fig\.[S2](https://arxiv.org/html/2610.02439#A11.F2)B\)\. Lastly, 76 of 80 networks fit within the error range expected when the model is true \(Fig\.[S2](https://arxiv.org/html/2610.02439#A11.F2)C\)\.

## 6Discussion

Here we reiterate that comparisons with LGGM\-X and WMGM serve as reference points for transfer and fit\. As such, we do not claim that our formulation is intended to supersede deep or mechanistic models in general\. Rather, it is a complementary method suited to settings requiring both interpretable parameters and amortized inference, such as in scientific applications when N=1 and networks are of sizes unseen during training\. We also do not claim that MWNOT is the optimal inference method for this framework\. Our approach has several limitations\. For one, the model assumes the same connection rule repeats across scales\. Thus, networks whose wiring rule varies irregularly across scales fall outside this assumption\. In addition, the inverse operator is not permutation\-equivariant\. We resolve this with a deterministic degree\-sorting canonicalization, which is always available without any metadata \(Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)\)\. Lastly, the model needs improvement in the spectral metrics\. Despite these limitations, the results presented here indicate that mechanistic interpretability and amortized inference can be effectively combined within a unified framework to improve our understanding of complex networks\. Future work includes making the generative rule more expressive and comparing other inverse\-operator architectures to identify which yields the most accurate inference\.

## AI Use Statement

In this work, we used generative AI tools to design and provide feedback on research methodology and experiments, implement methods, and interpret results\. We have not used generative AI tools to help develop theoretical models or conceptual frameworks, formulate mathematical claims, provide critical ingredients for proving mathematical claims, assist in the writing of proofs, propose or refine hypotheses, or clean and reformat datasets, and generating synthetic datasets, assisting with translation, and supporting qualitative and thematic data analysis are not applicable to this work\. Additionally, we used generative AI tools to suggest experimental parameters, create and edit software code, create scientific figures, draft parts of the paper, edit the paper to improve readability, identify relevant literature, and format references\. We have reviewed all AI\-assisted work\. All AI\-generated code, analyses and text edits were reviewed by the authors\. We take responsibility for the final content of this work, including text, claims or artifacts produced with the aid of generative AI\.

## Reproducibility Statement

We describe the full experimental protocol in the Appendix, including preprocessing, training settings, hyperparameters, and evaluation procedures\. Anonymized source code to reproduce all experiments is provided in the supplementary materials, with instructions for running each experiment\.

## References

- Achlioptaset al\.\(2005\)D\. Achlioptas, A\. Clauset, D\. Kempe, and C\. MooreOn the bias of traceroute sampling: or, power\-law degree distributions in regular graphs\.InProceedings of the thirty\-seventh annual ACM symposium on Theory of computing,pp\. 694–703\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Airoldiet al\.\(2013\)E\. M\. Airoldi, T\. B\. Costa, and S\. H\. ChanStochastic blockmodel approximation of a graphon: theory and consistent estimation\.Advances in Neural Information Processing Systems26\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px1.p1.1),[§3\.1\.1](https://arxiv.org/html/2610.02439#S3.SS1.SSS1.p1.1),[§3\.1](https://arxiv.org/html/2610.02439#S3.SS1.p1.1)\.
- Akarcaet al\.\(2021\)D\. Akarca, P\. E\. Vértes, E\. T\. Bullmore, and D\. E\. AstleA generative network model of neurodevelopmental diversity in structural brain organization\.Nature communications12\(1\),pp\. 4216\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Aldous \(1981\)D\. J\. AldousRepresentations for partially exchangeable arrays of random variables\.Journal of Multivariate Analysis11\(4\),pp\. 581–598\.Cited by:[§3\.1\.1](https://arxiv.org/html/2610.02439#S3.SS1.SSS1.p1.1),[§3\.1](https://arxiv.org/html/2610.02439#S3.SS1.p2.1)\.
- Alpert \(1993\)B\. K\. AlpertA class of bases inL2\{L\}^\{2\}for the sparse representation of integral operators\.SIAM journal on Mathematical Analysis24\(1\),pp\. 246–262\.Cited by:[§4](https://arxiv.org/html/2610.02439#S4.p2.1)\.
- Azizpouret al\.\(2024\)A\. Azizpour, N\. Zilberstein, and S\. SegarraScalable implicit graphon learning\.arXiv preprint arXiv:2410\.17464\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px3.p1.1)\.
- Barabási and Albert \(1999\)A\. Barabási and R\. AlbertEmergence of scaling in random networks\.Science286\(5439\),pp\. 509–512\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p2.1),[§1](https://arxiv.org/html/2610.02439#S1.p3.1)\.
- Batty and Longley \(1994\)M\. Batty and P\. A\. LongleyFractal cities: a geometry of form and function\.Academic press\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Bleiet al\.\(2017\)D\. M\. Blei, A\. Kucukelbir, and J\. D\. McAuliffeVariational inference: a review for statisticians\.Journal of the American Statistical Association112\(518\),pp\. 859–877\.External Links:[Document](https://dx.doi.org/10.1080/01621459.2017.1285773)Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Borgset al\.\(2008\)C\. Borgs, J\. T\. Chayes, L\. Lovász, V\. T\. Sós, and K\. VesztergombiConvergent sequences of dense graphs i: subgraph frequencies, metric properties and testing\.Advances in Mathematics219\(6\),pp\. 1801–1851\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px1.p1.1),[§3\.1\.1](https://arxiv.org/html/2610.02439#S3.SS1.SSS1.p1.1)\.
- Brandeset al\.\(2008\)U\. Brandes, D\. Delling, M\. Gaertler, R\. Gorke, M\. Hoefer, Z\. Nikoloski, and D\. WagnerOn modularity clustering\.IEEE transactions on knowledge and data engineering20\(2\),pp\. 172–188\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Broido and Clauset \(2019\)A\. D\. Broido and A\. ClausetScale\-free networks are rare\.Nature communications10\(1\),pp\. 1017\.Cited by:[§5\.1](https://arxiv.org/html/2610.02439#S5.SS1.p2.1)\.
- Casaliet al\.\(2013\)A\. G\. Casali, O\. Gosseries, M\. Rosanova, M\. Boly, S\. Sarasso, K\. R\. Casali, S\. Casarotto, M\. Bruno, S\. Laureys, G\. Tononi, and M\. MassiminiA theoretically based index of consciousness independent of sensory processing and behavior\.Science Translational Medicine5\(198\),pp\. 198ra105–198ra105\.Cited by:[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p5.1)\.
- Chan and Airoldi \(2014\)S\. Chan and E\. AiroldiA consistent histogram estimator for exchangeable graph models\.InInternational Conference on Machine Learning,pp\. 208–216\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px1.p1.1)\.
- Chenet al\.\(2019\)S\. Chen, D\. Liao, C\. Chen, T\. Wang, and K\. ChenConstruction and analysis of protein\-protein interaction network of heroin use disorder\.Scientific reports9\(1\),pp\. 4980\.Cited by:[footnote 2](https://arxiv.org/html/2610.02439#footnote2)\.
- Chennuet al\.\(2016\)S\. Chennu, S\. O’Connor, R\. Adapa, D\. K\. Menon, and T\. A\. BekinschteinBrain connectivity dissociates responsiveness from drug exposure during propofol\-induced transitions of consciousness\.PLoS computational biology12\(1\),pp\. e1004669\.Cited by:[Appendix K](https://arxiv.org/html/2610.02439#A11.SS0.SSS0.Px1.p1.1),[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p4.1)\.
- Daudinet al\.\(2008\)J\. Daudin, F\. Picard, and S\. RobinA mixture model for random graphs\.Statistics and computing18\(2\),pp\. 173–183\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Dorkenwaldet al\.\(2024\)S\. Dorkenwald, A\. Matsliah, A\. R\. Sterling, P\. Schlegel, S\. Yu, C\. E\. McKellar, A\. Lin, M\. Costa, K\. Eichler, Y\. Yin,et al\.Neuronal wiring diagram of an adult brain\.Nature634\(8032\),pp\. 124–138\.Cited by:[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p1.1)\.
- Erdős and Rényi \(1960\)P\. Erdős and A\. RényiOn the evolution of random graphs\.Publ\. Math\. Inst\. Hung\. Acad\. Sci5\(1\),pp\. 17–60\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p2.1),[§1](https://arxiv.org/html/2610.02439#S1.p3.1)\.
- Guptaet al\.\(2021\)G\. Gupta, X\. Xiao, and P\. BogdanMultiwavelet\-based operator learning for differential equations\.Advances in neural information processing systems34,pp\. 24048–24062\.Cited by:[§D\.1](https://arxiv.org/html/2610.02439#A4.SS1.p3.1),[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px2.p1.1),[§4\.1](https://arxiv.org/html/2610.02439#S4.SS1.p1.1),[§4](https://arxiv.org/html/2610.02439#S4.p2.1)\.
- Hajek \(1988\)B\. HajekCooling schedules for optimal annealing\.Mathematics of operations research13\(2\),pp\. 311–329\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Haoet al\.\(2023\)Z\. Hao, Z\. Wang, H\. Su, C\. Ying, Y\. Dong, S\. Liu, Z\. Cheng, J\. Song, and J\. ZhuGnot: a general neural operator transformer for operator learning\.InInternational conference on machine learning,pp\. 12556–12569\.Cited by:[§D\.1](https://arxiv.org/html/2610.02439#A4.SS1.p5.1)\.
- Kipf and Welling \(2016\)T\. N\. Kipf and M\. WellingVariational graph auto\-encoders\.InNeurIPS Workshop on Bayesian Deep Learning,Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p2.1)\.
- Kolda and Bader \(2009\)T\. G\. Kolda and B\. W\. BaderTensor decompositions and applications\.SIAM review51\(3\),pp\. 455–500\.Cited by:[§C\.2\.1](https://arxiv.org/html/2610.02439#A3.SS2.SSS1.p3.2)\.
- Kovachkiet al\.\(2023\)N\. Kovachki, Z\. Li, B\. Liu, K\. Azizzadenesheli, K\. Bhattacharya, A\. Stuart, and A\. AnandkumarNeural operator: learning maps between function spaces with applications to pdes\.Journal of Machine Learning Research24\(89\),pp\. 1–97\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px2.p1.1)\.
- La Roccaet al\.\(2018\)D\. La Rocca, N\. Zilber, P\. Abry, V\. van Wassenhove, and P\. CiuciuSelf\-similarity and multifractality in human brain activity: a wavelet\-based analysis of scale\-free brain dynamics\.Journal of neuroscience methods309,pp\. 175–187\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Lee and Jung \(2006\)C\. Lee and S\. JungStatistical self\-similar properties of complex networks\.Physical Review E—Statistical, Nonlinear, and Soft Matter Physics73\(6\),pp\. 066102\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Leeet al\.\(2006\)S\. H\. Lee, P\. Kim, and H\. JeongStatistical properties of sampled networks\.Physical Review E—Statistical, Nonlinear, and Soft Matter Physics73\(1\),pp\. 016102\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Leskovecet al\.\(2010\)J\. Leskovec, D\. Chakrabarti, J\. Kleinberg, C\. Faloutsos, and Z\. GhahramaniKronecker graphs: an approach to modeling networks\.\.Journal of Machine Learning Research11\(2\)\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p3.1),[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Liet al\.\(2018\)Y\. Li, O\. Vinyals, C\. Dyer, R\. Pascanu, and P\. BattagliaLearning deep generative models of graphs\.arXiv preprint arXiv:1803\.03324\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p2.1),[§1](https://arxiv.org/html/2610.02439#S1.p5.1)\.
- Liet al\.\(2021\)Z\. Li, N\. Kovachki, K\. Azizzadenesheli, B\. Liu, K\. Bhattacharya, A\. Stuart, and A\. AnandkumarFourier neural operator for parametric partial differential equations\.InInternational Conference on Learning Representations,Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px2.p1.1)\.
- Liaoet al\.\(2019\)R\. Liao, Y\. Li, Y\. Song, S\. Wang, W\. Hamilton, D\. K\. Duvenaud, R\. Urtasun, and R\. ZemelEfficient graph generation with graph recurrent attention networks\.Advances in neural information processing systems32\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p5.1)\.
- Liuet al\.\(2017\)J\. Liu, J\. Wang, Z\. Yu, and X\. XieFractal and multifractal analyses of bipartite networks\.Scientific reports7\(1\),pp\. 45588\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Long and Chen \(2021\)Y\. Long and Y\. ChenMultifractal scaling analyses of urban street network structure: the cases of twelve megacities in china\.PLoS One16\(2\),pp\. e0246925\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Lovász and Szegedy \(2006\)L\. Lovász and B\. SzegedyLimits of dense graph sequences\.Journal of Combinatorial Theory, Series B96\(6\),pp\. 933–957\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px1.p1.1),[§3\.1\.1](https://arxiv.org/html/2610.02439#S3.SS1.SSS1.p1.1)\.
- Luet al\.\(2021\)L\. Lu, P\. Jin, G\. Pang, Z\. Zhang, and G\. E\. KarniadakisLearning nonlinear operators via deeponet based on the universal approximation theorem of operators\.Nature Machine Intelligence3\(3\),pp\. 218–229\.Cited by:[§D\.1](https://arxiv.org/html/2610.02439#A4.SS1.p5.1),[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px2.p1.1),[§4\.1](https://arxiv.org/html/2610.02439#S4.SS1.p2.1)\.
- Luppiet al\.\(2021\)A\. I\. Luppi, M\. M\. Craig, I\. Pappas, P\. Finoia, G\. B\. Williams, J\. Allanson, J\. D\. Pickard, A\. M\. Owen, L\. Naci, D\. K\. Menon, and E\. A\. StamatakisBrain network integration dynamics are associated with loss and recovery of consciousness induced by sevoflurane\.Human Brain Mapping42\(9\),pp\. 2802–2822\.Cited by:[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p5.1)\.
- Molinaroet al\.\(2023\)R\. Molinaro, Y\. Yang, B\. Engquist, and S\. MishraNeural inverse operators for solving pde inverse problems\.arXiv preprint arXiv:2301\.11167\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px2.p1.1),[§3\.1\.5](https://arxiv.org/html/2610.02439#S3.SS1.SSS5.p1.1)\.
- Niuet al\.\(2020\)C\. Niu, Y\. Song, J\. Song, S\. Zhao, A\. Grover, and S\. ErmonPermutation invariant graph generation via score\-based generative modeling\.InInternational conference on artificial intelligence and statistics,pp\. 4474–4484\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p5.1)\.
- Pallaet al\.\(2010\)G\. Palla, L\. Lovász, and T\. VicsekMultifractal network generator\.Proceedings of the National Academy of Sciences107\(17\),pp\. 7640–7645\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p3.1),[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Purdonet al\.\(2013\)P\. L\. Purdon, E\. T\. Pierce, E\. A\. Mukamel, M\. J\. Prerau, J\. L\. Walsh, K\. F\. K\. Wong, A\. F\. Salazar\-Gomez, P\. G\. Harrell, A\. L\. Sampson, A\. Cimenser,et al\.Electroencephalogram signatures of loss and recovery of consciousness from propofol\.Proceedings of the National Academy of Sciences110\(12\),pp\. E1142–E1151\.Cited by:[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p4.1)\.
- Shihet al\.\(2025\)B\. Shih, A\. Peyvan, Z\. Zhang, and G\. E\. KarniadakisTransformers as neural operators for solutions of differential equations with finite regularity\.Computer Methods in Applied Mechanics and Engineering434,pp\. 117560\.Cited by:[§D\.1](https://arxiv.org/html/2610.02439#A4.SS1.p5.1),[§4\.1](https://arxiv.org/html/2610.02439#S4.SS1.p2.1),[§4](https://arxiv.org/html/2610.02439#S4.p2.1)\.
- Simonovsky and Komodakis \(2018\)M\. Simonovsky and N\. KomodakisGraphvae: towards generation of small graphs using variational autoencoders\.InInternational conference on artificial neural networks,pp\. 412–422\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p5.1)\.
- Stumpfet al\.\(2005\)M\. P\. Stumpf, C\. Wiuf, and R\. M\. MaySubnets of scale\-free networks are not scale\-free: sampling properties of networks\.Proceedings of the National Academy of Sciences102\(12\),pp\. 4221–4224\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Vignacet al\.\(2022\)C\. Vignac, I\. Krawczuk, A\. Siraudin, B\. Wang, V\. Cevher, and P\. FrossardDigress: discrete denoising diffusion for graph generation\.arXiv preprint arXiv:2209\.14734\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p5.1),[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p2.1)\.
- Vincket al\.\(2011\)M\. Vinck, R\. Oostenveld, M\. Van Wingerden, F\. Battaglia, and C\. M\. PennartzAn improved index of phase\-synchronization for electrophysiological data in the presence of volume\-conduction, noise and sample\-size bias\.Neuroimage55\(4\),pp\. 1548–1565\.Cited by:[Appendix K](https://arxiv.org/html/2610.02439#A11.SS0.SSS0.Px2.p1.1)\.
- Wanget al\.\(2025a\)Y\. Wang, R\. A\. Rossi, N\. Park, H\. Chen, N\. K\. Ahmed, P\. Trivedi, F\. Dernoncourt, D\. Koutra, and T\. DerrA large\-scale training paradigm for graph generative models\.InThe Thirteenth International Conference on Learning Representations,Cited by:[Appendix H](https://arxiv.org/html/2610.02439#A8.SS0.SSS0.Px1.p1.1),[Appendix H](https://arxiv.org/html/2610.02439#A8.p1.1),[§5\.3\.1](https://arxiv.org/html/2610.02439#S5.SS3.SSS1.p1.1),[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p2.1)\.
- Wanget al\.\(2025b\)Z\. Wang, Z\. Liu, T\. Ma, J\. Li, Z\. Zhang, X\. Fu, Y\. Li, Z\. Yuan, W\. Song, Y\. Ma,et al\.Graph foundation models: a comprehensive survey\.arXiv preprint arXiv:2505\.15116\.Cited by:[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p1.1)\.
- Watts and Strogatz \(1998\)D\. J\. Watts and S\. H\. StrogatzCollective dynamics of ‘small\-world’ networks\.Nature393\(6684\),pp\. 440–442\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p3.1)\.
- Wilsenachet al\.\(2025\)J\. B\. Wilsenach, C\. M\. Deane, G\. Reinert, and K\. WarnabyGraph models of brain state in deep anesthesia reveal sink state dynamics of reduced spatiotemporal complexity\.Network Neuroscience9\(4\),pp\. 1176–1198\.Cited by:[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p5.1)\.
- Wu \(1983\)C\. F\. J\. WuOn the convergence properties of the em algorithm\.The Annals of Statistics11\(1\),pp\. 95–103\.External Links:[Document](https://dx.doi.org/10.1214/aos/1176346060)Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p4.1)\.
- Xiaet al\.\(2023\)X\. Xia, G\. Mishne, and Y\. WangImplicit graphon neural representation\.InInternational Conference on Artificial Intelligence and Statistics,pp\. 10619–10634\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px3.p1.1)\.
- Xiaoet al\.\(2023\)X\. Xiao, D\. Cao, R\. Yang, G\. Gupta, G\. Liu, C\. Yin, R\. Balan, and P\. BogdanCoupled multiwavelet neural operator learning for coupled partial differential equations\.arXiv preprint arXiv:2303\.02304\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px2.p1.1),[§4](https://arxiv.org/html/2610.02439#S4.p2.1)\.
- Xuet al\.\(2021\)H\. Xu, P\. Zhao, J\. Huang, and D\. LuoLearning graphon autoencoders for generative graph modeling\.arXiv preprint arXiv:2105\.14244\.Cited by:[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px3.p1.1)\.
- Xue and Bogdan \(2017\)Y\. Xue and P\. BogdanReliable multi\-fractal characterization of weighted complex networks: algorithms and implications\.Scientific reports7\(1\),pp\. 7487\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1),[§2](https://arxiv.org/html/2610.02439#S2.SS0.SSS0.Px1.p1.1),[§3\.1\.3](https://arxiv.org/html/2610.02439#S3.SS1.SSS3.p1.1),[§3\.1](https://arxiv.org/html/2610.02439#S3.SS1.p1.1)\.
- Yang and Bogdan \(2020\)R\. Yang and P\. BogdanControlling the multifractal generating measures of complex networks\.Scientific reports10\(1\),pp\. 5541\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Yanget al\.\(2021\)R\. Yang, F\. Sala, and P\. BogdanHidden network generating rules from partially observed complex networks\.Communications Physics4\(1\),pp\. 199\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1),[§1](https://arxiv.org/html/2610.02439#S1.p3.1),[§1](https://arxiv.org/html/2610.02439#S1.p4.1),[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p3.1)\.
- Yehet al\.\(2018\)F\. Yeh, S\. Panesar, D\. Fernandes, A\. Meola, M\. Yoshino, J\. C\. Fernandez\-Miranda, J\. M\. Vettel, and T\. VerstynenPopulation\-averaged atlas of the macroscale human structural connectome and its network topology\.Neuroimage178,pp\. 57–68\.Cited by:[footnote 2](https://arxiv.org/html/2610.02439#footnote2)\.
- Yinet al\.\(2020\)C\. Yin, X\. Xiao, V\. Balaban, M\. E\. Kandel, Y\. J\. Lee, G\. Popescu, and P\. BogdanNetwork science characteristics of brain\-derived neuronal cultures deciphered from quantitative phase imaging data\.Scientific reports10\(1\),pp\. 15078\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p1.1)\.
- Youet al\.\(2018\)J\. You, R\. Ying, X\. Ren, W\. Hamilton, and J\. LeskovecGraphrnn: generating realistic graphs with deep auto\-regressive models\.InInternational conference on machine learning,pp\. 5708–5717\.Cited by:[§1](https://arxiv.org/html/2610.02439#S1.p5.1),[§5\.3\.2](https://arxiv.org/html/2610.02439#S5.SS3.SSS2.p2.1)\.

Supplemental Material

## Appendix AModel Selection

##### SelectingMM\.

MMis the number of base communities\. The base graphon splits the latent space\[0,1\]\[0,1\]intoMMintervals and each recursion step splits every interval intoMMsub\-intervals\. The depth\-KKgraphon therefore partitions the latent space intoMKM^\{K\}blocks, and each block must contain enough nodes for its parameters to be reflected in the observed graph\. For our smallest training graphs at the largest depth \(N=512N=512,K=6K=6; training pairsK=6K=6only withN≥512N\\geq 512\),M=2M=2yields26=642^\{6\}=64blocks, whereasM=4M=4would yield46=40964^\{6\}=4096blocks, far exceeding the number of nodes and leaving most blocks empty\. Supporting largerMMwould therefore require exponentially larger training graphs\.

##### Expressivity\.

FixingM=2M=2does not restrict the model to two communities\. Expanding two recursive steps,l\(1\)⊗l\(1\)l^\{\(1\)\}\\otimes l^\{\(1\)\}andp\(1\)⊗p\(1\)p^\{\(1\)\}\\otimes p^\{\(1\)\}form a validM=4M=4base rule, so anM=2M=2model at depthKKis equivalent to the subset of theM=4M=4model class at depthK/2K/2whose base rule is Kronecker\-factorizable\. What an unrestrictedM=4M=4model adds is capacity for non\-factorizable branching rules\. Periodic multirule generators \(Appendix[L](https://arxiv.org/html/2610.02439#A12)\) offer intermediate expressivity without requiring larger training graphs, and are a natural direction for future work\.

## Appendix BClosed\-Form Network Statistics under the Multifractal Step Graphon

A consequence of the Kronecker structure ofWKW\_\{K\}is that several network statistics admit closed\-form expressions that factorize across scales\. We derive two such expressions here: the expected edge density and the expected degree of a node assigned to a given block\.

### B\.1Expected edge density

The expected edge density ofWKW\_\{K\}is

ρ⁡\(WK\)=∫01∫01WK​\(x,y\)​𝑑x​𝑑y=∑q,hlq\(K\)​lh\(K\)​pq​h\(K\)\.\\rho\(W\_\{K\}\)=\\int\_\{0\}^\{1\}\\int\_\{0\}^\{1\}W\_\{K\}\(x,y\)\\,dx\\,dy=\\sum\_\{q,h\}l^\{\(K\)\}\_\{q\}\\,l^\{\(K\)\}\_\{h\}\\,p^\{\(K\)\}\_\{qh\}\.\(S1\)Becausel\(K\)l^\{\(K\)\}andp\(K\)p^\{\(K\)\}each factorize asKK\-fold Kronecker products, the sum decomposes into a product ofKKidentical scale\-1 sums:

ρ⁡\(WK\)=∏k=1K\(∑i,jli\(1\)​lj\(1\)​pi​j\(1\)\)=ρ​\(W1\)K\.\\rho\(W\_\{K\}\)=\\prod\_\{k=1\}^\{K\}\\left\(\\sum\_\{i,j\}l^\{\(1\)\}\_\{i\}\\,l^\{\(1\)\}\_\{j\}\\,p^\{\(1\)\}\_\{ij\}\\right\)=\\rho\(W\_\{1\}\)^\{K\}\.\(S2\)Edge density at depthKKis therefore theKK\-th power of the depth\-1 density, making explicit the recursive nature of the construction\.

### B\.2Expected degree of a node in blockqq

For a graph ofNNnodes with latent coordinates drawn iid fromUniform​\(0,1\)\\text\{Uniform\}\(0,1\), the expected out\-degree of a node assigned to blockq=\(q1,…,qK\)q=\(q\_\{1\},\\ldots,q\_\{K\}\)is

𝔼⁡\[dqout\]=N​∑hlh\(K\)​pq​h\(K\)=N​∏k=1K\(∑jlj\(1\)​pqk​j\(1\)\)\.\\mathbb\{E\}\[d^\{\\text\{out\}\}\_\{q\}\]=N\\sum\_\{h\}l^\{\(K\)\}\_\{h\}\\,p^\{\(K\)\}\_\{qh\}=N\\prod\_\{k=1\}^\{K\}\\left\(\\sum\_\{j\}l^\{\(1\)\}\_\{j\}\\,p^\{\(1\)\}\_\{q\_\{k\}j\}\\right\)\.\(S3\)The expression for in\-degree is analogous, withpqk​j\(1\)p^\{\(1\)\}\_\{q\_\{k\}j\}replaced bypj​qk\(1\)p^\{\(1\)\}\_\{jq\_\{k\}\}\. Each scale\-kkfactor depends only onqkq\_\{k\}, so the expected degree of a node decomposes into independent contributions from each level of the recursion\.

### B\.3Multifractal spectrum

The multifractal structure ofWKW\_\{K\}admits a closed\-form characterization through the partition function\. Define the scale\-1 moment sum

Z⁡\(q\)=∑i,jli\(1\)​lj\(1\)​\(pi​j\(1\)\)q\.Z\(q\)=\\sum\_\{i,j\}l^\{\(1\)\}\_\{i\}\\,l^\{\(1\)\}\_\{j\}\\,\\left\(p^\{\(1\)\}\_\{ij\}\\right\)^\{q\}\.\(S4\)Because the multiplicative structure ofWKW\_\{K\}makes the depth\-KKpartition function factorize asZK​\(q\)=Z​\(q\)KZ\_\{K\}\(q\)=Z\(q\)^\{K\}, the mass exponentτ⁡\(q\)\\tau\(q\)takes the closed form

τ⁡\(q\)=−log⁡Z⁡\(q\)log⁡M\.\\tau\(q\)=\-\\frac\{\\log Z\(q\)\}\{\\log M\}\.\(S5\)The Lipschitz–Hölder exponentα⁡\(q\)\\alpha\(q\)and multifractal spectrumf⁡\(α\)f\(\\alpha\)are then obtained via the Legendre transform:

α⁡\(q\)=d​τd​q,f⁡\(α⁡\(q\)\)=q​α​\(q\)−τ⁡\(q\)\.\\alpha\(q\)=\\frac\{d\\tau\}\{dq\},\\qquad f\(\\alpha\(q\)\)=q\\,\\alpha\(q\)\-\\tau\(q\)\.\(S6\)The full multifractal spectrum is therefore determined byθ\(1\)\\theta^\{\(1\)\}alone, independent ofKKin the limit\. This closed form makes precise the claim in Section[3\.1\.3](https://arxiv.org/html/2610.02439#S3.SS1.SSS3)thatθ\(1\)\\theta^\{\(1\)\}controls whereWKW\_\{K\}lies on the monofractal–multifractal continuum: monofractal behavior corresponds toZ⁡\(q\)=cqZ\(q\)=c^\{q\}for some constantcc, in which caseτ⁡\(q\)\\tau\(q\)is linear andf⁡\(α\)f\(\\alpha\)collapses to a single point; heterogeneity inp\(1\)p^\{\(1\)\}produces curvature inτ⁡\(q\)\\tau\(q\)and a non\-degenerate spectrum\.

### B\.4Implication

These factorizations show that macroscopic network statistics inherit the recursive structure of the generator: depth\-KKobservables are powers or products of depth\-1 quantities, and the multifractal spectrum is fully characterized byθ\(1\)\\theta^\{\(1\)\}alone\. This provides analytical tools for relatingθ\(1\)\\theta^\{\(1\)\}to observable network properties, and confirms that the recursive construction is not merely a notational device but produces statistics that scale predictably with depth\.

## Appendix CProof of Proposition 1 \(Identifiability\)

We establish identifiability through three steps: \(i\) the step graphonWKW\_\{K\}determines its block structure uniquely up to permutation, \(ii\) the Kronecker factorization of the block parameters is unique, and \(iii\) the combined result yields identifiability of the base parameters up to community label permutation\.

### C\.1Step 1: A step graphon determines its block parameters up to permutation

A step graphonWKW\_\{K\}is a piecewise\-constant function on\[0,1\]2\[0,1\]^\{2\}defined by a partition of\[0,1\]\[0,1\]intoMKM^\{K\}contiguous intervals with lengthsl\(K\)∈int⁡\(ΔMK−1\)l^\{\(K\)\}\\in\\mathrm\{int\}\(\\Delta^\{M^\{K\}\-1\}\)and a block interaction matrixp\(K\)∈\(0,1\)MK×MKp^\{\(K\)\}\\in\(0,1\)^\{M^\{K\}\\times M^\{K\}\}\. We show that these parameters are determined byWKW\_\{K\}up to a simultaneous permutation of the blocks\.

Since all entries ofl\(K\)l^\{\(K\)\}are strictly positive, the partition consists ofMKM^\{K\}intervals each with positive Lebesgue measure\. The functionWKW\_\{K\}takes the constant valuepq​h\(K\)p^\{\(K\)\}\_\{qh\}on the rectangleIq×IhI\_\{q\}\\times I\_\{h\}\.

Suppose two parameterizations\(l1\(K\),p1\(K\)\)\(l^\{\(K\)\}\_\{1\},p^\{\(K\)\}\_\{1\}\)and\(l2\(K\),p2\(K\)\)\(l^\{\(K\)\}\_\{2\},p^\{\(K\)\}\_\{2\}\)yield the same graphonWKW\_\{K\}almost everywhere\. Because the graphon is piecewise constant with finitely many pieces, equality a\.e\. implies equality on every rectangle of positive measure\.

The partition boundaries are determined by the discontinuities ofWKW\_\{K\}\. For a fixedyyin the interior of some intervalIhI\_\{h\}, the functionWK​\(⋅,y\)W\_\{K\}\(\\cdot,y\)is piecewise constant inxx, taking the valuepq​h\(K\)p^\{\(K\)\}\_\{qh\}whenx∈Iqx\\in I\_\{q\}\. The discontinuities of this function occur at the boundaries between intervalsIqI\_\{q\}andIq\+1I\_\{q\+1\}wheneverpq​h\(K\)≠p\(q\+1\)​h\(K\)p^\{\(K\)\}\_\{qh\}\\neq p^\{\(K\)\}\_\{\(q\+1\)h\}\. By varyingyyacross all intervals, a boundary betweenIqI\_\{q\}andIq\+1I\_\{q\+1\}is undetectable only ifpq​h\(K\)=p\(q\+1\)​h\(K\)p^\{\(K\)\}\_\{qh\}=p^\{\(K\)\}\_\{\(q\+1\)h\}for everyhh— that is, blocksqqandq\+1q\+1have identical interaction profiles with all other blocks\.

In this case, blocksqqandq\+1q\+1are indistinguishable: they have the same size \(since both parameterizations must agree on the set of interval widths\) and the same interaction profile\. Swapping their labels produces the same graphon\. This is precisely a permutation of the block indices\.

Because the intervals are contiguous \(ordered left\-to\-right on\[0,1\]\[0,1\]\), the set of partition boundaries is a property of the function itself, not of the parameterization\. Both parameterizations must use the same boundary points and therefore the same set of interval widths\. The only freedom is in which interval is assigned which label\. Thus, there exists a permutationπ∈SMK\\pi\\in S\_\{M^\{K\}\}such that:

l1,q\(K\)=l2,π⁡\(q\)\(K\),p1,q​h\(K\)=p2,π⁡\(q\)​π​\(h\)\(K\)l^\{\(K\)\}\_\{1,q\}=l^\{\(K\)\}\_\{2,\\pi\(q\)\},\\quad p^\{\(K\)\}\_\{1,qh\}=p^\{\(K\)\}\_\{2,\\pi\(q\)\\pi\(h\)\}\(S7\)for allq,h∈\{1,…,MK\}q,h\\in\\\{1,\\dots,M^\{K\}\\\}\. ∎

### C\.2Step 2: Uniqueness of the Kronecker factorization

We now show that if two base parameterizations\(l1,p1\)\(l\_\{1\},p\_\{1\}\)and\(l2,p2\)\(l\_\{2\},p\_\{2\}\)yield the same expanded parameters underKK\-fold Kronecker product \(up to a permutationπ∈SMK\\pi\\in S\_\{M^\{K\}\}\), then they are related by a base\-level permutationσ∈SM\\sigma\\in S\_\{M\}\.

#### C\.2\.1Step 2a: Factorization of the length vector

By construction,l\(K\)=l\(1\)⊗⋯⊗l\(1\)l^\{\(K\)\}=l^\{\(1\)\}\\otimes\\cdots\\otimes l^\{\(1\)\}\(KKtimes\)\. The entries ofl\(K\)l^\{\(K\)\}are indexed by multi\-indicesq=\(q1,…,qK\)∈\{1,…,M\}Kq=\(q\_\{1\},\\dots,q\_\{K\}\)\\in\\\{1,\\dots,M\\\}^\{K\}and take the form:

lq\(K\)=∏k=1Klqk\(1\)l^\{\(K\)\}\_\{q\}=\\prod\_\{k=1\}^\{K\}l^\{\(1\)\}\_\{q\_\{k\}\}\(S8\)These are the entries of a rank\-1 symmetric tensorl\(1\)⊗⋯⊗l\(1\)∈ℝM×⋯×Ml^\{\(1\)\}\\otimes\\cdots\\otimes l^\{\(1\)\}\\in\\mathbb\{R\}^\{M\\times\\cdots\\times M\}\. We show thatl\(1\)l^\{\(1\)\}is recoverable from this tensor up to permutation of its entries\.

Consider theM×MK−1M\\times M^\{K\-1\}matricization ofl\(K\)l^\{\(K\)\}, denotedL\(1\)L\_\{\(1\)\}, obtained by unfolding along the first mode:

L\(1\)=l\(1\)⋅\(l\(1\)⊗⋯⊗l\(1\)\)⊤L\_\{\(1\)\}=l^\{\(1\)\}\\cdot\\left\(l^\{\(1\)\}\\otimes\\cdots\\otimes l^\{\(1\)\}\\right\)^\{\\top\}\(S9\)Since all entries ofl\(1\)l^\{\(1\)\}are strictly positive,L\(1\)L\_\{\(1\)\}has rank 1 and its column space is spanned byl\(1\)l^\{\(1\)\}\. The simplex constraint∑ili\(1\)=1\\sum\_\{i\}l^\{\(1\)\}\_\{i\}=1fixes the scale uniquely\. The only remaining ambiguity is permutation of the entries ofl\(1\)l^\{\(1\)\}\.

To see why the permutation must act consistently across all modes: a permutationπ∈SMK\\pi\\in S\_\{M^\{K\}\}of the multi\-indices that maps one Kronecker product to another must act independently on each coordinate of the multi\-index\. Ifπ\\pimaps\(q1,…,qK\)\(q\_\{1\},\\dots,q\_\{K\}\)such that the first coordinateq1q\_\{1\}is permuted byσ\\sigmain one multi\-index but left unchanged in another, the multiplicative structurelq\(K\)=∏klqk\(1\)l^\{\(K\)\}\_\{q\}=\\prod\_\{k\}l^\{\(1\)\}\_\{q\_\{k\}\}would be violated\. Formally, the logarithms of the entries form an additive rank\-1 tensor:

log⁡lq\(K\)=∑k=1Klog⁡lqk\(1\)\\log l^\{\(K\)\}\_\{q\}=\\sum\_\{k=1\}^\{K\}\\log l^\{\(1\)\}\_\{q\_\{k\}\}\(S10\)A permutation preserving this additive structure must act as the same base permutationσ\\sigmaon each coordinate independently\([Kolda and Bader, 2009](https://arxiv.org/html/2610.02439#bib.bib22)\)\.

Therefore, ifl1\(K\)=π⁡\(l2\(K\)\)l\_\{1\}^\{\(K\)\}=\\pi\(l\_\{2\}^\{\(K\)\}\), thenπ\\piis induced by a base permutationσ∈SM\\sigma\\in S\_\{M\}acting asσ⁡\(q1,…,qK\)=\(σ⁡\(q1\),…,σ⁡\(qK\)\)\\sigma\(q\_\{1\},\\dots,q\_\{K\}\)=\(\\sigma\(q\_\{1\}\),\\dots,\\sigma\(q\_\{K\}\)\), andl1=σ⁡\(l2\)l\_\{1\}=\\sigma\(l\_\{2\}\)\.

#### C\.2\.2Step 2b: Factorization of the interaction matrix

Similarly,p\(K\)=p\(1\)⊗⋯⊗p\(1\)p^\{\(K\)\}=p^\{\(1\)\}\\otimes\\cdots\\otimes p^\{\(1\)\}\(KKtimes\)\. The entries satisfy:

pq​h\(K\)=∏k=1Kpqk​hk\(1\)p^\{\(K\)\}\_\{qh\}=\\prod\_\{k=1\}^\{K\}p^\{\(1\)\}\_\{q\_\{k\}h\_\{k\}\}\(S11\)Letvec⁡\(p\(1\)\)∈\(0,1\)M2\\mathrm\{vec\}\(p^\{\(1\)\}\)\\in\(0,1\)^\{M^\{2\}\}denote the vectorization ofp\(1\)p^\{\(1\)\}\. Then:

vec\(p\(K\)\)=vec\(p\(1\)\)⊗⋯⊗vec\(p\(1\)\)\\mathrm\{vec\}\(p^\{\(K\)\}\)=\\mathrm\{vec\}\(p^\{\(1\)\}\)\\otimes\\cdots\\otimes\\mathrm\{vec\}\(p^\{\(1\)\}\)\(S12\)By the same rank\-1 tensor uniqueness argument as in Step 2a,vec⁡\(p\(1\)\)\\mathrm\{vec\}\(p^\{\(1\)\}\)is recoverable up to permutation\. However, this permutation is not arbitrary overM2M^\{2\}entries — it must be consistent with the base\-level permutationσ\\sigmaestablished in Step 2a\.

The permutationπ\\piacts simultaneously on bothl\(K\)l^\{\(K\)\}andp\(K\)p^\{\(K\)\}, because both are defined over the same block indices\. Blockqqhas sizelq\(K\)l^\{\(K\)\}\_\{q\}and interaction valuespq​h\(K\)p^\{\(K\)\}\_\{qh\}\. Any permutation that preserves both must be factorizable asσ\\sigmaacting independently on each level of the multi\-index\. The permutationσ\\sigmafrom Step 2a induces the action:

p1,i​j=p2,σ⁡\(i\)​σ​\(j\)p\_\{1,ij\}=p\_\{2,\\sigma\(i\)\\sigma\(j\)\}\(S13\)Thereforep1=σ⁡\(p2\)p\_\{1\}=\\sigma\(p\_\{2\}\)under the same permutationσ\\sigma\. ∎

### C\.3Step 3: Combining the results

From Step 1,WK\(1\)=WK\(2\)W\_\{K\}^\{\(1\)\}=W\_\{K\}^\{\(2\)\}a\.e\. implies that\(l1\(K\),p1\(K\)\)\(l\_\{1\}^\{\(K\)\},p\_\{1\}^\{\(K\)\}\)and\(l2\(K\),p2\(K\)\)\(l\_\{2\}^\{\(K\)\},p\_\{2\}^\{\(K\)\}\)agree up to a block permutationπ∈SMK\\pi\\in S\_\{M^\{K\}\}\. From Step 2, the Kronecker structure forcesπ\\pito be induced by a base\-level permutationσ∈SM\\sigma\\in S\_\{M\}, yieldingl1=σ⁡\(l2\)l\_\{1\}=\\sigma\(l\_\{2\}\)andp1=σ⁡\(p2\)p\_\{1\}=\\sigma\(p\_\{2\}\)\.

This establishes that the multifractal step graphonWKW\_\{K\}uniquely determines the base parametersθ\(1\)=\(l\(1\),p\(1\)\)\\theta^\{\(1\)\}=\(l^\{\(1\)\},p^\{\(1\)\}\)up to community label permutation\.

### C\.4Discussion of Assumptions

Strict positivity ofl\(1\)l^\{\(1\)\}\.The requirementl\(1\)∈int⁡\(ΔM−1\)l^\{\(1\)\}\\in\\mathrm\{int\}\(\\Delta^\{M\-1\}\)ensures that every community has positive measure and contributes to the graphon\. If someli\(1\)=0l^\{\(1\)\}\_\{i\}=0, the corresponding community is absent and the effective model has fewer thanMMcommunities, reducing to an identifiable model of smaller dimension\.

Strict positivity ofp\(1\)p^\{\(1\)\}\.The requirementpi​j\(1\)∈\(0,1\)p^\{\(1\)\}\_\{ij\}\\in\(0,1\)ensures that the vectorized interaction matrix has no zero entries, which is needed for the rank\-1 tensor factorization to be unique\. Ifpi​j\(1\)=0p^\{\(1\)\}\_\{ij\}=0for some\(i,j\)\(i,j\), multiple factorizations may produce the same zero entry\. In practice, exact zeros are measure\-zero in the parameter space and do not arise in the sampling procedure described in Appendix[F](https://arxiv.org/html/2610.02439#A6)\.

K≥2K\\geq 2\.The result requires at least two Kronecker expansions\. ForK=1K=1, the graphon is a standard step graphon withMMblocks, and identifiability reduces to the classical step graphon case \(covered by Step 1 alone\)\.

Note this proposition does not claim identifiability over arbitrary graphon representatives under all measure\-preserving transformations, nor does it imply statistical stability from finite noisy graph samples\.

## Appendix DMultiwavelet Neural Operator Transformer

### D\.1Architecture

Given an observed adjacency matrixA∈ℝN×NA\\in\\mathbb\{R\}^\{N\\times N\}, the model processes the graph topology through the following operations\.

Channel lift\.A learned linear map lifts each entry of the inputAAtok2k^\{2\}channels, reshaping it intox\(0\)∈ℝN×N×k2x^\{\(0\)\}\\in\\mathbb\{R\}^\{N\\times N\\times k^\{2\}\}\. The Legendre polynomial basis enters through the fixed multiwavelet filters\.

Multiwavelet decomposition\.A 2D discrete multiwavelet transform is applied overLLspatial scales using a non\-standard wavelet representation\([Gupta et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib10)\)\. The multiwavelet basis is constructed via Gram–Schmidt orthogonalization of shifted Legendre polynomials, ensuring that high\-frequency details are strictly orthogonal to the coarse approximation\. At each scaleℓ∈\{1,…,L\}\\ell\\in\\\{1,\\dots,L\\\}, the feature grid is decomposed by a coarse filterℋ\\mathcal\{H\}and a composite detail filter𝒢\\mathcal\{G\}:

s\(ℓ\)\\displaystyle s^\{\(\\ell\)\}=x\(ℓ−1\)∗ℋ\\displaystyle=x^\{\(\\ell\-1\)\}\*\\mathcal\{H\}\(S14\)d\(ℓ\)\\displaystyle d^\{\(\\ell\)\}=x\(ℓ−1\)∗𝒢\\displaystyle=x^\{\(\\ell\-1\)\}\*\\mathcal\{G\}\(S15\)The coarse coefficientss\(ℓ\)s^\{\(\\ell\)\}capture macroscopic connectivity, whiled\(ℓ\)d^\{\(\\ell\)\}concatenates horizontal \(LH\), vertical \(HL\), and diagonal \(HH\) high\-frequency variations\. The outputx\(ℓ\)=s\(ℓ\)x^\{\(\\ell\)\}=s^\{\(\\ell\)\}serves as input for the next coarser scale\. This hierarchical decomposition mirrors the multiplicative construction of the multifractal step graphon\. WhenNNis not a power of two, the canonicalised input \(Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)\) is cropped to the largest dyadic size≤N\\leq N, which drops the lowest\-degree nodes\.

Multiscale tokenization\.The wavelet coefficients \(detailsd\(ℓ\)d^\{\(\\ell\)\}forℓ≤L\\ell\\leq Land the final coarses\(L\)s^\{\(L\)\}\) are flattened into sequences and concatenated with positional grids and scale\-indicator labels to form multiscale tokensT\(ℓ\)T^\{\(\\ell\)\}\.

Cross\-attention for operator approximation\.The model parameterizes the inverse operator via cross\-attention in a branch\-trunk architecture\([Lu et al\., 2021](https://arxiv.org/html/2610.02439#bib.bib16);[Shih et al\., 2025](https://arxiv.org/html/2610.02439#bib.bib12);[Hao et al\., 2023](https://arxiv.org/html/2610.02439#bib.bib15)\)\. The multiscale tokens encode the input function \(branch\), while learnable query embeddings evaluate the target parameters \(trunk\)\.

LetQl,Qp∈ℝM×DQ\_\{l\},Q\_\{p\}\\in\\mathbb\{R\}^\{M\\times D\}denote query matrices forl\(1\)l^\{\(1\)\}andp\(1\)p^\{\(1\)\}\. At each scaleℓ\\ell, tokensT\(ℓ\)T^\{\(\\ell\)\}are projected into KeysK\(ℓ\)∈ℝNℓ×DK^\{\(\\ell\)\}\\in\\mathbb\{R\}^\{N\_\{\\ell\}\\times D\}and ValuesV\(ℓ\)∈ℝNℓ×SV^\{\(\\ell\)\}\\in\\mathbb\{R\}^\{N\_\{\\ell\}\\times S\}\. The cross\-attention aggregates these into fixed\-size latent representations:

zs\(ℓ\)​\(m\)=∑d=1D\(1Nℓ​∑j=1Nℓkd​\(wj\(ℓ\)\)T​vs​\(wj\(ℓ\)\)\)​qd​\(m\)z\_\{s\}^\{\(\\ell\)\}\(m\)=\\sum\_\{d=1\}^\{D\}\\left\(\\frac\{1\}\{N\_\{\\ell\}\}\\sum\_\{j=1\}^\{N\_\{\\ell\}\}k\_\{d\}\\left\(w\_\{j\}^\{\(\\ell\)\}\\right\)^\{T\}v\_\{s\}\\left\(w\_\{j\}^\{\(\\ell\)\}\\right\)\\right\)q\_\{d\}\(m\)\(S16\)wherem∈\{1,…,M\}m\\in\\\{1,\\dots,M\\\}indexes the target community andwj\(ℓ\)w\_\{j\}^\{\(\\ell\)\}is thejj\-th wavelet token withinT\(ℓ\)T^\{\(\\ell\)\}\. The inner summation averages overNℓN\_\{\\ell\}grid locations, forming a data\-dependent basis from the wavelet coefficients\. The query embeddingqd​\(m\)q\_\{d\}\(m\)extracts the relevant features, yielding attention outputsZl\(ℓ\),Zp\(ℓ\)∈ℝM×SZ\_\{l\}^\{\(\\ell\)\},Z\_\{p\}^\{\(\\ell\)\}\\in\\mathbb\{R\}^\{M\\times S\}\.

Decoding\.Scale\-specific representations are concatenated into global embeddingsZp=\[Zp\(1\),…,Zp\(L\)\]Z\_\{p\}=\[Z\_\{p\}^\{\(1\)\},\\dots,Z\_\{p\}^\{\(L\)\}\]andZl=\[Zl\(1\),…,Zl\(L\)\]Z\_\{l\}=\[Z\_\{l\}^\{\(1\)\},\\dots,Z\_\{l\}^\{\(L\)\}\]\. MLPs decode these to predict the base parameters:

p^\(1\)\\displaystyle\\hat\{p\}^\{\(1\)\}=MLPp​\(Zp\)∈\[0,1\]M×M\\displaystyle=\\text\{MLP\}\_\{p\}\(Z\_\{p\}\)\\in\[0,1\]^\{M\\times M\}\(S17\)l^\(1\)\\displaystyle\\hat\{l\}^\{\(1\)\}=MLPl​\(Zl\)∈ΔM−1\\displaystyle=\\text\{MLP\}\_\{l\}\(Z\_\{l\}\)\\in\\Delta^\{M\-1\}\(S18\)

## Appendix EEmpirical Validation of the Multiwavelet Basis

We empirically verify that the multiwavelet basis used in MWNOT is responsive to the recursive multiscale construction of the multifractal step graphon to support the inductive bias of using this basis for multifractal objects\.

##### Method\.

For each recursion depthK∈\{1,…,6\}K\\in\\\{1,\\ldots,6\\\}, we construct the continuous graphonWKW\_\{K\}on anN×NN\\times Ngrid \(N=1024N=1024\) using fixed base parametersl\(1\)=\[0\.4,0\.6\]l^\{\(1\)\}=\[0\.4,0\.6\]andp\(1\)=\[0\.90\.10\.10\.9\]p^\{\(1\)\}=\\begin\{bmatrix\}0\.9&0\.1\\\\ 0\.1&0\.9\\end\{bmatrix\}for demonstration\. We projectWKW\_\{K\}onto the constant component of the Legendre multiwavelet basis and recursively apply the 2D multiwavelet transformLLtimes\. At each scaleℓ\\ell, we record the squared\-norm energy of the detail coefficientsE\(ℓ\)=∑\|d\(ℓ\)\|2E^\{\(\\ell\)\}=\\sum\|d^\{\(\\ell\)\}\|^\{2\}, and after the final iteration, we record the energy of the coarse approximationEcoarse\(L\)=∑\|s\(L\)\|2E^\{\(L\)\}\_\{\\text\{coarse\}\}=\\sum\|s^\{\(L\)\}\|^\{2\}\. From these we compute two summary quantities:

ρfine​\(K\)\\displaystyle\\rho\_\{\\text\{fine\}\}\(K\)=E\(1\)\+E\(2\)Ecoarse\(L\)\+∑ℓ=1LE\(ℓ\),\\displaystyle=\\frac\{E^\{\(1\)\}\+E^\{\(2\)\}\}\{E^\{\(L\)\}\_\{\\text\{coarse\}\}\+\\sum\_\{\\ell=1\}^\{L\}E^\{\(\\ell\)\}\},\(S19\)ρdetail​\(K\)\\displaystyle\\rho\_\{\\text\{detail\}\}\(K\)=∑ℓ=1LE\(ℓ\)Ecoarse\(L\)\+∑ℓ=1LE\(ℓ\)\.\\displaystyle=\\frac\{\\sum\_\{\\ell=1\}^\{L\}E^\{\(\\ell\)\}\}\{E^\{\(L\)\}\_\{\\text\{coarse\}\}\+\\sum\_\{\\ell=1\}^\{L\}E^\{\(\\ell\)\}\}\.\(S20\)

##### Result\.

Figure[S1](https://arxiv.org/html/2610.02439#A5.F1)reports both quantities as a function ofKK\. Bothρfine\\rho\_\{\\text\{fine\}\}andρdetail\\rho\_\{\\text\{detail\}\}grow approximately exponentially withKK\(linear on a semi\-log scale\), spanning roughly two orders of magnitude over the recursion depths used in training\.

##### Interpretation\.

Each recursion level introduces finer\-scale block boundaries toWKW\_\{K\}, which the multiwavelet basis encodes through its detail coefficients\. Deeper recursion therefore produces an exponentially larger share of energy in the detail components\. MWNOT’s encoder and cross\-attention operate on these detail coefficients; the architecture’s representational capacity is thus allocated to the components where the multifractal structure ofWKW\_\{K\}resides\. We interpret this as an empirical match between the basis and the construction, supporting the multiwavelet basis as a principled choice for representing this object class\.

Figure S1:Multiwavelet energy response to recursion depth\.\(a\) Fine\-scale detail shareρfine\\rho\_\{\\text\{fine\}\}\. \(b\) Total detail energy shareρdetail\\rho\_\{\\text\{detail\}\}\. Both quantities grow approximately exponentially withKK, consistent with the multiplicative scaling structure of the multifractal step graphon\.
### E\.1Input Ordering

Function\-space methods operate on inputs aligned to spatial axes\. More specifically, the multiwavelet basis decomposes the adjacency matrix as a discretized function on\[0,1\]2\[0,1\]^\{2\}, which benefits from inputs whose block structure is aligned with the spatial axes\. As such, function\-space formulations trade architectural permutation equivariance for a parameterization that allows for generalization to unseen graph sizes at inference\.

This consideration concerns only the finite\-sample to graphon inverse map learned by the architecture\. It is distinct from the population\-level identifiability ofWk↦θ\(1\)W\_\{k\}\\mapsto\\theta^\{\(1\)\}established in Proposition[1](https://arxiv.org/html/2610.02439#Thmproposition1), which holds at the level of graphons, and from the exchangeability of finite graphs sampled from a graphon under the Aldous framework, which guarantees that*generation*of finite samples under our model remains permutation\-invariant\.

Canonicalization\.We resolve the ordering dependence with a single deterministic, parameter\-free rule applied to every input, synthetic or empirical: nodes are sorted in descending order of degree \(weighted strength for weighted inputs\)\. Ties are broken by a stable sort in the zero\-shot benchmark and arbitrarily elsewhere\. The rule uses no metadata, is available for any graph, and makes the model input a function of the graph alone: two labelings of the same graph yield the same sorted adjacency matrix up to the order of equal\-degree nodes, so the inferredθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}is invariant to node relabeling up to that residual tie ambiguity\. Degree sorting also recovers an approximately block\-aligned form for the wavelet basis, since expected degree is monotone in community membership under the model \(Appendix[B](https://arxiv.org/html/2610.02439#A2)\)\.

Training\.Synthetic graphs are produced by multinomial community assignment followed by geometric edge sampling \(Appendix[F](https://arxiv.org/html/2610.02439#A6)\)\. This randomizes which node sits in which community, producing adjacency matrices that resemble static noise rather than block\-structured signals; the canonicalization above is applied to every generated training graph before augmentation \(Appendix[F](https://arxiv.org/html/2610.02439#A6)\)\.

Inference\.The same canonicalization is applied to every empirical network before inference, in the zero\-shot benchmark \(Appendix[H](https://arxiv.org/html/2610.02439#A8)\), the single\-instance reconstructions \(Appendix[I](https://arxiv.org/html/2610.02439#A9)\), and the EEG analysis\. Dataset\-native orderings, whether anatomical, temporal or spatial, are never used; all reported numbers are produced under the metadata\-free ordering\. The MMD statistics of the benchmark are permutation\-invariant, so canonicalizing the model input leaves the reference side of every comparison unchanged\.

## Appendix FTraining Details

### F\.1Data Generation

We generate synthetic training pairs\(A,\{p\(1\),l\(1\)\}\)\(A,\\\{p^\{\(1\)\},l^\{\(1\)\}\\\}\)by constructingKK\-level parametersl\(K\)l^\{\(K\)\}andp\(K\)p^\{\(K\)\}via Kronecker products of randomly sampled base measures\. The base parameters are sampled as follows\. For the community distributionl\(1\)∈ΔM−1l^\{\(1\)\}\\in\\Delta^\{M\-1\}, each entry is initialized to a minimum value of0\.10\.1to ensure that every community has non\-negligible measure, and the remaining probability mass\(1−0\.1​M\)\(1\-0\.1M\)is distributed by drawingMMuniform random values and normalizing them to sum to the residual\. The interaction matrixp\(1\)∈\(0,1\)M×Mp^\{\(1\)\}\\in\(0,1\)^\{M\\times M\}is sampled entry\-wise fromUniform​\(0\.1,1\.0\)\\text\{Uniform\}\(0\.1,1\.0\), ensuring strict positivity as required by the identifiability result \(Proposition 1\)\. Parameter sets with\|l0\(1\)−l1\(1\)\|≤0\.01\|l^\{\(1\)\}\_\{0\}\-l^\{\(1\)\}\_\{1\}\|\\leq 0\.01are rejected\. The targets are ordered by expected degree, so community 0 is the higher\-degree community\.

Nodes are assigned to communities from a multinomial distribution parameterized byl\(K\)l^\{\(K\)\}, and edge weights are drawn asei​j∼Geometric​\(1−pq​h\(K\)\)e\_\{ij\}\\sim\\text\{Geometric\}\(1\-p^\{\(K\)\}\_\{qh\}\)\. Non\-zero edge weights are log\-binned to form the input adjacency matrix:

Ai​j=⌊log2⁡\(ei​j\)⌋\+1forei​j\>0A\_\{ij\}=\\lfloor\\log\_\{2\}\(e\_\{ij\}\)\\rfloor\+1\\quad\\text\{for\}\\quad e\_\{ij\}\>0\(S21\)Inputs are canonicalised as in Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)before log\-binning\.

##### Augmentation\.

Each training graph is augmented with probability 0\.5 by node dropout, edge dropout and spurious edges, each at a rate drawn fromUniform⁡\(0,0\.3\)\\mathrm\{Uniform\}\(0,0\.3\)\. The edge masks are symmetric with zero diagonal\. Augmented graphs are binarized and are not re\-sorted after augmentation\.

### F\.2Loss Function

The model is trained with a composite loss\. The interaction matrix is optimized via binary cross\-entropy with logits \(ℒBCE\\mathcal\{L\}\_\{\\text\{BCE\}\}\) on the raw MLP outputsZ^p\\hat\{Z\}\_\{p\}, while the node distributionl^\(1\)∈ΔM−1\\hat\{l\}^\{\(1\)\}\\in\\Delta^\{M\-1\}is optimized via KL divergence \(ℒKL\\mathcal\{L\}\_\{\\text\{KL\}\}\):

ℒ\(θ\(1\),θ^\(1\)\)=ℒBCE\(p\(1\),Z^p\)\+ℒKL\(l\(1\)∥l^\(1\)\)\\mathcal\{L\}\(\\theta^\{\(1\)\},\\hat\{\\theta\}^\{\(1\)\}\)=\\mathcal\{L\}\_\{\\text\{BCE\}\}\(p^\{\(1\)\},\\hat\{Z\}\_\{p\}\)\+\\mathcal\{L\}\_\{\\text\{KL\}\}\(l^\{\(1\)\}\\parallel\\hat\{l\}^\{\(1\)\}\)\(S22\)

### F\.3Optimization and Hyperparameters

The model is trained for 2000 epochs using AdamW with a base learning rate of10−310^\{\-3\}, beginning with a 10\-epoch linear warmup from10−510^\{\-5\}to10−310^\{\-3\}, followed by cosine annealing to10−510^\{\-5\}over the remaining epochs\. The batch size is 64 \(32 per GPU on two GPUs\)\.

Each epoch draws 2000 new graphs from the generator\. Training graphs haveN∈\{128,512,1024\}N\\in\\\{128,512,1024\\\}andK∈\{2,4,6\}K\\in\\\{2,4,6\\\}\. Validation uses 100 fixed graphs atK=4K=4andN=512N=512, and we keep the epoch with the lowest validation MAE ofp\(1\)p^\{\(1\)\}plus that ofl\(1\)l^\{\(1\)\}\(epoch 1816\)\. The MWNOT architecture uses the following hyperparameters: query/key dimensionD=128D=128, latent feature/value dimensionS=128S=128,L=4L=4wavelet decomposition scales, and wavelet orderk=6k=6\(corresponding to order\-6 Legendre polynomial bases\)\. The number of base communities isM=2M=2throughout\. Recovery error \(Figure[4](https://arxiv.org/html/2610.02439#S5.F4)\) is computed on 100 new, unaugmented graphs per\(N,K\)\(N,K\)\(random seed 42 for the whole run\), with the depth sweep atN=1024N=1024, after matching predicted to true communities by the best of theM\!M\!permutations\.

### F\.4Compute resources\.

Training was performed on two NVIDIA V100 GPUs and completed in approximately 8 hours\. Inference is fast enough to run on a single GPU or modern CPU: per\-graph inference time is reported in Table[2](https://arxiv.org/html/2610.02439#S5.T2)\(“Time \(ms\)” column\) as one forward pass on the cropped input, averaged over 10 repeats after one warm\-up pass, excluding normalisation, sorting and log\-binning\.

## Appendix GExpressivity Sweep Details

For both sweeps, we setM=2M=2,K=4K=4, andN=1024N=1024, averaging over 10 realizations per setting\. The heterogeneity sweep uses𝐥\(1\)=\[0\.25,0\.75\]\\mathbf\{l\}^\{\(1\)\}=\[0\.25,0\.75\]and interpolates𝐩\(1\)\\mathbf\{p\}^\{\(1\)\}from a uniform matrix\[\[0\.5,0\.5\],\[0\.5,0\.5\]\]\[\[0\.5,0\.5\],\[0\.5,0\.5\]\]to a heterogeneous matrix\[\[0\.99,0\.05\],\[0\.05,0\.01\]\]\[\[0\.99,0\.05\],\[0\.05,0\.01\]\]aswhetw\_\{\\mathrm\{het\}\}varies from 0 to 1 over 30 steps\. Panel \(d\) of Figure[3](https://arxiv.org/html/2610.02439#S4.F3)is drawn from a separate sweep, in which𝐩\(1\)\\mathbf\{p\}^\{\(1\)\}is interpolated in 12 steps towards each of three targets \(left\-skewed\[\[0\.99,0\.2\],\[0\.2,0\.01\]\]\[\[0\.99,0\.2\],\[0\.2,0\.01\]\], symmetric\[\[0\.99,0\.01\],\[0\.01,0\.99\]\]\[\[0\.99,0\.01\],\[0\.01,0\.99\]\]and right\-skewed\[\[0\.01,0\.6\],\[0\.6,0\.99\]\]\[\[0\.01,0\.6\],\[0\.6,0\.99\]\]\), with one realization per step andq∈\[−3\.5,3\.5\]q\\in\[\-3\.5,3\.5\]\. The structure sweep uses𝐥\(1\)=\[0\.4,0\.6\]\\mathbf\{l\}^\{\(1\)\}=\[0\.4,0\.6\]and interpolates between a disassortative matrix\[\[0\.01,0\.99\],\[0\.99,0\.01\]\]\[\[0\.01,0\.99\],\[0\.99,0\.01\]\]atwstruct=−1w\_\{\\mathrm\{struct\}\}=\-1, the uniform matrix atwstruct=0w\_\{\\mathrm\{struct\}\}=0, and an assortative matrix\[\[0\.99,0\.01\],\[0\.01,0\.99\]\]\[\[0\.99,0\.01\],\[0\.01,0\.99\]\]atwstruct=\+1w\_\{\\mathrm\{struct\}\}=\+1over 21 steps\. In panels \(b\) and \(c\), degree assortativity is computed on the directed weighted graph, and Louvain modularity and the normalised\-Laplacian spectral gap on the weighted undirected graph, with isolated nodes removed;𝐩\(1\)\\mathbf\{p\}^\{\(1\)\}is clamped to\[0\.01,0\.99\]\[0\.01,0\.99\]\. The multifractal spectra in these sweeps useqqgrids of 51 points on\[−10,10\]\[\-10,10\]or 40 points on\[−8,8\]\[\-8,8\], instead of the grid of Appendix[J](https://arxiv.org/html/2610.02439#A10)\. Panel \(a\)’s degree\-distribution regimes use distinct hand\-tuned\(𝐥\(1\),𝐩\(1\),K\)\(\\mathbf\{l\}^\{\(1\)\},\\mathbf\{p\}^\{\(1\)\},K\)configurations to target each canonical shape \(full parameters in the released code\)\.

The Lipschitz–Hölder exponentα\\alphameasures local singularity strength: smallα\\alphacorresponds to hubs, largeα\\alphato peripheral nodes\.f⁡\(α\)f\(\\alpha\)gives the fractal dimension of the set of nodes at eachα\\alpha\. The widthΔ​α=αmax−αmin\\Delta\\alpha=\\alpha\_\{\\max\}\-\\alpha\_\{\\min\}summarizes multiscale heterogeneity\. Full definitions are in Appendix[J](https://arxiv.org/html/2610.02439#A10)\.

## Appendix HZero\-shot benchmark: evaluation protocol and full results

This appendix states the protocol behind Table[1](https://arxiv.org/html/2610.02439#S5.T1)and reports the full table with seed\-to\-seed variability \(Table[S1](https://arxiv.org/html/2610.02439#A8.T1)\)\. We follow the evaluation protocol of the LGGM benchmark\([Wang et al\., 2025a](https://arxiv.org/html/2610.02439#bib.bib38)\), such as the test splits, metrics, and released evaluation code\. The paragraphs below state how our model is applied within that protocol\.

##### Data\.

We use the twelve domains and the held\-out test graphs of the LGGM benchmark\([Wang et al\., 2025a](https://arxiv.org/html/2610.02439#bib.bib38)\)processed exactly as the benchmark’s evaluation pipeline \(undirected, self\-loops removed\)\. For ECO, as in the released evaluation code, the reference combines the training, validation and test splits \(251 graphs, six unique\)\.

##### Inference\.

The inverse operator receives the canonicalised \(Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)\) binary adjacency matrix, truncated to the largest dyadic size2⌊log2⁡N⌋2^\{\\lfloor\\log\_\{2\}N\\rfloor\}\(graphs with fewer than2L=162^\{L\}=16nodes are zero\-padded to 16 nodes\) as required by the multiwavelet decomposition\. One forward pass returnsθ^\(1\)=\(l^\(1\),p^\(1\)\)\\hat\{\\theta\}^\{\(1\)\}=\(\\hat\{l\}^\{\(1\)\},\\hat\{p\}^\{\(1\)\}\)\. No per\-graph optimisation is performed andθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}is used as returned\.

##### Recursion depth\.

The depthKKused to realise a graph fromθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}is set per reference graph by the closed\-form density relation of Appendix[B](https://arxiv.org/html/2610.02439#A2),ρemp=ρ​\(W1\)K\\rho\_\{\\mathrm\{emp\}\}=\\rho\(W\_\{1\}\)^\{K\}withρ⁡\(W1\)=∑i​jl^i​l^j​p^i​j\\rho\(W\_\{1\}\)=\\sum\_\{ij\}\\hat\{l\}\_\{i\}\\hat\{l\}\_\{j\}\\hat\{p\}\_\{ij\}, solved asK⋆=round⁡\(log⁡ρemp/log⁡ρ⁡\(W1\)\)K^\{\\star\}=\\mathrm\{round\}\\big\(\\log\\rho\_\{\\mathrm\{emp\}\}/\\log\\rho\(W\_\{1\}\)\\big\)and clamped to\[2,6\]\[2,6\]\. Hereρemp\\rho\_\{\\mathrm\{emp\}\}is the number of non\-zero entries divided byN2N^\{2\}, diagonal included, on the full uncropped network; rounding is half\-to\-even, andK=4K=4is used ifρemp\\rho\_\{\\mathrm\{emp\}\}orρ⁡\(W1\)\\rho\(W\_\{1\}\)lies outside\(0,1\)\(0,1\)\.K⋆K^\{\\star\}depends only onθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}and the observed edge density\. It is fixed before any graph is generated\. Therefore, no evaluation metric is minimised overKK\.

##### Generation and metrics\.

For each reference graph and each seed, one graph of the same size is sampled fromθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}at depthK⋆K^\{\\star\}\. These graphs are directed and may contain self\-loops\. Following the benchmark’s protocol, they are binarized, symmetrised by union \(A\+A⊤\>0A\+A^\{\\top\}\>0\) and stripped of self\-loops before calculating the reconstruction metrics\. We repeat the sampling for ten seeds \(42–51\) and report mean and standard deviation\. DiGress and LGGM\-X entries are the values published with the benchmark\. Ours\-Mono is the ablation of Ours, trained withp\(1\)p^\{\(1\)\}fixed at0\.50\.5and evaluated under the identical protocol\.

Table S1:Full version of Table[1](https://arxiv.org/html/2610.02439#S5.T1)with seed\-to\-seed standard deviations\. Lower is better\.DEG, CC, Spec, Orb: MMD of degree, clustering coefficient, Laplacian eigenvalues and orbit counts\. DiGress and LGGM\-X are published values\. Ours and Ours\-Mono are mean±\\pmSD over 10 seeds; the SD of the ALL row is the mean of the twelve per\-domain SDs\. Bold marks the best method where the mean±\\pmSD intervals do not overlap\. On ECO the orbit\-count kernel saturates, so the orbit MMD is a constant set by the reference graphs for any distinct, dense generated graphs; it is not used to mark the best method\.

## Appendix ISingle\-instance reconstruction: protocol

This appendix states the protocol behind Table[2](https://arxiv.org/html/2610.02439#S5.T2)\. Its metrics are defined in Appendix[J](https://arxiv.org/html/2610.02439#A10)\. Each dataset is a single observed network\. The comparison asks the reconstruction fidelity of networks generated from the zero\-shot estimateθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}, with a per\-instance expectation–maximisation fit \(WMGM\) as a reference for what direct optimisation on the same network achieves\.

##### Data\.

We use six networks, each observed once\. Four are directed \(C\. elegans, Drosophila medulla, Allen mouse and P2P\-Gnutella\)\. Two are undirected \(Power Grid and Minnesota Roads\)\. Only C\. elegans and Allen mouse carry edge weights\. The C\. elegans file is a383×575383\\times 575connection table whose rows are the 383 neurons\. We keep only the neuron\-to\-neuron connections, the leading383×383383\\times 383block\. The Allen mouse, Drosophila and C\. elegans networks contain self\-loops \(163, 105 and 8\), which are kept\. They enter the edge density, the node strengths used by the multifractal MSE and the degrees, where each self\-loop counts once\. Before use, we divide every weight by the smallest non\-zero weight in that network, so the weakest edge has weight 1\. We then order the nodes as described in Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)\. This full\-size network is the reference that every generated network is compared against\. Both methods receive this network after log\-binning \(Eq\. \([S21](https://arxiv.org/html/2610.02439#A6.E21)\)\)\. WMGM is fitted on the full network\. Our operator receives the crop to the largest dyadic size for compatibility with multiwavelet transform, as described in Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)\.

##### Recursion depth\.

As in Appendix[H](https://arxiv.org/html/2610.02439#A8),K⋆K^\{\\star\}is set by the closed\-form density relation fromθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}and the observed edge density, and is fixed before generation\. For WMGM, which is fitted separately at each depthK∈\{2,…,6\}K\\in\\\{2,\\dots,6\\\}, the reported depth is the one at which the fitted parameters are self\-consistent with the same relation, i\.e\. theKKminimising\|K⋆​\(θ^K\(1\)\)−K\|\|K^\{\\star\}\(\\hat\{\\theta\}^\{\(1\)\}\_\{K\}\)\-K\|\. No evaluation metric enters the choice for either method\.

##### Generation and evaluation\.

For each of ten seeds \(42–51\), one weighted network of the original size is sampled from each method’s parameters at its depth\. Samples are directed and may contain self\-loops\. The degree and spectral distances are computed on the binarized sample with self\-loops kept\. The sample is symmetrised by union for the spectral distance and, for the undirected references \(Power Grid and Minnesota Roads\), also for the degree distance\. The multifractal MSE uses the out\-strengths of the directed weighted sample\. Each sample is compared with the reference network using multifractal MSE, normalised degree Wasserstein distance, and spectral Wasserstein distance of the normalised Laplacian\. We report mean and standard deviation over seeds\.

##### WMGM fitting\.

WMGM is run on the same canonicalised, log\-binned matrix at full size withM=2M=2, up to 1000 EM iterations and an Adam step size of10−310^\{\-3\}, once per depth and seed\. Fitting stops early when the change in the evidence lower bound falls below10−410^\{\-4\}, and the step size is halved when the objective plateaus \(patience 50, cooldown 20, minimum10−610^\{\-6\}\)\. The parameters are initialized as follows:l\(1\)l^\{\(1\)\}as uniform,p\(1\)p^\{\(1\)\}fromUniform⁡\(0,1\)\\mathrm\{Uniform\}\(0,1\), and the community responsibilities fromDirichlet⁡\(1\)\\mathrm\{Dirichlet\}\(1\)\.p\(1\)p^\{\(1\)\}is clamped to\[10−6,1−10−6\]\[10^\{\-6\},1\-10^\{\-6\}\]\. Ties in the depth rule are resolved to the smallestKK\. Its reported time is the mean wall\-clock time per fit over the five depths and ten seeds, measured on the same GPU as our inference\.

Table S2:Full version of Table[2](https://arxiv.org/html/2610.02439#S5.T2)with seed\-to\-seed standard deviations\.MFrac: multifractal MSE; DegW: degree Wasserstein; SpecW: spectral Wasserstein\. Values are mean±\\pmSD over 10 seeds\. Bold marks the lower value where the two methods’ mean±\\pmSD intervals do not overlap\.

## Appendix JDefinitions of Evaluation Metrics for Single Observation Datasets

### J\.1Multifractal Spectrum \(f⁡\(α\)f\(\\alpha\)\) Estimation

To estimate the continuous multifractal spectrum of the networks, we apply the Chhabra\-Jensen method\. LetAAbe the adjacency matrix of a network withNNactive nodes, and letdi=∑jAi​jd\_\{i\}=\\sum\_\{j\}A\_\{ij\}be the degree of nodeii\. We first define the degree probability distributionPi=di/∑diP\_\{i\}=d\_\{i\}/\\sum d\_\{i\}\.

To isolate different scaling regimes, we define a one\-parameter family of normalized measures,μi​\(q\)\\mu\_\{i\}\(q\), parameterized by the momentq∈ℝq\\in\\mathbb\{R\}:

μi​\(q\)=Piq∑j=1NPjq\\mu\_\{i\}\(q\)=\\frac\{P\_\{i\}^\{q\}\}\{\\sum\_\{j=1\}^\{N\}P\_\{j\}^\{q\}\}\(S23\)By varyingqq, we heavily weight the massive hubs \(q\>0q\>0\) or the low\-degree periphery \(q<0q<0\)\. The Lipschitz\-Hölder exponentα⁡\(q\)\\alpha\(q\)and its corresponding fractal dimensionf⁡\(α⁡\(q\)\)f\(\\alpha\(q\)\)are then computed parametrically as:

α\(q\)=−1log⁡N∑i=1Nμi\(q\)logPi\\alpha\(q\)=\-\\frac\{1\}\{\\log N\}\\sum\_\{i=1\}^\{N\}\\mu\_\{i\}\(q\)\\log P\_\{i\}\(S24\)f\(α\(q\)\)=−1log⁡N∑i=1Nμi\(q\)logμi\(q\)f\(\\alpha\(q\)\)=\-\\frac\{1\}\{\\log N\}\\sum\_\{i=1\}^\{N\}\\mu\_\{i\}\(q\)\\log\\mu\_\{i\}\(q\)\(S25\)TheMultifractal MSEis calculated as the mean squared error between the empirical and synthetic continuous curves,fe​m​p​\(α\)f\_\{emp\}\(\\alpha\)andfs​y​n​\(α\)f\_\{syn\}\(\\alpha\), interpolated over a sharedα\\alphasupport space\. We use 30 values ofqqevenly spaced in\[−4,4\]\[\-4,4\]and interpolate both curves linearly on 50 points spanning the overlap of theirα\\alpharanges\.

### J\.2Degree Wasserstein Distance

To evaluate the absolute local hub structure, we utilize the 1\-Wasserstein distance \(Earth Mover’s Distance\) between the sets of empirical and synthetic degrees\. LetU⁡\(x\)U\(x\)andV⁡\(x\)V\(x\)be the cumulative distribution functions of the empirical degreesde​m​pd^\{emp\}and synthetic degreesds​y​nd^\{syn\}, respectively\. The unscaled degree distance is defined as:

W1​\(de​m​p,ds​y​n\)=∫−∞∞\|U⁡\(x\)−V⁡\(x\)\|​𝑑xW\_\{1\}\(d^\{emp\},d^\{syn\}\)=\\int\_\{\-\\infty\}^\{\\infty\}\|U\(x\)\-V\(x\)\|\\,dx\(S26\)For cross\-dataset comparisons independent of absolute scale \(as presented in the main text\), we compute the Normalized Degree Wasserstein by transforming the degrees into probability mass functions:d~i=di/∑dj\\tilde\{d\}\_\{i\}=d\_\{i\}/\\sum d\_\{j\}prior to calculating the Wasserstein integral\.

### J\.3Spectral Wasserstein Distance

To evaluate the global community structure, bottlenecking, and information flow of the networks, we analyze the eigenspectrum of the Normalized Laplacian matrix\. The Normalized Laplacian is defined as:

ℒ=I−D−1/2AD−1/2\\mathcal\{L\}=I\-D^\{\-1/2\}AD^\{\-1/2\}\(S27\)whereIIis the identity matrix andDDis the diagonal degree matrix\. Because the eigenvalues ofℒ\\mathcal\{L\}are strictly bounded within the interval\[0,2\]\[0,2\]regardless of network size or absolute density, they provide a perfectly scaled invariant measure of global topology\.

Letλe​m​p\\lambda^\{emp\}andλs​y​n\\lambda^\{syn\}denote the sorted sets of eigenvalues derived from the empirical and synthetic Normalized Laplacians\. The Spectral Wasserstein distance is calculated as the 1\-Wasserstein distance between these two discrete distributions:

W1​\(λe​m​p,λs​y​n\)=infγ∈Γ⁡\(λe​m​p,λs​y​n\)∫∫⁡\|x−y\|​𝑑γ​\(x,y\)W\_\{1\}\(\\lambda^\{emp\},\\lambda^\{syn\}\)=\\inf\_\{\\gamma\\in\\Gamma\(\\lambda^\{emp\},\\lambda^\{syn\}\)\}\\int\\int\|x\-y\|\\,d\\gamma\(x,y\)\(S28\)whereΓ\\Gammadenotes the collection of all joint probability measures with marginals equivalent to the empirical and synthetic eigenvalue distributions\.

The degree and spectral distances are computed on binarized graphs\. For directed networks, the degree distance uses out\-degrees and the spectral distance uses the union\-symmetrised graph\. Standard deviations over seeds \(Tables[S1](https://arxiv.org/html/2610.02439#A8.T1)and[S2](https://arxiv.org/html/2610.02439#A9.T2)\) are population SDs; Table[3](https://arxiv.org/html/2610.02439#S5.T3)reports sample SDs across subjects\.

## Appendix KEEG Analysis

##### Data\.

We use the propofol sedation dataset of[Chennu et al\. \(2016\)](https://arxiv.org/html/2610.02439#bib.bib40), which includes 20 healthy volunteers who were recorded at rest with 91\-channel scalp EEG \(250 Hz\) across four sessions: before propofol administration \(baseline\), under mild sedation, under moderate sedation, and after recovery\.

##### Network construction\.

Each epoch is band\-pass filtered to the alpha band \(8–15 Hz, fourth\-order Butterworth, zero\-phase\) and Hilbert\-transformed\. Functional connectivity between every pair of channels is the weighted phase\-lag index \(wPLI\) computed per epoch and averaged over the epochs of the session\([Vinck et al\., 2011](https://arxiv.org/html/2610.02439#bib.bib59)\)\. Edges are retained by a surrogate test\. For each epoch we generate 1000 phase\-randomised surrogates of every channel \(Fourier amplitudes preserved, phases drawn uniformly\), compute their wPLI matrices, and average them over epochs\. An edge is kept, with its wPLI value as weight, if the observed session\-mean wPLI exceeds the 95th percentile of the surrogate distribution for that channel pair\.

Figure S2:Supplementary EEG analyses\.\(A\)Predictive validity\. ROC curves of the leave\-one\-subject\-out \(LOSO\) logistic regression\. The legend gives the AUC with its 95% bootstrap interval over subjects\. Dotted line represents chance\.\(B\)Fit by state\. Reconstruction errorDD\(multifractal MSE, mean over ten generated networks\) of every EEG network, grouped by sedation state \(n=20n=20per state\)\. Boxes represent median and interquartile range\. Whiskers represent1\.5×1\.5\\timesIQR\. Dashed line represents the pass threshold of \(C\)D=0\.0185D=0\.0185\. \(Kruskal–Wallisp=0\.40p=0\.40; Friedman over subjectsp=0\.42p=0\.42; two\-sided Mann–Whitney tests of moderate sedation against each other state, Holm\-corrected,p=0\.65p=0\.65for all three; baseline vs recoveryp=0\.86p=0\.86\)\.\(C\)Synthetic reference range\. Distribution oflog10⁡D\\log\_\{10\}Dfor 500 model\-generated networks whose density lies in the EEG range \(gray\) and for the 80 EEG networks \(blue\), for the multifractal MSE \(MFrac; the decision metric\), the normalised degree Wasserstein distance \(DegW\) and the spectral Wasserstein distance \(SpecW\)\. Dashed lines represent the 95th percentile of the model\-generated networks \(0\.01850\.0185,0\.00160\.0016and0\.0170\.017\)\. The number of EEG networks below the line pass are as follows: 76 of 80 on MFrac, 71 of 80 on DegW, and 27 of 80 on SpecW\.
##### Inference\.

Each network is processed as in Appendix[I](https://arxiv.org/html/2610.02439#A9)\. Weights are normalised so that the weakest edge has weight 1, nodes are canonicalised \(Appendix[E\.1](https://arxiv.org/html/2610.02439#A5.SS1)\), the matrix is cropped to 64 nodes, which removes the 27 lowest\-strength channels, and log\-binned\. One forward pass of the operator returnsθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}\. Fromθ^\(1\)=\(l^\(1\),p^\(1\)\)\\hat\{\\theta\}^\{\(1\)\}=\(\\hat\{l\}^\{\(1\)\},\\hat\{p\}^\{\(1\)\}\)we report two scalar summaries\. One is the assortativity indexp^00\+p^11−p^01−p^10\\hat\{p\}\_\{00\}\+\\hat\{p\}\_\{11\}\-\\hat\{p\}\_\{01\}\-\\hat\{p\}\_\{10\}, which is positive when within\-class connections dominate\. The other is multifractal widthΔ​α=log2⁡\(p^max/p^min\)\\Delta\\alpha=\\log\_\{2\}\(\\hat\{p\}\_\{\\max\}/\\hat\{p\}\_\{\\min\}\), which is the width of the spectrum of Appendix[B](https://arxiv.org/html/2610.02439#A2)forM=2M=2\. Classical metrics \(mean clustering, mean shortest path, greedy modularity, spectral gap of the normalised Laplacian\) are computed on the same networks before cropping and after binarisation \(Ai​j\>0A\_\{ij\}\>0\)\. The shortest path and spectral gap use the largest connected component\.

##### Group comparison\.

The four sessions are compared as groups of 20 networks\. For each measure we report the following: 1\) mean and SD per state 2\) Cohen’sddwith pooled SD for each state minus moderate sedation and for baseline minus recovery 3\) a two\-sided Mann–WhitneyUUtest per contrast with Holm correction over the four contrasts\.

##### Predictive validity\.

To test whether the inferred rule summaries generalize to new individuals, we train a logistic regression on the twoθ^\(1\)\\hat\{\\theta\}^\{\(1\)\}\-derived summaries to distinguish moderate sedation from each other state, holding out one subject at a time\. Features are standardised within each training fold, and the model uses the default L2 penalty \(C=1C=1\)\. Predictions are pooled across held\-out subjects and scored by AUC with 95% bootstrap intervals over the 20 subjects \(1000 resamples\)\. The AUC is0\.690\.69\[0\.550\.55,0\.830\.83\] for moderate sedation versus baseline,0\.720\.72\[0\.580\.58,0\.850\.85\] versus mild sedation and0\.750\.75\[0\.620\.62,0\.890\.89\] versus recovery\. All three intervals exclude chance \(0\.50\.5\)\.

##### Synthetic reference range\.

We ask whether the reconstruction error of each EEG network is within the range seen when the model is true\. We draw parameter sets from the training prior and a depthKKuniformly from\{2,…,6\}\\\{2,\\dots,6\\\}, and generate one 91\-node network from each\. Unlike the EEG networks, which are symmetric, loop\-free and wPLI\-weighted, these networks are directed, may contain self\-loops and have integer weights\. As in Table[2](https://arxiv.org/html/2610.02439#S5.T2), the degree and spectral distances are computed on binarized, union\-symmetrised graphs with self\-loops kept, and the multifractal MSE on the directed weighted samples\. We keep a network only if its edge density lies in the observed EEG range, until 500 networks are kept\. Each network, synthetic or EEG, is then processed as described above, given to MWNOT as inputs, and ten networks are generated from each inferred rule \(seeds 42 to 51\)\. The generated graphs are each compared with the input network using the multifractal MSE of Appendix[J](https://arxiv.org/html/2610.02439#A10)\. The mean over the ten comparisons is the errorDD\. The pass threshold is the 95th percentile ofDDover the 500 synthetic networks,0\.01850\.0185\. The model’s fit on a given EEG is considered acceptable if itsDDis below the threshold\.

Of the 80 EEG networks, 76 pass \(18 baseline, 19 mild, 20 moderate and 19 recovery\)\. The medianDDof the EEG networks,0\.001390\.00139, is close to the synthetic median,0\.001300\.00130\. The degree Wasserstein and spectral Wasserstein distances are reported for completeness but do not enter the decision\. Under the same rule, 71 networks pass on the degree distance and 27 on the spectral distance\. The spectral gap between real and synthetic networks is consistent with Table[2](https://arxiv.org/html/2610.02439#S5.T2), where spectral error is also the weakest term\.

##### Fit and sedation state\.

State differences in the inferred rule could reflect worse fit under moderate sedation rather than a real change\. To rule this out, we compared fit across states\. We useDDas the fit score\. Its median is0\.00160\.0016at baseline,0\.00140\.0014under mild sedation,0\.00050\.0005under moderate sedation and0\.00220\.0022at recovery\. A Kruskal–Wallis test across the four states givesp=0\.40p=0\.40\. A Friedman test over the 20 subjects givesp=0\.42p=0\.42\. Two\-sided Mann–Whitney tests of moderate sedation against each other state, Holm\-corrected over the four contrasts of Table[3](https://arxiv.org/html/2610.02439#S5.T3), givep=0\.65p=0\.65for all three\. Baseline versus recovery givesp=0\.86p=0\.86\. We found no evidence that fit differs across states\.

## Appendix LCompression of a periodic multirule generator

###### Proposition 2\(Compression of a periodic multirule generator\)\.

For each scalekk, let

l\[k\]∈int⁡\(ΔMk−1\),p\[k\]∈\(0,1\)Mk×Mk,l^\{\[k\]\}\\in\\mathrm\{int\}\(\\Delta^\{M\_\{k\}\-1\}\),\\qquad p^\{\[k\]\}\\in\(0,1\)^\{M\_\{k\}\\times M\_\{k\}\},and define the multirule recursive expansion

l\(K\)=⨂k=1Kl\[k\],p\(K\)=⨂k=1Kp\[k\]\.l^\{\(K\)\}=\\bigotimes\_\{k=1\}^\{K\}l^\{\[k\]\},\\qquad p^\{\(K\)\}=\\bigotimes\_\{k=1\}^\{K\}p^\{\[k\]\}\.Assume that the rule sequence is periodic with periodRR, so that

\(l\[k\+R\],p\[k\+R\]\)=\(l\[k\],p\[k\]\)for all​k,\(l^\{\[k\+R\]\},p^\{\[k\+R\]\}\)=\(l^\{\[k\]\},p^\{\[k\]\}\)\\qquad\\text\{for all \}k,and letK=m​RK=mR\. Define

l~\(1\):=⨂r=1Rl\[r\],p~\(1\):=⨂r=1Rp\[r\],M~:=∏r=1RMr\.\\tilde\{l\}^\{\(1\)\}:=\\bigotimes\_\{r=1\}^\{R\}l^\{\[r\]\},\\qquad\\tilde\{p\}^\{\(1\)\}:=\\bigotimes\_\{r=1\}^\{R\}p^\{\[r\]\},\\qquad\\tilde\{M\}:=\\prod\_\{r=1\}^\{R\}M\_\{r\}\.Then

l~\(1\)∈int⁡\(ΔM~−1\),p~\(1\)∈\(0,1\)M~×M~,\\tilde\{l\}^\{\(1\)\}\\in\\mathrm\{int\}\(\\Delta^\{\\tilde\{M\}\-1\}\),\\qquad\\tilde\{p\}^\{\(1\)\}\\in\(0,1\)^\{\\tilde\{M\}\\times\\tilde\{M\}\},and the original multirule generator is exactly representable as

l\(K\)=l~\(1\)⊗⋯⊗l~\(1\)⏟m​times,p\(K\)=p~\(1\)⊗⋯⊗p~\(1\)⏟m​times\.l^\{\(K\)\}=\\underbrace\{\\tilde\{l\}^\{\(1\)\}\\otimes\\cdots\\otimes\\tilde\{l\}^\{\(1\)\}\}\_\{m\\text\{ times\}\},\\qquad p^\{\(K\)\}=\\underbrace\{\\tilde\{p\}^\{\(1\)\}\\otimes\\cdots\\otimes\\tilde\{p\}^\{\(1\)\}\}\_\{m\\text\{ times\}\}\.Hence a periodic multirule recursive generator reduces exactly to a single\-basis generator with effective branching number

M~=∏r=1RMr\.\\tilde\{M\}=\\prod\_\{r=1\}^\{R\}M\_\{r\}\.

###### Proof\.

Since eachl\[r\]l^\{\[r\]\}is a strictly positive probability vector, their Kronecker productl~\(1\)\\tilde\{l\}^\{\(1\)\}is again a strictly positive probability vector, now of dimension∏r=1RMr\\prod\_\{r=1\}^\{R\}M\_\{r\}\. Likewise,p~\(1\)\\tilde\{p\}^\{\(1\)\}has dimensionM~×M~\\tilde\{M\}\\times\\tilde\{M\}and entries in\(0,1\)\(0,1\)\.

Because the rule sequence is periodic andK=m​RK=mR, the full expansion can be grouped intommidentical blocks:

l\(K\)=\(⨂r=1Rl\[r\]\)⊗⋯⊗\(⨂r=1Rl\[r\]\)=\(l~\(1\)\)⊗m\.l^\{\(K\)\}=\\Bigl\(\\bigotimes\_\{r=1\}^\{R\}l^\{\[r\]\}\\Bigr\)\\otimes\\cdots\\otimes\\Bigl\(\\bigotimes\_\{r=1\}^\{R\}l^\{\[r\]\}\\Bigr\)=\(\\tilde\{l\}^\{\(1\)\}\)^\{\\otimes m\}\.The same argument gives

p\(K\)=\(p~\(1\)\)⊗m\.p^\{\(K\)\}=\(\\tilde\{p\}^\{\(1\)\}\)^\{\\otimes m\}\.Therefore the multirule generator is exactly representable as anmm\-level recursion of the single effective basis\(l~\(1\),p~\(1\)\)\(\\tilde\{l\}^\{\(1\)\},\\tilde\{p\}^\{\(1\)\}\)\. ∎

Obstruction to finite single\-basis reduction:If a multirule recursive generator does not admit a decomposition into repeated identical finite blocks, then it cannot be represented exactly by a finite\-dimensional single basis repeated across scales\.

###### Proof\.

Any finite\-dimensional single basis repeated across scales induces a rule sequence obtained by repetition of one finite block\. The claim follows by contrapositive\. ∎

###### Theorem 1\(Minimal effective branching number\)\.

Consider a multirule recursive generator of depthKK,

l\(K\)=⨂k=1Kl\[k\],p\(K\)=⨂k=1Kp\[k\],l^\{\(K\)\}=\\bigotimes\_\{k=1\}^\{K\}l^\{\[k\]\},\\qquad p^\{\(K\)\}=\\bigotimes\_\{k=1\}^\{K\}p^\{\[k\]\},where

l\[k\]∈int⁡\(ΔMk−1\),p\[k\]∈\(0,1\)Mk×Mk\.l^\{\[k\]\}\\in\\mathrm\{int\}\(\\Delta^\{M\_\{k\}\-1\}\),\\qquad p^\{\[k\]\}\\in\(0,1\)^\{M\_\{k\}\\times M\_\{k\}\}\.Assume that the rule sequence admits an exact single\-basis reduction\. Letq∗q\_\{\\ast\}be the smallest positive divisor ofKKsuch that the sequence

\(\(l\[1\],p\[1\]\),…,\(l\[K\],p\[K\]\)\)\\bigl\(\(l^\{\[1\]\},p^\{\[1\]\}\),\\dots,\(l^\{\[K\]\},p^\{\[K\]\}\)\\bigr\)can be partitioned intoK/q∗K/q\_\{\\ast\}consecutive identical blocks of lengthq∗q\_\{\\ast\}, up to one consistent permutation of community labels\.

Then the smallest effective branching number among all exact single\-basis reductions is

M~min=∏k=1q∗Mk\.\\tilde\{M\}\_\{\\min\}=\\prod\_\{k=1\}^\{q\_\{\\ast\}\}M\_\{k\}\.Moreover, the corresponding minimal effective basis is

l~min\(1\)=⨂k=1q∗l\[k\],p~min\(1\)=⨂k=1q∗p\[k\],\\tilde\{l\}\_\{\\min\}^\{\(1\)\}=\\bigotimes\_\{k=1\}^\{q\_\{\\ast\}\}l^\{\[k\]\},\\qquad\\tilde\{p\}\_\{\\min\}^\{\(1\)\}=\\bigotimes\_\{k=1\}^\{q\_\{\\ast\}\}p^\{\[k\]\},up to permutation of theM~min\\tilde\{M\}\_\{\\min\}effective labels, and satisfies

l\(K\)=\(l~min\(1\)\)⊗K/q∗,p\(K\)=\(p~min\(1\)\)⊗K/q∗\.l^\{\(K\)\}=\\bigl\(\\tilde\{l\}\_\{\\min\}^\{\(1\)\}\\bigr\)^\{\\otimes K/q\_\{\\ast\}\},\\qquad p^\{\(K\)\}=\\bigl\(\\tilde\{p\}\_\{\\min\}^\{\(1\)\}\\bigr\)^\{\\otimes K/q\_\{\\ast\}\}\.

###### Proof\.

Letq∗q\_\{\\ast\}be the shortest repeating block length\. Grouping each block into a single effective refinement step gives

l~min\(1\)=⨂k=1q∗l\[k\],p~min\(1\)=⨂k=1q∗p\[k\],\\tilde\{l\}\_\{\\min\}^\{\(1\)\}=\\bigotimes\_\{k=1\}^\{q\_\{\\ast\}\}l^\{\[k\]\},\\qquad\\tilde\{p\}\_\{\\min\}^\{\(1\)\}=\\bigotimes\_\{k=1\}^\{q\_\{\\ast\}\}p^\{\[k\]\},which have effective dimension

M~min=∏k=1q∗Mk\.\\tilde\{M\}\_\{\\min\}=\\prod\_\{k=1\}^\{q\_\{\\ast\}\}M\_\{k\}\.Since the sequence consists ofK/q∗K/q\_\{\\ast\}identical consecutive blocks,

l\(K\)=\(l~min\(1\)\)⊗K/q∗,p\(K\)=\(p~min\(1\)\)⊗K/q∗\.l^\{\(K\)\}=\\bigl\(\\tilde\{l\}\_\{\\min\}^\{\(1\)\}\\bigr\)^\{\\otimes K/q\_\{\\ast\}\},\\qquad p^\{\(K\)\}=\\bigl\(\\tilde\{p\}\_\{\\min\}^\{\(1\)\}\\bigr\)^\{\\otimes K/q\_\{\\ast\}\}\.Thus an exact single\-basis reduction with branching numberM~min\\tilde\{M\}\_\{\\min\}exists\.

Now let

l\(K\)=\(l^\(1\)\)⊗J,p\(K\)=\(p^\(1\)\)⊗Jl^\{\(K\)\}=\(\\hat\{l\}^\{\(1\)\}\)^\{\\otimes J\},\\qquad p^\{\(K\)\}=\(\\hat\{p\}^\{\(1\)\}\)^\{\\otimes J\}be any other exact single\-basis reduction\. Then the original rule sequence must decompose intoJJidentical consecutive blocks, each of lengthq=K/Jq=K/J\. By minimality ofq∗q\_\{\\ast\}, one hasq≥q∗q\\geq q\_\{\\ast\}\. Therefore no exact reduction can use a shorter repeating block thanq∗q\_\{\\ast\}, and the construction above attains the smallest possible effective branching number\. Hence

M~min=∏k=1q∗Mk\.\\tilde\{M\}\_\{\\min\}=\\prod\_\{k=1\}^\{q\_\{\\ast\}\}M\_\{k\}\.∎