summaries: The authors present a variational-autoencoder-based model for learning and generating graph structures. In particular, they propose a multiresolution graph network (MGN) that encodes a given graph in a hierarchical manner, i.e., at different levels of resolution.  Training the nodes of the coarsened graphs as latents of a variational auto-encoder allows for sampling a new graph in a bottom-up manner, i.e., at increasingly higher resolution.  The main contribution of this work consists of developing the framework of MGNs and demonstrating that they can be trained as hierarchical variational autoencoders, which yields the multiresolution graph variational autoencoders (MG-VAEs). The MG-VAEs represent generative models that allow for generating graphs in a multiresolution and equivariant manner. The authors evaluate these models in a range of different settings, from unsupervised representation learning, over link prediction on citation graphs, to the generation of molecules or graph-based image generation.

summaries: This study presents a multi-scale graph VAE with hierarchical graph coarsening. The Gumbel-max trick is employed to generate learnable hard partition for clustering, and permutation equivariant tensor operations are used (Kondor et al., 2018) to construct the group equivariance network. The authos find that the proposed framework exhibits competitive performance in graph generation, molecular generation, molecular representation learning, link prediction and graph-based image generation. Given the applicability of this framework, this proposed framework may be interesting to the graph generation community, and this paper is well-written.

summaries: This submission considers a multiresolution graph auto-encoder framework which is equivariant with respect to node permutation. Contrary to prior work on multi-scale graph generation, in this work the hierarchical structure in the encoder is learned through differentiable graph coarsening and is encouraged to produce balanced partitions of the graph. In the decoder, the graph is generated hierarchically and at each level, local adjacency matrices need to be predicted, rather than directly predicting the adjacency matrix of the full graph.  The method is evaluated on molecular generation, community graph generation,  and citation network generation in the main paper. The supplementary materials contains additional experiments on link prediction, unsupervised and supervised molecular property predictions, and graph-based image generation.

summaries: This paper proposes Multiresolution Equivariant Graph Variational Autoen-coders (MGVAE) which can learn and generate graphs in a multiresolution and equivariant fashion. Higher order message passing is used to encode the graph while maintaining learning mutually exclusive clusters so that coarsening into a lower resolution, thus creating a hierarchy of latent distributions. The model also maintain an end-to-end permutation equivariant with respect to node ordering. Their experimental results show that MGVAE achieves competitive re- sults with several generative tasks and graph link prediction etc.

abstractText: In this paper, we propose Multiresolution Equivariant Graph Variational Autoencoders (MGVAE), the first hierarchical generative model to learn and generate graphs in a multiresolution and equivariant manner. At each resolution level, MGVAE employs higher order message passing to encode the graph while learning to partition it into mutually exclusive clusters and coarsening into a lower resolution that eventually creates a hierarchy of latent distributions. MGVAE then constructs a hierarchical generative model to variationally decode into a hierarchy of coarsened graphs. Importantly, our proposed framework is end-to-end permutation equivariant with respect to node ordering. MGVAE achieves competitive results with several generative tasks including general graph generation, molecular generation, unsupervised molecular representation learning to predict molecular properties, link prediction on citation graphs, and graph-based image generation.

authors:
Peter Battaglia
Razvan Pascanu
Matthew Lai
Danilo Jimenez Rezende
Koray Kavukcuoglu

title: Interaction networks for learning about objects, relations and physics

venue: In Advances in Neural Information Processing Systems,

authors:
David M. Blei
Andrew Y. Ng
Michael I. Jordan

title: Latent dirichlet allocation

venue: J. Mach. Learn. Res.,

authors:
Nicola De Cao
Thomas Kipf

title: Molgan: An implicit generative model for small molecular graphs, 2018

authors:
Xu Chen
Xiuyuan Cheng
Stephane Mallat

title: Unsupervised deep haar scattering on graphs

authors:
Xiuyuan Cheng
Xu Chen
Stéphane Mallat

title: Deep haar scattering networks

authors:
Kai-Yang Chiang
Joyce Jiyoung Whang
Inderjit S. Dhillon

title: Scalable clustering of signed networks using balance normalized cut

venue: In Proceedings of the 21st ACM International Conference on Information and Knowledge Management,

authors:
Taco S. Cohen
Max Welling

title: Group equivariant convolutional networks

venue: Proceedings of The 33rd International Conference on Machine Learning,

authors:
Ronald R. Coifman
Mauro Maggioni

title: Diffusion wavelets

venue: Applied and Computational Harmonic Analysis,

authors:
Hanjun Dai
Azade Nazi
Yujia Li
Bo Dai
Dale Schuurmans

title: Scalable deep generative modeling for sparse graphs

venue: In Proceedings of the 37th International Conference on Machine Learning,

authors:
Inderjit Dhillon
Yuqiang Guan
Brian Kulis

title: A fast kernel-based multilevel algorithm for graph clustering

venue: In Proceedings of the Eleventh ACM SIGKDD International Conference on Knowledge Discovery in Data Mining,

authors:
Inderjit S. Dhillon
Yuqiang Guan
Brian Kulis

title: Weighted graph cuts without eigenvectors a multilevel approach

venue: IEEE Transactions on Pattern Analysis and Machine Intelligence,

authors:
Adji B. Dieng
Francisco J.R. Ruiz
David M. Blei
Michalis K. Titsias

title: Prescribed generative adversarial networks, 2019

authors:
David K Duvenaud
Dougal Maclaurin
Jorge Iparraguirre
Rafael Bombarell
Timothy Hirzel
Alan Aspuru-Guzik
Ryan P Adams

title: Convolutional networks on graphs for learning molecular fingerprints

authors:
Vijay Prakash Dwivedi
Chaitanya K. Joshi
Thomas Laurent
Yoshua Bengio
Xavier Bresson

title: Benchmarking graph neural networks

authors:
Jack Edmonds
Richard M. Karp

title: Theoretical improvements in algorithmic efficiency for network flow problems

venue: Journal of the ACM (JACM),

authors:
Alex Fout
Jonathon Byrd
Basir Shariat
Asa Ben-Hur

title: Protein interface prediction using graph convolutional networks

venue: In Proceedings of the 31st International Conference on Neural Information Processing Systems,

authors:
J. Gilmer
S.S. Schoenholz
P.F. Riley
O. Vinyals
G.E. Dahl

title: Neural message passing for quantum chemistry

authors:
Justin Gilmer
Samuel S. Schoenholz
Patrick F. Riley
Oriol Vinyals
George E. Dahl

venue: In Proceedings of the 34th International Conference on Machine Learning,

authors:
Xavier Glorot
Yoshua Bengio

title: Understanding the difficulty of training deep feedforward neural networks

venue: In Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics,

authors:
Ian Goodfellow
Jean Pouget-Abadie
Mehdi Mirza
Bing Xu
David Warde-Farley
Sherjil Ozair
Aaron Courville
Yoshua Bengio

title: Generative adversarial nets

title: Statistical theory of extreme values and some practical applications: a series of lectures

venue: US Govt. Print. Office, Number

authors:
Rafael Gmez-Bombarelli
Jennifer N. Wei
David Duvenaud
Jos Miguel Hernndez-Lobato
Benjamn Snchez-Lengeling
Dennis Sheberla
Jorge Aguilera-Iparraguirre
Timothy D. Hirzel
Ryan P. Adams
Aln Aspuru-Guzik

title: Automatic chemical design using a data-driven continuous representation of molecules

venue: ACS Central Science,

authors:
David K. Hammond
Pierre Vandergheynst
Rmi Gribonval

title: Wavelets on graphs via spectral graph theory

authors:
Martin Heusel
Hubert Ramsauer
Thomas Unterthiner
Bernhard Nessler
Sepp Hochreiter

title: Gans trained by a two time-scale update rule converge to a local nash equilibrium

authors:
P. Hohenberg
W. Kohn

title: Inhomogeneous electron gas

venue: Phys. Rev.,

authors:
Truong Son Hy
Shubhendu Trivedi
Horace Pan
Brandon M. Anderson
Risi Kondor

title: Predicting molecular properties with covariant compositional networks

venue: The Journal of Chemical Physics,

authors:
John Ingraham
Debora Marks

title: Variational inference for sparse and undirected models

authors:
Eric Jang
Shixiang Gu
Ben Poole

title: Categorical reparameterization with gumbel-softmax

venue: In ICLR,

authors:
Wengong Jin
Regina Barzilay
Tommi Jaakkola

title: Junction tree variational autoencoder for molecular graph generation

venue: In Proceedings of the 35th International Conference on Machine Learning,

authors:
Steven Kearnes
Kevin McCloskey
Marc Berndl
Vijay Pande
Patrick Riley

title: Molecular graph convolutions: moving beyond fingerprints

