{'fraction_non_alphanumeric': 0.03384269017053075, 'fraction_numerical': 0.01063422111515616, 'mean_word_length': 4.824103780164598, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Semi-Implicit Graph Variational Auto-Encoders

Semi-implicit graph variational auto-encoder (SIG-VAE) is proposed to expand the flexibility of variational graph auto-encoders (VGAE) to model graph data. SIG-VAE employs a hierarchical variational framework to enable neighboring node sharing for better generative modeling of graph dependency structure, together with a Bernoulli-Poisson link decoder. Not only does this hierarchical construction provide a more flexible generative graph model to better capture real-world graph properties, but also does SIG-VAE naturally lead to semi-implicit hierarchical variational inference that allows faithful modeling of implicit posteriors of given graph data, which may exhibit heavy tails, multiple modes, skewness, and rich dependency structures. Compared to VGAE, the derived graph latent representations by SIG-VAE are more interpretable, due to more expressive generative model and more faithful inference enabled by the flexible semi-implicit construction. Extensive experiments with a variety of graph data show that SIG-VAE significantly outperforms state-of-the-art methods on several different graph analytic tasks.


Introduction
Analyzing graph data is an important machine learning task with a wide variety of applications. Transportation networks, social networks, gene co-expression networks and recommendation systems are a few example datasets that can be modeled as graphs [16,14,15,19], where each node represents an agent (e.g., road intersection, person, and gene) and the edges manifest the interactions between the agents. The main challenge for analyzing graph datasets for link prediction, clustering, or node classification, is how to deploy graph structural information in the model. Graph representation learning aims to summarize the graph structural information by a feature vector in a low-dimensional latent space, which can be used in downstream analytic tasks.
While the vast majority of existing methods assume that each node is embedded to a deterministic point in the latent space [3,1,32,37,13,17,8], modeling uncertainty is of crucial importance in many applications, including physics and biology. For example, when link prediction in Knowledge Graphs is used for driving expensive pharmaceutical experiments, it would be beneficial to know what is the confidence level of a model in its prediction. To address this, variational graph auto-encoder (VGAE) [23] embeds each node to a random variable in the latent space. Despite its popularity, 1) the Gaussian assumption imposed on the variational distribution restricts its variational inference flexibility when the true posterior distribution given a graph clearly violates the Gaussian assumption; 2) the adopted inner-product decoder restricts its generative model flexibility. While recent study tries to address the first problem by changing the prior distribution but does not show much practical success [9], the latter one is not well-studied yet to the best of our knowledge.
Inspired by recently developed semi-implicit variational inference (SIVI) [47] and normalizing flow (NF) [34,21,30], which offer the interesting combination of flexible posterior distribution and effective optimization, we propose a hierarchical variational graph framework for node embedding of graph structured data, notably increasing the expressiveness of the posterior distribution for each node in the latent space. SIVI enriches mean-field variational inference with a flexible (implicit) mixing distribution. NF transforms a simple Gaussian random variable through a sequence of invertible differentiable functions with tractable Jacobians. While NF restricts the mixing distribution in the hierarchy to have an explicit probability density function, SIVI does not impose such a constraint. Both SIVI and NF can model complex posterior distributions, which will help when the underlying true embedded node distribution exhibits heavy tails and/or multiple modes. We further argue that the graph structure cannot be fully exploited by the posterior distribution from the trivial combination of SIVI and/or NF with VGAE, if not integrating graph neighborhood information. On the other hand, it does not address the flexibility of the generative model as stated as the second VGAE problem.
To address the aforementioned issues, instead of explicitly choosing the posterior distribution family in previous works [23,9], our hierarchical variational framework adopts a stochastic generative node embedding model that can learn implicit posteriors while maintaining simple optimization. Specifically, we innovate a semi-implicit hierarchical construction to model the posterior distribution to best fit both the graph topology and node attributes given graphs. With SIVI, even if the posterior is not tractable, its density can be evaluated with Monte Carlo estimation, enabling efficient model inference on top of highly enhanced model flexibility/expressive power. Our semi-implicit graph variational auto-encoder (SIG-VAE) can well model heavy tails, skewness, multimodality, and other characteristics that are exhibited by the posterior but failed to be captured by existing VGAEs. Furthermore, a Bernoulli-Poisson link function [49] is adopted in the decoder of SIG-VAE to increase the flexibility of the generative model and better capture graph properties of real-world networks that are often sparse. SIG-VAE facilitates end-to-end learning for various graph analytic tasks evaluated in our experiments. For link prediction, SIG-VAE consistently outperforms state-of-the-art methods by a large margin. It also achieves state-of-the-art performances for node classification without careful tuning. We further show that the new decoder is able to generate sparse random graphs whose statistics closely resemble those of real-world graph data. These results clearly demonstrate the great practical values of SIG-VAE.

