Modeling interacting dynamic networks: II. Systematic study of the statistical properties of cross-links between two networks with preferred degrees

For the homogeneous population in the steady state, there is only one ρ^ss(k), since all links are connected to nodes in the same network. Of all the distributions introduced above for two interacting networks, the most appropriate one for this comparative study is the global distribution, $\rho _1^{\rm ss}$, which should be same as $\rho _2^{\rm ss}$ for this symmetric case. As we found, these quantities settle down quite rapidly, so that we were able to exploit the relatively large systems (L = 1000 and κ = 250) used in previous studies [1]. Starting with empty networks, we discard the first 1K/2K MCS for the single-/two-network model. Thereafter, we measure the quantities of interest every 100 MCS and compile the average over 10⁴ measurements. The resultant ρ^sss are shown in figure 1. It is clear that there is no discernible difference between all three distributions, as might be expected. For the reader's convenience, we recapitulate known properties of ρ^ss and how its form (neither Poisson nor Gaussian) can be understood  [1]. We see that ρ^ss ∝e^{−μ|k−κ|} is a double exponential, i.e., a Laplacian. Our data indicate μ = 1.08 ± 0.01, which can be explained using a crude mean-field approximation to estimate the rates for our node to gain or lose a link. Respectively, these are Θ(κ − k) + 1/2 and Θ(k − κ) + 1/2, where Θ is the Heavyside function. The Θ terms correspond to how a chosen node will act, while the 1/2 accounts for how its partner nodes will act. In the steady state, each node is 'content' on the average; thus, the probability to add/cut is just 1/2. Balancing the gain/lose probability currents, we arrive at

$$\begin{equation} \frac {\rho ^{\rm ss}(k+1)}{\rho ^{\rm ss}(k)}=\frac {\Theta (\kappa -k)+1/2}{\Theta (k+1-\kappa )+1/2} \end{equation} \tag{ 5 }$$

which leads to ρ^ss(k) ∝ 3^−|k−κ|, in excellent agreement with the data.

Zoom In Zoom Out Reset image size

Needless to say, the same argument can be advanced for the symmetric two-network system, arriving at a similar result. However, as presented below, drastic differences between the single network and the interacting networks emerge when we measure another quantity: X.

In a previous paper [1], we investigated this quantity briefly. In particular, we found that X(t) displays very large fluctuations and, for systems with L = 1000, it takes very long times (≫3 × 10⁶ MCS) for X to reach the limits of its range. To build reliable histograms, it would take even longer to collect enough data. Thus, we consider smaller systems for the remainder of this paper, namely, L = 100 and κ = 25. We also discard the initial 10⁷ MCS, to let the system reach steady state, before taking measurements every 100 MCS for the next 2 × 10⁷ MCS. With the resultant time trace of 2 × 10⁵ points, we compile a histogram, which leads to P^ss(X). We also construct a power spectrum. Specifically, we divide the entire time trace into 20 shorter (10⁴ ≡ T) traces and obtain Fourier transformations of each: $\tilde {X} (\omega )\equiv \sum _{t=1}^{T}{X(t){\rm e} ^{{\rm i} {\omega }t}}$, where ω = 2πm/T (m∈[0, T − 1]). The power spectrum, I(ω), is defined as the average $\langle |\tilde {X}(\omega )^{2}|\rangle$ over these 20 FTs.

The results are plotted (in red) in figures 2 and 3. The presence of a broad plateau in P^ss(X) motivates us to explore the dynamics of X(t), to see if it simply performs an unbiased random walk within the confines of two 'soft walls'. This conjecture is confirmed by the latter plot, in which we see that I(ω) indeed follows 1/ω² quite well, crossing over (for small ω) to a constant dictated by the limits of X, namely, 0 and L².

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In stark contrast, the data displayed in green are indicative of very different behavior. Let us emphasize what X is—in this homogeneous population of 200 nodes, all preferring degree 25. First, we randomly partition the system into two sets of 100 nodes, labeling one set 'blue' and the other 'red'. At any t, the total number of 'red–blue' links is defined as X. Clearly, if we focus on any one node, there can be up to 199 contacts, though κ will limit the average to about 25.5. Thanks to homogeneity and the randomness in the dynamics, we can expect half of these to be cross-links. Thus, we are not surprised by the peak of P^ss(X) being located at ∼100 × 25.5/2 = 1275. We can go further, to estimate the observed standard deviation (σ_X ∼ 25) by the following crude argument. Denoting by x the number of cross-links of any node, we already arrived at its mean, i.e., 〈x〉 = 12.75. If we assume that the distribution of x is a binomial distribution (with probability 1/2 that the node acts on a cross-link), then this standard deviation is, roughly, $\sigma _{x}=\sqrt {25.5/4}$. Invoking the central limit theorem for 100 nodes, we easily find $\sigma _{X}=\sqrt {100}\sigma _{x}\simeq 25$, in excellent agreement with observation. As for its dynamics, we expect X(t) to be governed by white noise, as confirmed by the green line in figure 3.