venue: Journal of Computer-Aided Molecular Design,

authors:
Diederik P. Kingma
Jimmy Ba

title: Adam: A method for stochastic optimization

venue: In 3rd International Conference on Learning Representations,

authors:
Diederik P. Kingma
Max Welling

title: Auto-Encoding Variational Bayes

venue: In 2nd International Conference on Learning Representations,

authors:
Thomas N. Kipf
Max Welling

title: Variational graph auto-encoders, 2016

authors:
Alexej Klushyn
Nutan Chen
Richard Kurle
Botond Cseke
Patrick van der Smagt

title: Learning hierarchical priors in vaes

authors:
Daphne Koller
Nir Friedman

title: Probabilistic Graphical Models: Principles and Techniques

authors:
Risi Kondor
Truong Son Hy
Horace Pan
Brandon M. Anderson
Shubhendu Trivedi

title: Covariant compositional networks for learning graphs, 2018

authors:
Nils M. Kriege
Pierre-Louis Giscard
Richard C. Wilson

title: On valid optimal assignment kernels and applications to graph classification

venue: In Proceedings of the 30th International Conference on Neural Information Processing Systems,

authors:
Yujia Li
Richard Zemel
Marc Brockschmidt
Daniel Tarlow

title: Gated graph sequence neural networks

venue: In Proceedings of ICLR’16,

authors:
Yujia Li
Oriol Vinyals
Chris Dyer
Razvan Pascanu
Peter W. Battaglia

title: Learning deep generative models of graphs. ICML’18

authors:
Renjie Liao
Yujia Li
Yang Song
Shenlong Wang
Will Hamilton
David K Duvenaud
Raquel Urtasun
Richard Zemel

title: Efficient graph generation with graph recurrent attention networks

authors:
Zinan Lin
Ashish Khetan
Giulia Fanti
Sewoong Oh

title: Pacgan: The power of two samples in generative adversarial networks

authors:
Jenny Liu
Aviral Kumar
Jimmy Ba
Jamie Kiros
Kevin Swersky

title: Graph normalizing flows

authors:
Qi Liu
Miltiadis Allamanis
Marc Brockschmidt
Alexander Gaunt

title: Constrained graph variational autoencoders for molecule design

authors:
Jianxin Ma
Peng Cui
Kun Kuang
Xin Wang
Wenwu Zhu

title: Disentangled graph convolutional networks

venue: In Proceedings of the 36th International Conference on Machine Learning,

authors:
Haggai Maron
Heli Ben-Hamu
Hadar Serviansky
Yaron Lipman

title: Provably powerful graph networks

authors:
Haggai Maron
Heli Ben-Hamu
Nadav Shamir
Yaron Lipman

title: Invariant and equivariant graph networks

venue: In International Conference on Learning Representations,

authors:
Haggai Maron
Ethan Fetaya
Nimrod Segol
Yaron Lipman

title: On the universality of invariant networks

authors:
Haggai Maron
Or Litany
Gal Chechik
Ethan Fetaya

title: On learning sets of symmetric elements

authors:
F. Monti
D. Boscaini
J. Masci
E. Rodola
J. Svoboda
M.M. Bronstein

title: Geometric deep learning on graphs and manifolds using mixture model cnns

venue: IEEE Conference on Computer Vision and Pattern Recognition (CVPR),

title: Chapter 19: Undirected graphical models (markov random fields)

venue: In Machine Learning: A Probabilistic Perspective,

authors:
M. Niepert
M. Ahmed
K. Kutzkov

title: Learning convolutional neural networks for graphs

venue: In Proceedings of the International Conference on Machine Learning,

authors:
Mathias Niepert
Mohamed Ahmed
Konstantin Kutzkov

venue: In Proceedings of The 33rd International Conference on Machine Learning,

authors:
B. Perozzi
R. Al-Rfou
S. Skiena

title: Deepwalk: Online learning of social representations

venue: Proceedings of the 20th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (KDD),

authors:
Alec Radford
Luke Metz
Soumith Chintala

title: Unsupervised representation learning with deep convolutional generative adversarial networks, 2016

authors:
Raghunathan Ramakrishnan
Pavlo O Dral
Matthias Rupp
O Anatole von Lilienfeld

title: Quantum chemistry structures and properties of 134 kilo molecules

venue: Scientific Data,

authors:
Rajesh Ranganath
Dustin Tran
David Blei

title: Hierarchical variational models

authors:
Lars Ruddigkeit
Ruud van Deursen
Lorenz C. Blum
Jean-Louis Reymond

title: Enumeration of 166 billion organic small molecules in the chemical universe database gdb-17

venue: Journal of Chemical Information and Modeling,

authors:
Havard Rue
Leonhard Held

title: Gaussian markov random fields: Theory and applications

venue: In Monographs on Statistics and Applied Probability,

authors:
Raif Rustamov
Leonidas J Guibas

title: Wavelets on graphs via deep learning

authors:
Franco Scarselli
Marco Gori
Ah Chung Tsoi
Markus Hagenbuchner
Gabriele Monfardini

title: The graph neural network model

venue: IEEE Transactions on Neural Networks,

title: pytorch-fid: FID Score for PyTorch

authors:
P. Sen
G.M. Namata
M. Bilgic
L. Getoor
B. Gallagher
T. Eliassi-Rad

title: Collective classification in network data

venue: AI Magazine,

authors:
Hadar Serviansky
Nimrod Segol
Jonathan Shlomi
Kyle Cranmer
Eilam Gross
Haggai Maron
Yaron Lipman

title: Set2graph: Learning graphs from sets

authors:
Nino Shervashidze
Pascal Schweitzer
Erik Jan van Leeuwen
Kurt Mehlhorn
Karsten M. Borgwardt

title: Weisfeiler-lehman graph kernels

venue: Journal of Machine Learning Research,

authors:
Chence Shi
Minkai Xu
Zhaocheng Zhu
Weinan Zhang
Ming Zhang
Jian Tang

title: Graphaf: a flow-based autoregressive model for molecular graph generation

authors:
Donghyuk Shin
Si Si
Inderjit S. Dhillon

title: Multi-scale link prediction

authors:
Si Si
Donghyuk Shin
Inderjit S Dhillon
Beresford N Parlett

title: Multi-scale spectral decomposition of massive graphs

authors:
Martin Simonovsky
Nikos Komodakis

title: Graphvae: Towards generation of small graphs using variational autoencoders

authors:
Akash Srivastava
Lazar Valkov
Chris Russell
Michael U. Gutmann
Charles Sutton

title: Veegan: Reducing mode collapse in gans using implicit variational learning

authors:
Teague Sterling
John J. Irwin

title: Zinc 15 ligand discovery for everyone

authors:
L. Tang
H. Liu

title: Leveraging social media networks for classification

venue: Data Mining and Knowledge Discovery,

authors:
Erik Henning Thiede
Truong Son Hy
Risi Kondor

title: The general theory of permutation equivarant neural networks and higher order graph variational encoders

authors:
Arash Vahdat
Jan Kautz

title: Nvae: A deep hierarchical variational autoencoder

authors:
Petar Velikovi
Guillem Cucurull
Arantxa Casanova
Adriana Romero
Pietro Li
Yoshua Bengio

title: Graph attention networks

authors:
Martin J. Wainwright
Michael I. Jordan

title: A variational principle for graphical models

venue: New Directions in Statistical Signal Processing,

authors:
Hao Wu
Jonas Köhler
Frank Noe

title: Stochastic normalizing flows

authors:
Bingbing Xu
Huawei Shen
Qi Cao
Yunqi Qiu
Xueqi Cheng

title: Graph wavelet neural network

authors:
Yiding Yang
Zunlei Feng
Mingli Song
Xinchao Wang

title: Factorizable graph convolutional networks

authors:
Zhitao Ying
Jiaxuan You
Christopher Morris
Xiang Ren
Will Hamilton
Jure Leskovec

title: Hierarchical graph representation learning with differentiable pooling

authors:
Jiaxuan You
Rex Ying
Xiang Ren
William Hamilton
Jure Leskovec

title: GraphRNN: Generating realistic graphs with deep auto-regressive models

authors:
Jiaxuan You
Rex Ying
Xiang Ren
William L. Hamilton
Jure Leskovec

title: Code for graphrnn: Generating realistic graphs with deep auto-regressive model. 2018b

authors:
Dawei Zhou
Lecheng Zheng
Jiejun Xu
Jingrui He

title: Misc-gan: A multi-scale generative model for graphs

venue: Frontiers in Big Data,

heading: 1 INTRODUCTION