Background
Variational graph auto-encoder (VGAE). Many node embedding methods derive deterministic latent representations [13,17,8]. By expanding the variational auto-encoder (VAE) notion to graphs, [23] propose to solve the following problem by embedding the nodes to Gaussian random vectors in the latent space.
Problem 1. Given a graph G = (V, E) with the adjacency matrix A and M -dimensional node attributes X ∈ R N ×M , find the probability distribution of the latent representation of nodes Z ∈ R N ×l , i.e., p(Z | X, A).
Finding the true posterior, p(Z | X, A), is often difficult and intractable. In [23], it is approximated by a Gaussian distribution, Here, µ i and σ i are l-dimensional mean and standard deviation vectors corresponding to node i, respectively. The parameters of q(Z | ψ), i.e., ψ = {ψ i } N i=1 , are modeled and learned using two graph convolutional neural networks (GCNs) [22]. More precisely, µ = GCN µ (X, A), log σ = GCN σ (X, A) and µ and σ are matrices of µ i 's and σ i 's, respectively. Given Z, the decoder in VGAE is a simple inner-product decoder as The parameters of the model are found by optimizing the well known evidence lower bound (ELBO) [20,5,6,41]: Despite promising results shown by VGAE, a well-known issue in variational inference is underestimating the variance of the posterior. The reason behind this is the mismatch between the representation power of the variational family to which q is restricted and the complexity of the true posterior, in addition to the use of KL divergence, which is asymmetric, to measure how different q is from the true posterior.
Semi-implicit variational inference (SIVI). To well characterize the posterior while maintaining simple optimization, semi-implicit variational inference (SIVI) has been proposed by [47], which is also related to the hierarchical variational inference [33] and auxiliary deep generative models [25]; see [47] for more details about their connections and differences. It has been shown that SIVI can capture complex posterior distributions like multimodal or skewed distributions, which can not be captured by a vanilla VI due to its restricted exponential family assumption over both the prior and posterior in the latent space. SIVI assumes that ψ, the parameters of the posterior, are drawn from an implicit distribution rather than being analytic. This hierarchical construction enables flexible mixture modeling and allows to have more complex posteriors while maintaining simple optimization for model inference. More specifically, Z ∼ q(Z | ψ) and ψ ∼ q φ (ψ) with φ denoting the distribution parameters to be inferred. Marginalizing ψ out leads to the random variables Z drawn from a distribution family H indexed by variational parameters φ, expressed as The importance of semi-implicit formulation is that while the original posterior q(Z | ψ) is explicit and analytic, the marginal distribution, h φ (Z) is often implicit. Note that, if q φ equals a delta function, then h φ is an explicit distribution. Unlike regular variational inference that assumes independent latent dimensions, semi-implicit does not impose such a constraint. This enables the semi-implicit variational distributions to model very complex multivariate distributions.
Normalizing flow (NF) [30] also enriches the posterior distribution families. Compared to SIVI, NF imposes explicit density functions for the mixing distributions in the hierarchy while SIVI only requires q φ to be reparameterizable. This makes SIVI more flexible, especially when using it for graph analytics as explained in the next section, since the SIVI posterior can be generated by transforming random noise using any flexible function, for example a neural network.
3 Semi-implicit graph variational auto-encoder (SIG-VAE) Figure 1: SIG-VAE diffuses the distributions of the neighboring nodes, which is more informative than sharing deterministic features, to infer each node's latent distribution.
To address Problem 1 while well characterizing the posterior with modeling flexibility in the VGAE framework, we propose SIG-VAE that takes the semiimplicit variational distribution following the hierarchical formulation in SIVI: by introducing the implicit prior distribution parametrized by ψ, which can be sampled from the reparametrizable q φ (ψ | X, A). Such a hierarchical semi-implicit construct not only leads to flexible mixture modeling of the posterior but also enables efficient model inference, for example, with φ being parameterized by deep neural networks. To have tractable posterior inference, we construct SIG-VAE using a hierarchy of multiple stochastic layers. Specifically, the first stochastic layer q(Z | X, A) is reparameterizable and has an analytic probability density function. The layers added after are reparameterizable and computationally efficient to sample from. More specifically, we adopt a hierarchical encoder in SIG-VAE that injects random noise at L different stochastic layers: , where GNN is any type of existing graph neural networks, such as graph convolutional neural network (GCN) [22], GCN with Chebyshev filters [10], GraphSAGE [17], jumping knowledge (JK) networks [45], and graph isomorphism network (GIN) [44]. Note that in the equations above µ and Σ are random variables and thus q(Z | X, A) is not necessarily Gaussian after marginalization; u is N -dimensional random noise drawn from a distribution q u ; and q u is chosen such that the samples drawn from it are the same type as X, for example if X is categorical, Bernoulli is a good choice for q u . By concatenating the random noise and node attributes, the output of GNNs are random variables rather than deterministic vectors. Their expressive power is inherited in SIG-VAE to go beyond Gaussian, exponential family, or von Mises-Fisher [9] posterior distributions for the derived latent representations.
In SIG-VAE, when inferring each node's latent posterior, we incorporate the distributions of the neighboring nodes, better capturing graph dependency structure than sharing deterministic features from GNNs. More specifically, the input to our model at stochastic layer u is CONCAT(X, u ) so that the outputs of the subsequent stochastic layers give mixing distributions by integrating information from neighboring nodes (Fig. 1). Through experiments, we show that this neighborhood sharing is key for SIG-VAE to achieve superior graph analysis performance.
We further argue that increasing the flexibility of variational inference is not enough to better model real-world graph data as the optimal solution of the generative model does not change. In SIG-VAE, the Bernoulli-Poisson link [49] is adopted for the decoder to further increase the expressiveness of the generative model. Potential extensions with other decoders can be integrated with SIG-VAE if needed. Let where R ∈ R l×l + is a diagonal matrix with diagonal elements r k .