By contrast, σ_X for the symmetric two-network system is over an order of magnitude larger: ∼300! (Of course, by symmetry, the average 〈X〉 is expected to be the similar, i.e., ∼1250.) The reason behind the failure of crude arguments for the symmetric, interacting two-network system is quite subtle. Exploring this non-trivial issue will be the topic of the third paper in the series [26]. There, we will study the opposite (indeed, extreme) limit of this system, thereby bringing the essence of this remarkable behavior into sharp focus. For the remainder of this paper, we will report on more systematic studies of X in typical asymmetric two-network models.

We begin by varying χ₁ = χ₂ = χ, with N_α and κ_α kept at 100 and 25, respectively. Recall that χ controls effectively the interaction between two networks, since individuals with larger χ are more likely to take action on cross-links. Thus, χ = 0 represent two independent networks, while only cross-links are present in χ = 1 networks. For general χ, however, it is not directly related to the value of X (but only to the likelihood for X to change). Therefore, our expectation is that changing χ would not affect 〈X〉 or σ_X. Figure 4 shows the simulation results for χ = 0.1, ..., 0.9. Our expectation is largely borne out, especially for σ_X. There appears to be a slight rising trend in 〈X〉, ∼15% over this entire range of χ. It is difficult to explain these variations in detail, though the typical values deviate little from that in the symmetric case (χ = 0.5; 〈X〉∼ 1272), as predicted above. We will next see that more interesting behavior appears when we vary the other two control parameters.

Zoom In Zoom Out Reset image size

The next pair we vary is κ₁ = κ₂ = κ, with fixed N_α = L = 100 and χ_α = χ = 0.5. Figure 5 shows the simulation results for κ = 25,50,75, ..., 175. Needless to say, if κ exceeds the total population size, (200 in simulations here), every link will be established quickly, and X will be a constant L² = 10⁴ and σ_X ≡ 0. Thus, our simulations sample essentially the entire range of meaningful κs here. Now, κ controls the preferred degree of a node, so it is not surprising that 〈X〉∝ κ. In particular, since ρ_α(k) is sharply peaked at κ_α, a simple-minded estimate is that each node has κ_αχ_α cross-links. Thus, we arrive at 〈X〉≃ κ_αN_α/2 (=50κ here), plotted as a solid line in figure 5 and in surprisingly good agreement with the data. We are unable to find a similar estimate for the behavior of σ_X. The most striking features there are: (i) σ_X varies linearly with κ and (ii) σ_X is symmetric around its peaks at 100. That σ_X vanishes at both end points is clear, but it is unclear how other features arise. A phenomenological formula which accounts for these features is σ_X ∝min(〈X〉, L² −〈X〉), but how it emerges from the underlying dynamics is unknown. In stark contrast, since this system is also 'symmetric', an argument similar to the above (for its counterpart in a homogeneous population) would provide $\sqrt {N\kappa /4}= \sqrt {25\kappa }$, which clearly fails to match the data. Understanding the properties of σ_X here remains a challenge.

Zoom In Zoom Out Reset image size

