From subcritical behavior to a correlation-induced transition in rumor models | Nature Communications
Let us first define our model. In alignment with the DK and MT models, here we also have the same set of states (ignorant, spreader, or stifler), which are modeled by associating to each node i three Bernoulli random variables, Xi, Yi, and Zi. If node i does not know the rumor, it is classified as an ignorant (Xi = 1, Yi = Zi = 0). If it knows the rumor and is spreading it, it is called a spreader (Yi = 1, Xi = Zi = 0). However, if it knows the rumor but does not spread it, it is classified as stifler (Zi = 1, Xi = Zi = 0). Note that Xi + Yi + Zi = 1, implying that, for a fixed node, only one variable will be one, while the others will be zero. The spreading evolves through the contact between nodes defined by an undirected network, which is codified by the adjacency matrix A, whose entries Aij are equal to one if there is an edge between nodes i and j and zero otherwise. Our process is defined in continuous-time as a collection of Poisson processes. If the contact is between a spreader and an ignorant, the second node will learn the rumor and become another spreader at rate λ. On the other hand, if the contact happens between a spreader and someone that already knows the rumor (spreader or stifler), then the spreader that initiated the contact will lose interest in the rumor, thus becoming a stifler at a rate α. We also introduce a forgetting mechanism, where each stifler spontaneously becomes ignorant at a rate δ. We assume that only stiflers can forget the rumor as spreaders are actively trying to transmit the information. Denoting the state of the node i as \({({X}_{i},{Y}_{i},{Z}_{i})}_{i}\), the above described local rules are expressed as
A graphical representation of these transitions is presented in Fig. 1. Please see the Supplementary Videos for a visual example of the temporal evolution of the Monte Carlo Simulations for the standard Maki–Thompson model near the critical point and in the supercritical regime.
In this figure, the curly arrows represent transitions driven by contact, while the normal arrow represents a spontaneous transition. To simplify the notation, the states of the nodes are represented in the center of the circles, which show the only Bernoulli random variable that is equal to one.
We remark that, if δ = 0, we recover the original MT model. Notice that δ can be interpreted either as a forgetting mechanism or, for evolving rumors as, individuals who do not know the rumor’s current state that is, someone that has not been informed since the rumor changed. In this case, we could simplify the process by treating an updated rumor as, in fact, a new one. Either way, the essence of the process is captured by the forgetting mechanism (with parameter δ).
In summary, our model combines the “reluctance to tell stale news”, the annihilation by contact, with the eventual forgetting of the rumor. The latter transition assumes that one can only forget a rumor if she/he is not spreading it. The removal mechanism is the central difference between our model and the SIRS model. Namely, in our modified model, the transition from spreaders to stiflers depends on a interaction-based removal mechanism, while in the disease case, this transition is spontaneous. Although these models share similar features, we observed that their differences are enough to generate a different behavior, which is especially clear in the subcritical regime. In this regime, the SIRS has an exponential decay, while our modified rumor model presents a power-law decay to the absorbing state.
From a Markovian point of view, the modified model has a single absorbing state, while the original MT model has infinitely many absorbing states. In other words, by including of the transition \({Z}_{i}\mathop{\to }\limits^{\delta }{X}_{i}\) only the rumor-free state in which all the individuals are ignorant (i.e., Xi = 1 for all i ∈ {1, 2, …N}) traps the dynamics. Note that, for the MT model any state in which Yi = 0 for all i ∈ {1, 2, …N} is absorbing. In this case, as the system size increases, the number of absorbing states also increases. However, from a practical point of view, single absorbing state dynamics are easier to deal with as we can employ the quasi-stationary algorithm42 (for more, see section “Monte Carlo simulations”) and the lifespan method43 (for more, see section “The lifespan method”), allowing for a more precise characterization of the process.
To analyze the model, we first need to define a phase transition and corresponding critical point in our context. As MT and SIR dynamics have infinitely many absorbing states, in this case, the critical point is defined as the parameter that separates two scaling regimes. Before this point, the final number of stiflers (or recovered in the SIR) when the process reaches an absorbing state does not scale with the system size and hence its fraction goes to zero in the thermodynamic limit. After the critical point, the number of stiflers (recovered) scales with the system size. On the other hand, for processes like SIS (Susceptible-Infected-Susceptible) and our modified MT model with δ > 0 where we have a single absorbing state (disease-free state), the critical point is defined as the point where above it we have an active state in the thermodynamic limit, and below it, the system goes to the absorbing state in the thermodynamic limit. Note that our model has infinitely many absorbing states if δ = 0 and a single absorbing state if δ > 0.
To characterize the MT model, we can simulate the process beginning with many different initial conditions and measure the final fraction of stiflers and the time necessary to reach the absorbing state. Note that the initial condition must correspond to a single spreader in an ignorant population to capture the transition. Thus, the initial spreader node changes in the different independent runs of our simulation. Formally, the quantities of interest can be defined as
where \(\inf (\cdot )\) is the infimum and \(\left\langle \cdot \right\rangle\) is the average among simulations. In this context, ρZ measures how far the rumor spreads, while τ measures how fast it spreads. If δ > 0, our modified model has a single absorbing state, and it thus does not make sense to characterize the behavior of the model based on the final state as in Eqs. (2) and (3). Instead, we measure how far the rumor spreads by estimating the density of spreaders, ρY, in the metastable state using the quasi-stationary algorithm42 (for more, see section “Monte Carlo simulations”). Further, the time to reach the absorbing state, τ (see Eq. (3)), diverges exponentially in the supercritical regime when δ > 0. We instead estimate the lifespan of “finite” realizations, τf, using the method proposed in ref. 43 (for more, see section “The lifespan method”).
We remark that both the MT model and our modified version are represented by a continuous-time Markov chain with 3N possible micro-states. Thus, it is theoretically possible to write the infinitesimal generator, leading to a linear system that solves these processes exactly. However, this is not possible in practice due to the prohibitive computational cost. Notice that this is the case of the SIR model or any three-state dynamics. Note also that the concept of phase transitions is only valid in the thermodynamic limit. However, as it is commonly done in complex systems17,18,21,24,26, we also use the term phase transition in finite systems as we observe two different scaling behaviors, as mentioned above. Despite the fact that all our simulations are performed in a finite system, for this theoretical reason, most of our analysis implicitly assumes the thermodynamic limit. The exception here is section “Modified rumor model: Asymptotic analysis”, which follows an individual-based approach.
Maki–Thompson model (δ = 0)
Our first main result is showing that the MT model has a phase transition. Figure 2 shows the average fraction of stiflers ρZ as a function of λ in the standard MT model with α = 1.0 for random regular networks with different sizes. In this network, all nodes have the same number of neighbors 〈k〉k, and the connections between them are made at random, avoiding both self and multiple connections. We observe that, as the system size increases, for small λ the fraction of stiflers decreases with system size, whereas for large λ it is larger than zero and weakly dependent on the system size (see Fig. 2 for an example). Note that, as we increase the system size, for a small λ regime, i.e., λ < λc, where λc is the critical point, the fraction of stiflers goes to zero. On the other hand, for the larger λ regime, i.e., λ > λc, the curves converge to the value of the fraction in the thermodynamic limit. This behavior suggests a phase transition in the MT model as a function of λ. Complementary to the phase diagram, in Fig. 3, we show the corresponding lifespan. For a continuous phase transition, the peak of lifespan diverges in the thermodynamic limit43,44. By comparing Figs. 2 and 3 we can see that the curves for the order parameter meet at the point in which the lifespan diverges. This behavior is compatible with a continuous phase transition43,44. Moreover, we also observe an unexpected subcritical behavior, where the lifespan increases as λ decreases. For the sake of comparison, in an SIS or SIR dynamics the lifespan would be an increasing function of λ in the subcritical regime. Although these results are in striking contradiction with the mean-field approximations17,18,21,25, we highlight that the mean-field approach neglects correlations between node states which are crucial in rumor models given the contact-driven stifling mechanism. Further, the first-order mean-field approximation is not accurate for rumor models18.
Results for α = 1 and different sizes on a random regular networks with 〈k〉k = 10.
Results for α = 1 and different sizes on random regular networks with 〈k〉k = 10. The dashed line follows τf ~ λ−1.
Modified rumor model: Critical behavior
To better understand the transition in rumor models, we move to our modified model and concentrate our efforts on the α ≫ δ regime, which includes the MT model. In this case, in the transient period, we can neglect the forgetting transitions and assume that there are only two competing mechanisms, the spreading and the stifling. Thus, for locally tree-like networks and near the absorbing state, we can estimate the expected number of newly informed nodes that result from the initial spread of the rumor to a node of degree k. Note that, at least one spreading event is required to reach the absorbing state from an initial state with one informed node due to the contact-driven stifling process. Formally, we have the following recurrent expression
where qk(i) is an iterative expression that calculates the expected number of newly informed individuals in a process where the seed has k neighbors and, among them i are spreaders. This quantity should be iterated k times and accounts for: (a) the probability of spreading the rumor to one of its k − i available neighbors (ignorants), \(\frac{(k-i)\lambda }{i\alpha +(k-i)\lambda }\); which is multiplied by (b) the probability that the process stops through a stifling event due to one of its i spreader neighbors, \(i\frac{(i+1)\alpha }{(i+1)\alpha +(k-i-1)\lambda }\), (c) or that the process continues, which is encoded in the term qk(i + 1). In the case of the MT model on an infinite tree, the rumor propagation is exactly described by a branching process45. Here, we assume that each node’s degree is sampled from a fixed degree distribution that is independent of the node and that there are no correlations between the degrees of neighboring nodes. By averaging over the degree distribution, the condition that establishes the transition between a phase where the rumor dies out with probability 1 and a phase where there is a non-zero probability of infinite propagation [45 Theorem 1, Chap. I] is
where 〈⋅〉k is the expectation on the degree distribution. We remark that Eq. (5) is not a closed expression, but can be numerically solved.
In other words, Eq. (4) estimates the expected number of newly informed individuals as a result from a single initial spreader event by weighting the probabilities that the rumor stops (accounted by item (a)) or that it continues the spreading (accounted by qk(i + 1)). Moreover, the transition can be obtained from the condition at which the expected number of newly informed nodes is larger than one, i.e., are more likely to spread the information than to stop the spreading. Note also that this approximation is based on the fact that only a single event occurs at a time, which is a property of continuous Markov chains.
In Fig. 4 we compare the solution of q(1) with Monte Carlo critical point estimations for random regular networks with 〈k〉k = 10 for both δ = 0 and δ = 1 (see section “Monte Carlo simulations” for the simulation details and the Supplementary Information, Sec. I, for the individual susceptibility curves and critical point estimations). In the same figure, we also present the naive estimation, considering a first-order approximation, given as \({\lambda }^{* }=\frac{\alpha }{({\langle k\rangle }_{k}-1)}\). Note that this expression accounts for the competition between the spreading processes, with the rate 2(〈k〉k − 1), and the annihilation, with the rate 2α. Moreover, λ* is the rate at which spreading and stifling have the same probability. So, from this approximation, one would expect that, if λ > λ*, the process can spread to a fraction of the population, while if λ < λ* it should be constrained to a finite fraction of nodes. We observe that q(1) seems to be a reasonable approximation as its precision increases with α. On the other hand, λ* might be a reasonable approximation only for small enough α. These solutions suggest that the process can not be reduced to a first-order approximation, and frustrated trials to get to the absorbing state are expected and should be accounted for. Here, a frustrated trial is defined as the scenario in which the system needs one event to be forced into the absorbing state, but a spreading event takes place instead, thus avoiding the absorbing state. A crucial difference between q(1) and the mean-field approaches that fail to predict a phase transition is that q(1) takes into account that the rate at which a spreader turns into a stifler increases with the number of spreader neighbors.
Results for random regular networks with 〈k〉k = 10 and δ = 1 and N = 106. The continuous line expresses the value of λc obtained as a solution of q(1) = 1, from Eq. (5). In contrast, the dashed line represents the naive approximation that accounts only for the probability that the next event is spreading or stifling. In the inset we present the comparison for the low α regime.
Results for random regular networks with 〈k〉k = 10 and δ = 1 and N = 106. The continuous line expresses the value of λc obtained as a solution of q(1) = 1, from Eq. (5). In contrast, the dashed line represents the naive approximation that accounts only for the probability that the next event is spreading or stifling. In the inset we present the comparison for the low α regime.
Typically, heterogeneity plays a major role on dynamical processes. To understand its effects in the phase transition, we analyzed the critical properties of uncorrelated random networks with a power-law degree distribution, P(k) ~ k−γ, generated by the algorithm proposed by Catanzaro et al.46, named uncorrelated configuration model. We remark that, here, the degree correlations can be measure by means of the conditional probability \(P(k^{\prime} | k)\) that a node with degree k is connected with a vertex with degree \(k^{\prime}\). For uncorrelated networks, this probability can be estimated as the probability that any edge points to a vertex with degree \(k^{\prime}\), leading to \({P}_{{{{{\rm{unc}}}}}}(k^{\prime} | k)=k^{\prime} P(k^{\prime} )/{\langle k\rangle }_{k}.\) These networks are important from a numerical point of view because it is possible to test the behavior of dynamical systems whose theoretical solution is obtained only in the absence of correlations. Using the quasi-stationary method (see section “Monte Carlo simulations”), we estimated the critical point and focused our analysis on its behavior as a function of the network size. These results are summarized in Fig. 5 (please see the Supplementary Information, Sec. I, for the individual susceptibility curves and critical point estimations). Importantly, we observe that in power-law networks the threshold seems to vanish for γ < 3.0, and converges to a non-null value for γ > 3.0. This behavior is at odds with the behavior of the SIS model, in which the power-law networks have a vanishing critical point for any value of γ18,47,48. Conversely, it is similar to the SIRS (Susceptible-Infected-Recovered-Susceptible) model49, contact process50, the generalized SIS model with weighted infection rates51, and also modified versions of the SIS model32. In all of these models, the phase transition can be associated with a collective phenomenon that involves the activation of the whole network, whereas the phase transition in the standard SIS model shows an unusual behavior related to mutual reinfection of hubs43,52,53,54.
We plot λc as a function of N and for different values of γ and α, considering δ = 1.
Modified rumor model: Subcritical behavior
Another important result regards the subcritical regime, where the rumor can last for a “very long” time in the small λ regime. To the best of our knowledge, this behavior remained unexplored until now. We briefly mentioned this phenomenon in the analysis of Fig. 3 for the MT model. For our modified model, these behaviors can be observed in Fig. 6, which shows the lifespan τf as a function of the spreading rate λ for an uncorrelated power-law network with P(k) ~ k−γ with γ = 2.25 and N = 106. In the main figure, we show the subcritical regime, whose lifespan follows τf ~ λ−η, while in the inset, we show the peak for τf, which suggests a continuous phase transition. To calculate the lifespan we used the method proposed by Boguñá and collaborators43 (for more, see section “The lifespan method”).
Results for δ = 1, α = 0.5 on an uncorrelated power-law network with P(k) ~ k−γ with γ = 2.25, N = 106. In the main panel, we show a wide range of λ, emphasizing the sub-critical behavior, while in the inset we show the peak that suggests a continuous phase transition. The blue curve (dot dashed line) follows τf ~ λ−1 and the orange curve (dashed line) follows τf ~ λ−0.88, obtained from a fitting of the lifespan obtained using Monte Carlo simulations (the gray curve).
Results for δ = 1, α = 0.5 on an uncorrelated power-law network with P(k) ~ k−γ with γ = 2.25, N = 106. In the main panel, we show a wide range of λ, emphasizing the sub-critical behavior, while in the inset we show the peak that suggests a continuous phase transition. The blue curve (dot dashed line) follows τf ~ λ−1 and the orange curve (dashed line) follows τf ~ λ−0.88, obtained from a fitting of the lifespan obtained using Monte Carlo simulations (the gray curve).
Complementarily, Fig. 7 shows the temporal behavior of the density of spreaders ρ, where we can see that the time to reach the absorbing state is larger for λ = 0.0001 than for λ = 0.01, at odds with exponential decay. Phenomenologically, the transition from the spreader to the stifler depends on at least two individuals who are not ignorant. In the λ ≪ δ regime, the forgetting mechanism reduces the number of stiflers faster than the creation of new spreaders, thus increasing the number of necessary events to reach the absorbing state.
Result for a network with N = 106, γ = 2.25, α = 0.5 and δ = 1 in the regime λ ≪ δ for values of λ near the critical point λc ≈ 0.015.
The subcritical regime can be characterized in terms of the time to reach the absorbing state, Tabs. We remark that τ in Eq. (3) should converge to Tabs. This quantity can be approximated under the assumption that there is a single spreader and the population is locally tree-like, where nodes have 〈k〉k neighbors. In this case, assuming that node i is a spreader, the shortest sequence of events that will drive the process to the absorbing state consists of: (1) node i spreads the information to node j, (2) either i or j turn into stifler, (3) the other node turns into a stifler, (4) either i or j turn into ignorant, and (5) the other node turns into an ignorant. Since our process is a continuous Markov chain, only one event occurs at a time, the respective probabilities P and expected times 〈τ〉 of these events are
Therefore, the probability and the expected time for reaching the absorbing state through this chain are
From these quantities, we can approximate the average time necessary to reach the absorbing state, \(\left\langle {T}_{{{{{\rm{abs}}}}}}\right\rangle\). In the λ ≪ δ and λ ≪ α regime, \(\left\langle {\tau }_{1}\right\rangle \gg \left\langle {\tau }_{\ell }\right\rangle\), for ℓ = 2, 3, 4, 5. Thus, the time spent in the frustrated trials can be approximated by \(\left\langle {\tau }_{1}\right\rangle\). Here, we call as frustrated trials any sub-chain of events that does not lead to the absorbing state. So, \(\left\langle {T}_{{{{{\rm{abs}}}}}}\right\rangle\) can be approximated by counting the number of times the process fails to reach the absorbing state plus the time it succeed, which is given by
which solves as
and for the regime λ ≪ 1 and λ ≪ α, δ, follows
Note that Eq. (13) is not defined for the standard MT model, i.e., δ = 0. In the regime λ ≪ α, we expect that T ~ λ−1 and only two nodes learn the rumor. It is noteworthy that, in the subcritical regime, the time to reach the absorbing state is dominated by the time the system has to wait before a spreading event happens. In order to reach the absorbing state, whenever we have a single spreader, it first needs to inform a neighbor and, only then, the process will be allowed to reach the absorbing state. Notice that, after the spreading event, the stifling events are much faster than the spreading ones. As the spreading events are the slowest ones in the subcritical regime, the rate of these processes dominates Eq. (13), as shown in Eq. (14). Also, frustrated trials to reach the absorbing state are possible. However, in this case, we are only repeating the above-described process. Figure 8 shows T and T* as a function of λ, where the subcritical behavior is dominated by the time spent on the spreading events. Although our assumptions do not cover heterogeneous networks, T has a similar behavior as the one observed in Fig. 6 for a power-law network, suggesting that our assumptions are reasonable and might be applicable in different substrates.
Results for δ = 1.0 and α = 0.5. The orange curve is the result of Monte Carlo simulations in a random regular network with the same parameters, 〈k〉k = 10 and N = 105.
Note that in our modified model, δ > 0, the rumor might still wander around in the network before die-out. For every frustrated attempt to reach the absorbing state, the spreader individual might be a different one. In this way, although the rumor is fated to disappear, during the transient period, many different individuals might be informed about the rumor before it dies out. This effect is a consequence of the “reluctance to tell stale news” (the annihilation by contact) with the forgetting mechanism. The first mechanism imposes that at least two spreaders are necessary to reach the absorbing state, while the second mechanism enforces a single absorbing state.
Modified rumor model: Asymptotic analysis
From the analytical viewpoint, the analysis of arbitrary networks is much more difficult, and even the proof of the existence of a critical point is very hard or maybe impossible in a first-order approximation. This argument might be of interest to our community as it points to interesting future directions. In a first-order approximation we assume that the probabilities are independent. Denoting \(\left\langle {X}_{i}\right\rangle\), \(\left\langle {Y}_{i}\right\rangle\) and \(\left\langle {Z}_{i}\right\rangle\) by xi, yi, and zi, respectively, we have
Here we follow an asymptotic analysis, considering that \({y}_{i}={y}_{i}^{(1)}{\epsilon }^{c}+O({\epsilon }^{2c})\), \({z}_{i}={z}_{i}^{(1)}{\epsilon }^{k}+O({\epsilon }^{2k})\) and xi ∈ O(1), where ϵ ≪ 1. Moreover, we have to consider a scaling of the parameters as \(\lambda =\tilde{\lambda }{\epsilon }^{m}\) and \(\alpha =\tilde{\alpha }{\epsilon }^{n}\) to be able to balance the equations. Thus, from Eq. (15) in the steady-state and neglecting the higher-order terms, we have
This relation imposes that \(k=n+c+\min (c,k)\) and k = m + c, establishing three different possible regimes: (i) 0 < c < k, where c = m − n and k = 2m − n, (ii) c = k, where m = 0 and n = −k, and (iii) 0 < k < c, where c = −n and k = m − n. Without loss of generality, we also assume δ = 1. In the first regime, we have that \({y}_{i}^{(1)}=\frac{\tilde{\lambda }}{\tilde{\alpha }}\) and hence \({y}_{i} \sim \frac{\lambda }{\alpha }\). For the second and third regimes we have to assume a regular network (i.e., all nodes have the same degree) to obtain a closed-form solution, resulting in \({y}_{i} \sim \frac{\lambda }{\alpha \left(\lambda {{{\Lambda }}}_{\max }+1\right)}\) for the second, and \({y}_{i} \sim \frac{1}{\alpha {{{\Lambda }}}_{\max }}\) for the third regime, where \({{{\Lambda }}}_{\max }\) is the leading eigenvalue of the adjacency matrix. We remark that \({{{\Lambda }}}_{\max } = k\) in regular networks. For the details and derivations of this analysis, see the “Methods” section. A similar asymptotic analysis was made for the SIRS model in the Supplementary Information, Sec. II, for the sake of comparison.
Note that, we can go from one regime to another by controlling the parameters λ and α. However, regardless of the regimes, yi is always positive and larger than zero for any positive non-null rates in the three possible regimes, thus implying that the phase transition is not captured in the first-order approximation. Nonetheless, the transition was observed in our numerical experiments. Thus, we can conclude that the assumptions on the first-order mean-field approach are not enough to capture the essence of rumor models and correlations should be included.