Inference
To derive the ELBO for model inference in SIG-VAE, we must take into account the fact that ψ has to be drawn from a distribution. Hence, the ELBO moves beyond the simple VGAE as where h φ is defined in (1). The marginal probability density function h φ (Z|X, A) is often intractable, so the Monte Carlo estimation of the ELBO, L, is prohibited. To address this issue, SIVI derives a lower bound for the ELBO and optimizes this lower bound instead of optimizing the ELBO itself, which is tractable and asymptotically equals to the ELBO. SIG-VAE requires q(Z | ψ) to be explicit, and also requires it to either be reparameterizable or the ELBO under q(Z | ψ) to be analytic, while q φ (ψ | X, A) is required to be reparameterizable but not necessarily explicit. This captures the idea that combining an explicit q(Z | ψ) with an implicit q φ (ψ | X, A) is as powerful as needed, but makes the computation tractable. Following [47], we can derive a lower bound for the ELBO as follows This can be proved based on the first theorem in [47], which shows Unlike L, a Monte Carlo estimation of L only requires q φ (Z | ψ) to have an analytic density functions and q φ (ψ | X, A) to be convenient to sample from.
Directly optimizing L without early stopping could lead to a point mass density as q φ (ψ | X, A). This degenerates SIG-VAE to vanilla VGAE. To avoid degeneracy, a regularizer term can be added to the L.
We define a regularized lower bound as It has been proved by [27] that L K is a monotonic lower bound of the ELBO, satisfying L K ≤ L K+1 ≤ L. Therefore, setting K to zero means that L 0 = L, and as K goes to infinity L converges to the exact ELBO, i.e., lim K→∞ L K = L. representations (middle and right) are heat maps in R 3 . We expect that the embedding of the Swiss roll graph with inner-product decoder to be a curved plane in R 3 , which is clearly captured better by SIG-VAE. and VGAE (red). SIG-VAE clearly infers more complex distributions that can be multi-modal, skewed, and with sharp and steep changes. This helps SIG-VAE to better represent the nodes in the latent space.

