Recurrent Dirichlet Belief Networks for interpretable Dynamic Relational Data Modelling

1 Introduction

Dynamic data is a common feature in many real-world applications, including relational data analysis [Mucha et al., 2010; Phan and Airoldi, 2015; Yang and Koeppl, 2018] for learning time-varying node interactions, and text modelling [Guo et al., 2018; Schein et al., 2019] for exploring topic evolution. Modelling dynamic data has become a vibrant research topic, with popular techniques ranging from non-Bayesian methods, such as Collaborative Filtering with Temporal Dynamics (SVD++) [Koren, 2009], to Bayesian deep probabilistic frameworks such as Deep Poisson-Gamma Dynamical Systems (DPGDS) [Guo et al., 2018]. The main advantage of Bayesian deep probabilistic frameworks is the flexible model design and the strong modelling performance. However, most of these frameworks are static so that they cannot account for the evolution of relationships over time. It would be highly beneficial if the frameworks can be extended to the dynamic setting to enjoy the modelling advantages.

∗Corresponding Author

The Dirichlet Belief Network (DirBN) [Zhao et al., 2018] has been proposed recently as a promising deep probabilistic framework for learning interpretable deep latent structures. To date, the DirBN has mainly been used in two applications: (1) topic structure learning [Zhao et al., 2018], where latent representations are used to model the word distribution for topics; and (2) relational models [Fan et al., 2019a], where latent representations model the nodes' membership distribution over communities. By constructing a deep architecture for latent distributions, the DirBN can model high-order dependence between topic-word distributions (in topic models) and nodes' membership distributions (in relational models).

In this work, we propose a Recurrent Dirichlet Belief Network (Recurrent-DBN) to explore the complex latent structures in dynamic relational data. In addition to constructing an interpretable deep architecture for the data within individual time steps, we also study the temporal dependence in the dynamic relational data through (layer-to-layer) connections crossing consecutive time steps. Consequently, our Recurrent-DBN can describe long-term temporal dependence (i.e., the dependence between the current variables and those in the previous several time steps), improving over the one-order Markov structures that usually describe the dependence between the current variables and those in the previous one time step only.

For model inference, we further develop an efficient Gibbs sampling algorithm. Besides upward propagating latent counts as done by DirBN, we also introduce a backward step to propagate the counts from the current time step to the previous time steps. Our experiments on real-world dynamic relational data show significant advantages of the Recurrent-DBN over the state-of-the-art models in tasks of interpretable latent structure discovery and link prediction. Similar to DirBN that can be considered as a self-contained module [Zhao et al., 2018], our Recurrent-DBN could be flexibly adapted to account for dynamic data other than evolving relational data, such as time-varying counts and dynamic drifting text data. We summarise this paper's main merits as follows:

• Model Recurrent structures are designed to model long term temporal dependence. Also, interpretable and organised latent structures are well explored;
• Inference An efficient Gibbs sampling method is devised that first upward-backward propagates latent counts and then downward-forward samples variable;

Results Significantly improved model performance in realworld dynamic relational models compared to the stateof-the-art, including better link prediction performance and enhanced interpretable latent structure visualisation.

2 Background Information of DirBN

We first give a brief review of the DirBN model. In general, the DirBN constructs a multi-stochastic layered architecture to represent interpretable latent distributions for objects. We describe it within the relational data setting for illustrative purposes. Given a binary observed linkage matrix R ∈ {0, 1}N×N for N nodes, where Rij denotes whether node i has a relation to node j, the DirBN constructs an L-layer and K-length community membership distribution πi ={π (l) i }Ll=1 for each node i. The generative process for the membership distributions {π (l) i }Ll=1 , as well as the observed matrix R, can be briefly described as:

where α1×K is a concentration parameter generating the membership distribution in the 1st layer, β(l-1) i′ i represents the information propagation coefficient from node i′ to node i in the (l - 1)th layer, c and d are the hyper-parameters generating these propagation coefficients, Xi is the latent count information for node i and M is the sum of these counts, and f(Xi,Xj) represents the probabilistic function mapping a pair of membership distributions to a linkage probability. A larger value of β(l-1) i′ i indicates higher influence of π (l-1) i′ on the generation of π (l) i . Therefore, β(l-1) i′ i is set to 0 if node i′ is not connected to node i in the observed data R.