Finally, we vary L = N₁ = N₂ while holding κ_α = κ = 25 and χ_α = χ = 0.5. Due to limited computing power, we only explored an order of magnitude around 100: L ∈ [50,500], the results of which are shown in figure 6. Again, the dependence of the mean, 〈X〉, is easy to understand, namely, κχL ≃ 12.5L (solid line). However, similar to the systems reported above, the behavior of σ_X is more intriguing. In an inset of figure 6(b), we show a log–log plot of σ_X-L and see that σ_X scales well as L^0.63. Such anomalous scaling hints at critical phenomena and deserves a thorough investigation. In particular, the incidence matrix associated with the two communities can be viewed as an L × L Ising model (in the lattice gas language, with 0,1 in the entries). Then, X maps into 2M − L², with M being the total magnetization, while the variance $\sigma _{X}^{2}$ corresponds to the Ising susceptibility (×L²). So, our findings here imply a decreasing 'susceptibility' (∼L^−0.74). This remarkable property can be argued qualitatively, as follows. Note that κχ controls the creation of a cross-link, so that any particular link occurs, roughly, with probability κχ/L. Thus, increasing L while holding κχ fixed corresponds to an increasingly stronger magnetic field, which in turn leads to a decreasing magnetic susceptibility. In the next paper of this series [26], such a mapping can be established analytically (for the two-population model with maximally opposite preferences), while a system with N₁ = N₂ can be interpreted as a peculiar critical point [24].

Zoom In Zoom Out Reset image size

In section 3.2.1, we explored the region around a special symmetric two-network system by varying one pair of parameters (i.e., along certain lines). Here we turn our attention to a more general case, in which the networks differ by two parameters, specifically, just the sizes and the preferences. To make comparisons with the previously studied systems, we restrict ourselves to systems with fixed sums: N₁ + N₂ = 200, κ₁ + κ₂ = 50, and χ₁ + χ₂ = 1. Moreover, from the results above (especially figure 4), we see that the effects of changing χ are minimal. Thus, we simply fix both χs to 0.5. As we vary N_α, over the entire range, we need to keep, say, only κ₁ ≤ κ₂ to access the whole subspace of interest. This choice allows us to refer to α = 1,2 as the introverts and extroverts, respectively. In figure 7, we provide results for 〈X〉 and σ_X as a function of N₂, for three different pairs of (κ₁, κ₂). We choose (5,45), (15,35), and (25,25) partly for convenience and partly for having κ₁ ratios at 1:3:5.

Zoom In Zoom Out Reset image size

Focusing on the last pair first (blue points), we note that both curves are symmetric around N₂ = 100. Since the simulations were performed for pairs of N₂ around 100, the observed symmetry is an indication of the level of our statistical errors. The values of 〈X〉 and σ_X at the center are entirely consistent with the findings shown above: approximately 1300 and 300 respectively. Away from the center, 〈X〉 appears to decrease slowly at first, but turns up when N₂ reaches to within ∼10% of the boundaries. The slow decrease is related to N₁N₂ = (200 − N₂)N₂, the maximum allowed value for X. Indeed, in this regime, the fraction

$$\begin{equation}f\equiv \langle X\rangle /N_{1}N_{2} \end{equation} \tag{ 6 }$$

hovers around 13%, i.e., ∼κχ divided by the average number of nodes in each community. On the other hand, the non-monotonic behavior can hardly be expected. In section 3.2.3, we will present an approximation scheme which can provide some insight into this remarkable phenomenon.

Turning to the two asymmetric cases, we are not surprised by the lack of symmetry in the curves. However, there are prominent and interesting features, the origins of which do not readily come to mind. First, the levels of 〈X〉 are generally reduced, despite both χs being held at 0.5 and the average preference remaining at 25. Specifically, for a wide range of N₂ ≳ 100, the fs are relatively constant, matching roughly the ratio 1:3:5. Thus, it appears that the introverts are controlling the level of cross-links. Why the extroverts play a lesser role may be argued as follows: when the number of extroverts is well above κ₂, they can be 'content' by maintaining more links to other extroverts, instead of adding cross-links. Of course, it would be highly desirable to formulate an analytic and more convincing approximation scheme. Second, in addition to the sharp upturn for N₂ ∼ 200, the peaks of 〈X〉 are shifted to smaller N₂, to approximately 20 and 50 for κ₁ = 5 and 15 respectively. These features can also be roughly argued. Given that the preferred degrees are generally less than the number of nodes, it is understandable that the introverts are more likely to be frustrated, by the eagerness of extroverts to make cross-links. As we decrease N₂, (especially to values below N₁) we can expect such frustration to decrease, leading to the rise in 〈X〉 observed. Arguably, an 'optimal' state might be characterized by a balance between the total number of cross-links the introverts prefer, N₁κ₁, and that of the extroverts, N₂κ₂. (Note that we can ignore χ for this rough argument, since it affects predominantly the rate of actions on cross-links.) This balance occurs at N₂ = 4κ₁ for our parameters, i.e., 20 and 60 for the runs with black and red points, respectively. While such arguments produce a rough understanding of the data, more quantitative improvements are clearly needed.