Other VGAE extensions
It is possible to enable VGAE model flexibility by other variational inference methods, for example using NF. However, as NF requires deterministic transforms whose Jacobians shall be easy to compute, it cannot be combined with VGAEs in a similar way as SIG-VAE does. We indeed have constructed another non-Gaussian VGAE, i.e. NF-based variational graph auto-encoder (NF-VGAE) as follows: where the posterior distribution q K (Z (K) |A, X) is obtained by successively transforming a Gaussian random variable Z (0) with distribution q 0 through a chain of K transformations f k . NF-VGAE is a two-step inference method that 1) starts with Gaussian random variables and then 2) transforms them through a series of invertible mappings. Note that (7) is different from the SIG-VAE construction (3). In NF-VGAE, the GNN output layers are deterministic without neighborhood distribution sharing due to the deterministic nature of the initial density parameters in q 0 . The flexibility of SIG-VAE directly working on the stochastic distribution parameters in (3)(4) allows neighborhood sharing to achieve better performance in graph analytic tasks. SIVI can be integrated with VGAE in a similar way as in NF-VGAE, without neighborhood sharing: the features from multiple layers of GNNs can be aggregated and then transformed via multiple fully connected layers after being concatenated by random noise to derive the posterior distribution for each node separately. We call this Naive SIG-VAE. More specifically, Naive SIG-VAE injects random noise at C different stochastic fully connected layers for each node independently: where T t , g µ , and g σ are all deterministic neural networks and i is the node index. Please note that the first layer of SIVI can be integrated with NF rather than simple Gaussian. We leave that for future study. While these two models are able to approximate more flexible and complex posterior, such trivial combinations may fail to fully exploit graph dependency structure.

Experiments
We test the performances of SIG-VAE on different graph analytic tasks: 1) interpretability of SIG-VAE compared to VGAE, 2) link prediction in various real-world graph datasets including graphs with node attributes and without node attributes, 3) graph generation, 4) node classification in the citation graphs with labels. In all of the experiments, GCN [22] is adopted for all the GNN modules in SIG-VAE, Naive SIG-VAE, and NF-VGAE, implemented in Tensorflow [26]. The PyGSP package [11] is used to generate synthetic graphs. Implementation details for all the experiments, together with graph data statistics, can be found in the supplement.  We first demonstrate the expressiveness of SIG-VAE by illustrating the approximated variational distributions of node latent representations. We show that SIG-VAE captures the graph structure better and has a more interpretable embedding than VGAE on a generated Swiss roll graph with 200 nodes and 1244 edges (Fig. 2). In order to provide a fair comparison, both models share an identical implementation with the innerproduct decoder and same number of parameters. We simply consider the identity matrix I N as node attributes and choose the latent space dimension to be three in this experiment. This graph has a simple plane like structure. As the inner-product decoder assumes that the information is embedded in the angle between latent vectors, we expect that the node embedding to map nodes of the Swiss roll graph into a curve in the latent space. As we can see in Fig. 2, SIG-VAE derives a clearly more interpretable planar latent structure than VGAE. We also show the posterior distributions of five randomly selected nodes from the graph in Fig. 3. As we can see, SIG-VAE is capable of inferring complex distributions. The inferred distributions can be multi-modal, skewed, non-symmetric, and with sharp and steep changes. These complex distributions help the model to get a more realistic embedding capturing the intrinsic graph structure. To explain why multi-modality may arise, we used Asynchronous Fluid [31] to visualize the Swiss Roll graph by highlighting detected communities with different colors in Figure 4. The three red (two orange) nodes are the nodes with multi-modal (skewed) distributions. These nodes with multi-modal posteriors reside between different communities. The supplement contains additional results and discussions with a torus graph, with similar observations.