It is difficult to directly implement efficient Gibbs sampling for the DirBN because the prior and posterior distributions of the membership distributions π (l) i are not conjugate. To address this issue, a strategy of first upward propagating latent counts and then downward sampling variables has been developed in [Zhao et al., 2018]. Given the count information Xi for node i, the DirBN upward propagates Xi to all the nodes in the (L - 1)th layer through a Chinese Restaurant Table (CRT) distribution. Each node in the (L - 1)th layer collects these propagated counts and uses their sum as its latent count X (L-1) i in the (L - 1)th layer. This procedure is repeated until the counts have been assigned to all layers.

Thus, conjugate constructions can be created for each variable and thereby used to construct efficient Gibbs samplers.

3 Recurrent-DBN for Dynamic Relational Data Modeling

To handle dynamic relational data, we attach an index t to variables to denote the corresponding time step. Thus, the observed dynamic relational data can be described as R ∈ {0, 1}N×N×T for N nodes at T time steps, where Rij,t denotes whether node i has relation to node j at the tth time step. Each matrix {R-,t}t ∈ {0, 1}N×N can be either asymmetric (directional) or symmetric (non-directional) and we do not consider self-linkages {Rii,t}i,t.

3.1 Recurrent-DBN for Latent Structure Generation

In the Recurrent-DBN, we assume the time-dependent membership distribution of a node i in the l-th layer at time step t, π (l) i,t , follows a Dirichlet distribution. Its generative process can be described as below, with the propagation of π (l) i,t illustrated in Fig. 1 (Left). For notation convenience, any parameters with index 0 are set to zero. It is noted that we have already used the observed data into the generative process.

Here, β(l-1) i′ i,t ∈ R+ is the information propagation coefficient from node i′ in the (l-1)-th layer to node i in the l-th layer at the same time t, γ (l) i′ i,t-1 ∈ R+ is the information propagation coefficient from node i′ at time t - 1 to node i at time t in the same layer l, c(-) c , c(-) u , dc are the corresponding hyperparameters and α1×K is the concentration parameter for the membership distribution in the first layer. The larger the value of these coefficients, the stronger the connections between the two corresponding latent representations (i.e., π (l-1) i′ ,t and π (l) i,t , π (l) i′ ,t-1 and π (l) i,t ).

We restrict the two nodes to have information propagated only if they are observed with positive relationship (step (a).i and (a).ii). This can reduce the computational cost of calculating β (l) -,t , γ (l) -,t from O(N 2 ) to the scale of the number of positive relationships. Also, it encourages connected nodes to have more similar membership distributions and larger dependencies between each other.

The concentration parameter ψ (l) i,t for generating π (l) i,t comprises two parts: the information propagated from all other nodes' latent representations in the (l - 1)-th layer at time t, ∑i′ β(l-1) i′ i,t π (l-1) i′ ,t , and those in the l-th layer at time (t - 1), ∑i′ γ (l) i′ i,t-1π (l) i′ ,t-1 . In other words, ψ (l) i,t is a linear sum of all the previous-layers' information at the same time step and all the previous-time steps' information in the same layer. When the coefficients β dominate over γ, the hierarchical structure plays a more important role. Otherwise, the temporal dependence has higher influence.

3.2 Application to Dynamic Relational Data

After generating the membership distributions {π (l) i,t }, we use the Bernoulli-Poisson link function [Dunson and Herring, 2005; Zhou, 2015; Fan et al., 2019a] to generate the relational data at each time step:

1. Λk1k2 ∼ Gamma(λ1, λ0), ∀k1, k2 2. Mi,t ∼ Poisson(M), ∀i, t 3. Xi,t ∼Multinomial(Mi,t;π (L) i,t ), ∀i, t; 4. For t = 1, . . . , T, i, j = 1, . . . ,N, (a) Cij,k1k2,t ∼ Poisson(Xi,k1,tΛk1k2 Xj,k2,t), ∀k1, k2 (b) Rij,t = 1(∑k1,k2 Cij,k1k2,t > 0),

where Λk1k2 is a community compatibility parameter such that a larger value of Λk1k2 indicates a larger possibility of generating the links between communities k1 and k2, λ1, λ0,M are hyper-parameters, Mi,t is a scaling parameter for generating the related counting information for node i at time t, and Cij,k1k2,t is a community-to-community latent integer for linkage Rij at time t.

