Digital Ecology: Coexistence and Domination among Interacting Networks

The overwhelming success of Web 2.0, within which online social networks are key actors, has induced a paradigm shift in the nature of human interactions. The user-driven character of Web 2.0 services has allowed researchers to quantify large-scale social patterns for the first time. However, the mechanisms that determine the fate of networks at the system level are still poorly understood. For instance, the simultaneous existence of multiple digital services naturally raises questions concerning which conditions these services can coexist under. Analogously to the case of population dynamics, the digital world forms a complex ecosystem of interacting networks. The fitness of each network depends on its capacity to attract and maintain users’ attention, which constitutes a limited resource. In this paper, we introduce an ecological theory of the digital world which exhibits stable coexistence of several networks as well as the dominance of an individual one, in contrast to the competitive exclusion principle. Interestingly, our theory also predicts that the most probable outcome is the coexistence of a moderate number of services, in agreement with empirical observations.


I. INTRODUCTION
The rapid growth of online social networks (OSNs), such as Twitter or Facebook, has led to over two billion active accounts in 2014 [1], and hence they can be said to cover over one quarter of the world population and 72% of online U.S. adults [2].Bridging the gap between social sciences and information and communication technologies, OSNs constitute a crucial building block in the development of innovative approaches to the challenges our current society faces.OSNs have already changed the nature of human interactions on a worldwide scale.In contrast to the large-scale social patterns of individuals [3][4][5][6][7], the mechanisms underlying the fate of OSNs at the system level are not at all well understood.
Real-world social networked systems exhibit a very high level of complexity [8][9][10][11][12].In a recent study, we were able to identify the main mechanisms responsible for the evolution of quasi-isolated OSNs [13].However, most OSNs operate on a worldwide scale and are in constant competition for users' attention with numerous other services; a fact that makes it extremely challenging to model them.This competitive environment leads to the extinction of some networks, while others persist.This phenomenon suggests an ecological perspective on the interaction of multiple OSNs, from which networks are considered to form a complex digital ecosystem of interacting species that compete for the same resource: users' networking time.
In standard ecology theory, Gause's law of competitive exclusion [14] states that under constant environmental conditions, two species in competition for the same resource cannot coexist.Competitive exclusion is predicted * Electronic address: kkl@ffn.ub.edu by many theoretical models [15].However, many observations of natural ecosystems seem to contradict Gause's law, as in the case of the famous plankton paradox [16].Attempts to solve such paradoxes include the assumption of different roles (competition-colonization trade off [17,18]), the increase of the dimension of the systems, the inclusion of further species properties, etc. (see [19] and references within).However, such models allow for an unlimited number of coexisting species, which thereby creates a new paradox.Indeed, real ecosystems usually consist of a moderate number of coexisting species.Here, we show that the coexistence of networks that are in competition for the same resource, namely our society's networking time, is possible.Furthermore, our work predicts that the most probable outcome is the coexistence of a moderate number of networks.
Recent work [20] showed that the competition between Facebook and its competitors such as MySpace in the mid 2000s led to the extinction of Facebook's competitors and its own prevalence.However, the current existence of a large number of OSNs [21] suggests that the coexistence of multiple networks is indeed possible.This could be explained by analogy with the competition-colonization trade-off mentioned earlier, if we assume that different networks compete for different peer groups and hence one network can persist in each of these groups.Although the existence of different peer groups is certainly the case in reality, our aim in this paper is to introduce a general and concise theory for competition between identical networks that are in competition for the same set of potential users that allows either the coexistence of any number of networks or the domination of a single network.
We show that the coexistence of competing networks can indeed be modeled by allowing for the interplay of two very common mechanisms: preferential attachment and diminishing returns.Preferential attachment [22][23][24][25][26][27][28][29][30] is a fundamental principle that can be applied to growing networks and which states that newborn nodes are most likely to connect to the more popular nodes; this leads to a rich-get-richer effect.The principle of diminishing returns-or diminishing marginal returns-is widely used in economic theories and refers to the negative curvature of production functions.For example, suppose that sowing 1 kilogram of seed in a certain place yields a crop of one ton.However, 2 kilograms of seed produces only 1.5 tons of crop; and 3 kilograms of seed produces 1.75 tons of crop.Thus, the marginal return per increment of seed diminishes with the increasing amount of seed used.
In this paper, we demonstrate the following three points.First, multiple networks can coexist in a certain parameter region due to the interplay of a rich-get-richer mechanism and diminishing returns in the dynamics of the evolution of the networks.Second, we are most likely to observe only a moderate number of coexisting services.Finally, third, the influence of the mass media controls the observed diversity in the digital ecosystem.