Accurate link prediction
We further conduct extensive experiments for link prediction with various real-world graph datasets. Our results show that SIG-VAE significantly outperforms well-known baselines and state-of-the-art methods in all benchmark datasets. We consider two types of datasets, i.e., datasets with node attributes and datasets without attributes. We preprocess and split the datasets as done in [23] with validation and test sets containing 5% and 10% of network links, respectively. We learn the model parameters for 3500 epochs with the learning rate 0.0005 and the validation set used for early stopping. The latent space dimension is set to 16. The hyperparameters of SIG-VAE, Naive SIG-VAE, and NF-VGAE are the same for all the datasets. For fair comparison, all methods have the similar number of parameters as the default VGAE. The supplement contains further implementation details. We measure the performance by average precision (AP) and area under the ROC curve (AUC) based on 10 runs on a test set of previously removed links in these graphs.
With node attributes. We consider three graph datasets with node attribbutes-Citeseer, Cora, and Pubmed [35]. The number of node attributes for these dataset are 3703, 1433 and 500 respectively. Other statistics of the datasets are summarized in the supplement Table 1. We compare the results of SIG-VAE, Naive SIG-VAE, and NF-VGAE with six state-of-the-art methods, including spectral clustering (SC), DeepWalk (DW) [32] , GAE [23], VGAE [23], S-VGAE [9] and SEAL [48]. The inner-product decoder is also adopted in SIG-VAE to clearly demonstrate the advantages of the semi-implicit hierarchical variational distribution for the encoder.
We use the same hyperparameters for the competing methods as stated in [48,23,9]. As we can see in Table 1, SIG-VAE shows significant improvement in terms of both AUC and AP over state-of-the-art methods. Note the standard deviation of SIG-VAE is also smaller compared to other methods, indicating stable semi-implicit variational inference. Compared to the baseline VGAE, more flexible posterior in three proposed methods SIGVAE (with both inner-product and Bernoulli-Poisson link decoders), Naive SIG-VAE, and NF-VGAE can clearly improve the link prediction accuracy. This suggests that the Gaussian assumption does not hold for these graph structured data. The performance improvement of SIG-VAE with inner-product decoder (IP) over Naive SIG-VAE and NF-VGAE clearly demonstrates the advantages of neighboring node sharing, especially in the smaller graphs. Even for the large graph Pubmed, on which VGAE performs similar to S-VGAE, our SIG-VAE still achieves the highest link prediction accuracy, showing the importance of all modeling components in the proposed method including non-Gaussian posterior, using neighborhood distribution, and the sparse Bernoulli-Poisson link decoder. Without node attributes. We further consider five graph datasets without node attributes-USAir, NS [28], Router [36], Power [42] and Yeast [40]. The data statistics are summarized in the supplement Table 1. We compare the performance of our models with seven competing state-of-the-art methods including matrix factorization (MF), stochastic block model (SBM) [2], node2vec (N2V) [13], LINE [37], spectral clustering (SC), VGAE [23], S-VGAE [9], and SEAL [48].
For baseline methods, we use the same hyperparameters as stated in Zhang et al. [48]. For datasets without node attributes, we use a two-stage learning process for SIG-VAE. First, the embedding of each node is learned in the 128-dimensional latent space while injecting 5-dimensional Bernoulli noise to the system. Then the learned embedding is taken as node features for the second stage to learn 16 dimensional embedding while injecting 64-dimensional noise to SIG-VAE. We follow the same procedure for Naive SIG-VAE and NF-VGAE.
As we can see in Table 5, SIG-VAE again shows the consistent superior performance compared to the competing methods, especially over the baseline VGAE, in both AUC and AP. It is interesting to note that, while the proposed Berhoulli-Poisson decoder works well for sparser graphs, especially NS and Router datasets, SIG-VAE with inner-product decoder shows superior performance for the USAir graph which is much denser. Compared to the baseline VGAE, both Naive SIG-VAE and NF-VGAE improve the results with a large margin in both AUC and AP, showing the benefits of more flexible posterior. Comparing SIG-VAE with two other flexible inference methods shows not only SIG-VAE is not restricted to the Gaussian assumption, which is not a good fit for link prediction with the inner-product decoder [9], but also it is able to model flexible posterior considering graph topology. The results for the link prediction of the Power graph clearly magnifies this fact as SIG-VAE improves the accuracy by 34% compared to VGAE. The supplement contains the results with standard deviation values over different runs, showing the stability again.
Ablation studies have also been run to evaluate SIG-VAE with inner-product decoder in link prediction for citation graphs without using node attributes. The [AUC, AP] are [91.14, 90.99] for Cora and [88.72, 88.24] for Citeseer, lower than the values from SIG-VAE with attributes in Table 1 but are still competitive against existing methods (even with node attributes), showing the ability of SIG-VAE of utilizing graph structure. While some of the methods, like SEAL, work well for graphs without node attributes and some of others, like VGAE, get good performance for graphs with node attributes, SIG-VAE consistently achieves superior performance in both types of datasets. This is due to the fact that SIG-VAE can learn implicit distributions for nodes, which are very powerful in capturing graph structure even without any node attributes.

Graph generation
To further demonstrate the flexibility of SIG-VAE as a generative model, we have used the inferred embedding representations to generate new graphs. For example, SIG-VAE infers network parameters for Cora whose density and average clustering coefficients are 0.00143 and 0.24, respectively. Using the inferred posterior and learned decoder, a new graph is generated with corresponding r k to see if its graph statistics are close to the original ones. Please note that we have shrunk inferred r k 's smaller than 0.01 to 0. The density and average clustering coefficients of this generated graph based on SIG-VAE are 0.00147 and 0.25, respectively, which are very close to the original graph. We also generate new graphs based on SIG-VAE with the inner-product decoder and VGAE. The density and average clustering coefficients of the generated graphs based on SIG-VAE (IP) and VGAE are same, i.e. 0.1178 and 0.49, respectively, showing the inner-product decoder may not be a good choice for sparse graphs. The supplement includes more examples.