text:
Understanding graphs in a multiscale and multiresolution perspective is essential for capturing the structure of molecules, social networks, or the World Wide Web. Graph neural networks (GNNs) utilizing various ways of generalizing the concept of convolution to graphs (Scarselli et al., 2009) (Niepert et al., 2016b) (Li et al., 2016) have been widely applied to many learning tasks, including modeling physical systems (Battaglia et al., 2016), finding molecular representations to estimate quantum chemical computation (Duvenaud et al., 2015) (Kearnes et al., 2016) (Gilmer et al., 2017b) (Hy et al., 2018), and protein interface prediction (Fout et al., 2017). One of the most popular types of GNNs is message passing neural nets (MPNNs) that are constructed based on the message passing scheme in which each node propagates and aggregates information, encoded by vectorized messages, to and from its local neighborhood. While this framework has been immensely successful in many applications, it lacks the ability to capture the multiscale and multiresolution structures that are present in complex graphs (Rustamov & Guibas, 2013) (Chen et al., 2014) (Cheng et al., 2015) (Xu et al., 2019).
Ying et al. (2018) proposed a multiresolution graph neural network that employs a differential pooling operator to coarsen the graph. While this approach is effective in some settings, it is based on soft assigment matrices, which means that (a) the sparsity of the graph is quickly lost in higher layers (b) the algorithm isn’t able to learn an actual hard clustering of the vertices. The latter is important in applications such as learning molecular graphs, where clusters should ideally be interpretable as concrete subunits of the graphs, e.g., functional groups.
In contrast, in this paper we propose an arhictecture called Multiresolution Graph Network (MGN), and its generative cousin, Multiresolution Graph Variational Autoencoder (MGVAE), which explicitly learn a multilevel hard clustering of the vertices, leading to a true multiresolution hierarchy. In the decoding stage, to “uncoarsen” the graph, MGVAE needs to generate local adjacency matrices, which is inherently a second order task with respect to the action of permutations on the vertices, hence MGVAE needs to leverage the recently developed framework of higher order permutation equivariant message passing networks (Hy et al., 2018; Maron et al., 2019b).
Learning to generate graphs with deep generative models is a difficult problem because graphs are combinatorial objects that typically have high order correlations between their discrete substructures (subgraphs) (You et al., 2018a) (Li et al., 2018) (Liao et al., 2019) (Liu et al., 2019) (Dai et al., 2020). Graph-based molecular generation (Gmez-Bombarelli et al., 2018) (Simonovsky &
Komodakis, 2018) (Cao & Kipf, 2018) (Jin et al., 2018) (Thiede et al., 2020) involves further challenges, including correctly recognizing chemical substructures, and importantly, ensuring that the generated molecular graphs are chemically valid. MGN allows us to extend the existing model of variational autoencoders (VAEs) with a hierarchy of latent distributions that can stochastically generate a graph in multiple resolution levels. Our experiments show that having a flexible clustering procedure from MGN enables MGVAE to detect, reconstruct and finally generate important graph substructures, especially chemical functional groups.

heading: 2 RELATED WORK

text:
There have been significant advances in understanding the invariance and equivariance properties of neural networks in general (Cohen & Welling, 2016a) (Cohen & Welling, 2016b), of graph neural networks (Maron et al., 2019b), of neural networks learning on sets (Zaheer et al., 2017) (Serviansky et al., 2020) (Maron et al., 2020), along with their expressive power on graphs (Maron et al., 2019c) (Maron et al., 2019a). Our work is in line with group equivariant networks operating on graphs and sets. Multiscale, multilevel, multiresolution and coarse-grained techniques have been widely applied to graphs and discrete domains such as diffusion wavelets (Coifman & Maggioni, 2006); spectral wavelets on graphs (Hammond et al., 2011); finding graph wavelets based on partitioning/clustering (Rustamov & Guibas, 2013); graph clustering and finding balanced cuts on large graphs (Dhillon et al., 2005) (Dhillon et al., 2007) (Chiang et al., 2012) (Si et al., 2014); and link prediction on social networks (Shin et al., 2012). Prior to our work, some authors such as (Zhou et al., 2019) proposed a multiscale generative model on graphs using GAN (Goodfellow et al., 2014), but the hierarchical structure was built by heuristics algorithm, not learnable and not flexible to new data that is also an existing limitation of the field. In general, our work exploits the powerful group equivariant networks to encode a graph and to learn to form balanced partitions via back-propagation in a datadriven manner without using any heuristics as in the existing works.
In the field of deep generative models, it is generally recognized that introducing a hierarchy of latents and adding stochasticity among latents leads to more powerful models capable of learning more complicated distributions (Blei et al., 2003) (Ranganath et al., 2016) (Ingraham & Marks, 2017) (Klushyn et al., 2019) (Wu et al., 2020) (Vahdat & Kautz, 2020). Our work combines the hierarchical variational autoencoder with learning to construct the hierarchy that results into a generative model able to generate graphs at many resolution levels.

heading: 3 MULTIRESOLUTION GRAPH NETWORK

heading: 3.1 CONSTRUCTION