A. From quasi-isolated online social networks to interacting networks
The fate of a single network within the digital ecosystem depends crucially on the form of the interactions between it and its competitors, and the fitness of each of them.Nevertheless, without precise knowledge of the evolution of a single network in the absence of competitors, little insight can be gained into the fundamental interaction mechanisms.A theory of interacting networks must therefore be built on such precise understanding of the evolution of individual networks in isolation.In a recent study [13], we were able to gauge precisely the fundamental mechanisms driving the evolution of isolated networks, which we briefly summarize in what follows.
The evolution of an OSN is coupled to the underlying social structure.The following four dynamical processes drive the evolution of the system: 1. Viral activation: a susceptible node can be virally activated and added to the OSN by contact with an active neighbor in the traditional off-line network.
Such events happen at rate λ for each active link.
2. Mass media effect: each susceptible individual becomes active spontaneously at rate µ and may thus be added to the OSN in response to the visibility of the OSN.The balance between the mass media influence, µ, and the viral effect, λ, can be estimated from the topological evolution of the corresponding empirical network.This estimation can be performed by making use of the network exhibiting a dynamical percolation transition.The critical point of the transition depends on the ratio between λ and µ.This is due to the complementary roles the respective effects play in the topological evolution.
Matching the system size at the critical point then yields a linear relationship between λ and µ (see [13] for further details).
The macroscopic state of the system is characterized by the density of active nodes, defined as the quotient of active nodes over the total number of nodes, ρ a (t), and the density of passive nodes, ρ p (t).The density of susceptible nodes can be evaluated as: ρ s (t) = 1−ρ a (t)− ρ p (t).For a detailed discussion of the model we refer the reader to Ref. [13].In the present context, we want to emphasize that the model exhibits a threshold λ c below which the entire network eventually becomes passive.
Suppose now that, instead of a single network, there are n l networks competing for the same set of potential users.Each user can be active or passive in several networks simultaneously, as represented in Fig. 1, such that the long-term evolution of the fraction of active users in each layer determines the fate of the system: either several networks coexist or only a single network prevails.The first key point in the generalization of the model introduced in [13] concerns the role of the viral parameter λ.This parameter is a proxy for users' engagement in online activities, such as inviting their friends to participate in the network, generating or forwarding content, etc.However, such activities require users to spend a given amount of their time on them and their time is, obviously, bounded.This implies that when users are simultaneously active in two or more different networks or services, they are forced to decide the amount of time they devote to each of them.We model this effect by assuming that the viral parameter for each layer is λ i = λω i , where ω i a set of normalized weights (that is, n l i=1 ω i = 1) that quantify users' engagement with each OSN.In this way, n l i=1 λ i = λ is a conserved quantity related to the physical and cognitive limitations of users.The second key point in our generalization concerns the dependence of the share, λ i , of the total amount of virality for individual networks on the state of activity of the whole system, which is defined by the vector: ρ a = (ρ a 1 , ρ a 2 , • • • , ρ a n l ) T .We assume that the weights ω i are functions of ρ a that obey the following two conditions: 1. Symmetry: All networks are considered intrinsically equal.Therefore, the weight functions must satisfy the symmetry conditions: for any i, j, and k.This implies that when the fraction of active users is the same in all of them, the viral parameters λ i must also be equal in each network and, therefore, ω i = 1/n l ∀i.