Node classification & graph clustering
We also have applied SIG-VAE for node classification on citation graphs with labels by modifying the loss function to include graph reconstruction and semi-supervised classification terms. Results are summarized in Table 3. Our model exhibits strong generalization properties, highlighted by its competitive performance compared to the state-of-the-art methods, despite not being trained specifically for this task (note GAT uses 64 hidden features, while the other methods including SIG-VAE use 16). To show the robustness of SIG-VAE to missing edges, we randomly removed 10, 20, 50 and 70 (%) edges while keeping node attributes. The mean accuracy of 10 run for Cora (2 layers [32,16]) are 79.5, 78.7, 75.3 and 60.6, respectively. The supplement contains additional results and discussion for graph clustering, again without specific model tuning. Combining the advantages of semi-implicit hierarchical variational distribution and VGAE with a Bernoulli-Poisson link decoder, SIG-VAE is developed to enrich the representation power of the posterior distribution of node embedding given graphs so that both the graph structural and node attribute information can be best captured in the latent space. By providing a surrogate evidence lower bound that is asymptotically exact, the optimization problem for SIG-VAE model inference is amenable via stochastic gradient descent, without compromising the flexbility of its variational distribution. Our experiments with different graph datasets have shown the promising capability of SIG-VAE in a range of graph analysis applications with interpretable latent representations, thanks to the hierarchical construction that diffuses the distributions of neighborhood nodes in given graphs.

Supplementary Material
In this supplement, we first provide the detailed review of the related literature as well as the connection to our proposed work. Dataset statistics, network setups, and implementation details of performance evaluation experiments for different graph analytic tasks are then presented with richer experimental results in addition to the ones discussed in the main text.

A Related works
Variational graph auto-encoders (VGAE), proposed by [23], embed each node to a random variable in the latent space. VGAE, by extending the use of VAEs to graph structured data, is shown to be capable of learning interpretable latent representations for undirected graphs and getting competitive results in the link prediction task. However, the Gaussian assumption imposed on the variational distribution restricts the model flexibility when the true posterior distribution given a graph clearly violates the assumption. It also suffers from underestimating the variance of the posterior, which is a well-known issue of vanilla VAEs.
To better model graph data using variational distributions in VGAEs, [9] proposes the hyperspherical VGAE (S-VGAE), in which, instead of the Gaussian assumption for the posterior, the von Mises-Fisher distribution has been deployed. This assumption is not well-suited for all classes of graphs. For example, it has been proven that graphs with hierarchical tree-like structure have hyperbolic latent structures [29] which clearly cannot be represented well in a hyperspherical space. While S-VGAE outperforms vanilla VGAE in some graphs including Cora and Citeseer in terms of link prediction accuracy, its performance will be degraded for more complex graphs such as Pubmed. On the other hand, changing the prior is not going to change the flexibility and optimal solution of the generative model, but will affect the tightness of the ELBO and hence how well the generative model parameters can be inferred. This shows the necessity to develop a variational graph auto-encoders that not only is capable of inferring more flexible posteriors to represent a broader range of graphs, but also is able to have more flexible decoder especially for the real-world sparse graphs.
In this paper, we propose to develop a hierarchical variational model to increase the expressiveness of the posterior distribution for each node in the latent space. While Naive SIG-VAE and NF-VGAE can be used as a variational node embedding to effectively expand the variational posterior distribution family, SIG-VAE allows flexible implicit posteriors as well as exploitation of the neighbor dependency while maintaining simple optimization. We have further adopted a Bernoulli-Poisson link decoder to improve the flexibility of the generative model which has not been addressed in the previous studies.

