Reciprocity, community detection, and link prediction in dynamic networks

Many real networks are dynamical, i.e., the pattern of interactions between their nodes vary over time, e.g., network of exchanged emails in a company. The abundance of such datasets and the development of optimal numerical methods have led to a growing number of studies in this field [1–4]. In addition, interactions between nodes can be reciprocated, e.g., the people whom one retweets and the number of times she retweets them vary over time; so do the papers that researchers cite in their manuscripts and papers that cite one's scientific output. This latter issue has received little attention in previous studies.

Among the main approaches to study these systems, latent variable models assume that the existence of an edge between any pair of nodes is independent of other nodes, and is conditional on some latent variables which incorporate the hidden structure of the network. These techniques mainly focus on community membership as the main relevant latent variable, e.g., in the case of citations, the people who cite each other's works, inadvertently form a community. The stochastic block model (SBM) [5–7] and its variants provide flexible network generative models [8, 9]. In this framework, nodes are initially partitioned into communities, then edges are created between nodes, based on their community membership. There are several variants of dynamical equivalents of stochastic block model (DSBM) [10–14] which capture transition of community membership over time, reflecting the evolution of edge formation. Peixoto and Rosvall [15], and Matias et al [16] develop a non-parametric temporal SBM. Gauvin et al [17] consider non-negative tensor factorization, where communities are static but the affinity matrix changes over time. Bovet et al [18] use flow of random walkers co-evolving in the dynamic network to define communities. Various methods have been used to address whether the community membership or connectivity parameters could change over time, see [19] for a review. For instance, one could assume that communities are fixed in time but the connectivity parameters across groups changes, as in [11, 17], or that communities change in time [10, 20–22].

In Zhang et al [12], the authors extend some of the popular methods of modeling network structure, e.g., SBM, to represent dynamic networks. The main idea behind their Markovian approach is to find transition rates of appearance and disappearance of edges over time. Based on these rates, they were able to calculate the average probability of edges over all time steps, hence, they estimate a steady state probability distribution for each network model, depending on its structural parameters. Although the approach followed in [12] is efficient and analytically grounded, it was developed for models that incorporate communities as the only latent variable.

Nevertheless, in directed real-world networks, community membership may not be the only factor influencing network structure. Reciprocity, i.e., the tendency of a pair of nodes to form edges on both directions, has been subject of many studies [23–25] as a crucial factor to determine the structure of networks, in particular in social networks. Bartolucci et al [26] assume local conditional independence between pairs of edges, i.e., dyads, and extend the SBM to account for the reciprocal patterns in directed dynamical networks. Furthermore, they established various specifications of the proposed model corresponding to different reciprocal assumptions.

Recently, a generative model (CRep) has been introduced that, in addition to community membership, includes reciprocity as latent variable that dictates formation of edges between the nodes [25]. In other words, the appearance of a directed edge from node i to j not only depends on the community that the nodes belong to, but also is affected by the existence of the edge from j to i. In the case of citation network, it is more likely for an author to cite those other who already cited her, implying overlapping research areas.

In this work, following the approach in [12], we extend CRep and propose a continuous-time Markov process model for dynamic networks (DynCRep). Observing the system at discrete points in time, at each time step the transition rates of appearance and disappearance of a directed edge between two nodes depends on the current community membership of the nodes, as well as on the existence of a reciprocated edge between them.

We validate the applicability of the proposed model and its inference approach by performing experiments on real and synthetic networks for community detection and link prediction. We apply the model to synthetic datasets and observe that DynCRep shows a reasonable performance in terms of link prediction. Moreover, we test the model performance on real-world datasets in the domain of social and online communication to reproduce reciprocity, with promising results.

In our model, the temporally evolving network is captured in snapshots taken at fixed intervals, from t = 0 to T + 1. A(t) represents the dynamic adjacency matrix of the network, where a non-zero value of A_ij (t) represents a weighted edge from i to node j at time t, and A_ij (t) = 0 denotes no interaction. We assume that the total number of nodes is fixed over time, i.e., new nodes do not enter into the network, and nodes do not leave it; instead, existing edges can appear and disappear. We focus on directed, and weighted networks.