Preferential attachment:
We assume that users are in general more likely to subscribe to and participate in more active networks.Hence, the weight of a given network i must be a monotonically increasing function of ρ a i .Following the same line of reasoning, we also assume that a network with zero activity is not functional, so that ω i (ρ a i = 0) = 0.
Finally, consistent with observations in [13], we assume a linear relation between µ i and λ i where ν denotes the relative strength of the viral effect with respect to the mass media effect (in [13], we found ν ≈ 4 ∼ 5).These conditions can be interpreted as coarse-grained preferential attachment in the bipartite graph consisting of users and networks.Users are in general more prone to connect to networks which exhibit higher activity and, once active in more than one network, they are also more inclined to engage with the most active one more often.Notice that we are introducing a feedback loop between the global dynamics of the system and the microscopic parameters λ i .We are thus assuming that users are, somehow, able to sense the global activity of the system.This can be achieved in practice as a combination of the amounts of information that users receive from: the network itself [31][32][33], global media, the traditional off-line social network, etc.Although preferential attachment induces a rich-get-richer mechanism, in what follows we show that the interplay of this mechanism with the dynamics of the networks leads to the emergence of stable coexistence of multiple networks across a certain parameter region.

B. Mean-field approximation
The effects of complex topologies on epidemic-like spreading processes are well understood nowadays and cannot be ignored.However, the dynamics of our model is rich and complex enough on its own to be analyzed in isolation.Therefore, in this section we perform a meanfield analysis which provides important insight into the emergence and stability of a state of coexistence of multiple networks.In particular, we replace the real social contact network by a fully mixed population with an average number of contacts per user k .Section II C contains numerical simulations of our dynamics using a real social network [13,34].We can confirm in advance that the general picture drawn in this section is also observed in the real system.
For one network, the system is described by the following mean-field equations [13] The nontrivial steady-state solution is ρ s = 0 and ρ a = 1 − δ/λ k , which is stable only when λ ≥ δ/ k ≡ λ 1 c .This defines the critical value of λ below which activity is not possible, even in a single network.In the rest of the paper, we assume that λ > λ 1 c so that, even if coexistence is not possible, at least one network is always able to survive.Likewise, we also fix the timescale of our model by setting δ = 1 from now on.
In the case of an arbitrary number of OSNs, the system is characterized by the fraction of active and passive users in each layer, ρ a i and ρ p i , and the fraction of individuals in the traditional off-line social network that are susceptible to subscription in network i: ρ s i .We assume that the densities of active/passive/susceptible nodes are not correlated between different OSNs.Thus, the evolution equations in the mean-field approximation for the i-th layer are where we have used ρ p i = 1 − ρ a i − ρ s i .Note that the coupling between different OSNs is encoded in the weights, ω i (ρ a ).The stationary solution of Eqs.(3) that corresponds to the complete coexistence of all the n l networks is given by for λ > λ n l c ≡ n l k .This again defines a critical threshold for λ below which coexistence is impossible.At the opposite extreme, the stationary solution for the prevalence of just one single network, j, is and ρ s j = 0, for λ > λ 1 c .It is easy to see that this last solution is always stable when λ > λ 1 c .However, the stability of the coexistence solution depends, in general, on the particular form of the weights ω i (ρ a ).A detailed analysis of the Jacobian matrix of the system of Eqs.(3) evaluated at the coexistence point Eq.(4) shows that this state is stable only when The emergence of stable coexistence can be understood as the interplay between preferential attachment and diminishing returns.Preferential attachment affords an advantage in terms of respective weight, ω i , for networks which already exhibit higher activity; inducing a rich-get-richer effect.However, this is damped by the intrinsic dynamics of the system, which exhibits diminishing returns in terms of activity with respect to an enhancement of the corresponding weight ω i .As long as the preferential attachment mechanism is not strong enough to overcome this damping effect, any perturbation in the density of active nodes near the coexistence point will eventually decline.Hence, the coexistence is stable.From a mathematical point of view, this is equivalent to showing that, at the coexistence point, the function φ(ρ a * i ) is proportional to the dynamical return of the system when network i is perturbed.In other words, if the activity of network i is externally increased by a small amount ∆ρ a i , after some relaxation time, the dynamics brings the perturbation to the new value ∆ρ a i = φ(ρ a * i )∆ρ a i .Coexistence is stable whenever the dynamical perturbation ∆ρ a i is smaller than the external one ∆ρ a i (see Appendix B for details).It is possible to see that φ(ρ a * i ) diverges at ρ a * i = 0 and is zero when ρ a * i = 1, and thus there is always a value of λ above which the inequality ( 6) is fulfilled (see Appendix C for details).
Interestingly, a series of states of partial coexistence exist between the complete coexistence state and the prevalence of a single network, such that only a number n c < n l of OSNs coexist simultaneously.The symmetries of the weights ω i (ρ a ) imply that any such case is exactly the same as the complete coexistence state if we replace n l by n c in Eqs. ( 4) and (6).Finally, we note that the stability of the partial or complete coexistence solutions is independent of the value of ν (see Appendix A).Therefore, we can discuss the stability in the limit ν → ∞, which reduces the dimensionality of the dynamical system.
The symmetry and preferential attachment conditions of the weights ω i (ρ a ) combined with the normalization condition imply that, without loss of generality, ω i (ρ a ) can be written as where ψ can be any monotonically increasing function bounded on [0, 1] with ψ(0) = 0. To gain further insight, we consider the following form of function ψ(ρ a i ) = [ρ a i ] σ .By adjusting a single parameter this form allows us to describe a system between a set of decoupled networks, when σ = 0, and very strongly coupled ones, when σ 0. In this particular case, the stability condition of the coexistence state of n c networks is given by This inequality defines a set of n l − 1 critical lines σ c (λ; n c ) in the plane (λ, σ) that separate phases with n c and n c − 1 maximally coexisting networks.This is illustrated in Fig. 2a for the case of n l = 5 competing networks.
However, the stability of the coexistence solution does not guarantee that it is reached from arbitrary initial conditions because, as we show above, there are several other stable fixed points, each with its own basin of attraction.This is illustrated in Fig. 3, where we show the vector field in the plane (ρ a 1 , ρ a 2 ) for the case of two competing networks in the limit ν → ∞.For any fixed value of λ > λ 2 c and σ > σ c (λ; 2), the coexistence solution is an unstable saddle point.This implies that one of the networks will eventually prevail, independently of the initial conditions (Fig. 3 top right).At the critical point σ = σ c (λ; 2), the system undergoes a subcritical pitchfork bifurcation with the appearance of two unstable saddle points moving away from the (now stable) coexistence solution as σ is decreased (Fig. 3 top left and bottom).The subcritical character of the bifurcation is akin to first-order phase transitions.Indeed, an infinitesimal increase in the value of σ near the critical point makes the system jump from stable coexistence to the domination of one of the networks.Decreasing the value of σ afterwards does not, however, bring the system back into the coexistence state, as this type of bifurcation implies a hysteresis effect, as shown in the inset of Fig. 3.The two saddle points that emerge below the critical line determine the basin of attraction of the coexistence solution.This basin (depicted in blue in the top left plot of Fig. 3) is very narrow for low densities of active nodes, as found at the beginning of the evolution.This makes the system sensitive to stochastic fluctuations; a small perturbation of the initial conditions may push the system into a state of domination of one network.