Lastly, we turn to the fluctuations in the cross-links, characterized here by only the standard deviation σ_X, see figure 7(b). As pointed out above, these are far larger than naively expected and explaining their presence remains a challenge. Here, we merely highlight what is in the figure: the asymmetric curves for κ₁≠κ₂, the variation by over an order of magnitude, and the 'calming influence of the more introverted'. The last of these comments refers to the almost constant σ_X for the κ₁ = 5 case. We end this subsection with another observation. Noting that the peaks of σ_X for the asymmetric κs are also displaced toward the peaks of 〈X〉, we plot a 'normalized standard deviation', σ_X/〈X〉, see inset in figure 7(b). It is interesting that the peaks here are now located much closer to the center (though the curves remain asymmetric). Perhaps a detailed analysis of this quantity will facilitate the formulation of a viable theory.

In this last subsection, we present an approximation scheme which provides reasonably good agreement with data in a special regime. It is clear that, in a typical point in parameter space, there are so many competing features in our model that the contributions of all these factors will be difficult to untangle. A remarkable phenomenon—the upturn of 〈X〉 as N₂ nears its upper bound in figure 7(a)—is observed. Moreover, it appears that 〈X〉 assumes the same value at N₂ = 190 for all three pairs of κs! It behooves us, therefore, to explore this regime in more detail. The result, shown in figure 8, hints at the existence of some underlying 'universal' properties. Focusing on this regime, we discover that a simplification emerges, allowing us to gain some insight into this universality.

Zoom In Zoom Out Reset image size

In the regime of interest, introverts are seriously outnumbered and highly frustrated. In other words, they find themselves with far more contacts than they prefer (i.e., k ≫ κ₁), so that, when selected to act, they will always cut a link. A more precise and general characterization of this regime is N₁κ₁χ₁ ≪ N₂κ₂χ₂. The result of this activity is that there are no (or extremely few) I-I links in the system, so that we can attempt an approximation scheme for $\rho _1^{\rm ss}(k)$, the degree distribution of an I.⁴ From $\rho _1^{\rm ss}$, the average 〈k〉 can be computed. Since every link is a cross-link, 〈X〉 is just N₁〈k〉.

The strategy here is the same as the one used for the single network: finding the steady-state ρ^ss by balancing the rates for increasing and decreasing k, namely,

$$\begin{equation} \rho ^{\rm ss}(k)W[k\rightarrow k+1]=\rho ^{\rm ss}(k+1)W[k+1\rightarrow k] \label {rhoW} \end{equation} \tag{ 7 }$$

