Directed Graph Generation with Heat Kernels
Authors: Marc T. Law
Karsten Kreis
Haggai Maron
Venue: Accepted by TMLR ()

Abstract: Existing work on graph generation has, so far, mainly focused on undirected graphs. In this paper we propose a denoising autoencoder-based generative model that exploits the global structure of directed graphs (also called digraphs) via their Laplacian dynamics and enables one-shot generation.  Our noising encoder uses closed-form expressions based on the heat equation to corrupt its digraph input with uniform noise. Our decoder reconstructs the corrupted representation by exploiting the global topological information of the graph included in its random walk Laplacian matrix. Our approach generalizes a special class of exponential kernels over discrete structures, called diffusion kernels or heat kernels, to the non-symmetric case via Reproducing Kernel Banach Spaces (RKBS). This connection with heat kernels provides us with a geometrically motivated algorithm related to Gaussian processes and dimensionality reduction techniques such as Laplacian eigenmaps. It also allows us to interpret and exploit the eigenproperties of the Laplacian matrix. 
We provide an experimental analysis of our approach on different types of synthetic datasets and show that our model is able to generate directed graphs that follow the distribution of the training dataset even if it is multimodal.

Published in Transactions on Machine Learning Research (01/2025)
Directed Graph Generation with Heat Kernels
Mare T. Law NVIDIA
marel@nvidiaacommmmmmmmm com
Karsten Kreis
kkreis com
NVIDIA
Haggai Maron
honaron @nvidia. com
NVIDIA
Technion
Reviewed OpenReview:
openreview. https: h net/ forum? id= 600 lubhte
Abstract
Existing work on graph generation has, so far, mainly focused on undirected graphs. In this paper we propose a denoising autoencoder-basee generative model that exploits the global structure of directed graphs (also called digraphs) via their Laplacian dynamics and enables one-shot generation. Our noising encoder uses closod-form expressions based on the heat equation corrupt its digraph input with uniform noise. Our decoder reconstructs the corrupted representation by exploiting the global topological information of the graph included in its random walk Laplacian matrix. Our approach generalizes a special class of exponential kernels over discrete structures, called diffusion kernels or heat kernels, to the non-symmetric case via Reproducing Kernel Banach Spaces (RKBS). This connection with heat kernels provides us with a goometrically motivated algorithm related to Gaussian processes and dimensionality reduction techniques such as Laplacian eigenmaps. allows us interpret and exploit the eigenproperties of the Laplacian matrix. We provide an experimental analysis of our approach different types of synthetic datasets and show that model is able to generate directed graphs that follow the distribution of the training dataset even if it. is multimodal.
also
1 Introduction
The representation of directed graphs (or digraphs) has recently attracted interest from the machine learning community (Clough & Evans 2017 Sim et al. 2021) they naturally doscribe causal relations (Bombelli
et al. 1987 spatiotemporal events using chronological order (Law & Lucas 2023 or some stochastic processes such as Markov chains (Norris 1998 Existing work in digraph representation has focused on discrininative tasks such as classification or link prediction. In this work, we consider the task of digraph generation. One possible application is the modeling of new causal systems that follow the same distribution as some given training set.
Most of the machine learning literature on graph generation focuses on undirected graphs (Liao al. 2019 Nuu et al. 2020; You et al, 2018). Their goal is to generate plausible graphs of the same nature as those from some given training dataset (e.g., molecules (Vignac et al. 2023)). ). Existing approaches can be divided mostly into two categories: aufo-regressive and one-shot. Auto-regressive approaches (Liao et al. 2019 You et al.) 2018) start by generating small graphs to which sets of nodes and their corresponding edges are iteratively added until the final graph reaches a certain criterion (e.g., size). On the other hand, one-shot single step. Ono-shot approaches
approaches generate all the nodes and edges of the generated graphs were shown to be more efficient than rego-segees ones due the lack of intermediate steps that can also lead to worse generative performance because of error accumulation at each step and the fact that one-shot.

---

Published in Transactions on Machine Learning Research (01/2025)
Adjacency matrix A
Edge perturbation
Hoat diffusion encoder
(Closed-form solution)
Noisy representation X(2)
Learued docoder (Neural netwcark))
Reconstruction
0
0
0
0
1
0
0
1
0
0
0
Data
0
0
Figure 1: Our framework can be viewed denoising autoencoder. Our heat diffusion encoder maps a perturbed adjacency matrix À E {0, 1]mxm to noisy node representation matrix X(T) E [0, 1]mxd that is given as input of decoder that reconstructs the edges. (n = 5, d = 7 in the figure)
approaches can directly work with the global structure of the graph via ils Laplacian matrix instead of arbitrary subgraphs (Martinkus et al. 2022 Among the one-shot approaches, Spectre (Martinkus et al.) 2022h learns the distribution of the most informative eigenvectors of the Lapacian matrices. Nonetheless, Spectre does not generalize to digraphs as it relies on properties of Laplacian matrices that are satisfied only for undirected graphs (see explanation in Section 5 DiGress (Vignac al. 2023) considers separate representations for nodes and edges which discrete noise added. DiGress is formulated
classification problem such that denoising decoder classifies the category or existence of edges and nodes. However, DiGress also requires spectral features from [Bainn et al. (2021) that are valid only for undirected graphs they rely on symmetric scalar products. In conclusion, none of the existing onc-shot approaches can be casily adapted to digraphs.
Contributions. We propose a one-shot approach that generates digraphs in a single step. The training of our framework illustrated in Fig. [ Unlike previous approaches, we exploit the cigenproperties of the Laplacian matrix that are valid even when the graph is directed, so we can effectively use its global structure during the generation process. To this end, we propose denoising autoencoder (Vincent et al. 2008 whose noising encoder is not learned by neural network but exploits closed-form expressions based the heat equation for digraphs (Veerman & Lyons 2020), which effectively encodes the global topological information of the graph into node features. We propose to introduce noise via a nonhomogeneous term that makes our node represontations tend to a stochastic matrix with all its clements equal. Our denoising decoder is a neural net work trained to reconstruct the original node representations and adjacency matrix of the graph. We explain below the main setup and building blocks of our approach.
1.1 Problem definition and overview of our approach
During training, we are given a set of digraphs (Gi,...; Gm represented by their adjacency matrices Aj...., Am- Our goal at inference time is to generate new digraphs that follow the same distribution.
Training. At training time, a training graph G, is represented by its adjacency matrix A, First, this matrix perturbed using data augmentation techniques such as edge perturbation (Ding et al. 2022 The resulting perturbed adjacency matrix A is used as input of our first main component, called Heat Diffusion Encoder, which generates some noisy representation, called X(T). Our encoder produces outputs that are close to uniform distribution in order to approximate maximum entropy (or lack of information). A detailed description of the heat diffusion encoder can be found in Section 3 Next, our second component, called Denoising Decoder, takes X(") as input and is trained to predict an adjacency matrix that is as similar to the original (unperturbed) matrix A, as possible. The whole training process is illustrated in Fig. 1
Sampling new graphs after training. At inference time, an adjacency matrix is randomly generated and given as input of our autoencoder to generate novel digraph. We discuss other continuous sampling strategies in Appendix G.2
2

---

Published in Transactions on Machine Learning Research (01/2025)
2 Preliminaries on the Heat Equation
This section introduces our notation and a self-contained introduction to heat diffusion on digraphs. These concepts will be used in Section 3 to define a closed-form expression for our heat diffusion encoder. We
refer the reader to Veerman & Lyons (2020) for more details on Laplacian matrices of digraphs, to Chung & Yau y Kondor & Lafferty| (2002); Belkin & Niyogi (2003) for heat diffusion kernel-based methods on undirected graphs, and to Appendix for the details of the equations of this section and the next section. One main difference with other types of data such as images is that graphs lie in discrete space whose topological information depends on the adjacency of its nodes. Moreover, graphs are often sparse in terms of adjacency. Due to the discrete nature of graphs, calculating similarities between nodes has been a challenge. To tackle this, heat kernel approaches (Chung & Yau. 1999 Kondor & Lafferty| 2002) express the notion of discrete local neighborhood of nodes in terms of a global similarity over the nodes of an undirected graph. More exactly, heat kernels are generalizations of Gaussian kernels where the Euclidean distance that is used in Gaussian kernels is replaced by more general distance that takes into account the neighborhood information of the Laplacian matrix (Kondor & Lafferty 2002). The information in node features is propagated along their neighbors via the heat equation. We propose to add nonhomogencous term that introduces noise along the heal diffusion process. We explain in Appendix D how our approach generalizes heat kernels to digraphs.
Notation. We denote the identity matrix by I, and the all-ones vector by 1. We consider a graph G = (V, E) defined by its set of n nodes V = {er}"= and its set of odges E S V V. Its adjacency matrix A E {0, 1 satisfies A. = 1 ifl (vi, Uz) E E and A, = 0 otherwise. I[ G is undirected, then we can simply consider that A is symmetric (i.e., ji = At). However, we consider the more general case where A is not constrained to be symmetrie. Although optional, we add self-loops by constraining A to satisfy Vi, Ai have D == diag(1 A). R = is defined that Vi, E, Ajs. In other words,
The in-degree
diagonal matrix D
Let us define the matrix S == AD-1. S column stochastic (i.c., Vi, Ei,E, E, = and Vi,j, S, 2 0). We consider in this section that we are given some matrix N == (() E Raxd where the i-th row of N is the initial d-dimensional feature representation of e (i.e., X(t) with & = 0). The matrix N could be arbitrarily defined or given. In practice, we train N jointly with our denoising decoder, and we explain its training process in Section E and Appendix B In this section, we consider that. N fixed.
Laplacian dynamics. We define the negative of the random walk Laplacian matrix as -D The matrix L can be viewed as matrix form of the discrete Laplace operator which approximates the continuous Maplace-Beltrann differentiable real-valued function f defined on some manifold, the Laplace-Beltrami as the divergence of the gradient g and provides us with an estimate of how far apart f maps nearby points. Since L is not symmetric in general, we can use both L or L as they have different left and right cigenvectors. We then denote A E (I, In the main paper, we consider only the case A = L,, which is called diffusion model in the graph community (Veerman & Lyons 2020). To avoid confusion with diffusion models in the machine learning literature, we call it heal diffussio.. We give in Appendix Dal the formulae to solve our problem when A L (called consensus model (DeGroot 1974
operator in differential geometry (Belkin & Niyogi 2003). Given a twice-
Beltrami operator is defined
Heat equation. The global information of the graph G is diffused via the following heat equation:
d x(t) AX(t)= Q(t) where X(t) is the representation of nodes at time t 2 0, V 0, dt
(1)
d is a heat source term that introduces noise in the node representations. If Vt 2 0, (t) o, then Q is called homogeneous. It is called nonhomogeneous otherwise. The solution of
N X(0), and Qtt) E Rnd
equation 1 can be found in standard textbooks such as Edwards et al. (2020) and is written in equation 2
For any formulation of Q equation
I is solved by the following equation when A is constant over time:
Vt 2 0, X = c' x()) + L' clt Q(s)dss = Z(t) + Ft)).
(2)
where e denotes the matrix exponential of the matrix tA. To simplify notation, we define Z(t) == eAx((0). and F(!) = Sh' Q(s)ds
3

---

Published in Transactions on Machine Learning Research (01/2025)
Interpretation of the heat equation. The above heat equation can be interpreted as continuous way of performing message passing in the graph G over continuous steps by exploiting the neighborhood information included in A. Nonetheless, our continuous message passing approach introduces noise due to the presence of the nonhomogencous term.
Each column of N (0) contains some initial signal of the nodes, and the information of the Laplacian matrix and noise are jointly diffused in those signals as increases by following the heat equation in equation 2 One difference with standard diffusion processes is the use of the global topological information of the graph
via e' over time t. The noise is introduced via the term F(t). If Q is homogeneous, then Vt 2 0, F() = 0
(i.c., there is no noise) and equation 2
reduces to Vt 2 0, (t) Z(t).
3 Heat Diffusion Encoder
As in standard denoising autoencoders (Vincent et al. 2008), we define a noising process that maps ((0) to an informative representation close to some analytically tractable distribution (from which we could sample at inference time). To this end, we consider that equation
2 is the output of our heal diffusion encoder that is given x()) and A as input, and it constructs some noisy representation X(T) where T > is an arbitrary time constant defined in Section 4 In Proposition [ we propose to define the heat source term Q that X(T) is similar to some constant matrix that call M and that has all its elements equal to the same value in order to approximate maximum entropy (or lack of information) in the node representations when t = T. In Section 4 we train a decoder that reconstructs the nodes and edges when given some noisy X(T).
Formulation of the heat source term. Our goal is to formulate the heat source term Q so that X(T) tends to some non-informative matrix M as T tends to too. To this end, we first notice that the matrix e 'A is column stochastic for all 2 0 (see explanation Appendix F.). In order to have Z(t) = e also column stochastic for all i 0, we add the constraint that N ((0) is column stochastic (see proof in Appendix C.3) Since each column of the resulting node representations is a probability distribution vector, define each column of M as the uniform probability distribution vector, which corresponds to the maximum entropy probability distribution vector. In other words, we define the column stochastic uniform noise matrix as M : 111T E { By definition, the matrix M is constant. Each column of M also corresponds to the expected value of a random variable following a flat Dirichlet distribution.
Proposition 1. To satisfy limm X(T) M where M
1 111T, we formulate Q and F as follows: n
Q(s) = aa - (0) = F(1) = (1 - atte (R (Ax(0)
(3)
where >>0 0 is a noise diffusivity rate byperpavamneter, B 2 0 is hyperparameter that can be tuned to control the Laplacian dynamics further, and R := e CTAM E Raxd is a constont matrix for some arbitrary time constant > 0 defined in Section 4 See Appendic C for details and proofs.
With the above proposition, X(t can be written only as a function of X(0) and A in equation ?
Vt 0, (t) = et (X(0) + (e-at- 1) + (1 e -t R)
(4)
If we set 3 3=0 0, then .ع = I by definition of the matrix exponential, and we get simpler formulation of equation 4
0 = x(,(() = e - + - e ea)e AR.
In this case, we have the following simple formulation when = T::
0= - = - Z(T) + (1 - -T M
(5)
which is column stochastic, and we call 1 -aT the noise ratio time T. Since we have limr- e-a2 of 0, one can verify that Proposition I
satisfies our initial goal because equation 5 implies lim7 X(T) M.
Backward formulation. In the following, we consider that 3 = -0. However, our approach can be generalized to any B 2 0. A is given, we can write X(t) as a function of X(! for all time step T 2 0:
0=
0, x(1) at (1 + T) + (1 e A
(6)
4

---

Published in Transactions on Machine Learning Research (01/2025)
However, we assume that the denoising decoder that we train in Section E does not have access to A at inference time. Therefore, the decoder has to reconstruct the set of nodes and/or edges when given only X(T).
Efficient forward process. As explained in Vignac et al. 2023 we can build an efficient generative model for the following reasons:
The noisy representation X(t) has a closed-form expression that depends only on X(0) and A (see equation 4h. We can then directly calculate X(t) for all 2 0 without adding noise iteratively.
The ground truth node representation (t) can be written in closod-form when given only A and cither (0) or ((). In practice, we generate a perturbed adjacency matrix A, from which we generate its Laplaciann e- JM that is given input of a node decoder that has to reconstruct the ground truth node representation Z(t) c'Ax(()) for some & E [0, T7. This data augmentation technique acts as a regularizer that reconstructs node representations similar to the training distribution when given a noisy input.
- Á and noisy representation X(T) c (())
The limit distribution limm-++00 X(T) = M does not depend on X()). It is also worth noting that all the elements of X(T) are in the interval [(1 ([1-e-aT) /n, (1+ 1)e - /n]. In practice, we choose 0 so that sufficient information of the graph is preserved in the noisy
appropriate values of T > 0 and node representation, and so that denoised edges can be recovered from it. Specifically, we diffuse toward a distribution that can be well approximated by an analytic distribution (e.g., we can sample from (lat) symmetric Dirichlet distribution (Kotz et al. 2004)) while preserving sufficient information about X(0) to perform denoising. Morcover, ((t) is column stochastic when t == 0 and 2 T, but ((t) might contain some negative elements when & E (0, T) due to the formulation of R. This is not a problem in practice since our goal is to reconstruct (t) which is column stochastic for all 2 0.
4 Denoising Decoder and Sampling
The previous section defines a noising heat diffusion encoder that does not require learning neural network and is given an initial node representation matrix N = =(0) and a Laplacian matrix to generate in closed form some noisy representation al some arbitrary time T > O. We now propose a multi-task learning formulation to train decoders that reconstruct the original input. We assume that the decoders are not directly given the Laplacian matrix at inference time.
Our first. task learns neural network (called node decoder) that predicts the denoised node representation Z(t) where & E [0, (we recall that Z(0) X(0).. Our second and last task jointly learns another neural network (called edge decoder) to predict edges. This is similar to the way link prediction tasks are solved, and we observe in practice that the learned representations hold information from the graph in the form of the Laplacian singular vectors. In Section inference time.
4.3
we propose a sampling strategy to generate digraphs at
4.1 Training of the denoising decoders and data augmentation
given digraphs {G, = E))m drawn from some distribution G and of different sizes. Our goal is to generate graphs that follow the same distribution. Each graph G is represented by its adjacency matrix A E (0,1)%.xm, where n, == Vil is the number of nodes. The Laplacian matrix of G, is - I A (diag(11At))
Setup. We consider the task where, during training,
Node representation. We explain here how define the node representations of the training graphs. Let USS note Nlma == max; n, the number of nodes of the largest graph in the training set. We define some matrix O E Rmax where d > 0 is an arbitrary hyperparameter. For each graph G;, we define its column stochastic initial node representation matrix N = -X) E [0,1]m;xd as the uppersubmatrix of O whose columns are n-normalized with the softmax operator, which corresponds to applying mask and renormalizing. To simplify the notation, we write N instead of N since we consider that all the graphs of same size share the same initial node representation. N is trained by training O (see details about their training in Appendix B
5

---

Published in Transactions on Machine Learning Research (01/2025)
Algorithm 1 Generation of digraphs at inference time
Input: Node representations O C ["mmaxd, hyperparameters T, > O, J.PC [0,1]
1: Sample a 9max- Define N E Rexd upper submatrix of O followed by softmax operation for é--normalized columns
2: Generate discrete adjacency matrix A E {0. 13mxa such that vi # J. Ag Bernoull((x) and Vi, Aa =
3: Apply data augmentation to obtain perturbed matrix A (e.g., A = A @ C s.t. Vi & j. Ci ~ Bernoull(pp)
4: Calculate the diagonal matrix D C R4*" such that D;; = E, Aza- Define á = ÁD
5: Define B as erèc or optionally as the rank-s approximation of eTá via truncated SVD.
6: Give the matrix e- BN + (1 (1-e-"))M as input of the edge decoder that returns an adjacency matrix
We arbitrarily define the values of >>> 0 and >> 0 so that the noise ratio defined as (1 - e-aT) is close to (see Fig. 2 Following equation 5 we define X²(T) == - N (1- aa)MM In practice, we apply
data augmentation on (T) during training explained below.
Edge perturbation (Ding et al., 2022). One can apply data augmentation via edge perturbation which can be interpreted as injecting noise by considering the perturbed adjacency matrix A A @ C instead of A', where @ is the logical XOR operator, and the zero-diagonal corruption matrix E {0, [0,1]" has its non-diagonal elements equal to 1 with probability & E [o, I], and 0 with probability (1 - p). Following Veličković al. (2019). : we set p 1/n; and we sample a new matrix C each time C is called. We call A' = I A' (diag( 1 A the Laplacian from Á We obtain the formulation (T) = e eT'N, N, and ('(T) = = 2 (T) + (1 - e M. p =0 0, we have -=)))=X+0)
Permutation of adjacency matrices. N is the same for all graphs of same size, To promote permutation invariance of our model, we can replace A' and A by P A and P A where P E {0, []"x,,, is a (randomly sampled) permutation matrix. This is equivalent to replacing e eT'' by ETA'PN, This can be seen as augmenting the training set with adjacency matrices of isomorphic digraphs. In Appendix G.3 we experimentally show that using this kind of data augmentation technique does not have a negative impact on the optimization of our loss function.
X is given as input of a node decoder yp and edge decoder ab during training as described below.
Node decoding task. Our node decoder lakes the noisy node representation {((T) as input, and its goal is to reconstruct some target node representation matrix T that does not contain noise. In practice, we formulate the training loss of our node decoder as Lnoce (i) == II)((' T where || || is the Frobenius norm, and we arbitrarily define T := (L) = e N where AA' is the ground truth Laplacian matrix. Since each row of X'(T') represents node of the graph, we ideally want our model to be equivariant to the order of the rows of (T).. For this reason, we formulate as an ceanntion-baser permutation-ines invariant neural network called Set Transformer (Lee et al. 2019). &p considers each row of (T) as the element of a set. of node representations, and it is robust to the order of the rows. Implementation details can be found in Appendix B
Edge decoding task. We call our edge decoder vè and we denote the p-th row of (X)( by (xX
Our edge decoder predicts whether or not there exists a directed edge between pairs of nodes. We formulate the term:
Eppg H (w ([x(((())) (T)), v6 (X+(()) A Apy)) where [+]] denotes concatenation, w is a learned multilayer perceptron (MLP), and И the cross-entropy loss. Our edge decoder and node decoder share a common backbone (see architecture details in Appendix B. It is worth noting that if the goal is generate undirected graphs, then the concatenation operation can be replaced by symmetric operation such the addition. The training loss that we minimize is:
Ees Lodge (1) +y Liode
where y 2 a regularization parameter. Since both T and (T) depend on N, we optimize equation 7 by training jointly sp, v and N via gradient descent in order to reconstruct the training graphs. See Appendix for implementation details.
6

---

Published in Transactions on Machine Learning Research (01/2025)
4.2 Other considerations
Class-conditional generation. To add class label information, we give input of both decoders the concatenation of a matrix Y' E (0,1)*xXCC to (T) where |C] is the number of categories, and each row of Y' is one-hot vector whose nonzero index corresponds to the category of the graph G,. This sampling strategy is known as conditional sampling (Zhu et al. 2022).
The rest of the method is similar.
Choice of the final step T. The matrix (T) is given as input of decoders to reconstruct G,. We ideally
want (T) e- =e (T)) to the matrix M. This similarity depends
+ (1 e1-e-ad")M M to be similar
both T and a, and (+()) tends to M as T or a tend to +00.
0.8
f( 1-exp--x:1
We provide a detailed discussion about the impact of T and a in Appendix F We found that setting T = 1 and choosing a large enough works well in practice (e.g., a = 2.3 implies - e R 0.9, which means that about 90% of the values of X'(T) are noise). However, the optimal value of both T and a can be determined via cross-validation depending the task. Fig, illlustrates the ratio of noise for different values of a as function of T. It is worth noting that we want (T) to be similar to M so that sampling a similar matrix at inference time is easy. On the other hand, we also want (T) to
0.6
0.4
a
a
0.2
0.5
0.25
0.1
0
2
6
Figure 2: The noise ratio - e aT in X(T) as a function of for different values of a.
preserve enough information so that our neural networks can reconstruct T and A (i.e., the node and edge information of the graph) from it.
Learning N. It is known in the heat kernel literature (Belkin & Niyogi 2003 Chung & Yau 1999 that the coarse structure of an undirected graph is included in the subset of eigenvectors of its Laplacian matrix that correspond to its smallest eigenvalues and can be used for dimensionality reduction. PageRank (Page 1998 uses the same motivation and exploits the leading cigenvector of a nonsymmetric stochastic matrix. Spectre (Martinkus et al., 2022) exploits this observation by generating a symmetric Laplacian matrix spanned by a set of leading eigenvectors with bounded real eigenvalues. Since our eigenvalues and eigenvectors are usually complex and not unitary (when the adjacency matrix is not symmetrie), we consider related linear algebra properties such as column spaces and singular vectors, which allow us to work with real values. We experimentally observe in Sention 6.1 that the leading singular vectors of the learned matrix (1) = e' N and of e tend to be strongly correlated, which suggests that our model learns N so that it maps to the leading left singular vectors of e' This observation is mainly qualitative. We also observed that using large number of columns d to represent. N helps in practice to recover edges.
A Approximations. We show in Appendix E that our approach corresponds to a non-symmetric heat kernel method induced in a Reproducing Kernel Banach Space (RKBS) with a column stochastic and non-symmetric kernel matrix K e-aTT + I-c n; 0 (see Appendix F), and it then contains the structure of the graph. To scale our method to large graphs, we also propose to replace e'A' by its rank-s approximation obtained with truncated Singular Value Decomposition (SVD) where I n4. We replace e by the product of two rectangular rank-s matrices, which greatly reduces memory if 3 < n4. Although the obtained matrix is not stochastic when < n4, we observe in Section 6.2 that the training graphs can be fully reconstructed with our model while saving memory usage. The SVD can be preprocessed offline before training. Gaussian kernels are specific heat kernels that use a squared Euclidean distance, and the relationship between heat kernels and Gaussian Processes (GPs) is known in the literature (Kondor & Lafferty 2002). Using a truncated SVD corresponds to a low-rank approximation of the kernel matrix w.r.t. the Frobenius norm, which is similar to one of the approximation methods for GPs mentioned in Williams & Rasmussen [2006 Chapter 8.1).
11T, K has the same eigenvectors
the Laplacian -A' when

---

Published in Transactions on Machine Learning Research (01/2025)
4.3 Sampling strategy for generation of novel digraphs
Our sampling algorithm to generate novel digraphs is given in Algorithm while our training algorithm is detailed in Algorithm 2 During sampling at inference lime, we do not have access to inpul graphs from the dataset. We then construct graphs that, when given to our heat diffusion encoder, diffuse toward noisy graphs similar to those encountered after difffusing graphs during training. In this way, our denoising decoders can successfully produce denoised graphs similar to those of the dataset. This is conceptually similar to variational autoencoders (Kingma & Welling 2014), where during sampling the encoding distribution is approximated by simple prior distribution. How can we then analytically construct suitable input graphs during inference time? One solution is to generate matrix with each column sampled from flat Dirichlet distribution (Kotz et al. 2004) and give it as input of the decoders to generate a digraph. This works well when the training graphs (of a given category) are all similar to each other. However, it was observed in Vignac et al. (2023) that this kind of continuous sampling tends to destroy the graph's sparsity and creates very noisy graphs in practice. When the distribution of the graphs is multimodal, we found that sampling discrete adjacency matrices and applying standard data augmentation techniques for graphs both during training and sampling allows our model to sample graphs from the different modes. We study continuous sampling strategies in Appendix G.2
We propose one discrete sampling algorithm in Algorithm D Let us note μ E (0, 1J the ratio of pairs of distinct nodes that are adjacent in the training set. We first generate an adjacency matrix A E {0, 1|nx= such that each of its non-diagonal elements is assigned the value 1 with probability j4, and 0 otherwise. Following the motivation of denoising autoencoders, our decoders are trained to construct an (unperturbed) sample similar to training samples when given some noisy input. In Appendix G.2 we observe that our discrete sampling strategy is competitive with other continuous sampling strategies in terms of performance.
5 Related Work
Graph generalive approaches can be divided into two calegories which are uulo-regressive models and one- shot models. Auto-regressive models (Liao et al. 2019 You et al. 2018) generate a succession of graphs G1, G2...,, GT such that Vi, G, C Gi+1 and return the last generated graph Gr. At each iteration, the graph G2 is given as input of a neural network that generates Gi+1 by adding new nodes and their edges. Most of these models are typically slower than one-shot approaches that generate all the nodes and edges of graph in a single step. Three main one-shot approaches in the machine learning literature are Top-n (Vignac & Frossard 2022), Spectre (Martinkus al. 2022) and DiGress (Vignac et al. 2023). Other one-shot methods such as (Kwon et al., 2020 Mercado et al. 2021) are dedicated to molecular graphs and do not generalize to other tasks. Although Top-n is onc-shot, it assumes symmetric similarity functions between nodes inappropriate for digraphs.
Spectre (Martinkus et al. 2022 considers the generation of undirected graphs via their normalized graph Laplacian matrix L, := I D 1/2 A22AD-1 1/2, which is symmetric positive somi-definite and admits an cigen- decomposition of the form U- AU where U-1 UT and both the diagonal matrix A and U are real-valued. They exploit the intuition that coarse structure of the graph lies Stiefel manifold that contains the eigenvectors of the & smallest eigenvalues of the Laplacian. Spectre Martinkus et al. (2022) then trains a neural network that generates an adjacency matrix by sampling the i smallest eigenvalues and their corresponding eigenvectors by exploiting their Stiefel manifold structure. The authors mention that their work can be extended to the random-walk Laplacian matrix L, := I AD-1 since its right eigenvectors are the same as D-1/2U (up to column-wise l2 normalization). and its left eigenvectors can be formulated in a similar way. However, when the graph is directed and A is not symmetric, U is complex and not unitary. The information of the Laplacian matrix then does not lie complex Stiefel manifold, and Spectre can then not easily be extended to digraphs.
DiGress [Vignac et al. 2023 is denoising diffusion model for graphs. Instead of using the discrete Laplacian operator as we propose, they represent their nodes as a function of time E N as follows: x(1) (o'I [=-0)) 11 X(( 1) where X(I) is row-stochastic and (0, 1). DiGress relies on spectral
8

---

Published in Transactions on Machine Learning Research (01/2025)
(left) Erdös-Rény: graphs (right) stochastic block model graphs
Figure 3: Correlations between the left singular vectors of e' and e' N (more results in Fig. 4
features from Beaini et al. 2021)) designed for undirected graphs by using a symmetric similarity function between nodes. It is then not appropriate for digraphs.
Diffusion processes have been used in basee generative models (Ho al.
2020
Sohl-Dickstein et al,,
2015 Song et al. 2021 to generate undirected graphs (Jo et al. 2022 Niu et al
2020
adjacency matrices exploit Graph Neural
in general (Yan et al. 2024). Nonetheless, both Jo et al. (2022] and Niu et al. Network architectures that require a symmetric adjacency matrix so these methods are not casily adaptable to digraphs. Moreover, both approaches perturb each data dimension independently (i.e., without accounting for the global structure of the graph via its Laplacian). On the other hand, SwinGNN (Yan et al. 2024] considers diffusion models to reconstruct arbitrary adjacency matrices that do not have to be symmotric. Generally, in such score-based generative diffusion models, the input converges towards an entirely uninformative pure noise distribution, and generation requires a slow and iterative denoising process. In contrast, our approach encodes the global structure of the graph into the node representations via the Laplacian dynamics and effectively encodes the graph topology in a noisy, yet informative distribution when t = T, from which the denoised nodes and edges can be elficiently predicted in one shot. In other words, in contrast to score-based diffusion models, where the role of the diffusion process is purely to gradually perturb the data and destroy information, in our model the primary role of the (Laplacian) diffussion process is to additionally encode the graph structure into the node representations. Heat diffusion for the purpose of information propagation in that fashion has, for instance, been used for learning on curved 3D surfaces (Sharp al. 2022 Our method is a different approach compared to the widely used scoro-based diffusion models and specifically dosigned to avoid their slow, ileralive synthesis process.
2020
Instead, our model can be seen as a denoising autoencoder (Vincent et al., 2008), since it corrupts the input data with optional edge perturbation and with a nonhomogoncons heat diffusion process. During sampling, we also approximale the encoding distribution at ( T by using an analytically tractable distribution, which is similar variational autoencoders (Kingma & Welling 2014), where simple prior distribution in latent space models the encoding distribution. Hence, our approach can be scen related to GraphVAE (Simonovsky & Komodakis 2018), which, however, uses learned encoder neural networks instead of a heat diffusion process to encode small graphs in a latent space. Moreover, GraphVAF only tackles undirected graph synthesis. We emphasize that connection by proposing a sampling method that directly samples from a flat Dirichlet. distribution in Appendix G.2
It worth noting that was shown in Vincent (2011) that score matching techniques (Hyvärinen & Dayan 2005) can be seen as training specific type of denoising autoencoder. Our approach could be adapted to a score- based generative model with iterative synthesis. This could be done by replacing our deterministic nonhomogeneous heat source term by Wiener process using Dirichlet diffusion score done for example in Avdeyev et al. (2023). However, the method would be computationally expensive for large number of nodes since (1) each node would be associated with a different beta distribution, (2) there exists no known closed-form solution for this kind of linear stochastic differential equation, (3) and each intermediate time step would require calculating some matrix exponential involving A in our case. Instcad, we propose a framework similar to Martinkus et al. (2022) to efficiently work with singular vectors of the Laplacian matrix and generate (directed) graphs in one-shot manner.

---

Published in Transactions on Machine Learning Research (01/2025)
6 Experiments
We evaluate the generative power of our method that we call Digraph Generation with Diffusion Kernels (DGDK), focusing on digraphs only, for which our method was designed. We use Adam (Kingma & Ba 2014 as optimizer. DGDK works better as the number of columns d of N is larger. Experimental details can be found in Appendix Due to lack of space, we add other experiments in Appendix G
6.1 Column space of the learned representations
We adapt the intuition of Spectre [Martinkus et al. 2022) to digraphs. To this end, we experimentally observe in this subsection that our learned representations are strongly correlated with the leading singular vectors of e! In our first qualitative experiment, we apply conditional sampling (Zhu et. al. 2022) and consider two categories of digraphs: Erdös Erdös-Rényi (Erdös et al. 1960p with p 0.6 and stochastic block
model (Holland et al.) 1983) with 3 main blocks that are connected to cach other by following the transition probability matrix:
0.9
0.35
0.2
II
0.2
0.75
0.2
(8)
0.2
0.25
0.8)
The above transition matrix II between blocks cannot correspond to undirected graphs since undirected graphs would require that II symmetric. Each category contains m 100 training graphs with n 15 nodes cach, and we set the number of columns of N to d = 150. The size of the first two blocks is m 1 each and the size of the last block is n 2nn. In this setup, our training graphs are directed and their adjacency matrices are not symmetric. As explained in Section standard baselines such as Spectre cannot be exploited.
If we define T = N [0, 1j (i.e., we consider = 1), the goal of the node decoder in Section 4.1 is to reconstruct T whose column space is by definition included in the column space of e (i.c., spanned by its columns). T' is in general not a square matrix (i.e., n * d) and thus does not possess eigenvectors. We then consider its SVD e V and write the SVD of e' U2A2V The singular values of both A and A2 are ordered in descending order. Figure 3 illustrates the absolute values of the following matrix product U U1. Since both U1 and U, have their columns fo-normalized, each element of U.J U1 is the cosine between two singular vectors. cosine of O indicates orthogonality, hence independence, whereas higher absolute values indicate (cosine) correlations. The top left part of each plot corresponds to the leading singular vectors whereas the bottom right corresponds to the singular vectors with lower singular values. As one can see, there is strong cosine correlation between the leading singular vectors of e and of e' This suggests that N is learned to preserve the most informative singular vectors of c'A We exploit this observation in Section 6.2 by working with a low-rank approximation of !A'
6.2 Using approximations of the target matrix
We evaluate the generative power of DGDK in the class-conditional generation task. We use rank-s approximation of TA via a truncated SVD, and we formulate our target node matrix T' = Z'( where T= and a = 2.3. We consider the case where the number of nodes is Vi, n. = 21, and s = 15. The two categories follow the same propertics as in Section 6.1 and contain 3,000 non-isomorphic training graphs per category. During training, each mini-batch contains 10 graphs per category. Only the number of nodes per graph different. DGDK manages to reconstruct all the edges of the training graphs when given noisy representations (T)... This shows the effectiveness of our low-rank approximation approach.
We sample 10,000 test graphs per category by using Algorithm I with caas-conniiional generation (i.e., we provide the desired category as input). None of the generated test graphs are isomorphic to one another nor to the training graphs. This means that our model obtains a uniqueness and novelty scores of 100%. We also report in Table [ standard evaluation metrics based on the squared Maximum Mean Discrepancy (MMD) (O'Bray al. 2022) between the training set and the test set. These evaluation metrics mcasure the distance between the training and test distributions w.r.t. some graph properties. We adapt the descriptor
functions in O'Bray et al. (2022) to directed graphs.
10

---

Published in Transactions on Machine Learning Research (01/2025)
Table Squared MMD distances.
Dataset
Erdös-Rényi (p = 0.6)
Stochastic block model (3 blocks)
MMD metric
Degree
Clustering
Spectrum
Degree
Clustering
Spectrum
DGDK (ours)
1.1 10
1.0 10 3
1.3 10 6
1.2 10
2.0 10 4
6.2 10-6 G
SwinGNN Yan et 2024)
1.3 x 10-4
2.4 10-2
4.2 x 10-3
1.2 10-4
3.7 10-2
5.1 x 10-3
GRAN (Lioo et al. 2019)
1.5 x 10-4
5.6 x 10-2
1.9 x 10-4
1.3 x 10-4
8.3 x 10-2
3.6 x 10-4
Degree distribution histogram. Given graph G = (V, F),, we create a n-dimensional histogram by evaluating the in-degree deg(v) for v E V. The i-th position of the resulting histogram is the number of nodes with in-degrce i. We fi-normalize the histogram so that it sums to
Clustering coefficient. The local clustering coefficient for a directed graph of a node vs is formulated: C, = Ee E [0, 1] where N, {uy : (0), E E or (1y,0s) E Ej. It measures to what extent oi forms a clique. The different values Ci. Cn are binned into b-dimensional histogram. We set
= 100.
Laplacian spectrum. In the directed case, the eigenvalues of = L + I are complex but their absolute value is upper bounded by since S column stochastic. We bin their absolute values into a. 100-dimensional histogram.
We report in Table 1 the scores of DGDK, GRAN (Liao et al. 2019 (the state-of-the-at auto-regressive baseline that can be extended to digraphs), and SwinGNN (Yan et al. 2024) (a diffusion approach that can be applied to arbitrary adjacency matrices). It is worth noting that we train different GRAN and SwinGNN model for each category instead of using a class-conditional generation approach. Our MMD scores are close to 0, which means that the training and test graphs follow similar distributions. Moreover, DGDK outperforms both baselines in the clustering and spectrum evaluation metrics. This suggests that learning the global structure of the graph via its Laplacian in a single shot is beneficial for generation (for example compared to GRAN that sequentially considers multiple local subproblems).
6.3 Additional experiments
we consider an experiment similar to the one in Section 6.2 with digraphs containing different numbers of nodes (from 180 to 200). We compare different sampling strategies and show that strategies that exploit the learned matrix N tend to perform better. Tn Appendix G.4 we consider case where the training distribution contains multiple modes, and DGDK generates samples similar to the different modes.
We present additional experiments in Appendix G In Appendix G.2
7 Discussion
We discuss here different aspects of our approach that explained in detail in the appendix.
7.1
Scalability.
One main advantage of our approach is that use efficient closed-form solutions for the encoder that can be preprocessed offline to obtain efficient training times. As explained in Appendix B in practice GeForce RTX 3090 with 24 GB or VRAM, we manage to train 10,000 iterations of the datasct introduced in 6.2 in about hour when the batch size is 20 (10 training graphs per category). We can also train in hour 10,000 iterations of the dataset introduced in Section G.2] with larger graphs containing up to 200 nodes and a batch size of 2 (1 training graph per category). The main limitation scale to larger graphs is memory. Indeed, our low-rank approximation technique presented in Section 4.2 allows us to save memory that grows linearly in the size of the training set.. Nonetheless, larger graphs require larger VRAM.
NVIDIA
Section
11

---

Published in Transactions on Machine Learning Research (01/2025)
7.2 Connection with Heat Kernels and Reproducing Kernel Banach Spaces (RKBS)
The understanding of heat kernels and RKBS is not essential to run our method. Nonetheless, it provides connection with existing work that considers graphs and kernel methods (Kondor & Lafferty 2002 Belkin & Niyogi 2003). In particular, this connection is used in Belkin & Niyogi (2003) to create low-dimensional representations of the nodes, and in Kondor Lafferty (2002) to perform classification by using the kernel matrix as a similarity measure between nodes. In both cases, their graphs are undirected and their kernel matrix is symmetric, which allows them exploit Reproducing Kernel Hilbert Spaces (RKHS).
In our case, propose to add nonhomogeneous term to introduce noise, and show that the resulting non-symmetric matrix K still corresponds to non-symmetric kernel matrix induced by a RKBS. This motivates our low-rank approximation, inspired by kernel methods, to reduce used memory that grows linearly in the number of nodes instead of quadratically.
It is worth noting that since the kernel matrix K non-symmetric, the standard framework of RKHS that relies on (symmetric) inner products cannot be used. Our RKBS interpretation introduces a elegant alternative that is useful both theoretically and in practice.
8 Conclusion and Limitations
We have proposed a one-shot generative model that samples digraphs and is similar in essence to denoising autoencoders. Our encoder exploits closed-form expressions to add noise to digraph, and our decoder is trained to recover the global structure of the graph via its Laplacian dynamics. We show how our framework generalizes heat kernels and able to simulate the training distribution. We also propose a low-rank approximation of the heat kernel matrix to possibly scale to large graphs.
Limitations
Although our approach is scalable and produces digraphs of different sizes, one limitation is that our decoder is deterministic and its output is determined by its input. Nonetheless, different sampling strategies that we study thoroughly in Appendix G.2 can be used generate a large set of different types of input. This allows our model to generate diverse set of outputs.
Another limitation of our work is the lack of experiments on real-world datasets. This limitation is due to the fact that we did not find an appropriate real-world dataset to try our approach. For instance, even the genome dataset introduced in Marbach et al. (2012) contains three different graphs, each for different type of gene. Our goal is to work with a large number of graphs that follow the same type of distribution. On the other hand, in order to work with the dataset in Marbach et al. (2012), we would need to generate multiple subgraphs from the three large graphs, which is different from the task we are interested in.
Broader Impact Statement
Our work is not focused on applications. There are many potential societal consequences of our work, none which we feel must be specifically highlighted in the paper.
Acknowledgments
We thank Nicholas Sharp for helpful discussions about this project and the anonymous reviewers for their feedback.
References
Pavel Avdeyev, Chenlai Shi, Yuhao Tan, Kseniia Dudnyk, and Jian Zhou. Dirichlet diffusion score model for biological sequence generation. International Conference on Machine Learning, 2023.
Dominique Beaini, Saro Passaro, Vincent Lêtourneau, Will Hamilton, Gabriele Corso, and Pietro Liò. Directional graph networks. Tn International Conference Machine Learning, pp. 748 758. PMLR, 2021.
12

---

Published in Transactions on Machine Learning Research (01/2025)
Mikhail Belkin and Partha Niyogi. Laplacian eigenmaps for dimensionality reduction and data representation. Neural compulation, 15(6):1373-1333 2003.
Luca Bombelli, Joohan Lee, David Meyer, and Rafael D Sorkin. Space-time as a causal set. Physical review letters, 59)52,, 1987.
Fan Chung and S-T Yau. Coverings, heat kernels and spanning trees. the electronic journal of combinatorics, pp. R12-R12, 1999.
James R Clough and Tim S Evans. Embedding graphs in lorentzian spacetime. PloS one, 12(11):o0187311, 2017.
Morris H DeGroot. Reaching a consensus. Journal of the American Stalistical association, 69(445):118-215 1974.
Kaize Ding, Zhe Xu, Hanghang Tong, and Huan Liu. Data augmentation for deep graph learning: A survey. SIGKDD Explorations Newsletter, 24(2):61-77, dee 2022. ISSN 1931-0145. doi: 10.1145/3555637.3375660 3355644.. URL https://doi.orgg:
011114/5555677357755555 3575646)
C Henry Edwards, David Penney, and David T Calvis. Differential Equations and Linear Algebra, eBook. Pearson Higher Ed, 2020.
Paul Erdös, Alfréd Rényi, et al. On the evolution of random graphs Publ. Math. Inst. Hung. Acad. Sci, 5(1): 17 60, 1960.
Jonathan Ho, Ajay Jain, and Picter Abbeol. Denoising diffusion probabilistic models. In Advances in Neural Information Processing Systems, 2020.
Paul W Holland, Kathryn Blackmond Laskey, and Samuel Leinhardt. Stochastic blockmodels: First steps. Social networks, 5(2):109-137, 1983.
Aapo Hyvärinen and Peter Dayan. Estimation of non-normalized statistical models by score matching. Journal of Machine Learning Research, 6(4), 2005.
Jaehyeong Jo, Seul Lee, and Sung Ju Hwang. Score-based generative modeling of graphs via the system of
stochastic differential equations. In International Conference on Machine Learning, 2022.
Diederik P. Kingma and Jimmy Ba. Adam: A method for stochastic optimization. arXiv preprint arXiv:14 14N. 2014.
Dicderik Kingma and Max Welling. Auto-cncoding variational bayes. International Conference on Learning Representations, 2014.
Risi Imre Kondor and John D. Lafferty. Diffusion kernels on graphs and other discrete input spaces. In Proceedings of the Nineteenth International Conference Machine Learning, ICML '02, pp. 315 22-22,, San Francisco, CA, USA, 2002. Morgan Kaufmann Publishers Inc. ISBN 1558608737.
Samuel Kotz, Narayanaswamy Balakrishnan, and Norman Johnson. Continuous multivariate distributions, Volume 1: Models and applications, volume 1. John Wiley & Sons, 2004.
Youngchun Kwon, Dongseon Lee, Youn-Suk Choi, Kyoham Shin, and Seokho Kang. Compressed graph representation for scalable molecular graph generation. Journal of Cheminformatics, 12(1):1-8, 8, 2020.
Mare T. Law and James Lucas. Spacetime representation learning. In The Eleventh International Conference on Learning Representations, 2023. URL
BaapadiggV___hhh__aa M_rhYajc
Juho Lee, Yoonho Lec, Jungtack Kim, Adam Kosiorek, Scungjin Choi, and Yee Whye Tch. Set transformer: A framework for attention-based permutation-invarinn neural networks. In International conference on PMT.R, 2019.
machine learning, pp.
13

---

Published in Transactions on Machine Learning Research (01/2025)
Renjie Liao, Yujia Li, Yang Song, Shenlong Wang, Charlic Nash, William L. Hamilton, David Duvenaud, Raquel Urtasun, and Richard Zemel. Efficient graph generation with graph recurrent attention networks. In NewrIPS, 2019.
Daniel Marbach, James C Costello, Robert Küffner, Nicole M Vega, Robert J Prill, Diogo M Camacho, Kyle R Allison, Manolis Kellis, James J Collins, and Gustavo Stolovitzky. Wisdom of crowds for robust gene network inference. Nature methods, 9(8):796 804, 2012.
Karolis Martinkus, Andreas Loukas, Nathanaël Perraudin, and Roger Wattenhofer. Spectre: Spectral conditioning helps to overcome the expressivity limits of one-shot graph generators. In International Conference Machine Learning, pp. 15179. PMLR, 2022.
Rocío Mercado, Tobias Rastemo, Edvard Lindelöf, Günter Klambauer, Ola Engkvist, Hongming Chen, and Esben Jannik Bjerrum. Graph networks for molecular design. Machine Learning: Science and Technology, 2(2):025023, 2021.
James Mercer. Xvi. functions of positive and negative type, and their connection the theory of integral equations. Philosophical transactions of the royal society of London. Series A, containing papers of a mathematical or physical character, 200(441-555):4111-44 1909.
Chenhao NIu, Yang Song, Jiaming Song, Shengjia Zhao, Aditya Grover, and Stefano Ermon. Permutation invariant graph generation via score-based generative modeling, In International Conference on Artificial Intelligence and Statistics, pp. 4474-4484. PMLR, 2020.
James R. Norris. Markov chains. Number 2. Cambridge university press, 1998.
Leslie O'Bray, Max Horn, Bastian Rieck, and Karsten Borgwardt. Evaluation metrics for graph generative models: Problems, pitfalls, and practical solutions. International Conference on Learning Representations, 2022. URT. hattps://openrevieww wtt.:tt/oooorn...ddttttt/.......
Larry Page. The pagerank citation ranking: Bringing order to the web. Technicol Report, Slanford Universily, 1998.
Caroline E. Seely. Non-symmetric kernels of positive type. Annals of Mathemalics, pp. 172-176, 1919.
Nicholas Sharp, Souhaib Attaíki, Keenan Crane, and Maks Ovsjanikov. Diffusionnet: Discretization agnostic learning on surfaces. ACM Trans. Graph., 41(3), mar 2022.
Aaron Sim, Maciej Wiatrak, Angus Brayne, Páidí Creed, and Saee Paliwal. Directed graph embeddings in pscudo-ricmannian manifolds. International Conference on Machine Learning (ICML), 2021.
Martin Simonovsky and Nikos Komodakis. Graphvae: Towards generation of small graphs using variational autocncoders. In Vera Kirková, Yannis Manolopoulos, Barbara Hammer, Lazaros Iliadis, and Ilis Maglogiannis (eds.), Artificial Newral Networks and Machine Learning ICANN 2018, pp. 412-422, Cham, 2018. Springer International Publishing.
Alex J Smola and Bernhard Schölkopf. Learning with kernels, volume 4. Citeseer, 1998.
Jascha Sohl-Dickstein, Eric Weiss, Niru Maheswaranathan, and Surya Ganguli. Deep unsupervised learning using nonequilibrium thermodynamics. In International Conference on Machine Learning, 2015.
Guohui Song, Haizhang Zhang, and Fred J Hickernell. Reproducing kernel banach spaces with the 11 norm. arXiv preprint arXix:1101.4988 2011.
Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations. In International Conference Learning Representations, 2021.
JJP Veerman and Robert Lyons. A primer on laplacian dynamics in directed graphs. arXiv preprint 02X6:22002.02655, 2020.
14

---

Published in Transactions on Machine Learning Research (01/2025)
Petar Veličković, William Fedus, William L. Hamilton, Pictro Liò, Yoshua Bengio, and R. Devon Hjelm. Deep graph infomax. In International Conference on Learning Represenlations, 2019. URL https //openreview.netee http://orrumm/ddd-111ile..
Clement Vignac and Pascal Frossard. Top-n: Equivariant set and graph generation without exchangeability.
In International Conference Learning Representations, 2022. URL https://oppareezee.nntt//ooo
id=-Gk IPJWvk
Clement Vignac, Igor Krawczuk, Antoine Siraudin, Bohan Wang, Volkan Cevher, and Pascal Frossard.
Digress: Discrete denoising diffusion for graph generation. In The Eleventh International Conference Learning Representations, 2023. URL http:ssss//ppeeeeeeeeeeeee/eoees
Pascal Vincent. A connection between score matching and denoising autoencoders. Neural computation, 23 (7):1661-1674, 2011.
Pascal Vincent, Hugo Larochelle, Yoshua Bengio, and Pierre-Amis Antoine Manzagol. Extracting and composing robust features with denoising autoencoders. In Proceedings of the 25th international conference on Machine learning, pp. 1096-1103, 2008,
Christopher KI Williams and Carl Edward Rasmussen. Gaussion processes for machine learning, volume 2. MIT press Cambridge, MA, 2006.
Qi Yan, Zhengyang Liang, Yang Song, Renjie Liao, and Lele Wang. SwinGNN: Rethinking permutation
invariance in diffusion models for graph generation. Transactions Machine Learning Research, 2024.
ISSN 2835-8856. URL https://openreviex. http:::/ooeeeervee..eeeeerrrrrrr net::/orrmmm1dddaaiiiippp...
Jiaxuan You, Rex Ying, Xiang Ren, William Hamilton, and Jure Leskovec. Graphrmm: Generating realistic graphs with deep auto-regressive models. In International conference on machine learning, pp. 5708-5717. PMLR, 2018.
Haizhang Zhang, Yuesheng Xu, and Jun Zhang. Reproducing kernel banach spaces for machine learning. Journal of Machine Learning Research, 10(12), 2009.
Yanqiao Zhu, Yuanqi Du, Yinkai Wang, Yichen Xu, Jieyu Zhang, Qiang Liu, and Shu WT. A survey on deep graph generation: Methods and applications. arXiv preprint arxi::2003 06714, 2022.
15

---

Published in Transactions on Machine Learning Research (01/2025)
Algorithm Training algorithm (for each mini-batch)
Input: Node representations
R"muux", Xd E
hyperparameters T
0,
0,
0,
Bernoulli factor
E
[0, !]
Initialize mini-batch as mini-batch
Graph C; with Adjacency matrix A' and Laplacian matrix A' in the mini- batch do Promote permutation invariance is true then
Generate random permutation matrix E (0, 1)9cxns.
A t P A'P.
A' -
end if
If Data augmentation matrix is true then Generate perturbation matrix C E (0, 1japon s.t. Vi, Ca 0 and Vi # J. Cit Bernoulli(p)
A A C.
Á A [diag(1 A9) (diag(1AA)))
else
Á A
end
Define N E Rexd as upper submatrix followed by softmax operation for éi-normalized columns
Define B Tá optionally as the rank-s approximation of eTA via truncated SVD.
(T) < B + e=TT)
Define E as e! or optionally as the rank-s approximation of e'' via truncated SVD. T
EEN
L mint- batch + mini-batch Ledgo (t) (node(f) where 0 is regularization parameter.
end for Optimize Emini-batch with Adam Kingma Ba 2014
A
Summary
Our appendix structured as follows:
Section
B
provides experimental details including the architecture of the neural networks.
Section
provides the details of our equations for our heat diffusion encoder.
Section
D
provides the necessary details to solve the consensus problem.
Section
explains the connection of our method with heat kernels.
Section
studies the impact of and a on the eigenvalues of the different matrices that are involved
our model.
Section
presents additional experimental results.
Experimental details
Setup. We ran all our experiments on a single desktop with a NVIDIA GeForce RTX 3090 GPU (with 24 GB of VRAM) and 64 GB of RAM. We coded our project in Pytorch. We use double precision format to define our tensors. This important in the current version of Pytorch to obtain an accurate version of the matrix exponential.
Our training algorithm is illustrated in Algorithu 2 and follows the different data augmentation techniques mentioned in Section 4
Node representation matrix N. We define Raax E N as the maximum number of nodes in graph of the training set, and d as the number of columns of the initial node representation O E Rm xd. In practice, we set 150 and we initialize each element of O by sampling from the normal distribution parameterized by a mean of and standard deviation of Other hyperparameter values could be used. For a given graph G, of n; nodes, define N Ra as the upper submatrix of O, and it is n--orrallieed by using a column-wise softmax operator. Equation 7 is minimized via standard gradient descent. by training jointly v. yp and N
16

---

Published in Transactions on Machine Learning Research (01/2025)
(and hence O). In practice, we find that the larger d, the better the performance. However, approach is limited by the amount of available memory.
In Section 6.2 we train the model for 60,000 iterations. We recall that each mini-batch contains 10 training graphs per category, hence we have 20 graphs per mini-batch. We make this choice because this is the maximum amount we manage to load in our GPU VRAM. The training algorithm takes about one hour for 10,000 iterations, so about 6 hours in total. We use regularization parameter of y = 100, and a step size/learning rate of 0.0001. Following O'Bray et al. (0022), calculate MMD evaluation metric by using radial basis function (RBF) kernel with 10. We report performance for other values of a in Appendix G.2
Backbone. The common backbone of the node and edge decoder is:
a linear layer Rd + Rd where d' = 1000.
set (ransformer (Lee al. 2019) with four attention blocks, each of dimensionality 1000, with one head, 20 inducing points and 20 seed vectors. Each row-wise feedforward layer of attention block contains linear layers d' x d' with ReLU activation function. We choose one head so that the global structure of the Laplacian matrix is not separately processed by the different heads.
a multilayer perceptron (MLP) with one linear layer Rdd Rdd (with d" = 600), followed by four linear layers Rd Rd", and one linear layer Rdd Rdd" where d" 200.
The head of the node decoder is an MLP with an initial layer Rd Ra followed by one linear layer Rat - Ra""" and one linear layer Ra" - Rd. If the input is of size Reixd this returns matrix of same size. We also use column-wise normalization with softmax the output.
Edge decoder. The output of the backbone described above are node representations of size d" = 200 for each node. We concatenate them (as described in Section 4.1) to obtain representations of pairs of nodes of size 2d"" = 400. They are given as input of an MLP with one linear layer R2d" - Rdd" (with d'n 600), followed by 5 linear layers Rdt" - "d followed by a linear layer Te" + RR. that is used for cross-entropy loss (one element used for the absence of edge, and the other element is used for the existence of edge).
Alternatively, one could output real value for each pair of nodes and use binary cross entropy loss.
We use ReLU
activation function between all the linear layers.
In all our experiments, we give the same weighting for the positive and negative edges (i.e., presence or absence of edge).
C Details of Equations for the Heat Diffusion Encoder
C.1 Nonhomogeneous heat source term
We now give the details of the equations in Section
We recall the formulation of Q in Proposition
Q(s) = ae dse (R - xA)())
(9)
which implies the following formulation of F:
(s)ds 's I'd c Q()dss
F(I)
(10)
((0))ds e/ (' aa- - x(0))ds
ac(t- - as R
(11)
= A -ae-as oAx(0))
= (1 -at a R eAx
(12)
17

---

Published in Transactions on Machine Learning Research (01/2025)
C.2 Node representation over time
We assume in this subsection that L 0. equation 4 is written:
x(1) e A (x) + (e 1)e 8 (() + (1-e* (l-e-aa)R)
(13)
For any nonnegative time step 2 0, we can write X(t + r) as function of ((t) and vice versa.
X(t (t)+ ((+r)) A (X(0) + (e (c-a(t+r) 1) 1)e (()) (0) + (1 e - R
=e = d (((0) + (e - 1)c (()) + (1 - c + (E - - -e - -R))
e (((t) + (e (e - ceatte-altti - d'a (R x
From the equation above, we find:
X(t) x(t)=e-<((tt)) - -eraa - -ead---attt - C (R - e
(14)
When 3 0, we have:
X(T) = -all e T (() + - -all' aa''
Z(T) + (1 E -
(15)
Let us assume that , = 0, cquation can be written as:
(() = i (e (()) + +1 (1- (a))))))
(16)
which implics
= e0 x()( + e1- e.
(17)
Equation 17 implies the following formulation of x() as a function of X(1 + T):
((0) ( + T) +1 (1 c R
(18)
From equation 14 and by setting 3 = 0, we obtain:
X(t) =
+ T) e -e-l-ealltty e (R-
(19)
X + +
(e -e eA +
(/ e'Ax(0)
(20)
By using equation 18 the last term of equation
20 can be rewritten:
(e-at- e- a((+)))) e (0)) (e (e-at - c aa(t e' (olltt))g ( +11 - ca(l (a(l-T)) ()
(21)
- ((x((())(((00t)-)))
eA
(22)
Equation 20 is then rewritten:
X(t) = e E x(++T)+ + ((-
- cat + (e-at_ -
-aaa-e-aa - (1 (((tt)))))) cA
(23)
= _ - (c<a((+)) - att 46 at - - calltt) --((tt)+1) + e'
(24)
X(I) =M (X(++)) (1- e R
(25)
C.3
Stochasticity of the node representation matrix
In Section we mention that if the matrices e'A and N are both column stochastic for all / 2 0 VVerman & Lyons. 2020 then Z(1) = A is also column stochastic for all /. 12 0.
This casily verified. Let us assume that two matrices B and C are column stochastic. They then have nonnegative elements, and they satisfy 1 B = 1 and _ = 1T. The matrix (BC) column stochastic because it has nonnegative elements and satisfies 1 (BC) = 1 1 =
18

---

Published in Transactions on Machine Learning Research (01/2025)
D Consensus
We now give the formulac and constraints for the consensus model (i.e., when A = =[ DDGroott 1974-
In this model, the matrix é'á is row stochastic for all & 2 0. We then constrain both N = (0) and the matrix
M := 1 E [ jadd to be row stochastic, this implies Z(t) row stochastic for all 2 0. In this case, we have vt 2 0, A M. We then also define R := M. In the consensus model, X(t) Z(t) (1 - -- is
row stochastic for all t 2 0.
A detailed comparison between the diffusion and consensus models is discussed in Veerman & Lyons 2020 Section 6). The consensus model would be appropriate in contexts where each row of (() is one-hot vector corresponding to the category of the node. However, it might suffer from degenerate cases where some rows of (() do not depend on the same rows at their initial time I due to the properties of the left eigenvectors of L.
E Non-symmetric Heat Kernels
We explain how the heat kernel framework in Kondor & Lafferty| 2002) is special case of our encoder (see Section 3) when the kernel function is symmetric and the source term Q is homogencous. This connection motivates the study of the columm space of the Laplacian matrix that is the foundation of many heat kernel methods (Blkiin & Niyogi 2003 Kondor & Lafferty| 2002).
Heat kernels (Kondor & Lafferty, 2002). We first explain how our approach can be seen as a non- symmetric heat kernel when the term Q is homogeneous. We recall that function K X' X' R defined on some nonempty set X' is called kernel function if it satisfics Mercer's theorem (Mereer 1909 In other words, it satisfies fx fx K (p.g)) (F)) (4)dpdg (())d(d 2 0 for every function (w) of integrable square, or EPEX Eyex K(p, 9).ff 2 0 for all sets of real coefficients {fp} in the discrete case, which is the case we are interested in.
Kernel functions are usually defined to be symmetric (i.e., Vp, q, K(p.g)=K(u.)) to define a Reproducing Kernel Hilbert Space (Smola & Schölkopf 1998), and the symmetry of the kernel function between pairs of nodes is reasonable when the graph is undirected. However, kernels are not necessarily symmetric (Seely 1919
and one can define Reproducing Kernel Banach Space (RKBS) on X' 2 V equipped with the & norm (Song let al., 2011 Zhang et al. (2009) by considering the kernel matrix K E 10,1 defined such that K = e with K pg = (p, Let us define [0, 1] ;= : + [0, L]} the Banach space of functions mapping X' into [o. 1]. We also define the lincar operators
(,,,, g(.) := -- 3,K(, x/) and the bilincar form (f.g) == E. Σx= Og xz), where n E N, m E N,a; E R, E R. Fn. 2i...., 2n. E X are arbitrary. The resulting RKBS generalizes the heat kernels in Kondor & Lafferty (2022) to directed graphs although it restricts the functions to map into [0, 1] instead of R.
When Q is nonhomogeneous, we have X(T) = e aT Z(T) ((-T))((-e-) (1-e aT = M since 1 '11" TNNN- = M when N is column stochastic. This leads the kernel matrix: n
aT TT + (1 - (/111/11 N
K k aT TA a- aT
(26)
As in Laplacian eigenmap approaches for dimensionality reduction (Ching & Yau 1999). K can be thought of as an operator on functions defined on nodes of the graph and we obtain the node representation matrix KN [o, 1jxxd The fact that KN column stochastic allows us to upper bound the & norm of its columns
by 1 (it is equal to 1), and then satisfy the properties of RKBS on X' with the e norm (Song et al. 2011).
We now explain how our framework falls into the framework of Song et al. 2011 Zhang et al. 2009
Let us consider the case V = X, which implies |X| = n. We consider the non-symmetric function K as Kpp = K(p,g) where the kernel matrix K= --=-ad e-aT TA + (1-e-oT)/n11 E (0, 1|xxx is column stochastic.
Since K(p.q) is nonnegative for all p E X and q E X, a sufficient condition to satisfy EpEx Eyex ((p,))( q)fpff 2 0 is to constrain (fplpex to be a set of nonnegative coefficients. In our ex- periments, set {pp to be set of nonnegative coefficients that sum to 1 (i.e., column stochastic).
19

---

Published in Transactions on Machine Learning Research (01/2025)
The above explanation did not require the notion of RKBS. We can nonetheless use Proposition 5 of Zhang et al. (009)) which is stated as follows: If the input space X is a finile sel, then any noutrivial function K X x X' is the reproducing kernel of some RKBS on X.
In our case, K nontrivial because K is full rank. Our kernel matrix (that is full rank, and column stochastic hence with nonnegative elements) naturally satisfies the first three requirements of Song et al. (2011). also satisfies the relaxation of their fourth requirement in their Section 6. We recall their requirements:
(A1) for all sequences {xp E {1,.....m) n}} S X' of pairwise distinct sampling points, the matrix K : (p. E R is non singular.
(A1) is satisfied because K is full rank in our case.
(A2) is bounded, namely, [sss t) & M for some positive constant M and all s, {l,...,
(A2) is satisfied because K column stochastic, so Vs,t, (K(, [)] This is satisfied if M.
(A3) for all pairwise distinct p X, N and having its --rrr finite, Ej=ee ez for all
E X implies = 0.
(A3) is satisfied because Vj. (xz,, > 0 and Vj # (i,K(r), 2 0 in our case. Therefore, (A3) can be satisfied only if c =0 0.
the relaxation of (A4) in their Section 6 can be formulated as follows: let us write as follows:
K 11,,-1), K kii(n11),n Kn,1((n-1)
K
(27)
K
where K1:(n (1:((=1).1:(a-11 [0, (0=1)_)))) and 11;(n-1).n [0, are submatrices of K. The relaxation of (A4) is satisfied if there exists some 8n such that: [[( "oipda((n))""K(nt-----
S 8n is satisfied. Since Ki,nn-).nn 1).n has its values in [0, 1], we have to be able to bound the values of (K (i)
We can
bound them by exploiting the (inverse of the) eigenvalues of that depend on T and a (see equation 29
since we know that Vr, [A, + +1| 1.
It is worth noting that this section has proven that it is possible to formulate RKBS represent our node similarities. However, the kernel matrix K is given input of our algorithm via the adjacency matrix A. One could define some node representation space that would induce the matrix K by using the theory of RKBS instead of considering that K is given as input of the algorithm.
F Difference between T and a
To understand the difference of impact between T and we need to study the eigenvalues of the kernel matrix described in Appendix K aa- + (1- - ==="" e 10,1
It worth noting that S, = and TA all have the same set of right and left eigenvectors. The only difference is their set of eigenvalues. Since is column stochastic, is diagonalizable, its spectral radius and it has at least one cigenvalue oqual with as left cigenvoctor. The number of eigenvalues of S that are equal to 1 is the number of reaches of the graph (Veerman & Lyons 2020p. A reach of a directed graph is maximal unilaterally connected set (see Veerman & Lyons (2020) for details). By definition of A S II, A and S have the same eigenvectors and those that correspond to the cigenvalue of S, correspond to the eigenvalue 0 of A, and to the eigenvalue e = 1 of e' for all / E IR The matrix e is column stochastic for all /. 2 0, because its spectral radius then 1, it has at least one eigenvalue equal to with 1 as left. eigenvector, and it has the same eigenvectors as the column stochastic matrix S.
Let us note A U the eigendecomposttio of A. The eigendecomposition of S S U(A IJU-
and the eigendecomposstion of e is elA U Let consider that the first of U-1 is 71
20

---

Published in Transactions on Machine Learning Research (01/2025)
where Y # 0 appropriate factor (i.e., the first row of U-1 is collinear to 1 17), and let note:
A1 0
0 2
A
(28)
An-1 0 An
We know that X since S is column stochastic. Morcover, both 1 1 and TA have 1 left eigenvector with corresponding eigenvalue equal to 1, so K = e- (1 with corresponding eigenvalue equal to 1. The eigendecomposition of K is then k= V$V-1 where V + U in general, but the first row of V-1 is collinear to 1 The diagonal matrix & is written:
a -e-aT)
also has 1 as left eigenvector
0 e T(22 a2
(29)
e T(n-11a)
e (TAn-al)
It is worth noting that all the nonzero eigenvalues A, of A have negative real part by definition of S (i.e., since the spectral radius of 1). If a = 0, then Vr, Ar 0 = ((Ar-aa) = 1. If a > 0 and T > 0, then the real part of Ap - a is negative for all r, which implies le (Ar-a)) < | If and T > 0, we also have for all r, limy- leT [ 0 and lim
T(r (e
=0
From equation 29 the main difference between T and the eigenvalues inside the exponential, whereas 0 only has impact the real part of the eigenvalue inside the exponential.
that T acts as multiplicative factor on
G Additional Experimental Results
G.1 Study of the column space of the learned representations
Fig illustrates additional qualitative results from the experiments in Scction 6.1 showing cosine values between the singular vectors (ordered by magnitude of their singular values) of e e and e'A'N. To improve visualization, we set all the absolute values lower than 0.3 to 0. One can see that most cosine correlations appear along the diagonal. This shows that the singular vectors of the different matrices are correlated. Although some absolute values are high in the bottom right corner, their corresponding singular values are much smaller compared to those in the top left corner, the overall importance of their correlation is weaker.
It is worth noting that we conditional sampling in this experiment so our model has to jointly learn representations that are relevant for both categories.
G.2 Approximation techniques for larger graphs
In this subsection, consider task similar to the described in Section 6.2 The goal in Section 6.2
to show that digraphs could be entirely reconstructed with our approach by using low-rank approximations of the Laplacian for small graphs. In other words, if we give some training adjacency matrix A E {0, 1)xn to our noising encoder, then our decoder is able to reconstruct all the elements of A. For simplicity, we did not consider data augmentation techniques in Section 6.2
We now describe our experimental setup to generate larger graphs of different sizes. As in Section 6.2 consider the class-conditional digraph generation with two categories. We use the following categories
21

---

Published in Transactions on Machine Learning Research (01/2025)
(left) Erdös-Rényi graphs
(right) stochastic block model graphs
Figure 4: Correlations between the left singular vectors of e'A' and e N.
Erdös-Rényi with == 0.4, and a stochastic block model with 5 blocks and the following transition matrix between the blocks:
0.9
0.2
0.4
0.2
0.4
0.3
0.9
0.15
0.5
0.45
II
0.4
0.95
0.05
0.4
(30)
0
0.3
0.4
0.75
0.45
0.1
0.4
0.4
0.15
0.7
All the graphs contain n; nodes where n, E {180, 181,
The first block contains 40 nodes, the second block contains 20 nodes, the third and fourth blocks contain 35 nodes each, and the last block contains from 50 to 70 nodes. We set μ = 0.4, p = 1/n, and 50 for the rank-s approximation of TT' e [0,1]xxm, in this experiment.
200).
At inference time, sample digraphs that contain ny nodes where n4 E (180, 181,.., 2003. Quantitative results reported in Table and Table 3 Once again, DGDK outperforms GRAN and SwinGNN in evaluation metrics that take into account global properties of the graph. Data augmentation slightly improves performance.
In both Table and Table 3 we report the MMD scores for different values of the variance parameter in the RBF kernel: -² 100, 10 and 1. We recall that we use の² 100 in Table We call Discrete DGDK the sampling strategy described in Algorithm where we sample discrete adjacency matrices E (0, 10,1) described in line 2 of Algorithm 1 Continuous DGDK corresponds to the same sampling strategy except.
22

---

Published in Transactions on Machine Learning Research (01/2025)
Table 2: Squared MMD distances over 5 random initializations (average standard deviation) for the
Rénói distribution :000.4.
MMD metric
Degree
Clustering
Spectrum
Sampling from flat. Dirichlet distribution
0.00091 + 0.00043
0.0079 + 0.0026
0.00193
0.0004
Continuous DGDK (with data augmentation)
0.00074 1 0.00004
0.0068 1 0.0005
0.00085
0.0001
Discrete DGDK (with data augmentation)
0.00073 1 0.00002
0.0067 1 0.0006
0.00085
0.0001
II
Discrete DGDK (without data augmentation)
0.00082 + 0.00003
0.0069 = 0.0004
0.00092 + 0.0003
6
SwinGNN
0.00091 + 0.00084
0.0091 = 0.0012
0.00110 + 0.0008
GRAN
0.00118 + 0.001 24
0.0138 = 0.0035
0.00143 + 0.0013
Sampling from flat Dirichlet distribution Contimuous DGDK (with data augmentation)
0.0093 0.0065
0.077 = 0.0298
0.0202
0.0074 1 0.0005
0.067 1 0.0049
0.0085 0.001
Discrete DGDK (with data augmentation)
0.0072 1. 0.0002
0.065 1 0.0009
0.0085 0.001
II
Discrete DGDK (without data augmentation)
0.0083 + 0.0004
0.071 F 0,0014
0.0123 1 0.002
SwinGNN
0.0090 + 0.0035
0.079 # 0.0021
0.0113 + 0.006
GRAN
0.0121 + 0.0103
0.143 H 0.0241
0.0149 + 0.012
Sampling from flat Dirichlet distribution Contimmous DGDK (with data augmentation)
0,103 0.098
0.57 0.14
0.184 = 0.0340
0.076 + 0.009
0.57 0.08
0.084 = 0.0122
Discrete DGDK (with data augmentation)
0.071 1 0.002
0.54 1 0.02
0.083 1 0.0092
II
Discrete DGDK (without data augmentation)
0.084 1 0.004
0.60 0.02
0.103 = 0.0104
SwinGNN GRAN
0.102 + 0.064
0.64 = 0.10
0.101 = 0.0126
0.142 + 0.105
1.39 = 0.22
0.145 = 0.1156
Table 3: Squared MMD distances over 5 random initializations (average + standard deviation) for the stochastic block model (5 blocks).
MMD metric
Degree
Clustering
Spectrum
Sampling from flat Dirichlet distribution
0.00031
0.0002
0.0069 + 0.0026
0,00232
0.0004
100
Continuous DGDK (with data augmentation)
0.00017
0.0001
0.0057
0.0002
0.00041
0.0003
Discrete DGDK (with data augmentatiinn)
0.00015
0.0001
0.0039
0.0023
0.00038
0.0003
II
Discrete DGDK (without data augmentation)
0.00031
0.0001
0.0046 1 0.0025
0.00041
0.0004
SwinGNN
0.00046
0.0009
0.0245 = 0.0094
0.00831
0.0103
GRAN
0.00053 0.0008
0.0654 = 0.0057
0.02472 + 0.0144
Sampling from flat Dirichlet distribution Continuous DGDK (with data augmentation)
0.0021 0.001
0.067 = 0,020
0.0187 0.002
0.0018 0.001
0.059 0.006
0.0040 0.003
Discrete DGDK (with data augmentationn]
0.0015 0.001
0.056 1 0.001
0.0039
II
Discrete DGDK (without data augmentation)
0.059 T 0.004
0.0043
SwinGNN
0.0041
0.264 1 0.062
0.0793
GRAN
0.0056 0.004
0.664 + 0.091
0.2629 0.124
Sampling from flat Dirichlet distribution Continuous DGDK (with data augmentation)
0.0255 0.0133
0.59112 + 0.342
0.172 + 0.031
0.0186 0.0132
0.58531 + 0.162
0.041 = 0.036
Discrete DGDK (with data augmentation)
0.0134 0.0093
0.48165 = 0.012
0.038 = 0.032
II
Discrete DGDK (without data augmentation)
0.0331 0.0103
0.56420 0.025
0.048 1 0.042
SwinGNN GRAN
0.0424 0.0251
0.68073
0.351
0.388 1 0.724
0.0533 + 0.0421
1.56743 1 0.931
1.253 + 1.123
that sample non-diagonal elements of A uniformly in the continuous interval [0, 1] instead of {0, 1}, while keeping the constraint Vi, Au = The two methods are competitive with baselines.
On the other hand, we also report scores when we sample each column of the input directly from a flat Dirichlet distribution. This strategy does not exploit the learned matrix N and is outperformed by our other sampling strategies although it still outperforms the GRAN deadline.
G.3 Data augmentation by using permutations of Laplacian matrices
In this subsoction, we study the impact of the data augmentation technique adding adjacency matrices of isomorphic digraphs as explained in Section 1.2 Figure illustrates the loss value of equation 7 as a function of the number of iterations with and without this data augmentation technique. In this setup, use 7 100.
23

---

Published in Transactions on Machine Learning Research (01/2025)
without data augmentation adding permulations of adjacency matrices with data augmentation ng permulations edjacency matrices
1069
1000
2000
4000
5003
2000
8003
5950
10002
number lterations
Figure 5: Loss values obtained when optimizing equation 7 with or without data augmentation by adding permutations of training adjacency matrices.
We use the experimental setup described in Section 6.2 with 3,000 training graphs per category, 10 graphs per category per batch. Each epoch corresponds to 300 iterations. Both loss curves follow similar patterns, although the one corresponding to data augmentation has slightly higher loss value.
One reason of the low impact of this data augmentation technique is that the matrix is jointly learned with nonlincar decoders o and & that are robust to this kind of transformation.
G.4 Additional experiment: multimodal categories
In Section 6.2 the categories are unimodal and it is assumed that they can be conditioned over at test time. We now show that if the distribution of the training set is multimodal and the modes are given to the model, model is able to sample graphs from the different modes.
We consider the following training set containing four modes of same size:
(1) a single block stochastic model with probability 0.28,
(2) two (disjoint and connected) components of same size with p = 0.48 cach,
(3) three components same size with p = 0.78 each,
(4) and four components of same size with p 0.97 each.
This corresponds to an average edge ratio of 0.24.
Fig. 7 illustrates some graphs generated with Am | when = and a = 1, they follow the training distribution. We emphasize that blocks are disjoints only in this subsection for visualization purpose, not in the previous subsections where some edges connect different blocks.
Impact of 7. The regularization parameter 7 in equation 7 acts as tradeoff between learning the edge decoder and the node decoder. If then the node decoder is not trained. Our goal is to predict the adjacency matrix from random matrix in [0, _]nxd given input of the decoders. Only the edge decoder is useful for this task. We found that the loss function does not converge during training when y 0. This suggests that the data augmentation technique that adds edge perturbation to the adjacency matrix is important, at least for this experiment. positive value of y is necessary for the loss to decrease during training and learn mcaningful decoders. This implics that the node decoder learns the global structure of the graph by reconstructing node representations that depend on the ground truth Laplacian. This is beneficial to the edge decoder as both decoders share common backbone.
24

---

Published in Transactions on Machine Learning Research (01/2025)
Figure 6: Digraphs generated when a 0.
Impact of a. Using 0 results in generated samples with only block (see Fig. 6)- As increases, the number of components per sample increases (sce Fig. 7 This is because the sampled adjacency matrices in Algorithm generated with Bernoulli distribution correspond to one component. As increases, the information of the sampled adjacency matrix gets partially lost and the decoder is able to reproduce graphs similar to the training distribution, and it samples graphs uniformly from all the modes. We provide some ablation study on the impact the noise diffusivity rate hyperparameter a 0.
When a 0, more than 97% of the generated graphs contain only one block. Some generated graphs with 0 are illustrated in Figure 0
When = 1 (i.e., 1 x 0.63), 19% of the generated graphs contain a single block, 26% contain 2 blocks, 30% contain 3 blocks, 20% contain 4 blocks, 2% contain 5 blocks, 2% contain 6 blocks, 1% contain 7 blocks. The distribution is similar to the training set that is uniformly distributed (i.e., 25% for each mode). Some generated graphs when a are illustrated Figure
When a 2.3 (i.c., 1 a-T 0.9), 17% of the generated graphs contain single block, 14% contain blocks, 20% contain blocks, 20% contain 4 blocks, 8% contain blocks, 5% contain blocks, 8% contain blocks, and the remaining 8% contains up to 15 blocks. Some graphs generated when 2.3 are illustrated in Figure Many generated graphs contain nodes that are isolated.
We report in Table the MMD scores that compare the training set with the generated graphs for different values of a. The performance shows that a has to be chosen carefully so that noise is introduced, but not in excess. It worth noting that GRAN is not appropriate in this experiment since has to be given the maxinum number of nodes in the graph. Our model is given number of nodes and has to generate one multiple (disconnected) graphs that follow the training distribution.
25

---

Published in Transactions on Machine Learning Research (01/2025)
Table 4: Squared MMD distances for the experiments on the multimodal dataset in Section G.4 for different values of
Dataset
Multimodal dataset
MMD metric
Degree
Clustering
Spectrum
2.6 x 10 10-3 3
1.7 x 10
3.3 10 4
0.9 x 10-3
6.3 x 10
2.3 x 10
2.3
1.6 x 10-3
1.0 x 10-3
5.9 10-4
Figure 7: Digraphs generated when = 1.0.
26

---

Published in Transactions on Machine Learning Research (01/2025)
Figure 8:
Digraphs generated when 2.3.
27


Reviews:
[{"review_id": "8IDBXXx2LE", "reviewer": "Reviewer_VUaS", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Pros:\n\n- The idea of using heat diffusion to encode the graph is interesting in the context of graph generation. There are a few benefits with such a design: 1) one can explicitly control the diffusion speed via the designed diffusion equation; 2) the closed-form solution to the heat diffusion allows fast computation.\n\n- The connection with heat kernel methods is nice, as the proposed encoder incorporates them as special cases. This part is currently presented in the Appendix. It would be great to move some parts in the main paper.\n\n- The experiments on small-scale graph datasets show promising results compared to the auto-regressive models.\n\nCons and Suggestions:\n\n- Although the proposed method is closely related to denoising auto-encoder, it lacks a principled probabilistic model. In particular, during generation, one has to first sample an adjacency matrix from the Erdős–Rényi model and then inject noise into it by randomly corrupting the adjacency matrix. This part is somewhat ad-hoc since the training distribution of adjacency matrices is typically very different from the Erdős–Rényi model, where edges are independent. The authors argue that this procedure is similar to using a prior in a variational auto-encoder (VAE). However, the proposed model does not have a probabilistic formulation like VAE, e.g., there is no ELBO to constrain the discrepancy between the prior and the encoder distribution. I hope the authors could discuss this part and provide more convincing arguments.\n\n- Reproducing Kernel Banach Spaces (RKBS) has been mentioned in the abstract as a way to generalize a special class of exponential kernels to the non-symmetric case. However, the authors do not explicitly discuss it in the full paper. Some discussion about the classic Heat kernel (Kondor & Lafferty, 2002) involves RKBS and is only provided in the appendix. It is unclear what the exact technical contributions the authors have made regarding this part. It would be great to highlight the exact technical contribution w.r.t. RKBS in the main paper.\n\n- It has been mentioned in Proposition 1 and some other places that Q is the so-called heat source term that introduces noise in heat diffusion. However, in the proposed definition of Q in Proposition 1, no variables are random. The description is confusing.\n\n- In principle, one can add any “noise” to the adjacency matrix within this framework. The authors mention that they followed Veličković et al. (2019) to set the corrupting probability. However, the reference is about unsupervised graph representation learning, whereas this paper is about graph generation. It is unclear if the setup can be directly transferred. It would be great to provide an ablation study on the effect of different types of noise and the noise level.\n\n- Learning a set of input node representations N shared by all graphs makes sense for the purpose of generating adjacency matrices. However, it seems hard to utilize it if one wants to additionally generate node attributes that differ from graph to graph, e.g., atom types in molecule graphs. I hope the authors could discuss the generalization of this framework to attributed graphs.\n\n- The experiments need more baselines. The proposed method is currently only compared with the auto-regressive models. There are quite a few latest methods, including discrete diffusion [1] and continuous diffusion [2]. In particular, for continuous diffusion [2], the denoiser treats the adjacency matrix as an image. Therefore, it can handle both directed and undirected graphs in a unified manner. It would be great to add additional comparisons. \n\n- The experiments need some real-world datasets. Currently, only synthetic datasets are used, and the graph sizes are somewhat limited, e.g., the number of nodes in Section 6.1 is 15, and the number of nodes in Section 6.2 is 21. It would be great to include experiments on larger directed graphs from real-world problems. This would help strengthen the motivation of studying directed graph generation.\n\nMinor: \n\n- In proposition 1, it would be better to restate the goal explicitly so that the proposition becomes self-contained.\n\n- In the 1st paragraph of Section 3, “To this end, we consider that Eq. equation 2 is the output of our heat diffusion encoder…” should be “To this end, we consider that Eq. 2 is the output of our heat diffusion encoder…”\n\n[1] Vignac, C., Krawczuk, I., Siraudin, A., Wang, B., Cevher, V. and Frossard, P., 2023 DiGress: Discrete Denoising diffusion for graph generation. In International Conference on Learning Representations.\n\n[2] Yan, Q., Liang, Z., Song, Y., Liao, R. and Wang, L., 2024. Swingnn: Rethinking permutation invariance in diffusion models for graph generation. Transactions on Machine Learning Research.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "No", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper presents a denoising auto-encoder-based generative model for directed graphs. The key idea is to exploit a heat-diffusion-based encoder that maps the “noisy” adjacency matrix to noisy node representations. The directed graphs are then reconstructed by a decoder. Experiments on small graph datasets show promising results.", "claims_and_evidence_yn": "No", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "The requested changes are in my comments (cons) in Strengths And Weaknesses.", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "yKkJZbZPKg", "reviewer": "Reviewer_e7P4", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "**Strengths**\n- The paper is clearly written and easy to follow. I was able to reimplement their work with little effort\n- DGDK is clearly a novel and interesting method\n- The work is extremely well motivated\n\n**Weaknesses**\n- Scalability. It appears that these experiments only used 21 vertices and a handful of edges.\n- Some of the notation is a bit hard to follow in spots. I might have missed it but it was challenging to realize that $M$ is a constant matrix, for example.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This work presents a novel approach for generating graphs in a physically inspired way. Their physical model of choice is the driven heat equation, and they develop a technique that can generate a general graph (either directed or undirected) that is similar to the training set in one shot.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "Overall the paper is well written and self contained. I have a few comments:\n\n- It would be nice if the author could comment on the scalability of this technique, or approximations that might be applicable. Additionally, it would be nice to have a runtime table in the text\n- One of the core properties of directed graphs are the presence of strongly connected components -- It would be nice to see an analysis of whether DGDK preserves these clusters. I understand that the MMD is low, but it is possible to have a low MMD *and* destroy the community structure.\n- Please comment on the importance of your link predictor? Does using a symmetric edge-encoding appreciably change the efficacy of your decoding procedure?\n- Could you present an ablation on $\\gamma$?\n- The noise parameter, $\\rho$, was set to $1/n_i$. How does the handle higher noise regimes? Or lower ones?\n- Would a directed node encoder (e.g. DirGNN or MagNET) improve this model?", "broader_impact_concerns": "None", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "q0YMf33002", "reviewer": "Reviewer_i1zX", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n- The method focuses on directed graph generation, which is an important problem, particularly as the majority of the methods focused solely on undirected graphs\n- The proposed approach is original and novel, particularly the application of heat equation for modeling noise\n- The empirical setup & evaluation is sound and show the benefit of the proposed method\n\nWeaknesses:\n- Perhaps the biggest weakness of the paper is that the evaluation is performed solely on the digraphs and synthetic datasets, thus not providing sufficient evidence to the community whether it could be applied in practice to real-world datasets. For example on gene network inference [1], where having direction is particularly important.\n- In general I have found sections 2 and 3 slightly hard to read, with often the mathematical notation not backed by motivation. It would be worth spending more time on explaining the reasoning behind the methods and move some of the items to the appendix.\n\n\n[1] Marbach et al. Wisdom of crowds for robust gene network inference. Nature Methods, 9(8):796–804, 2012. ISSN 1548-7091.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "The authors propose a method for one-shot directed graph generation with denoising auto encoder-based generative model using the heat equation with the graph Laplacian. The input graph is perturbed using augmentation techniques which is then modeled with an encoder which generates noise based on the heat equation, thus not relying on a neural network to add noise. The output of the encoder is then used to reconstruct the original input with diffusion leveraging random walk Laplacian. The framework is evaluated on synthethic datasets.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "- The paper would be significantly strengthened by adding evaluation on at least one real-world dataset", "broader_impact_concerns": "No concerns.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: Overall, all three reviewers find the work interesting and support acceptance of the paper.

The paper was improved during the discussion paper, and most of the reviewers' concerns were addressed. Reviewer VUaS asked for comparisons to ref. [2], and the authors agreed to include it in the next version of the paper. I encourage the authors to include those comparisons in the final version.
Decision: Accept (poster)