C. Real-world topology
The analysis presented in the previous section is based on two strong and unrealistic assumptions: the fully mixed hypothesis of the underlying off-line social network and the absence of fluctuations in the densities of active users.The first assumption has a strong impact on the value of the critical threshold λ 1 c and the fraction of active users in a single network when λ > λ 1 c .Fluctuations have an important impact mainly at the beginning of the evolution, when the number of active users is small.Such fluctuations can induce the system to change stochastically from one basin of attraction to another, leading the system to different steady stateseither coexistence or domination-even if it starts from the same initial configuration with identical parameters.To understand the effects of these assumptions within a real scenario, we performed large-scale numerical simulations of our model on a real social network, the Slovakian friendship-oriented OSN Pokec [13] in 2012.The size of this network (1.2 × 10 6 users) represents 25% of the population of Slovakia but demographic analysis shows that it covers a much larger fraction of the population susceptible to ever participate in OSNs.This makes Pokec a very good proxy of the underlying social structure.
We first study the coexistence space in the plane (σ, λ) in the case n l = 5.To do so, for each value of λ and σ we first set the system to the coexistence solution ρ a * .We then apply a small positive perturbation to one of the networks ρ a * 1 → ρ a * 1 + δρ a 1 .The evolution of the system after this perturbation can be used to determine the stability of the coexistence state (simulation details can be found in Appendix D 2).The results are shown in Fig. 2b.Even though the position of the critical point of a single network λ 1 c of the real Pokec network is extremely different from the mean-field prediction, the critical lines as a function of the ratio λ/λ 1 c follow a linear trend, as in the mean-field prediction.Interestingly, the slopes of these lines (although they are different from those in the mean-field case) scale with n l in the same way as in the mean-field case (see the inset in Fig. 2b).However, the stability of the coexistence solution per se does not guarantee that coexistence is reached from any initial configuration.This is particularly relevant when the evolution starts from empty networks, as fluctuations in the number of active users at the beginning of the evolution can induce the system to jump from one basin of attraction to another.Therefore, to determine the effective coexistence space in the plane (σ, λ), we evaluate the probability that a state of coexistence of a certain number of networks is reached when starting from empty networks.In the case of two competing networks, we define the effective critical line σ ef f c (λ; 2) as the line below which the probability of the two networks reaching coexistence is greater than 1/2.
Figure 2c shows the results of this program for two competing networks and ν = 4.The effective critical line follows the critical line in Fig. 2b for low values of λ and saturates at a constant value when λ/λ 1 c 1.This result can be understood in terms of the shape of the basin of attraction of the coexistence solution near the origin.Indeed, only in this region are fluctuations important enough to make the system change from one basin to the other.As an illustration, in the inset of Fig. 2c we show such a basin for n l = 2 and different values of λ in the mean-field approximation.As can be observed, the shape of the basin in the neighborhood of ρ a 1,2 ∼ 0 is almost independent of the value of λ, which explains why the probability of reaching the coexistence state saturates at a constant value.This saturation effect is similarly observed for systems of more networks, where the effective critical lines of higher coexistence states successively saturate at lower values; that is , which narrows the effective coexistence region in the plane (λ, σ) for large numbers of networks.This is particularly relevant because, although our theory allows for the coexistence of an arbitrarily large number of networks, the stochastic nature of the dynamics, coupled with the narrow form of the basin of attraction at low densities of active users, makes such coexistence highly improbable.Therefore, our model predicts-even without knowledge of the exact empirical parameters-a moderate number of coexisting networks in a large fraction of the parameter space.
The results shown in Fig. 2c are obtained for a fixed value of the parameter ν.While this parameter has no influence on the stability of the coexistence solution, and thus no effect on the results shown in Fig. 2b, it has a strong influence on the probability of reaching coexistence.Indeed, when ν is finite, the last term in Eq. ( 3) acts, at the beginning of the evolution, as a temporal boost that increases the fraction of active users in each network.This mechanism drives the system closer to the coexistence state where its attractor is broader.Figure 4 shows the simulation results of the probability of reaching coexistence as a function of ν for two competing networks.For small values of ν, the initial boost is large and the system almost always ends up in the coexistence state.For larger values of ν, the probability decreases significantly.We conclude that a higher boost-hence a smaller value of ν-favors the effective reachability of the coexistence state; whereas a small boost reduces that probability dramatically.Since ν is related to the influence of mass media, these results show that mass media influence plays a crucial role in the diversity of the digital ecosystem.
The temporal evolution of the process also shows interesting patterns.Figure 5 shows typical realizations of the process below and above the effective critical line in the case of two competing networks.It should be noted that in both cases, during the first stage of the evolution, the two networks acquire a very similar number of active users, making the forecasting of which network will eventually prevail very difficult.In a second stage, the symmetry is broken and one of the networks starts dominating, while the activity of the other declines.This pattern of "rise and fall" has been observed in many real OSNs [35].In our model, however, such behavior is a consequence of the non-linear coupling between the networks, without the need to introduce an exogenous mechanism to explain it [20].Meanwhile, the effective critical lines shown in Fig. 2c separate regions in a probabilistic way.This implies that in the vicinity of these lines, it is possible to find realizations that, with the same parameters and initial conditions, have opposite fates.This is illustrated in Fig. 6 where we show two different realizations of three competing networks.In the first column of Fig. 6, we show one such realization where two out of three networks coexist and, in the second column, a realization where only one of the three networks prevails.

