High prevalence regimes in the pair quenched mean-field theory for the susceptible-infected-susceptible model on networks
HHigh prevalence regimes in the pair quenched mean-field theory for thesusceptible-infected-susceptible model on networks
Diogo H. Silva, Francisco A. Rodrigues, and Silvio C. Ferreira
1, 3 Departamento de F´ısica, Universidade Federal de Vi¸cosa, 36570-900 Vi¸cosa, Minas Gerais, Brazil Instituto de Ciˆencias Matem´aticas e de Computa¸c˜ao,Universidade de S˜ao Paulo, S˜ao Carlos, SP, Brazil. National Institute of Science and Technology for Complex Systems, 22290-180, Rio de Janeiro, Brazil
Reckoning of pairwise dynamical correlations significantly improves the accuracy of mean-fieldtheories and plays an important role in the investigation of dynamical processes on complex networks.In this work, we perform a nonperturbative numerical analysis of the quenched mean-field theory(QMF) and the inclusion of dynamical correlations by means of the pair quenched mean-field (PQMF)theory for the susceptible-infected-susceptible (SIS) model on synthetic and real networks. We showthat the PQMF considerably outperforms the standard QMF on synthetic networks of distinct levelsof heterogeneity and degree correlations, providing extremely accurate predictions when the system isnot too close to the epidemic threshold while the QMF theory deviates substantially from simulationsfor networks with a degree exponent γ > .
5. The scenario for real networks is more complicated,still with PQMF significantly outperforming the QMF theory. However, despite of high accuracy formost investigated networks, in a few cases PQMF deviations from simulations are not negligible. Wefound correlations between accuracy and average shortest path while other basic networks metricsseem to be uncorrelated with the theory accuracy. Our results show the viability of the PQMFtheory to investigate the high prevalence regimes of recurrent-state epidemic processes on networks,a regime of high applicability.
I. INTRODUCTION
Pairwise approximation constitutes a valuable tool re-currently used for understanding dynamical processes ongraphs (networks or lattices) and, particularly, epidemicspreading on the top of complex networks [1–3]. Thisapproach outperforms ordinary mean-field approxima-tions extending dynamical equations from one-site to apair level [4]. Extensions to higher orders methods using n -cluster approximations [5] can lead to more accuratetheories at the cost of increasing the theoretical complex-ity. While being of limited application for low dimensionalsystems near to critical phase transitions [4], pair approx-imations can be remarkable improvements with respectto the one-site theory if either we are not too close to thetransition [6] or if the system dimension is large such asthe case of random graphs [7, 8].For dynamical processes on the top of complex net-works, heterogeneities play a central role [9, 10] thathas to be taken into account to reproduce the mostfundamental results [11, 12]. Particularly, the interplaybetween structural heterogeneity and dynamical corre-lations has been investigated using heterogeneous pair-approximations [7, 8, 13–15]. We consider the susceptible-infected-susceptible (SIS) model [9] whose dependenceon heterogeneities serves as reference for many other dy-namical processes [9, 16]. In the SIS model individualsare represented by vertices of a network and can be in ei-ther susceptible or infected states. Infected vertices healsspontaneously with rate µ and infect their susceptibleneighbors with rate λ per contact. A central aspect ofspreading phenomena is the epidemic threshold λ c abovewhich an extensive fraction of the population is infectedor, in other words, the epidemic prevalence is finite. Heterogeneities of networks can change drastically thebehavior of the threshold. If the network possesses aheavy-tailed degree distribution in the form of a power-law P ( k ) ∼ k − γ the epidemic threshold of SIS modelis zero in the thermodynamical limit when the networksize goes to infinity [17–19]. This involves a very spe-cial type of transition from an active and fluctuating toan absorbing state [20, 21] which can be knocked outwith small modifications of the SIS dynamics [22, 23].The SIS transition on uncorrelated power-law networkscan be of two types [20–22]: if the degree exponent issmall ( γ < .
5) the activation is triggered by a denselyconnected core identified by the maximal index of a k -core decomposition [20]. If the degree exponent is large( γ > .
5) then the activation is ruled by the hubs. Thelatter involves long-term epidemic activity on star sub-graphs, composed of a single hub (the center) and its k hub nearest-neighbors (the leaves), through a feedbackmechanism where the hub infects the leaves which in turnreinfect the hub. This activity must last for sufficientlylong times to permit the mutual activation of hubs whichare not directly connected (long-range) [17, 18, 21].Heterogeneities can be included in mean-field approxi-mations in different forms [10]. Two widely used approx-imations are the heterogeneous mean-field (HMF) [11]and quenched mean-field (QMF) [12, 24] theories. Theformer consists of a cross-graining where only the degreeof the nodes and the statistical degree correlations areincluded in the dynamical equations for the probabilitythat a node is infected [11, 25], and neglects the dynami-cal correlations. The latter includes the full connectivitystructure of the networks but still neglects dynamicalcorrelations assuming that the states of nearest-neighborsare independent. Due to the aforementioned nature of the a r X i v : . [ phy s i c s . s o c - ph ] J un SIS activation mechanisms, the explicit inclusion of thenetwork connectivity structure, as in the QMF approach,is imperative to construct mean-field theories of the SISmodel since heterogeneous mixing, which correspondsto an annealed network [20], will destroy localization ef-fects as for example the self-sustained activity in a starsubgraph.Despite of the detailed microscopic description of theQMF theory, neglecting dynamical correlations in SISmodel can lead to modest accuracy with significant devi-ations from simulations [26–28] if epidemic involves, forexample, activation localized in the hubs that spreadsto the rest of network [17, 20]. Indeed, the thresholdpredicted by the QMF theory for the SIS model, givenby the inverse of the largest eigenvalue of the adjacencymatrix [12, 19] (see Sec. II for details), involves a local-ized phase on random power-law networks with degreeexponent γ > . etal . [31] claimed that it has limitations to compute highepidemic incidence regimes and proposed that a micro-scopic Markov chain approach (MMCA) [35], which is adiscrete time version of the QMF theory, could be usedinstead. However, a nonperturbative approach is possiblethough numerical integration of both QMF and PQMFdynamical equations. Since large discrepancies betweendiscrete and continuous time approaches can be presentin the SIS dynamics [36], a nonperturbative analysis ofQMF and PQMF theories is necessary. We develop non-perturbative analyses of both QMF and PQMF theoriesfor SIS on networks using numerical integration of thecorresponding dynamical equations. We consider bothlarge synthetic networks generated with the Weber-Portomodel [37] and real networks with different levels of degreecorrelation. We address regimes not asymptotically closeto the epidemic threshold since the mean-field theories failin predicting very low densities of infected vertices [28, 38].However, numerical analyses supports that this asymp-totic scaling is confined into a small interval near to λ = λ c → γ < . γ > .
5) and correlations (assor-tative, disassortative, or uncorrelated) investigated whileQMF theory presents significant deviations for γ > . II. MEAN-FIELD THEORIES
We investigate the SIS model on a connected, undi-rected, and unweighted network with i = 1 , . . . , N verticeswhose structure is encoded in the adjacency matrix A ij defined by A ij = 1 if i and j are connected and A ij = 0otherwise. The healing rate is fixed to µ = 1 without lossof generality.The probability that a vertex i is infected, representedby ρ i , evolves as [27] dρ i dt = − ρ i + λ (cid:88) j φ ij A ij (1)where φ ij is the probability that a vertex i is susceptibleand its nearest-neighbor j is infected. Equation (1) is ex-act but not closed. A closed system is obtained taking theone-site approximation φ ij ≈ ρ i (1 − ρ j ) that correspondto the QMF theory [24, 29] dρ i dt = − ρ i + λ (1 − ρ i ) N (cid:88) j =1 A ij ρ j . (2)The QMF epidemic threshold is given by λ QMFc Λ (1) = 1where Λ (1) is the largest eigenvalue (LEV) of the adjacencymatrix A ij .The PQMF theory includes dynamical correlations con-sidering the evolution of φ ij which, in its complete butnot closed form, depends on the triplets; see Ref. [27].The PQMF theory consists of approximating the triplets[ A i , B j , C l ] in which i and l are both connected to j by[ A i , B j , C l ] ≈ [ A i , B j ][ B j , C l ][ B j ] . (3)Here A , B , and C are the states of the vertices that, inthe SIS case, can be either infected or susceptible. Theapproximation given by Eq. (3) considers that ( i, j, l ) doesnot form a triangle, i.e., i and l are connected to i but notto each other. Actually, the effects of clustering have beenrecently investigated [34] and it was shown that even innetworks with high cluster coefficient, plenty of triangles,the approximation given by Eq. (3) performs very wellfor the steady state. The final PQMF equation for φ ij becomes dφ ij dt = − (2 + λ ) φ ij + ρ j + λ (cid:88) l ω ij φ jl − ρ j ( A jl − δ il ) − λ (cid:88) l φ ij φ il − ρ i ( A il − δ lj ) , (4)in which ω ij = 1 − φ ij − ρ i . Equations (1) and (4) form aclosed system of N + M equations where M = (cid:80) j k i isthe number of edges of the network. More details of thederivation are given in Ref. [27].The epidemic threshold within the PQMF frameworkis given by the transcendent equation λ c Ω (1) ( λ c ) = 1 , [28]where Ω (1) is the largest eigenvalue of the weighted adja-cency matrix B ij given by B ij = 2 + λ λ + 2 A ij λ k i λ +2 . (5)See [28] for the derivation of Eq. (5).Very close to the threshold, the epidemic prevalence ρ , defined as ρ = N (cid:80) i ρ i , approaches zero following apower-law in the form ρ (cid:39) a ( λ − λ c ) β where β is a criticalexponent [4] and a is prefactor that may depend on thenetwork size. Either in QMF [29, 39] and PQMF [28]theories it can be shown that β QMF = β PQMF = 1 while a QMF or a PQMF can be expressed in terms of the prin-cipal eigenvector (PEV) { v (1) i } of either A ij or B ij ( λ c ),respectively, as [28, 29] a = (cid:80) Ni =1 v (1) i N (cid:80) Ni =1 (cid:104) v (1) i (cid:105) . (6)These results are straightforwardly derived when the net-work presents a spectral gap Λ (1) (cid:29) Λ (2) , where Λ (2) isthe second LEV of the adjacency matrix; see e.g. [28, 29].However, it was shown that β QMF = 1 is always true [40]and the same is expected for the PQMF theory since pairapproximations should not change the universality classpredicted by the one-vertex theory [5].The mean-field exponent β = 1 does not match therigorous results obtained by Mountford et al . [38] in ther-modynamical limit N → ∞ where ρ ∼ λ − γ if 2 < γ < / λ γ − (ln λ ) γ − if 5 / < γ < λ γ − (ln λ ) γ − if γ > , (7)according to which β > γ >
2. For largenetworks with γ < /
2, where the epidemic thresholdis very accurately reproduced by the QMF theory, thenumerical integration performed in Ref. [28] confirms thedeviation from the exact result for λ approaching λ QMFc = (1) while stochastic simulations are in agreement withthe rigorous results. However, the simulations show apre-asymptotic behavior fully consistent with β QMF = 1. × -4 × -4 × -4 × -4 × -4 λ−1/Λ (1) × -4 × -4 × -4 × -4 × -4 ρ QMF (N=10 )QMF (N=10 )N=10 N=10 FIG. 1. Epidemic prevalence around the transition threshold λ c = 1 / Λ (1) for UCM networks with γ = 2 . k min = 3, and k max = 2 √ N . Simulations are represented by solid lines withsymbols and numerical integration of the QMF equation (2)is given by the dashed lines. Horizontal arrows indicate theinterval where curves depart from linearity. In a linear scale , the region that departs from linearityis squeezed around λ = λ c → a QMF also does: we found a QMF = 0 . N = 10 and a QMF = 0 . N = 10 . Finally, we see that QMFis not able to capture quantitatively the amplitude of thelinear region observed in simulations reinforcing the needof nonperturbative analyses of the PQMF theory.A closed solution for Eqs (1) and (4) can be derivedfor the particular case of homogeneous networks where P ( k ) = δ k,m for which ρ i = ρ and φ ij = φ . The expressionfor stationary epidemic prevalence is ¯ ρ = λ − λ c m − + λ − λ c , λ c = 1 m − . (8) III. RESULTS
We numerically integrated QMF and PQMF equationsusing a fourth-order Runge-Kutta method with time step δt = 10 − to 10 − . Initial conditions consistent withthe exact closure equations relating pairwise and singlevertex probabilities such as [ S i , I j ] + [ I i , I j ] = [ I j ] must bechosen and the steady state is insensitive to a particularchoice. We performed stochastic simulations of the SIS Obviously, this region will not shrunk in a logarithm scale andthe asymptotic scaling is the theoretical one given by Eq. (7). The solution is same derived for the contact process in, e.g.,Ref. [41] replacing the infection rate λ CP = mλ SIS , where λ CP isthe infection rate for the contact process and λ SIS for SIS. dynamics on networks using an optimized Gillespie algo-rithm [42]; see Appendix C. The absorbing states, whichare rigorously the unique real stationary state in finitesize networks, were circumvented using quasistationary(QS) simulations [43]; see Appendix B.
A. Synthetic networks
A simple metrics to quantify the correlations is theaverage degree of the nearest-neighbors of vertices with agiven degree k [9, 44], represented by κ nn ( k ). The func-tional form of κ nn ( k ) reveals correlation patterns of thenetwork. If κ nn is an increasing function of k , the networkpresents assortative correlations where vertices of similardegree tend to be connected. Conversely, if κ nn decreaseswith k the network has disassortative correlations wherevertices of high degree tend to be connected with verticesof low degree. Finally, if κ nn does not depend of thedegree, the networks is said uncorrelated or neutral andassumes de value κ nn = (cid:104) k (cid:105) / (cid:104) k (cid:105) [9].We investigated networks with distributions P ( k ) ∝ k − γ and k = k min , . . . , k max where k min = 3. For γ < k max = 2 √ N that permits to build networks with-out degree correlation with the uncorrelated configurationmodel (UCM) [45]. The factor 2 helps to accelerate theconvergence to the asymptotic limit where both N → ∞ and k max → ∞ [28]. For γ >
3, a rigid cutoff givenby P ( k max ) N = 1 [46] was used to suppress multiple(localized) transitions [47] and facilitating threshold de-termination. Degree correlations were included using thebenchmark model proposed by Weber and Porto [37], here-after called Weber-Porto configuration model (WPCM);see Appendix A. The dependence κ nn ( k ) ∝ k α was in-vestigated where α < α = 0 and α > γ = 2 . γ , Figs 2(b) and (c). The accuracy of thetheories at low prevalence is better for assortative andworse for disassortative networks when compared withthe neutral case. Another interesting dependence on theassortativity can be observed in these curves. For lowprevalence, assortative and disassortative networks pos-sess, respectively, higher and lower densities if comparedwith the uncorrelated networks. At high prevalences,the converse is observed where disassortative networkspresent higher densities than the assortative and neutralnetworks. The same behavior is observed for all values of γ , indicating that it is related to the degree correlations.The behavior at low densities can be explained in terms λ ρ α=−0.2α=0α=0.2 QMFPQMF (a) λ ρ α=−0.2α=0α=0.2 QMFPQMF0.2 0.22 0.24 0.26 0.28 (b) λ ρ α=−0.2α=0α=0.2 QMFPQMF (c)
FIG. 2. Epidemic prevalence as a function of the infection ratefor WPCM networks with N = 10 vertices, degree exponent(a) γ = 2 .
3, (b) 2.8 and (c) 3.5 and different levels of degreecorrelations. Symbols represent stochastic simulation whilesolid and dashed lines the numerical integration of the QMFand PQMF theories, respectively. Bottom and top insets showzoom of low and high prevalence, respectively. of reduced capacity to transmit infection when hubs aresurrounded by low degree vertices in the disassortativecase rather than being directly connected more likely inthe assortative case. We cannot provide simple argumentsfor inverted dependence on the assortativity degree athigh densities, but it is very precisely reproduced by thePQMF theory. We also simulated the SIS on larger net-works with N = 10 vertices and the level of accuracy ofthe mean-field theories is similar.We also investigate the role of the heterogeneity compar-ing the PQMF theory with a pair homogeneous mean-field λ ρ γ=2.3 (PQMF) γ=2.3 (PHMF) γ=2.8 (PQMF) γ=2.8 (PHMF) γ=3.5 (PQMF) γ=3.5 (PHMF) FIG. 3. Comparison of pair quenched (dashed lines) andhomogeneous (solid lines) mean-field theories for uncorrelatednetworks with different levels of heterogeneity. The networksize is N = 10 . theory (PHMF) using Eq. (8) with m replaced by theaverage degree (cid:104) k (cid:105) of the network [48, 49]. The densityof infected vertices obtained in both pairwise approachesare shown in Fig. 3. Beyond the expected discrepancy fordescribing the low prevalence regimes, since one theorypredicts a finite while the other a vanishing threshold,the regime of high epidemic prevalence is affected by in-clusion of heterogeneity. As one could expect, the moreheterogeneous networks present the larger discrepanciesbetween homogeneous and heterogeneous theories. B. Real networks
Real networks usually present some degree of correla-tion and, in many case, the patterns can be quite complexexhibiting both assortative and disassortative correlationsfor distinct ranges of degree [9, 50]. Therefore the com-parison between mean-field theories and simulations arenecessary in order to determine in which extent the accu-racy observed in synthetic networks holds in the real-worldcounterparts. We selected some networks with differentlevels of heterogeneity, sizes, and correlations recentlyused in the investigation of epidemic processes [28, 51, 52].For detailed information about the original references forall the networks see Refs. [53, 54].Figure 4 presents the prevalence as a function of theinfection rate for 12 real networks. We remark that dataasymptotically close to the epidemic threshold are knownto mismatch simulations [28] and are beyond the scopeof the present work. In some cases, QMF and PQMFare indistinguishable from each other and agree almostperfectly with simulations in the scale presented in thesefigures. In other cases, QMF theory deviates considerablyfrom simulations while PQMF remains accurate. In orderto quantify the differences we define the relative deviationof densities obtained in simulations ( ρ sim ) and the QMF theory ( ρ QMF ) as η QMF = (cid:82) λ λ [ ρ QMF ( λ ) − ρ sim ( λ )] dλ (cid:82) λ λ ρ sim ( λ ) dλ , (9)where λ and λ are the initial and final infection ratesin the simulations presented in Fig. 4, and an equivalentdefinition of η PQMF for the PQMF theory. The rela-tive deviations are given in Table I. The positivity of η shows that the mean-field theories overestimate the den-sity obtained in simulations as expected since dynamicalcorrelations, pruned in mean-field theories, reduce thespreading capacity of the epidemic process. As can beseen, we have η PQMF (cid:28) η QMF such that PQMF is alwaysmore precise than QMF. Three networks present signifi-cant deviations of the PQMF theory from the simulations,namely, for Amazon customer, electronic circuit S838, andair traffic networks with 24%, 7% and 6.6% of deviation,respectively.
C. Accuracy versus structural properties
The gain of PQMF theory with respect to QMF in realnetworks is expressive but it still deviates from simulationsin some cases, as shown in Table I. One central questionis to determine when either QMF or PQMF performanceis satisfactory enough. Near to the transition point, whenthe prevalence is very small, a relation between accuracyand localization of the PEV of the weighted adjacency ma-trices obtained in the linearization has been proposed [28].This is justified by the fact that a leading contribution forthe probability that a vertex i is infected in the mean-fieldtheories is proportional to the corresponding PEV of A ij or B ij ( λ c ) for QMF and PQMF theories, respectively. Forsake of completeness of Ref. [28], in which the accuracyon epidemic threshold was discussed thoroughly, Fig. 5shows the steady state density calculated slightly abovethe epidemic threshold of the PQMF theory against thenetwork size. The PQMF theory is much more accuratethan QMF but also start to deviate from simulations asthe networks size increases. In both cases the accuracy islarger for networks with less localized PEV as quantifiedby the inverse participation ratio (IPR) [29] defined as Y (1) = N (cid:88) i =1 (cid:104) v (1) i (cid:105) , (10)where v (1) i is the normalized PEV defined in Sec. II . Aslarger as the IPR, more localized is the PEV. Insets ofFig. 5 shows the IPR for both QMF and PQMF theorieswhere we see that the latter is much less localized thanthe former, but still increases towards a finite value asthe network size increases indicating localization asymp-totically.However, the nonperturbative theory accounts forthe contributions of the complete basis of eigenvectors, ρ AS CaidaN=26475 Amazon 3/2/2003N=262111 AstroPhys-1993-2003N=17903 CiteSeerN=365154 ρ JungN=6120 S 838N=512 SlashdotN=51083 URV emailN=1133 λ ρ λ λ λ CoraN=23166 AirTrafficN=1226 Gnutella 8/24/2002 N=26498 cElgansNeuralN=297
FIG. 4. Epidemic prevalence on real networks. Symbols represent stochastic simulations while solid and dashed lines representthe numerical integration of the QMF and PQMF equations, respectively. In each panel, the usual name and size of networksare given.Network η QMF η PQMF Y (1)QMF Y (1)PQMF AS Caida 0.032 0.0022 0.0240 0.0139Amazon 3/2/03 0.75 0.24 0.106 0.0114AstroPhys 93-03 0.035 0.013 0.0045 0.0043CiteSeer 0.067 0.010 0.0177 0.0109Jung 0.0079 0.0021 0.0478 0.0335S 838 0.28 0.070 0.179 0.0340Slashdot 0.032 0.0011 0.144 0.0347URV email 0.022 0.0028 0.0096 0.0087Cora 0.152 0.037 0.0100 0.0090Air Traffic 0.38 0.066 0.0191 0.0154Gnutella 8/24/02 0.17 0.0066 0.214 0.0800cElegans Neural 0.057 0.0078 0.0189 0.0175TABLE I. Relative deviations and inverse participation ratiosfor QMF and PQMF theories applied to real networks. whether it is of A ij or B ij ( λ c ), and this comparison isnot justifiable anymore. Indeed, as shown in Table I,the accuracy of the QMF theory can be high even whenthe PEV is localized as, for example, in the case of theSlashdot network. In other cases, as Air Traffic and Coranetworks, the PEV localization corresponding to QMFand PQMF theories are similar but the performance ofthe latter is much better. We performed a statistical analysis of the correlations between ln η and ln Y (1) andfound no statistically significant p -values of p QMF = 0 . p PQMF = 0 .
46. It is worth to note that the statisticalanalyses considering the linear data present even smallerstatistical significance.We also checked (logarithm) statistical correlations of η with other basic network metrics, namely, the hetero-geneity coefficient ε = (cid:104) k (cid:105) / (cid:104) k (cid:105) , the modularity coefficient Q [55], the average clustering coefficient (cid:104) c (cid:105) , and averageshortest distances (cid:104) (cid:96) (cid:105) . We found statistical significancewith p < .
02 only with ε and (cid:104) l (cid:105) . Actually, ε and (cid:104) l (cid:105) are correlated since more heterogeneous networks tendto have shortest average distances due to the shortcutsintroduced by hubs [50]. The correlation between η and (cid:104) (cid:96) (cid:105) actually is not very surprising since one intuitivelyexpects that the shorter are distances more the mean-fieldhypothesis of neglecting long-range correlations becomesaccurate. One interesting feature is that the approxima-tion given by Eq. (3) in the PQMF theory discards thepossibility of triangles [27], in which the neighbors i and l of j are also connected. So, one could expect a worseperformance in networks with high clustering coefficientbut no statistical correlation with this metrics was found( p PQMF = 0 . N -3 -2 -1 ρ ∗ simulationsQMFPQMF 10 N (a) (b) FIG. 5. Finite size scaling of the steady-state density evaluatedat λ = 2 λ PQMFc . WPCM networks with degree exponent γ = 2 . α = − .
2) and (b)neutral ( α = 0) degree correlations are considered. Insets showthe corresponding IPR calculated for PEV of the correspondingmean-field theory. IV. CONCLUSIONS
Theoretical understanding of dynamical processes onnetworks constitutes a powerful tool for protectionagainst threats such as disease dissemination, misinforma-tion propagations, transportation infrastructure overload,among many other examples. Reliable theoretical approxi-mations are usually required to consider the heterogeneousstructure of the contact networks and the dynamical cor-relations, in which the states of neighboring individualsare statistically correlated. These features are explicitlyincluded in the PQMF theory [27]. However, this the-ory has been applied mainly to analyze the behavior ofepidemic processes in the neighborhood of the transitionfrom an endemic to a disease-free state trough perturba-tive analyses where the epidemic prevalence is very small.In this work, we contribute to fill this gap performing adetailed nonperturbative numerical analysis of the SISmodel on synthetic and real networks within a wide rangeof heterogeneities and assortativities.For synthetic networks generated with the Weber-Porto [37] configuration model, we report that the PQMFtheory predicts with great accuracy the regime of highprevalence observed in stochastic simulations in networkswith power-law degree distributions for all values of de-gree exponents investigated ( γ = 2 .
3, 2.8, and 3.5) anddegree correlations (disassortative, neutral and assorta-tive). In the case of large γ > .
5, where hubs tend tobe separated apart as the network size increases, we ob-served that the PQMF theory significantly outperformsthe simpler QMF theory where heterogeneity is fully con-sidered but dynamical correlations are neglected, beingthe discrepancy between theories larger for large γ . Thehigh accuracy of the PQMF theory at high prevalencescontrasts with its bad performance for asymptotically low densities where the theory is known to deviate fromexactly known critical behavior [38] where ρ ∼ λ β with β > β MF = 1 [28]. Weargue, however, that this mismatch is constrained to aregion very close to λ = λ c → + such that the regime ofnot too low density can still be accurately described bythe PQMF theory.In a set of real networks, where much more complexstructures and correlations can be present, we observedthat PQMF always outperforms (sometimes very signif-icantly) the QMF theory, but may still presents nonnegligible deviations from simulations in some cases; seeTable I. Differently from the low prevalence regime wherethe accuracy of mean-field theories is correlated withspectral properties of Jacobian matrices, only trivial sta-tistical correlations with simple network metrics could beidentified and the problem of predicting when the non-perturbative analysis is sufficiently accurate given certainnetwork properties remains open.Finally, we expect that our work will stimulate theapplication of nonperturbative approaches through thenumerical integration of continuous-time equations to ad-dress other fundamental problems of dynamical processeson networks. Appendix A: Weber-Porto configuration model
The WPCM networks are generated as follows. Thedegree of each vertex is drawn according to the degreedistribution P ( k ) such that each node has k unconnectedhalf-edges. Two half-edges are chosen and connected withprobability P link ( q (cid:48) , q ) = f ( q (cid:48) , q ) f max , (A1)where q and q (cid:48) are the respective degrees of the chosenvertices and f max is the maximum value of f ( q, q (cid:48) ) = 1 + ( κ nn ( q ) − (cid:104) k (cid:105) e )( κ nn ( q (cid:48) ) − (cid:104) k (cid:105) e ) (cid:104) kκ nn (cid:105) e − (cid:104) k (cid:105) , (A2)where (cid:104) A ( k ) (cid:105) e = (cid:80) k A ( k ) P e ( k ) where P e ( k ) = kP ( k ) / (cid:104) k (cid:105) is the probability that an edge ends on a vertex of degree k . Self and multiple connections are forbidden. In theabsence of degree correlations, we have κ nn = (cid:104) k (cid:105) e , im-plying f ( q, q (cid:48) ) = 1 and P link = 1. See Ref. [37] for moredetails. Appendix B: Quasi-stationary method
We applied the standard QS method [42, 43, 56] wherethe dynamics returns to a previously visited active config-uration with at least one infected vertex every time thesystem falls into the absorbing state where all vertices aresimultaneously susceptible. The method is implementedby building and constantly updating a list with M = 100active configurations. Every time the systems falls intothe absorbing state one of the M configurations is chosenwith equal chance to replace the absorbing state. The listis updated with probability 10 − by unit of time and theupdate consists of replacing a randomly selected config-uration of the list by the present state of the dynamics.The QS averages are computed during an averaging timevarying from t av = 10 µ − to 10 µ − after a relaxationtime t rlx = 10 , the longer times for the lower densitieswhere fluctuations are more relevant. Appendix C: Stochastic simulation of the SIS model
Simulations of SIS model were performed using the op-timized Gillespie algorithm described in [42]. The numberof infected vertices N inf and the total number of edges em-anating from them N SI are computed and kept updatedalong the simulations. In each time step, with probability q = µN inf µN inf + λN SI (C1)one infected vertex is chosen with equal chance and healed.With the complementary probability1 − q = λN SI µN inf + λN SI , (C2) one infected vertex i is chosen with probability propor-tional to its degree. One of nearest-neighbors of i , repre-sented by j , is chosen with equal chance. If j is susceptible,it becomes infected and, otherwise, no change of state isimplemented. The time is incremented by δt = − ln uµN inf + λN SI (C3)where u is a pseudo random number uniformly distributedin the interval (0 , ACKNOWLEDGMENTS
This work was partially supported by the Brazilianagencies CNPq and FAPEMIG. This study was financedin part by the Coordena¸c˜ao de Aperfei¸coamento de Pessoalde N´ıvel Superior - Brasil (CAPES) - Finance Code 001. [1] I. Kiss, J. Miller, and P. Simon,
Mathematics of Epi-demics on Networks: From Exact to Approximate Models (Springer International Publishing, Cham, Switzerland,2017).[2] W. Wang, M. Tang, H. Eugene Stanley, and L. A. Braun-stein, “Unification of theoretical approaches for epidemicspreading on complex networks,” Reports Prog. Phys. ,036603 (2017).[3] G. F. de Arruda, F. A. Rodrigues, and Y. Moreno, “Fun-damentals of spreading processes in single and multilayercomplex networks,” Phys. Rep. , 1 (2018).[4] J. Marro and R. Dickman, Nonequilibrium Phase Transi-tions in Lattice Models , Al´ea-Saclay (Cambridge Univer-sity Press, 2005).[5] D. Ben-Avraham and J. K¨ohler, “Mean-field ( n , m )-cluster approximation for lattice models,” Phys. Rev. A , 8358 (1992).[6] J. Joo and J. L. Lebowitz, “Pair approximation of thestochastic susceptible-infected-recovered-susceptible epi-demic model on the hypercubic lattice,” Phys. Rev. E ,036114 (2004).[7] J. P. Gleeson, “High-accuracy approximation of binary-state dynamics on networks,” Phys. Rev. Lett. ,068701 (2011).[8] J. P. Gleeson, “Binary-state dynamics on complex net-works: Pair approximation and beyond,” Phys. Rev. X ,021004 (2013).[9] A. Barrat, M. Barthlemy, and A. Vespignani, DynamicalProcesses on Complex Networks (Cambridge UniversityPress, New York, NY, USA, 2008). [10] R. Pastor-Satorras, C. Castellano, P. Van Mieghem, andA. Vespignani, “Epidemic processes in complex networks,”Rev. Mod. Phys. , 925 (2015).[11] R. Pastor-Satorras and A. Vespignani, “Epidemic Spread-ing in Scale-Free Networks,” Phys. Rev. Lett. , 3200(2001).[12] Yang Wang, D. Chakrabarti, Chenxi Wang, andC. Faloutsos, “Epidemic spreading in real networks: aneigenvalue viewpoint,” in (IEEE Comput. Soc, 2003) pp.25–34.[13] K. T. D. Eames and M. J. Keeling, “Modeling dynamicand network heterogeneities in the spread of sexuallytransmitted diseases,” Proc. Natl. Acad. Sci. , 13330(2002).[14] M. J. Keeling, “The effects of local spatial structure onepidemiological invasions,” Proc. R. Soc. London. Ser. BBiol. Sci. , 859 (1999).[15] A. S. Mata, R. S. Ferreira, and S. C. Ferreira, “Hetero-geneous pair-approximation for the contact process oncomplex networks,” New J. Phys. , 053006 (2014).[16] F. A. Rodrigues, T. K. D. Peron, P. Ji, and J. Kurths,“The Kuramoto model in complex networks,” Phys. Rep. , 1 (2016).[17] S. Chatterjee and R. Durrett, “Contact processes on ran-dom graphs with power law degree distributions havecritical value 0,” Ann. Probab. , 2332 (2009).[18] T. Mountford, J.-C. Mourrat, D. Valesin, and Q. Yao,“Exponential extinction time of the contact process on fi-nite graphs,” Stoch. Process. their Appl. , 1974 (2016). [19] C. Castellano and R. Pastor-Satorras, “Thresholds forEpidemic Spreading in Networks,” Phys. Rev. Lett. ,218701 (2010).[20] C. Castellano and R. Pastor-Satorras, “Competing acti-vation mechanisms in epidemics on networks,” Sci. Rep. , 24 (2012).[21] M. Bogu˜n´a, C. Castellano, and R. Pastor-Satorras, “Na-ture of the epidemic threshold for the susceptible-infected-susceptible dynamics in networks,” Phys. Rev. Lett. ,068701 (2013).[22] S. C. Ferreira, R. S. Sander, and R. Pastor-Satorras,“Collective versus hub activation of epidemic phases onnetworks,” Phys. Rev. E , 032314 (2016).[23] W. Cota, A. S. Mata, and S. C. Ferreira, “Robustness andfragility of the susceptible-infected-susceptible epidemicmodels on complex networks,” Phys. Rev. E , 012310(2018).[24] D. Chakrabarti, Y. Wang, C. Wang, J. Leskovec, andC. Faloutsos, “Epidemic thresholds in real networks,”ACM Trans. Inf. Syst. Secur. , 1 (2008).[25] M. Bogu˜n´a and R. Pastor-Satorras, “Epidemic spreadingin correlated complex networks,” Phys. Rev. E , 047104(2002).[26] S. C. Ferreira, C. Castellano, and R. Pastor-Satorras, “Epidemic thresholds of the susceptible-infected-susceptible model on networks: A comparison of numericaland theoretical results,” Phys. Rev. E , 041125 (2012).[27] A. S. Mata and S. C. Ferreira, “Pair quenched mean-fieldtheory for the susceptible-infected-susceptible model oncomplex networks,” EPL (Europhysics Lett. , 48003(2013).[28] D. H. Silva, S. C. Ferreira, W. Cota, R. Pastor-Satorras,and C. Castellano, “Spectral properties and the accu-racy of mean-field approaches for epidemics on correlatedpower-law networks,” Phys. Rev. Res. , 033024 (2019).[29] A. V. Goltsev, S. N. Dorogovtsev, J. G. Oliveira, andJ. F. F. Mendes, “Localization and Spreading of Diseasesin Complex Networks,” Phys. Rev. Lett. , 128702(2012).[30] G. St-Onge, J.-G. Young, E. Laurence, C. Murphy, andL. J. Dub´e, “Phase transition of the susceptible-infected-susceptible dynamics on time-varying configuration modelnetworks,” Phys. Rev. E , 022305 (2018).[31] J. T. Matamalas, A. Arenas, and S. G´omez, “Effectiveapproach to epidemic containment using link equationsin complex networks,” Sci. Adv. , eaau4212 (2018).[32] Q. Wu and T. Hadzibeganovic, “Pair quenched mean-fieldapproach to epidemic spreading in multiplex networks,”Appl. Math. Model. , 244 (2018).[33] Q. Wu, R. Zhou, and T. Hadzibeganovic, “Conditionalquenched mean-field approach for recurrent-state epi-demic dynamics in complex networks,” Phys. A Stat.Mech. its Appl. , 71 (2019).[34] Q. Wu and T. Hadzibeganovic, “An individual-basedmodeling framework for infectious disease spreading inclustered complex networks,” Appl. Math. Model. , 1(2020).[35] S. G´omez, A. Arenas, J. Borge-Holthoefer, S. Meloni,and Y. Moreno, “Discrete-time markov chain approachto contact-based disease spreading in complex networks,”EPL , 38009 (2010). [36] P. G. Fennell, S. Melnik, and J. P. Gleeson, “Limitationsof discrete-time approaches to continuous-time contagiondynamics,” Phys. Rev. E , 052125 (2016).[37] S. Weber and M. Porto, “Generation of arbitrarily two-point-correlated random networks,” Phys. Rev. E ,046111 (2007).[38] T. Mountford, D. Valesin, and Q. Yao, “Metastabledensities for the contact process on power law randomgraphs,” Electron. J. Probab. (2013), 10.1214/EJP.v18-2512.[39] P. Van Mieghem, “The N-intertwined SIS epidemic net-work model,” Computing , 147 (2011).[40] P. Van Mieghem, “The viral conductance of a network,”Comput. Commun. , 1494 (2012).[41] R. S. Ferreira and S. C. Ferreira, “Critical behavior ofthe contact process on small-world networks,” Eur. Phys.J. B , 462 (2013).[42] W. Cota and S. C. Ferreira, “Optimized Gillespie algo-rithms for the simulation of Markovian epidemic processeson large and heterogeneous networks,” Comput. Phys.Commun. , 303 (2017).[43] R. S. Sander, G. S. Costa, and S. C. Ferreira, “Samplingmethods for the quasistationary regime of epidemic pro-cesses on regular and complex networks,” Phys. Rev. E , 042308 (2016).[44] R. Pastor-Satorras, A. V´azquez, and A. Vespignani, “Dy-namical and Correlation Properties of the Internet,” Phys.Rev. Lett. , 258701 (2001).[45] M. Catanzaro, M. Bogu˜n´a, and R. Pastor-Satorras, “Gen-eration of uncorrelated random scale-free networks,” Phys.Rev. E , 027103 (2005).[46] S. N. Dorogovtsev, a. V. Goltsev, and J. F. F. Mendes,“Critical phenomena in complex networks,” Rev. Mod.Phys. , 1275 (2008).[47] A. S. Mata and S. C. Ferreira, “Multiple transitions ofthe susceptible-infected-susceptible epidemic model oncomplex networks,” Phys. Rev. E , 012816 (2015).[48] R. Juh´asz, G. ´Odor, C. Castellano, and M. A. Mu˜noz,“Rare-region effects in the contact process on networks,”Phys. Rev. E , 066125 (2012).[49] M. M. de Oliveira, S. G. Alves, and S. C. Ferreira, “Dy-namical correlations and pairwise theory for the symbioticcontact process on networks,” Phys. Rev. E , 052302(2019).[50] A.-L. Barab´asi and M. P´osfai, Network science (Cam-bridge University Press, Cambridge, 2016).[51] C. Castellano and R. Pastor-Satorras, “Relating Topolog-ical Determinants of Complex Networks to Their SpectralProperties: Structural and Dynamical Effects,” Phys. Rev.X , 041024 (2017).[52] R. Pastor-Satorras and C. Castellano, “Eigenvector Local-ization in Real Networks and Its Implications for EpidemicSpreading,” J. Stat. Phys. , 1110 (2018).[53] F. Radicchi, “Predicting percolation thresholds in net-works,” Phys. Rev. E , 010801 (2015).[54] F. Radicchi and C. Castellano, “Breaking of the site-bondpercolation universality in networks,” Nat. Commun. ,10196 (2015).[55] M. Newman, Networks: An Introduction (Oxford Univer-sity Press, Oxford New York, 2010).[56] M. M. de Oliveira and R. Dickman, “How to simulate thequasistationary state,” Phys. Rev. E71