where W[k→k'] specifies the probability for a node with degree k to become k'. To find the appropriate Ws here, we rely on the following argument. Focus on a particular introvert, i, with existing degree k. For the regime of interest, we can assume k > κ₁, so that the only way for it to gain a link is for an extrovert (not already connected to i) choosing to add a cross-link to it. Several factors contribute to the rate for such a process:

• (N₂ − k)/N, the probability for an extrovert not connected to i to be selected
• $w_2^{+}$, the probability that this extrovert will add a link
• χ₂, the probability for it to add a cross-link (1/2 here), and
• p, the probability that this cross-link is added to i.

Now, in this regime, we believe the extroverts should be mostly 'content' (κ₂ ≪ N). This belief is supported by the observed degree distributions, a typical case shown as red points in figure 9(a). Thus, we will approximate $w_{2}^{+}$ by 1/2. As for the last quantity, p⁻¹ should be the number of introverts unconnected to our extrovert. We may estimate it by N₁(1 − f), where f is the fraction of cross-links that are present. To be slightly more accurate, especially crucial for small N₁, we recognize that our extrovert is already connected to i, so we propose

$$\begin{equation} \frac {1}{p}\simeq 1+(N_{1}-1) (1-f) . \end{equation} \tag{ 8 }$$

Putting the factors together, we have

$$\begin{equation}W[k\rightarrow k+1]=\frac {N_{2}-k}{N}w_{2}^{+}\chi _{2}p\simeq \frac {(N_{2}-k) /4N}{1+(N_{1}-1) (1-f) }. \end{equation} \tag{ 9 }$$

On the other hand, the contributions to W[k + 1→k] come from two processes. One is node i being selected to cut a cross-link. The probability of selecting this particular node i is 1/N, and with rate 1 it will cut a link. But this action will be taken on a cross-link only with probability χ₁ (1/2 here). In the other process, an extrovert connected to i is chosen (probability k/N) and cuts one of its cross-links (probability $w_{2}^{-}\chi _{2}\simeq 1/4$). Following the argument above for p⁻¹, we have

$$\begin{equation} \frac {1}{p}\simeq 1+(N_{1}-1) f. \end{equation} \tag{ 10 }$$

Of course, here p corresponds to the probability that the extrovert cuts the cross-link to i. Therefore, we reach

$$\begin{equation}W[k+1\rightarrow k]=\frac {1}{2N}+\frac {k+1}{N}w_{2}^{-}\chi _{2}p\simeq \frac {1}{2N}+\frac {(k+1) /4N}{1+(N_{1}-1) f}. \end{equation} \tag{ 11 }$$

With explicit expressions for the W's, we exploit equation (7) to derive a recursion relation for the steady-state degree distribution:

$$\begin{equation} \rho _{1}^{\rm ss}(k+1) =\rho _{1}^{\rm ss}(k) R(k) \label {DD-rr} \end{equation} \tag{ 12 }$$

where

$$\begin{equation}R\left (k\right ) =\left \{ \frac {N_{2}-k}{1+(N_{1}-1) (1-f) }\right \} \left \{ 2+\frac {k+1}{1+(N_{1}-1) f}\right \} ^{-1}. \label {R} \end{equation} \tag{ 13 }$$

Thus, $\rho _{1}^{\rm ss}(k)$ is explicitly

$$\begin{equation} \rho _{1}^{\rm ss}(k) =\rho _{1}^{\rm ss}(0) \prod _{\ell =0}^{k-1}R\left (\ell \right ) \label {rho-ans} \end{equation} \tag{ 14 }$$

where $\rho _{1}^{\rm ss}(0)$ can be fixed by the normalization $\Sigma _{0}^{N_{2}}\rho _{1}^{\rm ss}(k) =1$.

Zoom In Zoom Out Reset image size

So far, f is not a known quantity. But it can determined self-consistently, through the equations

$$\begin{equation}f=\frac {\left \langle X\right \rangle }{N_{1}N_{2}}=\frac {\left \langle k\right \rangle }{N_{2}}=\frac {1}{N_{2}}\sum _{0}^{N_{2}}k\rho _{1}^{\rm ss}(k). \label {f} \end{equation} \tag{ 15 }$$

The result is a prediction (i.e., no fitting parameters) for $\rho _{1}^{\rm ss}$ and, hence, 〈X〉. Plotted as black lines in figures 9 and 8, these predictions agree remarkably well—provided we remain in the N₁ ≪ N₂ regime. Clearly, the deviations of 〈X〉 from the theoretical curve in figure 8 reflect the limits of this regime.

Finally, we return the issue of 'universality'. Note that results equations (13)–(15) are independent of the κs. Of course, their validity relies on the assumption that k₁ ≫ κ₁, which seems reasonable in these cases where all the introverts are highly frustrated. In figure 9(b), we show data for a more extreme system (N₁, N₂ = 5,195) in which this notion of universality is indeed well borne out. An alternative perspective is that, in this regime, the N₁ × N₁ block of the full N × N adjacency matrix is frozen at zero, so that much of what we wish to compute can be gleaned from the smaller N₁ × N₂ incidence matrix. As long as we restrict ourselves to considering 'zero N₁ × N₁ blocks', the role of κ₁ is entirely marginal.