III. DISCUSSION
OSNs constantly compete to attract and retain users' attention.From this point of view, OSNs and other digital services can be understood as forming a complex digital ecosystem of interacting species that compete for the same resource: our networking time.In this work, we have introduced a very general and concise theory of such an ecosystem.Akin to standard ecological theories of competing species, the fitness of OSNs increases with their performance following a preferential attachment (or rich-get-richer) mechanism.However, unlike the case of standard ecology, the total fitness of the system is a conserved quantity, which induces diminishing returns in the fitness of each network.Over a range of parameters, the combination of these two mechanisms leads to stable states of coexistence of many networks, in stark contrast to the competitive exclusion principle [14].
However, stable coexistence is only possible across a range of the parameter space, which is delimited by a critical line.At that critical line the system undergoes a subcritical pitchfork bifurcation akin to a first-order phase transition.Our model thus predicts that a minimal change or perturbation in the interactions between the different networks can have a catastrophic effect on the fate of the system.In any case, due to the stochastic nature of the dynamics and the multitude of fixed points, a stable coexistence solution is not always reached.The probability of reaching such a solution is an indicator of the diversity observed in the digital ecosystem.Interestingly, we find that over a large proportion of the parameter space the most probable outcome is the coexistence of a moderate number of digital services; in agreement with empirical observations.This number is, in general, greatly affected by the magnitude of the mass media influence.
The flexibility of our theory allows us to reproduce, with only three parameters, a large number of possible outcomes that have been observed empirically.Moreover, it can easily be modified to account for more complex situations, such as intrinsic differences between the networks or different launch times.This would allow an understanding to be gained of the extent to which higher intrinsic performance of one network can overcome the launch time advantage of another.Our model can also be extended to incorporate a description of the worldwide ecology of OSNs by incorporating different underlying societies that would represent different countries.It remains a task for future research to validate our assumptions regarding the coarse-grained coupling mechanism.
Appendix A: Jacobian matrix Here we analyze the Jacobian matrix of the dynamical system defined in Eq. ( 3) whose entries correspond to the following derivatives The Jacobian matrix (of dimension 2n l × 2n l ) can be written as where M η,θ represent n l × n l matrices with the following elements The matrix M s,s is a diagonal matrix with all its diag- s is also a diagonal matrix with all its diagonal elements equal to some value c and finally M s,a has all its elements equal to 0.
By using the Laplace expansion starting from the 2n l , 2n l entry and expanding row-wise, one finds after which means we have the eigenvalues Λ n l +1,...,2n l = d with degeneracy n l and the remaining eigenvalues are those of M a,a .Because d < 0, this means that the stability of the coexistence solution is exclusively determined by the dynamics in the limit ν → ∞, which reduces the dimensionality of the system from 2n l to n l and decouples the dynamics of ρ a i from ρ s i .See Fig. 3 for the dynamical properties.
The matrix M a,a has the form and its eigenvalues are Λ 1 = α + (n l − 1)β and Λ 2,...,n l = (α − β), which has degeneracy n l − 1.With ).From here, one obtains the result that Λ 1 is always negative and the stability is controlled by the eigenvalues Λ 2,...,n l , which are which have to be negative to satisfy stable coexistence.This leads to the condition n as described in Eq. ( 6).
Appendix B: Dynamical response to a small perturbation for n 1 To evaluate the stability of the system in the coexistence state, one has to perform an analysis of the Jacobian matrix as shown in Appendix A. In this section, we discuss the response of the system to a small perturbation, ∆ρ a i , of the activity of one network.In the limit n 1, we can neglect the effect of perturbing the i-th network on the remaining ones.The perturbation induces a shift in the corresponding weight according to Our initial perturbation triggers the dynamical response ∆ρ from the system given by where ρ a * i (ω i ) = 1 − 1/λ k ω i .With Eq. (B1), we obtain The coexistence solution is stable if the perturbation decreases; this means that the dynamical response ∆ρ a i has to be smaller than the initial perturbation ∆ρ a i .Mathematically, this leads to the condition which is equivalent to Eq. ( 6) in the limit n l 1.The left-hand side of Eq. (B4) is proportional the ratio between the dynamical response of the system and the initial perturbation.If this ratio is smaller than one, the initial perturbation will decrease and the coexistence state is stable.In the general case, one has to consider the eigenvalues of the Jacobian matrix as in Appendix A.
Appendix C: Existence of stable coexistence region Our assumptions of symmetry and normalization allow us to write where ψ(ρ a i ) is an arbitrary monotonically increasing function with ψ(0) = 0, which is bounded on the interval [0, 1].We have which we can plug into Eq.( 6) to obtain Since ψ(0) = 0 and ψ (0) = ψ(0), the left-hand side of Eq. (C3) diverges for ρ a * i → 0. Because ψ is bounded, we have lim Therefore, there is always a ρ a * i (and so a value of λ) for which the inequality (C3) is fulfilled.To simulate our model we use the Pokec network [13], a real OSN, as the underlying off-line network.We take advantage of the fact that the temporal events in our model are independent Poisson point processes, which allows us to use the Gillespie [36,37] algorithm (also known as the Doob-Gillespie [38] algorithm).For a single network, the algorithm works as follows: 1. Initialize the system and fix the rates corresponding to the respective events (here λ, µ, δ = 1) 3. Update the status of the system.So, if in step 2 a mass media activation was chosen, we randomly choose a susceptible node and change its status to active.For a deactivation, we randomly choose an active node which then becomes passive.In the case of viral activation, we randomly choose an edge connecting a susceptible and an active node and activate the susceptible node at the end of the link.Similarly, for viral reactivation we randomly choose an edge connecting a passive and an active node and the passive node at the end of this edge becomes active.We increase the time: t → t + τ .We iterate by returning to step 2 until the end of the simulation is reached.
Generalization of the algorithm to multiple layers is straightforward.One evaluates the probabilities of having a certain dynamical process in a certain layer; for example, mass media activation in layer i occurs with probability µ i N S,i , where N S,i denotes the number of susceptible nodes with respect to layer i (all the nodes which are in the underlying network but not in the i-th layer) and µ i is the corresponding rate in layer i. Accordingly, the probability of viral activation in layer i is given by λ i E SA,i , where E SA,i is the number of edges connecting active and susceptible nodes in layer i.One then chooses a dynamical processes in a certain layer in accordance with these probabilities.Finally, τ is given by the inverse of the sum over all these probabilities in all the layers.