Through the Multinomial distributions with π (L) i,t as event probabilities, Xi,t can be regarded as an estimator of π (L) i,t . Since the sumMi ∼ Poisson(M), according to the PoissonMultinomial equivalence, each Xi,k,t is equivalently distributed as Xi,k,t ∼ Poisson(Mπ(L) i,k,t ). Therefore, both the prior distribution for generating Xi,k,t and the likelihood based on Xi,k,t are Poisson distributions. We may form feaand marginalisation techniques [Fan et al., 2019a], which al-

sible categorical distribution on its posterior inference. This trick is inspired by the recent advances in data augmentation lows us to implement posterior sampling for Xi,k,t efficiently. The counts Xi,t lead to the generation of the K × K integer matrix Cij,t. Based on the Bernoulli-Poisson link function [Dunson and Herring, 2005; Zhou, 2015], the observed Rij,t is mapped to the latent Poisson count random variable matrix Cij,t. It is shown in [Fan et al., 2019a] that {Cij,k1k2,t}k1,k2 = 0 if Rij,t = 0. That is, only the non-zero links are involved during the inference for Cij,k1k2,t, which largely reduces the computational complexity, especially for large and sparse dynamic relational data.

Recurrent structure. Before describing the inference of Recurrent-DBN, we discuss the characteristic of the recurrent structure of ourmodel. Instead of using the one-orderMarkov property to describe the temporal dependence (assuming the state at time t depends on the states at time t-1 only), which is adopted by most probabilistic dynamic models, the deep structure of the Recurrent-DBN allows the latent variables at time t depend on those at time steps from t - 1 to t - L. For example, by using the law of total expectations, we can have the expectation of the latent count X-,t in a 2-layered Recurrent-DBN as (We use the notation - to denote a related parameter or variable hereafter):

In Eq. (2), E [X-,t|-] depends on both π-,t-1 and π-,t-2. This format can be extended straightforwardly to L-layers and involve more previous membership distributions. Such recurrent structures allow us to summarise and abstract those random variables, capturing both the hierarchical latent structures and the dynamic dependencies.

3.3 Inference

The joint distribution of the latent variables is expressed as:

While the DirBN only has upward-propagation for the latent counts and downward-sampling for the latent variables, for the Recurrent-DBN we develop an upward-backward propagation and forward-downward Gibbs sampling algorithm for count propagation and latent variable sampling. Posterior simulation for the Recurrent-DBN involves two key steps in each sampling iteration: (1) propagating the countsX upward and backward to the upper layers and previous time steps via a latent count variablem; (2) forward and downward sampling π,β,γ given the propagated latent countsm. Full updates for the other variables are similar to those in [Fan et al., 2019a].

Upward-Backward Propagating the Latent Counts

Figure 1 (right) illustrates the upward-backward propagation of counts X to the latent count variable m at each hidden layers. Generally speaking, for i, i′ = 1, . . . ,N, l = 1, . . . , L, t = 1, . . . , T, k = 1, . . . ,K, the latent variable ψ is generated as Eq. (1). m(l) i,k,t refers to the latent counts for the node i in layer l at time t for the k-th community. By integrating the m(l) i,k,t , the likelihood term of ψ (l) i,t can be calculated as:

where Γ(-) is a Gamma function.

By introducing the auxiliary variables q(l) i,t and y(l) i,k,t , the likelihood term of ψ (l) i,t can be further augmented as:

where the q(l) i,t and y(l) i,k,t can be generated as:

Consequently, y(l) i,k,t can be considered as the 'derived latent counts' for node i derived from the latent counts m(l) i,k,t . Each y(l) i,k,t can then be upward and backward distributed based on the probabilities of ψ(l) i,k,t as follows:

Here, the y(l) i,k,t is divided into two parts: one is delivered to each i′ at time t of layer l - 1 ((Z (l-1) i1,k,t , . . . , Z (l-1) iN,k,t )), and the other to each i′ at time t - 1 of layer l (A(l) i1,t-1,k , . . . , A(l) i1,t-1,k )). We denote them as Z (l-1) i-,k,t and A(l) i-,t-1,k respectively. The latent counts of lower layers and previous time steps can thus be calculated respectively as:

Let m(L) i,T =Xi,T , for t = T - 1, . . . , 2, i, i′ = 1, . . . ,N, the specification in terms of layer L is as follows,

To summarize, upward and backward propagation derives y (l) i,t from the latent countsm(l) i,t . Then, y (l) i,t is distributed to all i′ at time t of layer l-1 and time t-1 of layer l respectively as Z (l-1) i-,k,t andA(l) i-,t-1,k . Lastly, Z (l-1) -i,k,t andA(l) -i,t-1,k contribute to the generation of m(l-1) i,t and m(l) i,t-1 respectively. By repeating this process through layers and crossing time steps, we propagate the X to the m(l) upward and backward sequentially.

Forward-Downward Sampling Latent Variables

The generated ψ,q,m(l) , (Z,A) can enable to form closed Gibbs sampling algorithm for the following variables:

Sampling {π (l) i,t }i,t,l After obtaining the latent counts m(l) i,t for each layer and each time step, the posterior inference of π (l) i,t can be proceeded as:

Sampling {β(l) i′ i,t , γ (l) i′ i,t }i′ ,i,l,t The likelihood term of β(l) i′ i,t can be represented as:

The prior of β(l) i′ i,t , β(l) i′ i,t is Gam(γ (l) i , 1 c(l) ). Their posterior distribution is

4 RelatedWork

Several Bayesian deep probabilistic frameworks have been proposed to capture the temporal dependence in dynamic data [Gan et al., 2015; Gong, 2017; Henao et al., 2015]. The Deep Dynamic Sigmoid Belief Network [Gan et al., 2015] sequentially stacks models of sigmoid belief networks and uses the binary-valued hidden variables to depict the log-range dynamic dependence. The Deep Dynamic Poisson Factor Analysis (DDPFA) [Gong, 2017] incorporates the Recurrent Neural Networks (RNN) into the Poisson Factor Analysis (PFA) to depict the long-range dynamic dependence. However, in DDPFA, the parameters in RNN and the latent variables in PFA are optimized separately. Poisson Gamma Dynamic Systems (PGDS) [Schein et al., 2016] are developed to model the counting data through a "shallow" modelling strategy. Dynamic-PGDS (DPGDS) [Guo et al., 2018] is probably the closest work to our approach. Compared with DPGDS, our Recurrent-DBN differs in three aspects: (1) our RecurrentDBN generates normalized latent representations and thus provides more interpretable structures; (2) the count information is propagated in a different way; (3) our RecurrentDBN is devised in the setting of relational modelling, while DPGDS is for the topic modelling setting.

For modelling dynamic network data, many of the existing works are "shallow" probabilistic modelling. The dynamic Tensorial Mixed Membership Stochastic Block model (TMBM) [Tarr´es-Deulofeu et al., 2019] and the Fragmentation Coagulation Based MMSB (fcMMSB) [Yu and Fan, 2020] combine the notable mixed-membership stochastic block model with a dynamic setting. The Bayesian Poisson Tensor Factorization (BPTF) [Schein et al., 2015] and the Dependent Relational Gamma Process model (DRGPM) [Yang and Koeppl, 2018] are the representative works that use Poisson matrix factorization techniques to address dynamic counting data. There are also some models using the collaborative filtering techniques such as SVD++. Some methods are not developed for dynamic network data originally, but they have later been applied to the dynamic scenario, such as structurebased models like Common Neighbor (CN) [Newman, 2001], and network embedding models, including Scalable Multiplex Network Embedding (MNE) [Zhang et al., 2018] and DeepWalk [Perozzi et al., 2014]. It is noted that there is a recent trend in using the graphon theory [Lloyd et al., 2012; Orbanz and Roy, 2014] to model the network data [Fan et al., 2016; 2018b; 2018a; 2019b; 2020].

5 Experiments

We evaluate the performance of our proposed RecurrentDBN on five real-world data sets, by comparing with nine baseline methods: Mixed Membership Stochastic Block model (MMSB) [Airoldi et al., 2008], T-MBM, fcMMSB, BPTF, DRGPM, SVD++, CN, MNE and DeepWalk. Except MMSB, all of the other eight baseline models are implemented with the released code. For MMSB, we use Gibbs sampling for the inference of all variables.

5.1 Data Set and Experimental Setting