A matrix w(t) of dimensions K × K determines the evolving structure of the K communities over time and we refer to w(t) as the affinity matrix. Different assumptions about w(t) result in communities with different structures. For instance, in the case of diagonal entries being greater than off-diagonal ones, communities are assortative—that is, individuals are more inclined towards intra-community interactions than inter-community interactions. The K-dimensional vectors u_i (t) and v_i (t) denote the out-going and in-coming communities at time step t, respectively.

Here, we keep the community membership constant over time; hence, we drop the notion of time dependency. We develop the model in two different varieties: (1) the affinity matrix varies over time (w-DYN), i.e., the connectivity pattern between communities changes over time, for instance, a group of nodes which form a community at time step t could be peripheral nodes at another time step [11], and (2) the affinity matrix also remains static (w-STATIC).

Following the continuous-time Markov process approach in [12], we assume that networks evolve on the real-valued times; hence, the appearance and disappearance of the edges are continuous parameters. However, we observe the network at discrete time steps. At each time step, a Poisson distribution governs the existence of edges between nodes such that an edge between two nodes is formed at a rate ${\hat{\lambda }}_{ij}(t)$. This rate depends on both the community that nodes belong to, and the existence of the reciprocated tie at the previous time step:

$$\begin{align}\hfill {\hat{\lambda }}_{ij}(t)& ={\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t-1)\hfill \\ \hfill & \equiv \sum\limits _{k,q}{u}_{ik}{v}_{jq}\enspace {w}_{kq}(t)+\eta \enspace {A}_{ji}(t-1),\hfill \end{align} \tag{ 1 }$$

where η as a hyperparameter regulates the reciprocity effects, similarly as in [25]. The difference between equation (1) and the edge probability in [25] is that the dependency on the reciprocated tie is on the previous time step, while standard CRep considers only the same time t, being an approach valid for static networks. Furthermore, an edge could disappear with rate μ.

2.1. Dynamic CRep

The aim of this study is to infer the latent parameters of the model, namely, Θ ≡ {u, v, w, η, μ}, given the adjacency matrix observed at each time step. To this end, we perform this inference task by maximizing the log-likelihood. Given Θ, all the pairs of nodes are conditionally independent; as a result, the joint-probability of the node-pairs could be approximated by a factorized form. Here, we develop a Markov process, according to which, at every time step, the probability of edges depends only on the previous time step:

$$\begin{align}\hfill P(\left\{A(t)\right\}\vert {\Theta})& =P\left(\left\{A(t)\right\}\vert \left\{A(t-1)\right\},{\Theta}\right)\hfill \\ \hfill & =\prod\limits _{i,j}\left\{P\left({A}_{ij}(0)\vert {A}_{ji}(0),{\Theta}\right)\prod\limits _{t=1}^{T}\left\{P({A}_{ij}(t)\vert {A}_{ij}(t-1),{A}_{ji}(t-1),{\Theta})\right\}\right\}.\hfill \end{align} \tag{ 2 }$$

We further assume that at the initial time step the probability A_ij (0) of an edge between two nodes follows a Poisson distribution with mean ${\hat{\lambda }}_{ij}={\lambda }_{ij}(0)$, i.e., there is no reciprocated edge in the past:

$$\begin{equation}P({A}_{ij}(0)\vert {A}_{ji}(0),{\Theta})=\frac{{\mathrm{e}}^{{\lambda }_{ij}(0)}{\lambda }_{ij}{(0)}^{{A}_{ij}(0)}}{{A}_{ij}(0)!}.\end{equation} \tag{ 3 }$$

At each time-step, edges appear with rate ${\hat{\lambda }}_{ij}(t)$, and disappear with rate μ. We follow an approach similar to that of Zhang et al [12] and calculate the probability of the existence of edges by solving a master equation. Defining ${p}_{ij}^{k}(t)$ as the probability of having k edges, i.e., an edge with the weight equal to k, between nodes i, j at time t, this quantity satisfies the following master equation:

$$\begin{equation}\frac{\mathrm{d}{p}_{ij}^{k}(t)}{\mathrm{d}t}={\hat{\lambda }}_{ij}(t)\enspace {p}_{ij}^{k-1}(t)+(k+1)\mu {p}_{ij}^{k+1}(t)-\left({\hat{\lambda }}_{ij}(t)+k\mu \right){p}_{ij}^{k}(t).\end{equation} \tag{ 4 }$$

To solve this equation, we use a generating function approach [27], by defining $g(z,t)={\sum }_{k=0}^{\infty }{p}^{k}(t){z}^{k}$. The solution for the generating function,

$$\begin{equation}g(z,t)=f\left[(z-1){\text{e}}^{-\mu t}\right]{\mathrm{e}}^{\frac{(z-1){\hat{\lambda }}_{ij}(t)}{\mu }},\end{equation} \tag{ 5 }$$

could be expanded in terms of z to give us ${p}_{ij}^{t}$ (more details in section S1 (https://stacks.iop.org/JPCOMPLEX/03/015010/mmedia)). There are four possible transitions from time t − 1 to t: (1) there is no edge neither at time t − 1, nor at t; (2) the appearance of an edge from non-edge, (3) disappearance of an existing edge, and (4) an existing edge remains; with the following probabilities, respectively,

$$\begin{align}\hfill & {p}_{0\to 0}={\text{e}}^{-\beta ({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t))}\hfill \\ \hfill & {p}_{0\to 1}=\beta ({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t)){\text{e}}^{-\beta ({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t))}\hfill \\ \hfill & {p}_{1\to 0}=\beta \enspace {\text{e}}^{-\beta ({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t))}\hfill \\ \hfill & {p}_{1\to 1}=(1-\beta ){\text{e}}^{-\beta ({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t))},\hfill \end{align} \tag{ 6 }$$

where β = 1 − e^−μ . This leads to the following time-dependent, log-likelihood:

$$\begin{align}\hfill L(T,{\Theta})& =\mathrm{log}[P(\left\{A(t)\right\}\vert \left\{A(t-1)\right\},{\Theta})]\hfill \\ \hfill & =\sum\limits _{i,j}\left\{\mathrm{log}\left[{\mathrm{e}}^{-{\lambda }_{ij}(t)}{\lambda }_{ij}{(t)}^{{A}_{ij}(0)}\right]+\sum\limits _{t=1}^{T}\mathrm{log}\left[{\text{e}}^{-\beta \left({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t)\right)}\right.\right.\hfill \\ \hfill & \quad \left.\left.\times {\left[\beta \left({\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t)\right)\right]}^{\left(1-{A}_{ij}(t-1)\right){A}_{ij}(t)}{\beta }^{{A}_{ij}(t-1)\left(1-{A}_{ij}(t)\right)}\times {(1-\beta )}^{{A}_{ij}(t-1){A}_{ij}(t)}\right]\right\}\enspace .\hfill \end{align} \tag{ 7 }$$

We add parameters' regularization by assuming gamma-distributed priors for the membership vectors:

$$\begin{equation}P({u}_{ik};a,b)\propto {u}_{ik}^{a-1}\enspace {\mathrm{e}}^{-b{u}_{ik}},\end{equation} \tag{ 8 }$$

where a ⩾ 1, to ensure the maximization of the log-likelihood (the second derivative must be negative), similarly for the v_ik . This adds new terms to the log-likelihood:

$$\begin{equation}\mathcal{L}(T,{\Theta})=L(T,{\Theta})+(a-1)\sum\limits _{i,k}\enspace \mathrm{log}\enspace {u}_{ik}-b\sum\limits _{ik}{u}_{ik}+(a-1)\sum\limits _{i,k}\enspace \mathrm{log}\enspace {v}_{ik}-b\sum\limits _{ik}{v}_{ik}.\end{equation} \tag{ 9 }$$

In the experiments below we set the values of the hyper-priors to enforce sparsity, i.e., a = 1.5, b = 10.

Maximizing $\mathcal{L}(T,{\Theta})$ requires taking the derivative of equation (9) w.r.t. each parameter individually and setting them to zero. Because the summations in the logarithm render the calculations difficult, we employ a variational approximation using Jensen's inequality. Inference is then performed using the expectation-maximization algorithm (EM); details are provided in section S1A.

Hitherto, we have included all the dependencies on the reciprocated edge A_ji (t − 1) by considering the previous time step t − 1. However, the model still applies if we incorporate the reciprocated edge at the same time step, i.e., considering A_ji (t). This choice may depend on the application itself based on the expectations and insight of the practitioner from the reciprocity effects. Alternatively, one can choose between these two options with model-selection criteria. In our experiments on real data we deployed them both, and presented the version that performs best in cross-validation tasks (section S5A).

We continue with two specifications of the model with different assumptions on the temporal evolution of the affinity matrix. In the first approach, w-DYN, the affinity matrix is treated as a time-dependent variable; while the community membership vectors, u_i , v_i , are kept static over time. Notice that a similar scenario could be obtained by fixing w and changing u_i , v_i in time [11], our model can be easily adapted to accommodate this alternative interpretation. Our model assumes fixed number of communities K. As we consider a mixed-membership model, we have the flexibility of allowing nodes to belong to various communities and with various intensities, thus allowing to capture the likelihood of the data well by effectively changing how an entry u_ik or v_ik impacts the magnitude of λ_ij (t) via w(t) in the w-DYN scenario, while keeping K constant.

In the second scenario, w-STATIC, the affinity matrix is kept static as well. The purpose of considering these scenarios is to make the model flexible in dealing with various community structures (see sections S2 to S4 for more details on each scenarios). Notice that in the case of w-STATIC, although all the latent variables are fixed in time, the network can still evolve, as edges appear and disappear based on the parameters β and μ. This is also the case for the Markov model (without reciprocity) in [12].

For instance, the EM algorithm for w-STATIC yields:

$$\begin{equation}{u}_{ik}=\frac{a-1+\sum\limits _{j,q,t}\enspace {\rho }_{ij}^{(1)}(t)\enspace {\phi }_{ijkq}\enspace {\hat{A}}_{ij}(t)}{b+\sum\limits _{j,q}{v}_{jq}\enspace {w}_{kq}\enspace \left(1+\beta \enspace T\right)}\end{equation} \tag{ 10 }$$

$$\begin{equation}{v}_{jq}=\frac{a-1+\sum\limits _{i,k,t}\enspace {\rho }_{ij}^{(1)}(t)\enspace {\phi }_{ijkq}\enspace {\hat{A}}_{ij}(t)}{b+\sum\limits _{i,k}\enspace {u}_{ik}\enspace {w}_{kq}\enspace \left(1+\beta \enspace T\right)}\end{equation} \tag{ 11 }$$

$$\begin{equation}{w}_{kq}=\frac{\sum\limits _{i,j,t}{\rho }_{ij}^{(1)}(t){\phi }_{ijkq}{\hat{A}}_{ij}(t)}{\sum\limits _{i,j}\enspace {u}_{ik}\enspace {v}_{jq}\enspace \left(1+\beta \enspace T\right)}\end{equation} \tag{ 12 }$$

$$\begin{equation}\eta =\frac{\sum\limits _{i,j,t}{\rho }_{ij}^{(2)}(t){\hat{A}}_{ij}(t)}{\sum\limits _{i,j}\enspace \sum\limits _{t=1}^{T}\enspace \beta \enspace {A}_{ji}(t-1)},\end{equation} \tag{ 13 }$$

where we defined ${\hat{A}}_{ij}(t)={A}_{ij}(t)(1-{A}_{ij}(t-1))$ if t > 0, in which ${\hat{A}}_{ij}(0)={A}_{ij}(0)$ and we have the variational distributions

$$\begin{equation}{\rho }_{ij}^{(1)}(t)=\frac{{\lambda }_{ij}}{{\lambda }_{ij}+\eta \enspace {A}_{ji}(t-1)}\quad \end{equation} \tag{ 14 }$$

$$\begin{equation}{\rho }_{ij}^{(2)}(t)=\frac{\eta \enspace {A}_{ji}(t-1)}{{\lambda }_{ij}+\eta \enspace {A}_{ji}(t-1)}\quad \end{equation} \tag{ 15 }$$

$$\begin{equation}{\phi }_{ijkq}=\frac{{u}_{ik}{v}_{jq}{w}_{kq}}{\sum\limits _{k,q}{u}_{ik}{v}_{jq}{w}_{kq}}.\end{equation} \tag{ 16 }$$

The parameter β has no closed-form update:

$$\begin{align}\hfill & -\beta \left[T\sum\limits _{i,j}{\lambda }_{ij}+\sum\limits _{i,j,t=1}^{t=T}\left(\eta {A}_{ji}(t-1)\right)+\frac{1}{1-\beta }{A}_{ij}(t-1){A}_{ij}(t)\right]\hfill \\ \hfill & \qquad +\sum\limits _{i,j,t=1}^{t=T}\left[\hat{A}(t)+{A}_{ij}(t-1)(1-{A}_{ij}(t))\right]=0,\hfill \end{align} \tag{ 17 }$$

but this equation can be solved numerically using root-finding methods. The algorithm proceeds by randomly initializing the parameters u, v, w, η, β; then we estimate the variational distributions ρ⁽¹⁾, ρ⁽²⁾, and ϕ, using equations (14)–(16) (E-step), while keeping the parameters fixed. In the next step (M-step), we update the parameters, while keeping ρ⁽¹⁾, ρ⁽²⁾ and ϕ fixed. This procedure is repeated until the convergence of the likelihood in equation (9). An overview of the algorithm is described in algorithm 1.

Algorithm 1. DynCRep (w-DYN): EM algorithm.

Input: network $A(t)={\left\{{A}_{ij}(t)\right\}}_{i,j=1}^{N}$,

number of communities K.

Output: membership $u=\left[{u}_{ik}\right],\enspace v=\left[{v}_{ik}\right]$; network

affinity matrix $w(t)=\left[{w}_{kq}(t)\right]$; reciprocity

parameter η; edge disappearance rate β(t).

Initialize u, v, w(t), η, β(t) at random.

Repeat until $\mathcal{L}$ converges:

1. Calculate ρ1(t) and ϕ(t) (E-step):

$\begin{align*}\hfill & {\rho }_{ij}^{(1)}(t)=\frac{{\lambda }_{ij}(t)}{{\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t)}\enspace ,\quad {\rho }_{ij}^{(2)}(t)=\frac{\eta \enspace {A}_{ji}(t)}{{\lambda }_{ij}(t)+\eta \enspace {A}_{ji}(t)}\enspace ,\quad \hfill \\ \hfill & {\phi }_{ijkq}(t)=\frac{{u}_{ik}{v}_{jq}{w}_{kq}(t)}{\sum _{k,q}{u}_{ik}{v}_{jq}{w}_{kq}(t)}\enspace .\quad \hfill \end{align*}$

2. Update parameters Θ (M-step):

(i) For each node i and community k update memberships:

$\begin{align*}\hfill \quad {u}_{ik}=\frac{a-1+\sum _{j,q,t}\enspace {\rho }_{ij}^{(1)}(t)\enspace {\phi }_{ijkq}(t)\enspace {\hat{A}}_{ij}(t)}{b+\sum _{j,q}{v}_{jq}\enspace {\sum }_{t=0}^{T}\hat{\beta }(t)\enspace {w}_{kq}(t)}\\ \hfill \quad {v}_{ik}=\frac{a-1+\sum _{j,q,t}\enspace {\rho }_{ij}^{(1)}(t)\enspace {\phi }_{jiqk}(t)\enspace {\hat{A}}_{ij}(t)}{b+\sum _{j,q}{u}_{jq}\enspace {\sum }_{t=0}^{T}\hat{\beta }(t)\enspace {w}_{kq}(t)}\end{align*}$

(ii) For each pair (k, q) update affinity matrix:

$\quad {w}_{kq}(t)=\frac{\sum _{i,j}{\rho }_{ij}^{(1)}(t){\phi }_{ijkq}(t){\hat{A}}_{ij}(t)}{\sum _{i,j}{u}_{ik}\enspace {v}_{jq}\hat{\beta }(t)}$

(iii) Update reciprocity parameter:

$\quad \eta =\frac{\sum _{i,j,t}{\rho }_{ij}^{(2)}(t){\hat{A}}_{ij}(t)}{\sum _{i,j,t=1}\hat{\beta }(t)\enspace {A}_{ji}(t-1)}$

2.2. Applications

2.2.1. Synthetic networks: AUC

Having explained the nuts and bolts of our model, we now turn to its application on dynamic network data. We start by considering synthetic networks generated by section 2.1 with known community structure and reciprocity. We assess the ability of the model in predicting the network at future time steps using past observations. We look in particular at the impact of reciprocity in determining edges, by generating networks with varying $\eta \in \left\{0.05,0.2,0.5\right\}$, while keeping the other parameters fixed.

For the tests reported here we use N = 500, initial average degree ⟨k⟩ = 5, and β = 0.2. We generate K = 3 hard communities of equal size with assortative structure. Having fixed the parameters, we generate 20 samples of networks for each of the three values of η. For each network we generate an initial state followed by up to T = 6 further snapshots. The initial state is generated using only the community structure (no reciprocity) using equation (3). The successive snapshots are generated according to the instructions of section 2.1. In this study, to test the ability of our model in capturing the dynamical features, we generate the first three time snapshots (T = 1, 2, 3) with an assortative community structure and the rest of the snapshots (T = 4, 5, 6) with a disassortative community structure.

For each time step t ∈ [1, T], we hide the individual snapshot A(t) and fit the data using the previous snapshots A(0), ..., A(t − 1). We test whether a model that accounts for reciprocity is able to successfully predict the network's evolution. Success is measured using the area under the curve (AUC), i.e., the probability that a randomly selected edge has higher expected value than a randomly selected non-existing edge. A value of 1 means perfect reconstruction, while 0.5 is pure random chance. The expected value of an edge is computed using:

$$\begin{equation*} \mathbb{E}\left[{A}_{ij}(t)\right]=\begin{cases}\frac{{p}_{0\to 1}}{{p}_{0\to 1}+{p}_{0\to 0}}\quad \hfill & \quad \text{if}\enspace {A}_{ij}(t-1)=0\hfill \\ \frac{{p}_{1\to 1}}{{p}_{1\to 1}+{p}_{1\to 0}}\quad \hfill & \quad \text{if}\enspace {A}_{ij}(t-1)=1\hfill \end{cases}\end{equation*}$$

$$\begin{equation}=\begin{cases}\beta (t)({\lambda }_{ij}(t)+\eta {A}_{ji}(t-1))\quad \hfill & \quad \text{if}\enspace {A}_{ij}(t-1)=0\hfill \\ 1-\beta (t)\quad \hfill & \quad \text{if}\enspace {A}_{ij}(t-1)=1.\hfill \end{cases}\end{equation} \tag{ 18 }$$

Notice that while the expected value at time t uses explicitly only the network at the previous time step, all the parameters are inferred using the whole network history, i.e., the model is trained with $\left\{A(0),\dots ,A(t-1)\right\}$. We compare with a model that does not account for reciprocity, i.e., our model with η = 0 (DynCRep₀) [25].

Figure 1 shows the results of these tests. As we can see, the ability to predict future edges is greater for a model that accounts for reciprocity, and the performance gap increases for higher values of η. This gap is partially offset by increasing the number of snapshots, as both the models have access to more information to make their estimates. Remarkably, DynCRep has stronger performance also in the low-reciprocity regime, η = 0.05. This cannot be clearly seen by looking at figure 1, as the mean AUC of the two models are within the error bars due to random fluctuations of the network structure across samples. Instead, the stronger performance of DynCRep in the low-reciprocity regime is revealed by looking at the percentage of samples where DynCRep has higher AUC than DynCRep₀, on a trial-by-trial case (see table 1 for details). While w-STATIC, the static version of the algorithm, performs slightly better than its non-reciprocated version, with larger performance gap at later times, w-DYN, the algorithm with time-varying affinity matrix, outperforms its non-reciprocated equivalent at all time steps.

Zoom In Zoom Out Reset image size

Table 1. Edge prediction in synthetic networks. The stronger performance of DynCRep in the low-reciprocity regime, η = 0.05, is revealed by looking at the percentage of samples where DynCRep has higher AUC than DynCRep₀, on a trial-by-trial case, over 20 trials.

	w-DYN		w-STATIC
T	DynCRep	DynCRep0	DynCRep	DynCRep0
1	0.0	0.0	57.0	43.0
2	71.0	29.0	43.0	57.0
3	86.0	14.0	38.0	62.0
4	67.0	33.0	43.0	57.0
5	71.0	29.0	52.0	48.0
6	81.0	19.0	57.0	43.0

Although both variants of the algorithm give better performance than their non-reciprocated version, it could be seen from figure 1 that w-DYN is more robust in link prediction tasks as η increases, and as the planted evolving structure of the affinity matrix changes from assortative to disassortative over time (T = 4, 5, 6).

2.2.2. Real world data: reciprocity/AUC

To evaluate the capability of our proposed model in retrieving network features, we apply the model to real world datasets. In this case, we first apply the inference algorithm to each time snapshot of the dynamic real dataset and learn the network's latent variables, i.e., Θ. Then, we use these latent variables as the input for the generative model, section 2.2, to generate dynamic synthetic networks similar to the fitted real datasets. Thus, we can compare dynamic synthetic networks, here 5 samples, and the original network. In this paper, we study the performance of our model in reproducing reciprocity as a significant structural parameter of the network. We implement our algorithm on two social and communication datasets, namely, email Eu core network [28] and statistics citation networks [29] (see section S6C for details on data pre-processing).

EU email network

Email-Eu-core network (EU) is constructed from internal emails exchanged between members of a large European research institution. At each time step, there is a directed edge from i to j, if i sent an email to j. Reciprocity may play a role in that receiving incoming emails may, or not, trigger a response email, similarly to other types of social communication [30]. The recorded dataset spans over a period of 803 days. However, we studied the dynamics of the dataset by dividing it in both daily and monthly durations. In the first case, we divide the edges in daily intervals (EU-daily); then select the snapshots from 5 consecutive days, randomly. In the latter case, the intervals are monthly; we select the snapshots from the first recorded year (EU-monthly).

Figure 2 shows the performance of w-DYN and w-STATIC versions of DynCRep  in reproducing the reciprocity of the EU-daily network. As expected in email networks, the reciprocity is high in this case; hence, w-DYN and w-STATIC perform similarly in reproducing reciprocity. It is noticeable that the ability of reproducing reciprocity may change depending on how the network is built. For instance, if we consider the monthly time steps, EU-monthly network, we observe a different performance, see appendix S6B.

Zoom In Zoom Out Reset image size

Figure 3 indicates the captured AUCs, measuring performance in link prediction tasks. The AUC is calculated as described in section 2.2.1. We can notice the improvement over the time snapshots, and DynCRep tends to perform slightly better. Therefore by having access to the history of the dataset and accounting for reciprocity we can achieve better results in predicting future connections.

Zoom In Zoom Out Reset image size

It is worth mentioning that we performed the experiments for different values of the number of communities; however, the results do not show high sensitivity to this parameter. Therefore, we fixed K = 4 for the EU network, equivalent to the number of departments in the corresponding institute.

Statistics citation dataset

The second example of an empirical dataset is the citation networks for statisticians, which is based on the research papers published in four of the top journals in statistics from 2003 to the first half of 2012. We construct a network by selecting a sample of the data from 2003 to 2007 and dividing it into annual intervals. This way we will have a network of citations over 4 years, where nodes are authors and an edge from nodes i to j at time step T represents that i cites j's papers in that year. In this system, we may expect that reciprocity plays a role in that receiving a citation may trigger a citation back.

Despite the fact that the reciprocity in this dataset is much lower than EU-daily dataset, figure 4 shows that we are able to capture it competitively. In addition, although the two versions outperform each others at different time steps, they still behave similarly in reproducing the reciprocity. Moreover, in both empirical datasets, the best performance is obtained for the case that reciprocated edges presented at the same time step were used in the model.

Zoom In Zoom Out Reset image size

As it could be seen from figure 5, AUC values are always higher for DynCRep, showing that accounting for reciprocity improves link prediction tasks also for this dataset. It should be noted that, at each time step T we calculate AUC by having access to the edges up to time T − 1, then predicting edges at time T. Hence, the AUC cannot be calculated for the first time step. In this case we fix K = 3, the minimum number of communities with the highest performance, i.e., we perform five-fold cross validation [25] to calculate the value of AUC, then we choose K as the number of communities with the highest value for AUC.

Zoom In Zoom Out Reset image size

In this work, we study reciprocity in dynamic networks. In reality, many datasets, e.g., networks of friendship, of gene expression patterns or communication networks, describe interactions that evolve over time, thus making them unsuitable objects of analysis for aggregate methods. In addition, the interactions in these networks might not simply change over time, but their evolution could also be affected by their past reciprocated interactions; generally, such reciprocal interactions have received little attention as additional drivers of this dynamics.

To remedy this problem, we combine insights from previous works to incorporate reciprocity into a generative model approach with latent community structure. Specifically, we extend the assumptions formulated in [25] to situations where networks change in time. For this, we consider a Markovian transition matrix which governs the evolution of the parameters over time snapshots. Being a generative model, our approach can be used to build dynamic synthetic networks, with desired reciprocity and community structure. Its algorithmic implementation is based on an efficient EM algorithm, which can be applied to large systems. As we assume a chronological order in observing the reciprocated edges, we can estimate the joint probability distribution as a factorized distribution of time steps.

We consider two varieties of our model. In one case, community membership vectors remain static over time and only the affinity matrix contains temporal information. In the other case, the affinity matrix is treated as a static parameter, similarly as the community memberships; in both cases, reciprocity parameter and the rate of edge removal are kept static. These two scenarios enable us to thoroughly analyze the model and its performance in networks with different interaction patterns. For instance, in the case of a non-homogeneous community structure over time, the first version would be a more suitable approach, since it could capture the evolving community structures.

There are a number of directions in which this model could be extended. To capture more realistic properties of the real world datasets, we can generalize the model to the case of multilayered networks, where nodes can have more than one type of interaction. For instance, in a social network, an individual can have connections based on friendship, as well as her business affiliations.

In addition, considering a node related reciprocity parameter instead of a global reciprocity parameter could improve the applicability of the model. We have focused here on the case where edges change in time, but one can envisage situations where nodes appear and disappear as well. This would also be a natural model extension. Finally, we considered here reciprocity as main network structural property, but similar investigations can be performed for other properties involving more that one pair of nodes, as triadic closure or transitivity.

The authors thank the International Max Planck Research School for Intelligent Systems (IMPRS-IS) for supporting Martina Contisciani. All the authors were supported by the Cyber Valley Research Fund.

Any data that support the findings of this study are included within the article. An open source version of the code available online at https://github.com/hds-safdari/DynCRep.