Network activity
Figure 7: Steady state activity of our for a single network when using Pokec as the underlying off-line network as a function of the virality parameter λ (see [13]).Below a critical value λc, the network becomes entirely passive.

Empirical stability
Independent of the topological properties of the network, the activities for the steady state solution for an arbitrary number of layers is encoded in the activity curve of a single layer, which we show in Fig. 7.At the steady state of n c coexisting networks, each prevailing layer has the same share of the total virality λ i = ω i λ = λ nc , whereas the remaining ones have λ j = 0.The steady state activity of the i-th network is then given by the activity value of a single layer shown in Fig. 7 at λ = λ i .

Figure 1 :
Figure 1: Multiplex layout of two online social network layers.The bottom layer represents the underlying social structure and the remaining layers represent each OSN.

3 . 4 .
Deactivation: active users become spontaneously passive at rate δ and no longer trigger viral activations or reactivate other passive nodes.Viral reactivation: at rate λ, active users can reactivate their passive neighbors.The neighbor then becomes active and can trigger both viral activations and viral reactivations.

Figure 2 :
Figure 2: a) Regions of maximal possible coexistence in the mean-field approximation as a function of λ and σ for 5 networks evaluated from Eq. (8).b) Stability regions for the full stochastic model with a real underlying topology.For details see appendix D 2. The inset shows nc versus the inverse slope of linear fits to the respective lines.c)The most probable configuration reached from empty initial conditions for two networks.The dashed line corresponds to the empirical stability of the two networks.The insets (x and y axes each denote the activity from 0 to 1) show the basins of attraction in the mean-field approximation for σ = 0.8 and λ/λ 1 c = 4 (left), λ/λ 1 c = 6 (center), and λ/λ 1 c = 8 (right).