B Node embedding
Node embedding is to represent each node in a graph by a low-dimensional vector in a latent space. The geometric relations of vectors in the latent space reflect the probability of two corresponding nodes interacting with each other in the graph [18]. A good node embedding preserves node connectivity in graph as well as local neighborhood structures. More formally, node embedding can be formulated as follows.
Node embedding. Given a graph G = (V, E) where V is the set of nodes and E the set of edges, with the adjacency matrix A, X ∈ R N ×M denoting M -dimensional node attributes for N = |V| nodes, and a function s G : V × V → R measuring node similarity, find an encoder function, ENC : R N ×N + × R N ×M → R l , a decoder function, DEC : R l × R l → R + , and a latent representation of nodes Z ∈ R N ×l such that Z = ENC(A, X), where z i corresponds to the embedding representation of node v i ∈ V. Optimal parameters of ENC and DEC functions can be derived by finding the solutions to the following optimization problem where loss is a user-specified loss function based on the ultimate objective of network analysis.
Different node embedding methods vary in the choice of the loss function, s G , ENC, DEC and the optimization algorithm. For example, in graph factorization (GF) method [1], s G is defined based on the adjacency matrix, i.e., s G (v i , v j ) = A i,j ; loss is the mean squared error; and the inner-product decoder is adopted, i.e., DEC(z i , z j ) = z T i z j .

C Variational inference with normalizing flows
To increase the expressive power of a probabilistic model, a simple but powerful idea is to transform the corresponding random variables with complex deterministic and/or stochastic mappings. To construct flexible, arbitrarily complex and scalable approximate posterior distributions, normalizing flow (NF) transforms a simple random variable through a sequence of invertible differentiable functions with tractable Jacobians. More specifically, NF uses an invertible, smooth mapping f : R d → R d to transform a random variable z with distribution q(z) to the resulting random variable z = f (z) with the distribution: One may apply a chain of K transformations f k to obtain the density q K (z) from a random variable z 0 with distribution q 0 as: While normalizing flow helps to improve the model flexibility of the corresponding variational posterior, it requires the mapping to be deterministic and invertible, and the mixing distribution in the hierarchy to have an explicit density function. Removing these restrictions, there have been several recent attempts to define highly flexible variational posterior with implicit models. While an implicit variational distribution can be made highly flexible, it becomes necessary in each iteration to address the problem of density ratio estimation, which is often transformed into a problem related to learning generative adversarial networks [12]. SIVI addresses this issue by using an analytic conditional variational distribution which is not required to be reparameterizable. Table 4 provides the detailed statistics of the graph datasets used in our experiments.

E Experimental setups and hyperparameter tuning
Interpretable latent representations experiments. In these experiments, the code provided by [23] is used to derive the embedding for VGAE. The size of the first hidden layer of VGAE is 256 and the size of the output layer is 3. For SIG-VAE, two stochastic layers with sizes equal to [32,32] and an additional GCN layer of size 16 are used to model the µ. The dimension of injected standard Gaussian noises [ 1 , 2 ] are [32,32]. Covariance matrix Σ is deterministic and is inferred through two layers of GCNs with sizes equal to [32,16]. To remove the effect of decoder, we consider the inner-product decoder for this set of experiments. Link prediction with node attributes. For SIG-VAE, we use a stochastic layer with size equal to 32 and an additional GCN layer of size 16 is used to model µ. The dimension of injected Bernoulli noise for the stochastic layer is 64. For Naive SIG-VAE, we use two GCN layers with sizes equal to [32,16] followed by a fully connected layers with size 16 to infer µ. We inject 64-dimensional Bernoulli noise to the fully connected layer. We implement NF-VGAE by extending VGAE (two GCN layers with sizes equal to [32,16]) with invertible linear-time transformations of length 4 to keep its number of parameters close to the competing methods. We learn the model parameters for 3500 epochs with the learning rate 0.0005 and the validation set used for early stopping.
Link prediction without node attributes. For SIG-VAE, we use a stochastic layer with size equal to 32 and an additional GCN layer of size 16 is used to model µ. The dimensions of injected Bernoulli noise is 32. For Naive SIG-VAE, we use two GCN layers with sizes equal to [32,16] followed by a fully connected layer with sizes 16 to infer µ. We inject 32-dimensional Bernoulli noise to the fully connected layers. We learn the all model parameters for 2500 epochs with the learning rate 0.0005 and use the validation set for the early stopping. We use a two-stage learning process for SIG-VAE, Naive SIG-VAE, and NF-VGAE. First, the embedding of each node is learned in the 128-dimensional latent space while injecting 5-dimensional Bernoulli noise to the system in the case of SIG-VAE and Naive SIG-VGAE. Then we use the learned embedding as node features for the second stage to learn 16 dimensional embedding while injecting more noise to SIG-VAE. We follow the same procedure for Naive SIG-VAE too.
Graph generation. We have not specifically tuned the model but directly adopt the implementation setups for link prediction with and without node attributes.
Node classification and graph clustering. We use two GCN layers with sizes equal to [32,16] followed by a fully connected layer with sizes 16 to infer µ. We inject 64-dimensional Bernoulli noise to the GCN layers. Learning rate is set to be 0.0005.
Analysis of the complexity. For the analysis of the real-world graph dataset Cora on a single GeForce GTX 1080 GPU node, it took 24.5, 11.7 , and 9.5 seconds for SIG-VAE, NF-VGAE, and VGAE methods with 100 epochs, respectively. For the analysis of the small real-world graph dataset NS on a same GPU node, it took 7.23, 7.84, and 7.09 seconds for SIG-VAE, NF-VGAE, and VGAE methods with 100 epochs, respectively.