The real-world relational data sets used in this paper are: Coleman [Coleman, 1964], Mining Reality [Eagle and Pentland, 2006], Hypertext [Isella et al., 2011], Infectious [Isella et al., 2011] and Student Net [Fan et al., 2014]. The summarized statistics are detailed in Table 2. For the hyperparameters, we specify M ∼ Gamma(N, 1) for all data sets, {c(l) c , c(l) u }l, d, dc and Λk1,k2 are all given Gamma(1, 1) priors and L = 3. For MMSB, we set the membership distribution according to Dirichlet(11×K ).

5.2 Link Prediction

For link prediction, we randomly extract a proportion of 10% of relational data entries (either links or non-links) at each time step as the test set. The remaining 90% is used for training. The test relational data are not used to construct the information propagation matrix (i.e., we set {β(l) i′ i,t }l,t = 0 if Ri′ i is the testing data). We estimate the posterior mean of e-∑k1,k2 Xi,k1,tΛk1k2 Xj,k2,t as the linkage probability for each test data. These linkage probabilities are then used to calculate two evaluation metrics: the area under the curve of the receiver operating characteristic (AUC) and the precisionrecall (precision). Higher values of AUC and precision indicate better model performance.

The detail results are shown in Table 1. We report the average evaluation results for each model over 16 runs. Each run uses 3000 MCMC iterations with the first 1500 discarded as burn-in. Overall, Recurrent-DBN outperforms the baseline models for both metrics on almost all data sets. As might be expected, the value of AUC and precision increase with higher model complexity of Recurrent-DBN (i.e., larger values of K). For the other methods, fcMMSB is competitive with DRGPM and outperforms the other baselines. However, they all perform worse than Recurrent-DBN, especially for data sets with large numbers of N or T. We can see that the Recurrent-DBN has clear advantages in learning dynamic relational data, thanks to the deep hierarchical structure and recurrent long-term temporal dependence modelling.

5.3 Latent Variable Visualization

To gain further insights, we visualize the latent variables in Figure 2. It can be observed from the top part that: (1) for the same time step, the membership distributions change gradually with the increase of layers; (2) the membership distributions share some similarities for consecutive time steps and the similarities slowly shift along with the time. For example, the left bottom area of {π (3) i=1:30 } seems to have 3 different patterns: time steps t = 1, t = 2 ∼ 4, and t = 5 ∼ 10. The bottom part of Figure 2 visualizes propagation coefficients. It is reasonable to see the values of β in the first layer and γ in the first several time steps are small, since less information is propagated in these cases. The values become larger when more information is propagated. Also, the layer-wise propagation seems to have a larger influence than the cross-time propagation, with an average value of β/γ = 1.2 ∼ 1.4.

6 Conclusion

We have presented a probabilistic deep hierarchical structure named Recurrent Dirichlet Belief Networks (RecurrentDBN) for learning dynamic relational data. Through Recurrent-DBN, the evolution of the latent structure is characterized by both the cross-layer and the cross-time depen-

0.0

Figure 2: Top: visualizations of the membership distributions ({π (l) i=1:30 }3l=1) for the Infectious data set. Rows represent the nodes and columns represent the communities (with K = 10); Bottom: visualizations of average propagation coefficients β (l) t ,γ (l) t and their ratio. β (l) t ,γ (l) t are re-scaled for visualization convenience.

dencies. We also develop an upward-backward-forwarddownward information propagation to enable efficient Gibbs sampling for all variables. The experimental results on a variety of real data sets demonstrate the excellent predictive performance of our model, and the inferred latent structure provides a rich interpretation for both hierarchical and dynamic information propagation. Our Recurrent-DBN can be applied to tasks like dynamic topic models [Guo et al., 2018; Zhao et al., 2018]) and dynamic collaborative filtering. We keep these potential applications as the future work.

Acknowledgements

This work is partly supported by ARC Discovery Project DP180100966. Yaqiong Li is a recipient of UTS Research Excellence Scholarship. Xuhui Fan and Scott A. Sisson are supported by the Australian Research Council through the Australian Centre of Excellence in Mathematical and Statistical Frontiers (ACEMS, CE140100049). Bin Li is supported by Shanghai Municipal Science & Technology Commission (16JC1420401) and the Program for Professor of Special Appointment (Eastern Scholar) at Shanghai Institutions of Higher Learning.