Figure 3 :
Figure 3: Top: Mean-field approximation in the limit ν → ∞ (this reduces the system dimension from 4 to 2 and allows the diagram to be plotted, see appendix A).Left: Stable coexistence solution (λ/λ 1 c = 4, σ = 0.8).The basin of attraction for the coexistence solution is marked in blue.Right: Only the domination solution is stable (λ/λ 1 c = 4, σ = 1.2).Bottom: Bifurcation diagram for two OSN layers showing subcritical pitchfork bifurcation at σ = σc for λ/λ 1 c = 4.The inset shows the hysteresis induced by this type of bifurcation.

Figure 4 :
Figure 4: Probability of reaching the coexistence state for two networks for different values of ν, for λ/λ 1 c = 6 and σ = 0.70.The yellow area denotes one standard deviation (from top to bottom).

Figure 5 :
Figure 5: Evolution of the fraction of active users (top) and the total number of users (bottom) for two competing networks.The first column corresponds to the parameters λ/λ 1 = 5, σ = 0.5, and ν = 4 which lies in the coexistence region.The second column represents the parameters λ/λ 1 c = 5, σ = 0.75, and ν = 4, which lies in the dominance region.

Figure 6 :
Figure 6: Evolution of the fraction of active users (top) and the total number of users (bottom) for three competing networks.Both columns correspond to the same parameters λ/λ 1 c = 7.5, σ = 0.8, ν = 4, but are different realizations.

2 .•
(a) Evaluate the number of susceptible nodes (N S ), the number of active nodes (N A ) as well as the number of edges connecting susceptible and active nodes (E SA ) and the number of edges connecting active and passive nodes (E P A ). Evaluate the sum S = µN S + λ [E SA + E P A ] + δN A .(b) Generate random numbers to choose the next event.The probabilities for the events are the following: Mass media activation: µN S /S • Viral activation: λE SA /S • Viral reactivation: λE P A /S • Deactivation: δN A /S (c) Evaluate the corresponding time step τ .The corresponding time step is given by τ = S −1 .