text:
An undirected weighted graph G = (V,E ,A,Fv,Fe) with node set V and edge set E is represented by an adjacency matrix A ∈ N∣V ∣×∣V ∣, where Aij > 0 implies an edge between node vi and vj with weight Aij (e.g., Aij ∈ {0,1} in the case of unweighted graph); while node features are represented by a matrix Fv ∈ R∣V ∣×dv , and edge features are represented by a tensor Fe ∈ R∣V ∣×∣V ∣×de . The second-order tensor representation of edge features is necessary for our higher-order message passing networks described in the next section. Indeed, Fv can be encoded in the diagonal of Fe. Definition 1. A K-cluster partition of graph G is a partition of the set of nodes V into K mutually exclusive clusters V1, ..,VK . Each cluster corresponds to an induced subgraph Gk = G[Vk].
Definition 2. A coarsening of G is a graph G̃ of K nodes defined by a K-cluster partition in which node ṽk of G̃ corresponds to the induced subgraph Gk. The weighted adjacency matrix Ã ∈ NK×K of G̃ is
Ãkk′ = { 1 2 ∑vi,vj∈Vk Aij , if k = k ′,
∑vi∈Vk,vj∈Vk′ Aij , if k ≠ k ′,
where the diagonal of Ã denotes the number of edges inside each cluster, while the off-diagonal denotes the number of edges between two clusters.
Fig. 3.1 shows an example of Defs. 1 and 2: a 3-cluster partition of the Aspirin C9H8O4 molecular graph and its coarsening graph. Def. 3 defines the multiresolution of graph G in a bottom-up manner
in which the bottom level is the highest resolution (e.g., G itself) while the top level is the lowest resolution (e.g., G is coarsened into a single node).
Definition 3. An L-level coarsening of a graph G is a series of L graphs G(1), ..,G(L) in which
1. G(L) is G itself.
2. For 1 ≤ ` ≤ L − 1, G(`) is a coarsening graph of G(`+1) as defined in Def. 2. The number of nodes in G(`) is equal to the number of clusters in G(`+1).
3. The top level coarsening G(1) is a graph consisting of a single node, and the corresponding adjacency matrix A(1) ∈ N1×1.
Definition 4. An L-level Multiresolution Graph Network (MGN) of a graph G consists of L−1 tuples of five network components {(c(`),e(`)local,d (`) local,d (`) global,p (`))}L`=2. The `-th tuple encodes G (`) and transforms it into a lower resolution graph G(`−1) in the higher level. Each of these network components has a separate set of learnable parameters (θ(`)1 ,θ (`) 2 ,θ (`) 3 ,θ (`) 4 ,θ (`) 5 ). For simplicity, we collectively denote the learnable parameters as θ, and drop the superscript. The network components are defined as follows:
1. Clustering procedure c(G(`);θ), which partitions graph G(`) into K clusters V
(`) 1 , . . . ,V (`) K . Each cluster is an induced subgraph G (`) k of G (`) with adjacency matrix A
(`) k .
2. Local encoder elocal(G (`) k ;θ), which is a permutation equivariant (see Defs. 5, 6) graph
neural network that takes as input the subgraph G(`)k , and outputs a set of node latents Z
(`) k represented as a matrix of size ∣V (`) k ∣ × dz .
3. Local decoder dlocal(Z (`) k ;θ), which is a permutation equivariant neural network that tries
to reconstruct the subgraph adjacency matrixA(`)k for each cluster from the local encoder’s output latents.
4. (Optional) Global decoder dglobal(Z(`);θ), which is a permutation equivariant neural network that reconstructs the full adjacency matrixA(`) from all the node latents ofK clusters Z(`) =⊕kZ (`) k represented as a matrix of size ∣V (`)∣ × dz .
5. Pooling network p(Z(`)k ;θ), which is a permutation invariant (see Defs. 5, 6) neural network that takes the set of node latents Z(`)k and outputs a single cluster latent z̃ (`) k ∈ dz .
The coarsening graph G(`−1) has adjacency matrix A(`−1) built as in Def. 2, and the corresponding node features Z(`−1) =⊕k z̃ (`) k represented as a matrix of size K × dz .
Algorithmically, MGN works in a bottom-up manner as a tree-like hierarchy starting from the highest resolution graph G(L), going to the lowest resolution G(1) (see Fig. 3.1). Iteratively, at `-th level, MGN partitions the current graph into K clusters by running the clustering procedure c(`). Then, the local encoder e(`)local and local decoder d (`) global operate on each of the K subgraphs separately, and can be executed in parallel. This encoder/decoder pair is responsible for capturing the local structures. Finally, the pooling network p(`) shrinks each cluster into a single node of the next level.
Optionally, the global decoder d(`)global makes sure that the whole set of node latents Z (`) is able to capture the inter-connection between clusters.
In terms of time and space complexity, MGN is more efficient than existing methods in the field. The cost of global decoding the highest resolution graph is proportional to ∣V ∣2. For example, while the encoder can exploit the sparsity of the graph and has complexity O(∣E ∣), a simple dot-product decoder dglobal(Z) = sigmoid(ZZT ) has both time and space complexity ofO(∣V ∣2) which is infeasible for large graphs. In contrast, the cost of running K local dot-product decoders is O(∣V ∣2/K), which is approximately K times more efficient.

heading: 3.2 HIGHER ORDER MESSAGE PASSING

text:
In this paper we consider permutation symmetry, i.e., symmetry to the action of the symmetric group, Sn. An element σ ∈ Sn is a permutation of order n, or a bijective map from {1, . . . , n} to {1, . . . , n}. The action of Sn on an adjacency matrix A ∈ Rn×n and on a latent matrix Z ∈ Rn×dz are [σ ⋅A]i1,i2 = Aσ−1(i1),σ−1(i2), [σ ⋅Z]i,j = Zσ−1(i),j , σ ∈ Sn. Here, the adjacency matrix A is a second order tensor with a single feature channel, while the latent matrix Z is a first order tensor with dz feature channels. In general, the action of Sn on a k-th order tensor X ∈ Rn k ×d (the last index denotes the feature channels) is defined similarly as:
[σ ⋅X ]i1,..,ik,j = Xσ−1(i1),..,σ−1(ik),j , σ ∈ Sn.
Network components of MGN (as defined in Sec. 3.1) at each resolution level must be either equivariant, or invariant with respect to the permutation action on the node order of G(`). Formally, we define these properties in Def. 5. Definition 5. An Sn-equivariant (or permutation equivariant) function is a function f ∶Rn k ×d → Rn k′ ×d′ that satisfies f(σ ⋅ X ) = σ ⋅ f(X ) for all σ ∈ Sn and X ∈ Rn k ×d. Similarly, we say that f is Sn-invariant (or permutation invariant) if and only if f(σ ⋅X ) = f(X ).
Definition 6. An Sn-equivariant network is a function f ∶ Rn k ×d → Rn k′ ×d′ defined as a composition of Sn-equivariant linear functions f1, .., fT and Sn-equivariant nonlinearities γ1, .., γT :
f = γT ○ fT ○ .. ○ γ1 ○ f1.
On the another hand, an Sn-invariant network is a function f ∶ Rn k ×d → R defined as a composition of an Sn-equivariant network f ′ and an Sn-invariant function on top of it, e.g., f = f ′′ ○ f ′.
In order to build higher order equivariant networks, we revisit some basic tensor operations: the tensor product A⊗B and tensor contraction A↓x1,..,xp (details and definitions are Sec. A). It can be shown that these tensor operations respect permutation equivariance (Kondor et al., 2018). Based on these tensor contractions and Def. 5, we can construct the second-order Sn-equivariant networks as in Def. 6 (see Example 1): f = γ ○MT ○ .. ○ γ ○M1. The second-order networks are particularly essential for us to extend the original variational autoencoder (VAE) (Kingma & Welling,
2014) model that approximates the posterior distribution by an isotropic Gaussian distribution with a diagonal covariance matrix and uses a fixed prior distributionN (0,1). In constrast, we generalize by modeling the posterior byN (µ,Σ) in which Σ is a full covariance matrix, and we learn an adaptive parameterized priorN (µ̂, Σ̂) instead of a fixed one. Only the second-order encoders can output a permutation equivariant full covariance matrix, while lower-order networks such as MPNNs are unable to. See Sec. 4.2, B and C for details. Example 1. The second order message passing has the message H0 ∈ R∣V ∣×∣V ∣×(dv+de) initialized by promoting the node features Fv to a second order tensor (e.g., we treat node features as self-loop edge features), and concatenating with the edge features Fe. Iteratively,
Ht = γ(Mt), Mt =Wt[⊕ i,j (A⊗Ht−1)↓i,j],
where A ⊗ Ht−1 results in a fourth order tensor while ↓i,j contracts it down to a second order tensor along the i-th and j-th dimensions, ⊕ denotes concatenation along the feature channels, and Wt denotes a multilayer perceptron on the feature channels. We remark that the popular MPNNs (Gilmer et al., 2017b) is a lower-order one and a special case in which Mt = D−1AHt−1Wt−1 where Dii = ∑jAij is the diagonal matrix of node degrees. The messageHT of the last iteration is still second order, so we contract it down to the first order latent Z =⊕iHT ↓i.

heading: 3.3 LEARNING TO CLUSTER

text:
Definition 7. A clustering of n objects into k clusters is a mapping π ∶ {1, .., n}→ {1, .., k} in which π(i) = j if the i-th object is assigned to the j-th cluster. The inverse mapping π−1(j) = {i ∈ [1, n] ∶ π(i) = j} gives the set of all objects assigned to the j-th cluster. The clustering is represented by an assignment matrix Π ∈ {0,1}n×k such that Πi,π(i) = 1. Definition 8. The action of Sn on a clustering π of n objects into k clusters and its corresponding assignment matrix Π are
[σ ⋅ π](i) = π(σ−1(i)), [σ ⋅Π]i,j = Πσ−1(i),j , σ ∈ Sn. Definition 9. Let N be a neural network that takes as input a graph G of n nodes, and outputs a clustering π of k clusters. N is said to be equivariant if and only if N (σ ⋅ G) = σ ⋅N (G) for all σ ∈ Sn.
From Def. 9, intuitively the assignement matrix Π still represents the same clustering if we permute its rows. The learnable clustering procedure c(G(`);θ) is built as follows:
1. A graph neural network parameterized by θ encodes graph G(`) into a first order tensor of K feature channels p̃(`) ∈ R∣V (`) ∣×K . 2. The clustering assignment is determined by a row-wise maximum pooling operation:
π(`)(i) = arg max k∈[1,K] p̃ (`) i,k (1)
that is an equivariant clustering in the sense of Def. 9.
A composition of an equivariant function (e.g., graph net) and an equivariant function (e.g., maximum pooling given in Eq. 1) is still an equivariant function with respect to the node permutation. Thus, the learnable clustering procedure c(G(`);θ) is permutation equivariant.
In practice, in order to make the clustering procedure differentiable for backpropagation, we replace the maximum pooling in Eq. 1 by sampling from a categorical distribution. Let π(`)(i) be a categorical variable with class probabilities p(`)i,1 , .., p (`) i,K computed as softmax from p̃ (`) i,∶ . The Gumbel-max trick (Gumbel, 1954)(Maddison et al., 2014)(Jang et al., 2017) provides a simple and efficient way to draw samples π(`)(i):
Π (`) i = one-hot(arg max
k∈[1,K]
[gi,k + log p (`) i,k ]),
where gi,1, .., gi,K are i.i.d samples drawn from Gumbel(0,1). Given the clustering assignment matrix Π(`), the coarsened adjacency matrix A(`−1) (see Defs. 1 and 2) can be constructed as Π(`) T A(`)Π(`).
It is desirable to have a balanced K-cluster partition in which clusters V(`)1 , ..,V (`) K have similar sizes that are close to ∣V(`)∣/K. The local encoders tend to generalize better for same-size subgraphs. We want the distribution of nodes into clusters to be close to the uniform distribution. We enforce the clustering procedure to produce a balanced cut by minimizing the following Kullback–Leibler divergence:
DKL(P ∣∣Q) = K
∑ k=1
P (k) log P (k)
Q(k) where P = (
∣V (`) 1 ∣ ∣V(`)∣ , .., ∣V (`) K ∣ ∣V(`)∣ ), Q = ( 1 K , .., 1 K ). (2)
The whole construction of MGN is equivariant with respect to node permutations of G. In the case of molecular property prediction, we want MGN to learn to predict a real value y ∈ R for each graph G while learning to find a balanced cut in each resolution to construct a hierarchical structure of latents and coarsen graphs. The total loss function is
LMGN(G, y) = ∣∣f( L
⊕ `=1
R(Z(`))) − y∣∣ 2
2 +
L
∑ `=1
λ(`)DKL(P (`) ∣∣Q(`)), (3)
where f is a multilayer perceptron, ⊕ denotes the vector concatenation, R is a readout function that produces a permutation invariant vector of size d given the latent Z ∣V (`) ∣×d at the `-th resolution, λ(`) ∈ R is a hyperparamter, and DKL(P (`)∣∣Q(`)) is the balanced-cut loss as defined in Eq. 2.

heading: 4 HIERARCHICAL GENERATIVE MODEL

text: In this section, we introduce our hierarchical generative model for multiresolution graph generation based on variational principles.

heading: 4.1 BACKGROUND ON GRAPH VARIATIONAL AUTOENCODER

text:
Suppose that we have input data consisting of m graphs (data points) G = {G1, ..,Gm}. The standard variational autoencoders (VAEs), introduced by Kingma & Welling (2014) have the following generation process, in which each data graph Gi for i ∈ {1,2, ..,m} is generated independently:
1. Generate the latent variables Z = {Z1, ..,Zm}, where each Zi ∈ R∣Vi∣×dz is drawn i.i.d. from a prior distribution p0 (e.g., standard Normal distribution N (0,1)). 2. Generate the data graph Gi ∼ pθ(Gi∣Zi) from the model conditional distribution pθ.
We want to optimize θ to maximize the likelihood pθ(G) = ∫ pθ(Z)pθ(G∣Z)dZ . However, this requires computing the posterior distribution pθ(G∣Z) =∏mi=1 pθ(Gi∣Zi), which is usually intractable. Instead, VAEs apply the variational principle, proposed by Wainwright & Jordan (2005), to approximate the posterior distribution as qφ(Z ∣G) =∏mi=1 qφ(Zi∣Gi) via amortized inference and maximize the evidence lower bound (ELBO) that is a lower bound of the likelihood:
LELBO(φ, θ) = Eqφ(Z ∣G)[log pθ(G∣Z)] −DKL(qφ(Z ∣G)∣∣p0(Z))
= m
∑ i=1
[Eqφ(Zi∣Gi)[log pθ(Gi∣Zi)] −DKL(qφ(Zi∣Gi)∣∣p0(Zi))]. (4)
The probabilistic encoder qφ(Z ∣G), the approximation to the posterior of the generative model pθ(G,Z), is modeled using equivariant graph neural networks (see Example 1) as follows. Assume the prior over the latent variables to be the centered isotropic multivariate Gaussian pθ(Z) = N (Z; 0, I). We let qφ(Zi∣Gi) be a multivariate Gaussian with a diagonal covariance structure:
log qφ(Zi∣Gi) = logN (Zi;µi, σ 2 i I), (5)
where µi, σi ∈ R∣Vi∣×dz are the mean and standard deviation of the approximate posterior output by two equivariant graph encoders. We sample from the posterior qφ by using the reparameterization trick: Zi = µi + σi ⊙ , where ∼ N (0, I) and ⊙ is the element-wise product.
On the another hand, the probabilistic decoder pθ(Gi∣Zi) defines a conditional distribution over the entries of the adjacency matrix Ai: pθ(Gi∣Zi) = ∏(u,v)∈V2i pθ(Aiuv = 1∣Ziu,Ziv). For example, Kipf & Welling (2016) suggests a simple dot-product decoder that is trivially equivariant: pθ(Aiuv = 1∣Ziu,Ziv) = γ(Z T iuZiv), where γ denotes the sigmoid function.

heading: 4.2 MULTIRESOLUTION VAES

text:
Based on the construction of multiresolution graph network (see Sec. 3.1), the latent variables are partitioned into disjoint groups, Zi = {Z (1) i ,Z (2) i , ..,Z (L) i } where Z (`) i = {[Z (`) i ]k ∈ R∣[V (`) i ]k ∣×dz}k is the set of latents at the `-th resolution level in which the graph G (`) i is partitioned into a number of clusters [G(`)i ]k.
In the area of normalzing flows (NFs), Wu et al. (2020) has shown that stochasticity (e.g., a chain of stochastic sampling blocks) overcomes expressivity limitations of NFs. In general, our MGVAE is a stochastic version of the deterministic MGN such that stochastic sampling is applied at each resolution level in a bottom-up manner. The prior (Eq. 6) and the approximate posterior (Eq. 7) are represented by
p(Zi) = L
∏ `=1
p(Z (`) i ) =
L
∏ `=1 ∏ k
p([Z (`) i ]k), (6)
qφ(Zi∣Gi) = qφ(Z (L) i ∣G (L) i )
1
∏ `=L−1
qφ(Z (`) i ∣Z (`+1) i ,G (`) i ), (7)
in which each conditional in the approximate posterior are in the form of factorial Normal distributions, in particular
qφ(Z (`) i ∣Z (`+1) i ,G (`) i ) =∏
k
qφ([Z (`) i ]k∣Z (`+1) i , [G (`) i ]k),
where each probabilistic encoder qφ([Z (`) i ]k∣Z (`+1) i , [G (`) i ]k) operates on a subgraph [G (`) i ]k as follows:
• The pooling network p(`+1) shrinks the latent Z(`+1)i into the node features of G (`) i as in
the construction of MGN (see Def. 4).
• The local (deterministic) graph encoder d(`)local encodes each subgraph [G (`) i ]k into a mean
vector and a diagonal covariance matrix (see Eq. 5). A second order encoder can produce a positive semidefinite non-diagonal covariance matrix, that can be interpreted as a Gaussian Markov Random Fields (details in Sec. B). The new subgraph latent [Z(`)i ]k is sampled by the reparameterization trick.
The prior can be either the isotropic Gaussian N (0,1) as in standard VAEs, or be implemented as a parameterized Gaussian N (µ̂, Σ̂) where µ̂ and Σ̂ are learnable equivariant functions (details in Sec. C). The reparameterization trick for conventionalN (0,1) prior is the same as in Sec. 4.1, while the new one for the generalized and learnable priorN (µ̂, Σ̂) is given in Sec. B. On the another hand, the probabilistic decoder pθ(G (1) i , ..,G (L) i ∣Z (1) i , ..,Z (L) i ) defines a conditional distribution over all subgraph adjacencies at each resolution level:
pθ(G (1) i , ..,G (L) i ∣Z (1) i , ..,Z (L) i ) =∏
`
pθ(G (`) i ∣Z (`) i ) =∏
` ∏ k
pθ([A (`) i ]k ∣[Z (`) i ]k).
Extending from Eq. 4, we write our multiresolution variational lower bound LMGVAE(φ, θ) on log p(G) compactly as
LMGVAE(φ, θ) =∑ i ∑ ` [E qφ(Z (`) i ∣G (`) i ) [log pθ(G (`) i ∣Z (`) i )] −DKL(qφ(Z (`) i ∣G (`) i )∣∣p0(Z (`) i ))], (8)
where the first term denotes the reconstruction loss (e.g., ∣∣A(`)i − Â (`) i ∣∣ where A (L) i is Gi itself, A (`<L) i is the adjacency produced by MGN at level `, and Â (`) i are the reconstructed ones by the decoders); and the second term is indeed DKL(N (µ (`) i ,Σ (`) i )∣∣N (µ̂ (`), Σ̂(`))) where µ(`)i ∈ R ∣V (`) i ∣×d and Σ(`)i ∈ R ∣V (`) i ∣×∣V (`) i ∣×d are the mean and covariance tensors produced by the `-th encoder for graph Gi, while µ̂(`)and Σ̂(`) are learnable ones in an equivariant manner as in Sec. C. In general, the overall optimization is given as follows:
min φ,θ,{µ̂(`),Σ̂(`)}`
LMGVAE(φ, θ;{µ̂ (`), Σ̂(`)}`) +∑
i,`
λ(`)DKL(P (`) i ∣∣Q (`) i ), (9)
where φ denotes all learnable parameters of the encoders, θ denotes all learnable parameters of the decoders, andDKL(P (`) i ∣∣Q (`) i ) is the balanced-cut loss for graph Gi at level ` as defined in Sec. 3.3.

heading: 5 EXPERIMENTS

text: Many more experimental results and details are presented in the Sec. D of the Appendix.

heading: 5.1 MOLECULAR GRAPH GENERATION

text:
We examine the generative power of MGN and MGVAE in the challenging task of molecule generation, in which the graphs are highly structured. We demonstrate that MGVAE is the first hierarchical graph VAE model generating graphs in a permutation-equivariant manner that is competitive against autoregressive results. We train on two datasets that are standard in the field:
1. QM9 (Ruddigkeit et al., 2012) (Ramakrishnan et al., 2014): contains 134K organic molecules with up to nine atoms (C, H, O, N, and F) out of the GDB-17 universe of molecules. 2. ZINC (Sterling & Irwin, 2015): contains 250K purchasable drug-like chemical compounds with up to twenty-three heavy atoms.
We only use the graph features as the input, including the adjacency matrix, the one-hot vector of atom types (e.g., carbon, hydrogen, etc.) and the bond types (single bond, double bond, etc.) without any further domain knowledge from chemistry or physics. First, we train autoencoding task of reconstructing the adjacency matrix and node features. We use a learnable equivariant prior (see Sec. C) instead of the conventional N (0,1). Then, we generate 5,000 different samples from the prior, and decode each sample into a generated graph (see Fig. 3). We implement our graph construction (decoding) in two approaches:
1. All-at-once: We reconstruct the whole adjacency matrix by running the probabilistic decoder (see Sec. 4). MGVAE enables us to generate a graph at any given resolution level `. In this particular case, we select the highest resolution ` = L. This approach of decoding preserves equivariance, but is harder to train. On ZINC, we extract several chemical/atomic features from RDKit as the input for the encoders to reach a good convergence in training. 2. Autoregressive: The graph is constructed iteratively by adding one edge in each iteration, similarly to (Liu et al., 2018). But this approach does not respect permutation equivariance.
In our setting for small molecules, L = 3 and K = 2`−1 for the `-th level. We compare our methods with other graph-based generative models including GraphVAE (Simonovsky & Komodakis, 2018), CGVAE (Liu et al., 2018), MolGAN (Cao & Kipf, 2018), and JT-VAE (Jin et al., 2018). We evaluate the quality of generated molecules in three metrics: (i) validity, (ii) novelty and (iii) uniqueness as
the percentage of the generated molecules that are chemically valid, different from all molecules in the training set, and not redundant, respectively. Because of high complexity, we only train on a small random subset of examples while all other methods are trained on the full datasets. Our models are equivalent with the state-of-the-art, even with a limited training set (see Table 1). Admittedly, molecule generation is a somewhat subject task that can only be evaluated with objective numerical measures up to a certain point. Qualitatively, however the molecules that MGVAE generates are as good as the state of the art, in some cases better in terms of producing several high-quality drug-like molecules with complicated functional groups and structures. Many further samples generated by MGVAE and their analysis can be found in the Appendix.

heading: 5.2 GENERAL GRAPH GENERATION BY MGVAE

text:
We further examine the expressive power of hierarchical latent structure of MGVAE in the task of general graph generation. We choose two datasets from GraphRNN paper (You et al., 2018a):
1. Community-small: A synthetic dataset of 100 2-community graphs where 12 ≤ ∣V ∣ ≤ 20. 2. Ego-small: 200 3-hop ego networks extracted from the Citeseer network (Sen et al., 2008)
where 4 ≤ ∣V ∣ ≤ 18.
The datasets are generated by the scripts from the GraphRNN codebase (You et al., 2018b). We keep 80% of the data for training and the rest for testing. We evaluate our generated graphs by computing Maximum Mean Discrepancy (MMD) distance between the distributions of graph statistics on the test set and the generated set as proposed by (You et al., 2018a). The graph statistics are node degrees, clustering coefficients, and orbit counts. As suggested by (Liu et al., 2019), we execute 15 runs with different random seeds, and we generate 1,024 graphs for each run, then we average the results over 15 runs. We compare MGVAE against GraphVAE (Simonovsky & Komodakis, 2018), DeepGMG (Li et al., 2018), GraphRNN (You et al., 2018a), GNF (Liu et al., 2019), and GraphAF (Shi et al., 2020). The baselines are taken from GNF paper (Liu et al., 2019) and GraphAF paper (Shi et al., 2020). In our setting of (all-at-once) MGVAE, we implement only L = 2 levels of resolution and K = 2` clusters for each level. Our encoders have 10 layers of message passing. Instead of using a high order equivariant network as the global decoder for the bottom resolution, we only implement a simple fully connected network that maps the latent Z(L) ∈ R∣V ∣×dz into an adjacency matrix of size ∣V ∣ × ∣V ∣. For the ego dataset in particular, we implement the learnable equivariant prior as in Sec. B and Sec.C. Table 2 includes our quantitative results in comparison with other methods. MGVAE outperforms all competing methods. Figs. 10 11 show some generated examples and training examples on the 2-community and ego datasets.

heading: 6 CONCLUSION

text: We introduced MGVAE built upon MGN, the first generative model to learn and generate graphs in a multiresolution and equivariant manner. The key idea of MGVAE is learning to construct a series of coarsened graphs along with a hierarchy of latent distributions in the encoding process while learning to decode each latent into the corresponding coarsened graph at every resolution level. MGVAE achieves state-of-the-art results from link prediction to molecule and graph generation, suggesting that accounting for the multiscale structure of graphs is a promising way to make graph neural networks even more powerful.

heading: A BASIC TENSOR OPERATIONS

text:
In order to build higher order equivariant networks, we revisit some basic tensor operations: tensor product (see Def. 10) and tensor contraction (see Def. 11). It can be shown that these tensor operations respect permutation equivariance. Based on them, we build our second order message passing networks. Definition 10. The tensor product of A ∈ Rn a with B ∈ Rn b yields a tensor C = A ⊗ B ∈ Rn a+b where Ci1,i2,..,ia+b = Ai1,i2,..,iaBia+1,ia+2,..,ia+b Definition 11. The contraction of A ∈ Rn a
along the pair of dimensions {x, y} (assuming x < y) yields a (a − 2)-th order tensor
Ci1,..,ix−1,j,ix+1,..,iy−1,j,iy+1,..,ia = ∑ ix,iy Ai1,..,ia
where we assume that ix and iy have been removed from amongst the indices of C. Using Einstein notation, this can be written more compactly as
C{i1,i2,..,ia}∖{ix,iy} = Ai1,i2,..,iaδ ix,iy
where δ is the Kronecker delta. In general, the contraction of A along dimensions {x1, .., xp} yields a tensor C = A↓x1,..,xp ∈ R na−p where
A↓x1,..,xp =∑ ix1 ∑ ix2 ...∑ ixp Ai1,i2,..,ia
or compactly as A↓x1,..,xp = Ai1,i2,..,iaδ ix1 ,ix2 ,..,ixp .

heading: B MARKOV RANDOM FIELDS

text:
Undirected graphical models have been widely applied in the domains spatial or relational data, such as image analysis and spatial statistics. In general, k-th order graph encoders encode an undirected graph G = (V,E) into a k-th order latent z ∈ Rn k ×dz , with learnable parameters θ, can be represented as a parameterized Markov Random Field (MRF) or Markov network. Based on the Hammersley-Clifford theorem (Murphy, 2012) (Koller & Friedman, 2009), a positive distribution p(z) > 0 satisfies the conditional independent properties of an undirected graph G iff p can be represented as a product of potential functions ψ, one per maximal clique, i.e.,
p(z∣θ) = 1
Z(θ) ∏ c∈C
ψc(zc∣θc) (10)
where C is the set of all the (maximal) cliques of G, and Z(θ) is the partition function to ensure the overall distribution sums to 1, and given by
Z(θ) =∑ z ∏ c∈C ψc(zc∣θc)
Eq. 10 can be further written down as
p(z∣θ)∝∏ v∈V ψv(zv ∣θ) ∏ (s,t)∈E ψst(zst∣θ) ⋅ ⋅ ⋅ ∏ c=(i1,..,ik)∈Ck ψc(zc∣θ)
where ψv , ψst, and ψc are the first order, second order and k-th order outputs of the encoder, corresponding to every vertex in V , every edge in E and every clique of size k in Ck, respectively. However, factorizing a graph into set of maximal cliques has an exponential time complexity, since the problem of determining if there is a clique of size k in a graph is known as an NP-complete problem. Thus, the factorization based on Hammersley-Clifford theorem is intractable. The second order encoder relaxes the restriction of maximal clique into edges, that is called as pairwise MRF:
p(z∣θ)∝∏ s∼t ψst(zs, zt)
Our second order encoder inherits Gaussian MRF introduced by (Rue & Held, 2005) as pairwise MRF of the following form
p(z∣θ)∝∏ s∼t ψst(zs, zt)∏ t ψt(zt)
where ψst(zs, zt) = exp ( − 12zsΛstzt) is the edge potential, and ψt(zt) = exp ( − 1 2 Λttz 2 t + ηtzt) is the vertex potental. The joint distribution can be written in the information form of a multivariate Gaussian in which
Λ = Σ−1
η = Λµ
p(z∣θ)∝ exp(ηTz − 1
2 zTΛz) (11)
Sampling z from p(z∣θ) in Eq. 11 is the same as sampling from the multivariate GaussianN (µ,Σ). To ensure end-to-end equivariance, we set the latent layer to be two tensors µ ∈ Rn×dz and Σ ∈ Rn×n×dz that corresponds to dz multivariate Gaussians, whose first index, and second index are first order and second order equivariant with permutations. Computation of Σ is trickier than µ, simply because Σ must be invertible to be a covariance matrix. Thus, our second order encoder produces tensor L as the second order activation, and set Σ = LLT . The reparameterization trick from Kingma & Welling (2014) is changed to
z = µ +L , ∼ N (0,1)

heading: C EQUIVARIANT LEARNABLE PRIOR

text:
The original VAE published by Kingma & Welling (2014) limits each covariance matrix Σ to be diagonal and the prior to beN (0,1). Our second order encoder removes the diagonal restriction on the covariance matrix. Furthermore, we allow the prior N (µ̂, Σ̂) to be learnable in which µ̂ and Σ̂ are parameters optimized by back propagation in a data driven manner. Importantly, Σ̂ cannot be learned directly due to the invertibility restriction. Instead, similarly to the second order encoder, a matrix L̂ is optimized, and the prior covariance matrix is constructed by setting Σ̂ = L̂L̂ T . The
Kullback-Leibler divergence between the two distributions N (µ,Σ) and N (µ̂, Σ̂) is as follows:
DKL(N (µ,Σ)∣∣N (µ̂, Σ̂)) = 1
2 (tr(Σ̂
−1 Σ) + (µ̂ −µ)T Σ̂ −1 (µ̂ −µ) − n + ln(
det Σ̂ detΣ )) (12)
Even though Σ is invertible, but gradient computation through the KL-divergence loss can be numerical instable because of Cholesky decomposition procedure in matrix inversion. Thus, we add neglectable noise = 10−4 to the diagonal of both covariance matrices.
Importantly, during training, the KL-divergence loss breaks the permutation equivariance. Suppose the set of vertices are permuted by a permutation matrix P σ for σ ∈ Sn. Since µ and Σ are the first order and second order equivariant outputs of the encoder, they are changed to P σµ and P σΣP T σ accordingly. But
DKL(N (µ,Σ)∣∣N (µ̂, Σ̂)) ≠ DKL(N (P σµ,P σΣP T σ )∣∣N (µ̂, Σ̂))
To address the equivariance issue, we want to solve the following convex optimization problem that is our new equivariant loss function
min σ∈Sn
DKL(N (P σµ,P σΣP T σ )∣∣N (µ̂, Σ̂)) (13)
However, solving the optimization based on Eq. 13 is computationally expensive. One solution is to solve the minimum-cost maximum-matching in a bipartite graph (Hungarian matching) with the cost matrix Cij = ∣∣µi − µ̂j ∣∣ by O(n4) algorithm published by Edmonds & Karp (1972), that can be still improved further into O(n3). The Hungarian matching preserves equiariance, but is still computationally expensive. In practice, instead of finding a optimal permutation, we apply a free-matching scheme to find an assignment matrix Π such that: Πij∗ = 1 if and only if j∗ = arg minj ∣∣µi − µ̂j ∣∣, for each i ∈ [1, n]. The free-matching scheme preserves equivariance and can be done efficiently in a simple O(n2) algorithm that is also suitable for GPU computation.

heading: D EXPERIMENTS

heading: D.1 LINK PREDICTION ON CITATION GRAPHS

text:
We demonstrate the ability of the MGVAE models to learn meaningful latent embeddings on a link prediction task on popular citation network datasets Cora and Citeseer (Sen et al., 2008). At training time, 15% of the citation links (edges) were removed while all node features are kept, the models are trained on an incomplete graph Laplacian constructed from the remaining 85% of the edges. From previously removed edges, we sample the same number of pairs of unconnected nodes (non-edges). We form the validation and test sets that contain 5% and 10% of edges with an equal amount of non-edges, respectively.
We compare our model MGVAE against popular methods in the field:
1. Spectral clustering (SC) (Tang & Liu, 2011) 2. Deep walks (DW) (Perozzi et al., 2014) 3. Variational graph autoencoder (VGAE) (Kipf & Welling, 2016)
on the ability to correctly classify edges and non-edges using two metrics: area under the ROC curve (AUC) and average precision (AP). Numerical results of SC and DW are experimental settings are taken from (Kipf & Welling, 2016). We reran the implementation of VGAE as in (Kipf & Welling, 2016).
For MGVAE, we initialize weights by Glorot initialization (Glorot & Bengio, 2010). We repeat the experiments with 5 different random seeds and calculate the average AUC and AP along with their standard deviations. The number of message passing layers ranges from 1 to 4. The size of latent representation is 128. The number of coarsening levels is L ∈ {3,7}. In the `-th coarsening level, we partition the graph G(`) into 2` (for L = 7) or 4` (for L = 3) clusters. We train for 2,048 epochs using Adam optimization (Kingma & Ba, 2015) with a starting learning rate of 0.01. Hyperparameters optimization (e.g. number of layers, dimension of the latent representation, etc.) is done on the validation set. MGVAE outperforms all other methods (see Table 3).
We propose our learning to cluster algorithm to achieve the balancedK-cut at every resolution level. Besides, we also implement two fixed clustering algorithms:
1. Spectral: It is similar to the one implemented in (Rustamov & Guibas, 2013). • First, we embed each node i ∈ V into Rnmax as (ξ1(i)/λ1(i), .., ξnmax(i)/λnmax(i)),
where {λn, ξn}nmaxn=0 are the eigen-pairs of the graph LaplacianL = D −1(D−A) where
Dii = ∑jAij . We assume that λ0 ≤ .. ≤ λnmax . In this case, nmax = 10. • At the `-th resolution level, we apply the K-Means clustering algorithm based on the
above node embedding to partition graph G(`). 2. K-Means:
• First, we apply PCA to compress the sparse word frequency vectors (of size 1,433 on Cora and 3,703 on Citeseer) associating with each node into 10 dimensions.
• We use the compressed node embedding for the K-Means clustering.
Tables 4 and 5 show that our learning to cluster algorithm returns a much more balanced cut on the highest resolution level comparing to both Spectral and K-Means clusterings. For instance, we
have L = 7 resolution levels and we partition the `-th resolution into K = 2` clusters. Thus, on the bottom levels, we have 128 clusters. If we distribute nodes into clusters uniformly, the expected number of nodes in a cluster is 21.15 and 25.99 on Cora (2,708 nodes) and Citeseer (3,327 nodes), respectively. We measure the minimum, maximum, standard deviation of the numbers of nodes in 128 clusters. Furthermore, we measure the Kullback–Leibler divergence between the distribution of nodes into clusters and the uniform distribution. Our learning to cluster algorithm achieves low KL losses of 0.02 and 0.01 on Cora and Citeseer, respectively.

heading: D.2 MOLECULAR GRAPH GENERATION

text:
In this case, MGVAE and MGN are implemented withL = 3 resolution levels, and the `-th resolution graph is partitioned into K = 2`−1 clusters. On each resolution level, the local encoders and local decoders are second-order Sn-equivariant networks with up to 4 equivariant layers. The number of channels for each node latent dz is set to 256. We apply two approaches for graph decoding:
1. All-at-once: MGVAE reconstructs all resolution adjacencies by equivariant decoder networks. Furthermore, we apply learnable equivariant prior as in Sec. C. Our second order encoders are interpreted as Markov Random Fields (see Sec. B). This approach preserves permutation equivariance. In addition, we implement a correcting process: the decoder network of the highest resolution level returns a probability for each edge, we sort these probabilities in a descending order and gradually add the edges in that order to satisfy all chemical constraints. Furthermore, we investigate the expressive power of the second order Sn-equivariant decoder by replacing it by a multilayer perceptron (MLP) decoder with 2 hidden layers of size 512 and sigmoid nonlinearity. We find that the higher order decoder outperforms the MLP decoder given the same encoding architecture. Table 6 shows the comparison between the two decoding models.
2. Autoregressive: This decoding process is constructed in an autoregressive manner similarly to (Liu et al., 2018). First, we sample each vertex latent z independently. We randomly select a starting node v0, then we apply Breath First Search (BFS) to determine a particular node ordering from the node v0, however that breaks the permutation equivariance. Then iteratively we add/sample new edge to the existing graph Gt at the t-th iteration (given a randomly selected node v0 as the start graph G0) until completion. We apply second-order MGN with gated recurrent architecture to produce the probability of edge (u, v) where one vertex u is in the existing graph Gt and the another one is outside; and also the probability of its label. Intuitively, the decoding process is a sequential classification.
We randomly select 10,000 training examples for QM9; and 1,000 (autoregressive) and 10,000 (allat-once) training examples for ZINC. It is important to note that our training sets are much smaller comparing to other methods. For all of our generation experiments, we only use graph features as the input for the encoder such as one-hot atomic types and bond types. Since ZINC molecules are larger then QM9 ones, it is more difficult to train with the second order Sn-equivariant decoders (e.g., the number of bond/non-bond predictions or the number of entries in the adjacency matrices are proportional to squared number of nodes). Therefore, we input several chemical/atomic features

heading: Dataset Method Validity Novelty Uniqueness

text:
computed from RDKit for the all-at-once MGVAE on ZINC (see Table 7). We concatenate all these features into a vector of size 24 for each atom.
We train our models with Adam optimization method (Kingma & Ba, 2015) with the initial learning rate of 10−3. Figs. 4 and 5 show some selected examples out of 5,000 generated molecules on QM9 by all-at-once MGVAE, while Fig. 6 shows the molecules generated by autoregressive MGN. Qualitatively, both the decoding approaches capture similar molecular substructures (bond structures). Fig. 9 shows an example of interpolation on the latent space on ZINC with the all-at-once MGVAE. Fig. 7 shows some generated molecules on ZINC by the all-at-once MGVAE. Fig. 8 and table 8 show some generated molecules by the autoregressive MGN on ZINC dataset with high Quantitative Estimate of Drug-Likeness (QED) computed by RDKit and their SMILES strings. On ZINC, the average QED score of the generated molecules is 0.45 with standard deviation 0.21. On QM9, the QED score is 0.44 ± 0.07.

heading: D.3 GENERAL GRAPH GENERATION BY MGVAE

text: Figs. 10 11 show some generated examples and training examples on the 2-community and ego datasets, respectively.

heading: D.4 UNSUPERVISED MOLECULAR PROPERTIES PREDICTION ON QM9

text: Density Function Theory (DFT) is the most successful and widely used approach of modern quantum chemistry to compute the electronic structure of matter, and to calculate many properties of molecular systems with high accuracy (Hohenberg & Kohn, 1964). However, DFT is computationally expensive (Gilmer et al., 2017a), that leads to the use of machine learning to estimate the properties of compounds from their chemical structure rather than computing them explicitly with DFT (Hy et al., 2018). To demonstrate that MGVAE can learn a useful molecular representations and capture important molecular structures in an unsupervised and variational autoencoding manner, we extract the highest resolution latents (at ` = L) and use them as the molecular representations for the downstream tasks of predicting DFT’s molecular properties on QM9 including 13 learning

heading: Target Unit Mean STD Description

text:
targets. For the training, we normalize all learning targets to have mean 0 and standard deviation 1. The name, physical unit, and statistics of these learning targets are detailed in Table 9.
The implementation of MGVAE is the same as detailed in Sec. D.2. MGVAE is trained to reconstruct the highest resolution (input) adjacency, its coarsening adjacencies and the node atomic features. In this case, we do not use any chemical features: the node atomic features are just one-hot atomic types. After MGVAE is converged, to obtain the Sn-invariant molecular representation, we average the node latents at the L-th level into a vector of size 256. Finally, we apply a simple Multilayer Perceptron with 2 hidden layers of size 512, sigmoid nonlinearity and a linear layer on top to predict the molecular properties based on the extracted molecular representation. We compare the results in Mean Average Error (MAE) in the corresponding physical units with four methods on the same split of training and testing from (Hy et al., 2018):
1. Support Vector Machine on optimal-assignment Weisfeiler-Lehman (WL) graph kernel (Shervashidze et al., 2011) (Kriege et al., 2016)
2. Neural Graph Fingerprint (NGF) (Duvenaud et al., 2015) 3. PATCHY-SAN (PSCN) (Niepert et al., 2016a) 4. Second order Sn-equivariant Covariant Compositional Networks (CCN 2D) (Kondor et al.,
2018) (Hy et al., 2018)
Our unsupervised results show that MGVAE is able to learn a universal molecular representation in an unsupervised manner and outperforms WL in 12, NGF in 10, PSCN in 8, and CCN 2D in 8 out of 13 learning targets, respectively (see Table 10). There are other recent methods in the field that use several chemical and geometric information but comparing to them would be unfair.

heading: D.5 SUPERVISED MOLECULAR PROPERTIES PREDICTION ON ZINC

text:
To further demonstrate the comprehensiveness of MGN, we apply our model in a supervised regression task to predict the solubility (LogP) on the ZINC dataset. We use the same split of 10K/1K/1K
for training/validation/testing as in (Dwivedi et al., 2020). The implementation of MGN is almost the same as detailed in Sec. D.2, except we include the latents of all resolution levels into the prediction. In particular, in each resolution level, we average all the node latents into a vector of size 256; then we concatentate all these vectors into a long vector of size 256×L and apply a linear layer for the regression task. The baseline results are taken from (Yang et al., 2020) including:
1. Multilayer Perceptron (MLP), 2. Graph Convolution Networks (GCN), 3. Graph Attention Networks (GAT) (Velikovi et al., 2018), 4. MoNet (Monti et al., 2017), 5. Disentangled Graph Convolutional Networks (DisenGCN) (Ma et al., 2019), 6. Factorizable Graph Convolutional Networks (FactorGCN) (Yang et al., 2020), 7. GatedGCNE (Dwivedi et al., 2020) that uses additional edge information.
Our supervised result shows that MGN outperforms the state-of-the-art models in the field with a margin of 20% (see Table 11).

heading: D.6 GRAPH-BASED IMAGE GENERATION BY MGVAE

text:
In this additional experiment, we apply MGVAE into the task of image generation. Instead of matrix representation, an image I ∈ RH×W is represented by a grid graph of H ⋅W nodes in which each node represents a pixel, each edge is between two neighboring pixels, and each node feature is the corresponding pixel’s color (e.g., R1 in gray scale, and R3 in RGB scale). Fig. 12 demonstrates an exmaple of graph representation for images. Since images have natural spatial clustering, instead of learning to cluster, we implement a fixed clustering procedure as follows:
• For the `-th resolution level, we divide the grid graph of size H(`) ×W (`) into clusters of size h × w that results into a grid graph of size H
(`) h × W (`) w
, supposingly h and w are divisible by H(`) and W (`), respectively. Each resolution is associated with an image I(`)
that is a zoomed out version of I(`+1). • The global encoder e(`) is implemented with 10 layers of message passing that operates
on the whole H(`) ×W (`) grid graph. We sum up all the node latents into a single latent vector Z(`) ∈ Rdz . The global decoder d(`) is implemented by the convolutional neural network architecture of the generator of DCGAN model (Radford et al., 2016) to map Z(`) into an approximated image Î(`). The Sn-invariant pooler p(`) is a network operating on each small h × w grid graph to produce the corresponding node feature for the next level `+ 1. MGVAE is trained to reconstruct all resolution images. Fig. 13 shows an example of reconstruction at each resolution on a test image of MNIST (after the network converged).
We evaluate our MGVAE architecture on the MNIST dataset (LeCun et al.) with 60,000 training examples and 10,000 testing examples. The original image size is 28×28. We pad zero pixels to get the image size of 25 × 25 (e.g., H(5) =W (5) = 32). Each cluster is a small grid graph of size 2 × 2 (e.g., h = w = 2). Accordingly, the image sizes for all resolutions are 32 × 32, 16 × 16, 8 × 8, etc. In this case, the whole network architecture is a 2-dimensional quadtree. The latent size dz is selected as 256. We train our model for 256 epochs by Adam optimizer (Kingma & Ba, 2015) with the initial learning rate 10−3. In the testing process, for the `-th resolution, we sample a random vector of size dz from prior N (0,1) and use the decoder d(`) to decode the corresponding image. We generate 10,000 examples for each resolution. We compute the Frechet Inception Distance (FID) proposed by (Heusel et al., 2017) between the testing set and the generated set as the metric to evaluate the quality
of our generated examples. We use the FID implementation from (Seitzer, 2020). We compare our MGVAE against variants of Generative Adversarial Networks (GANs) (Goodfellow et al., 2014) including DCGAN (Radford et al., 2016), VEEGAN (Srivastava et al., 2017), PacGAN (Lin et al., 2018), and PresGAN (Dieng et al., 2019). Table 12 shows our quantitative results in comparison with other competing generative models. The baseline results are taken from Prescribed Generative Adversarial Networks paper (Dieng et al., 2019). MGVAE outperforms all the baselines for the highest resolution generation. Figs. 14 and 15 show some generated examples of the 32 × 32 and 16 × 16 resolution, respectively.