F.1 Interpretable latent representations
In addition to the results of the Swiss roll graph in the paper, we also compare the latent representations of SIG-VAE and VGAE for a torus graph with 256 nodes connected by 512 edges as illustrated in Figure 5. We consider the coordinates of each node in R 3 as node attributes for both methods in this experiment. We expect that the embedding of nodes to be symmetric since the graph itself is symmetric. We know that the inner-product decoder tries to embed a ring graph to a circle in space. Also, connected nodes should be in the same angle. Thus, the embedding of connected circles as in torus in R 3 should be some lines coming out of center while their altitude is changing periodically. As we can see in Figure 5, SIG-VAE demonstrates a better latent representation than VGAE. To gain more insights about the posterior distributions, we show the distributions inferred by SIG-VAE and VGAE for three nodes in Figure 6. The inferred distributions are indeed skewed and multi-modal, very different from Gaussian. Being able to capture complex non-Gaussian distributions helps the model to represent the graph structure in a more meaningful way. We expect that the embedding of the torus graph with the inner-product decoder to be multiple lines coming out of the center in R 3 , which is clearly better captured by SIG-VAE. : Latent representation distributions of three nodes in the torus graph using SIG-VAE (blue) and VGAE (red). SIG-VAE clearly infers more complex distributions that are multi-modal or skewed. This helps SIG-VAE to better represent the nodes in the latent space.

F.2 Link prediction
More complete link prediction results with the standard deviation values from different runs are presented here. As we can see in Tables 5 and 6, SIG-VAE shows the consistent superior performance compared to the competing methods, especially over the baseline VGAE, in terms of both AUC and AP. It is interesting to note that, while the proposed sparse decoder works well for the sparser graphs, especially NS and Router sparse datasets, SIG-VAE with the inner-product decoder shows superior performance for the USAir graph which is much denser. Compared to the baseline VGAE, both Naive SIG-VAE and NF-VGAE improve the results with a large margin in terms of both AUC and AP, showing the benefits of more flexible variational posterior. Comparing SIG-VAE with two other flexible inference methods shows that not only SIG-VAE is not restricted to the Gaussian assumption, which is not a good fit for link prediction with the inner-product  decoder [9], but also it is able to model flexible posterior considering graph topology. The results for the link prediction of the Power graph clearly magnifies this fact as SIG-VAE improves the accuracy by 34% compared to VGAE.

F.3 Graph generation
In addition to the results of Cora dataset in the paper, we also used the inferred embedding representations of different graph dataset with and without node attributes to generate new graphs. Results are summarized in Table 7. The SIG-VAE results are much closer to the real-world graph in terms of both graph density and average clustering for very sparse graphs. For the USAir dataset, which is much dense compare to othe graphs, the average clustering coefficient of SIG-VAE with inner-product decoder is closer to the read-world graph. This can be describe the better link prediction results of SIG-VAE for USAir dataset. On the other hand, the generated graph by SIG-VAE with the Bernoulli-Poisson link decoder is much sparser as its density is very closer to the read-world graph. This shows the benefit of the proposed decoder to improve the flexibility of the generative model.

F.4 Graph clustering
SIG-VAE can be applied in the other application including graph clustering. We first tried SIG-VAE for getting low-dimential feature space and then apply Gaussian mixture clustering (GMM) on citation graphs with labels including Cora and Citeseer and compare its results with VGAE. We consider same number of parameters and GCN layer for both model. Results are summarized in Table 8. We report the normalized mutual information (NMI) and unsupervised clustering accuracy (ACC) of 10 runs. The decoders for both methods are inner-product decoder.