Role of limiting dispersal on metacommunity stability and persistence

We use the well-known Rosenzweig-MacArthur model Rosenzweig and MacArthur (1963) to represent the prey-predator dynamics in each patch, where prey experiences logistic growth and predator exhibits a Holling Type II functional response. We consider two identical patches with diffusive coupling between them, whose governing equation is represented as

	$\displaystyle\frac{dV_{i}}{dt}$	$\displaystyle=rV_{i}\left(1-\frac{V_{i}}{k}\right)-\frac{aV_{i}}{V_{i}+b}H_{i}+d_{v}(\beta V_{j}-V_{i})\;,$		(1a)
	$\displaystyle\frac{dH_{i}}{dt}$	$\displaystyle=\frac{caV_{i}}{V_{i}+b}H_{i}-mH_{i}+d_{h}(\alpha H_{j}-H_{i})\;,$		(1b)

where $i,j=1,2$ denote the patch number and $i\neq j$. $V_{i}$ and $H_{i}$ are the prey and predator population densities, respectively. The parameter $r$ denotes the intrinsic growth rate, $k$ denotes the carrying capacity, $a$ and $c$ denote the maximum predation rate and predator efficiency, respectively. Half-saturation constant is denoted by $b$ and the mortality rate of the predator by $m$. These parameters determine the local dynamics of the prey and the predator in a given patch. The spatial interactions of the prey and predator between two patches are determined by the prey dispersal rate $d_{v}$, predator dispersal rate $d_{h}$ and the limiting factors $\alpha$ and $\beta$. Dispersal between the patches is controlled by the limiting factors, which can be decreased from $1$ to zero. For $\alpha=\beta=1$, the coupling is the usual diffusive coupling widely employed in the population ecology Goldwyn and Hastings (2008); Holland and Hastings (2008). $\alpha<1$ accounts for the limited predator dispersal by reducing the immigration of the predator from patch $j$ to patch $i$, while $\beta<1$ accounts for the limited prey dispersal by reducing the immigration of the prey from patch $j$ to patch $i$.

Note that the coupled Rosenzweig-MacArthur model with $\alpha=\beta=1$ has been widely employed to understand the intriguing collective dynamical behaviours under various coupling scenarios Dutta and Banerjee (2015); Arumugam et al. (2016); Arumugam and Dutta (2018); Arumugam et al. (2021). In particular, the coexistence of spatially synchronized state and death state, leading to chimeralike state, has been reported in nonlocally coupled Rosenzweig-MacArthur model Dutta and Banerjee (2015), the phenomenon of rhythmogenesis including amplitude and oscillation deaths has been reported under the environmental coupling Arumugam et al. (2016), the effect of trade-off between the mismatch in the timescale of the species, dispersal and external force on the collective dynamics has been reported in diffusively coupled Rosenzweig-MacArthur model Arumugam and Dutta (2018) and the effect of prey-predator dispersal on the metacommunity persistence in both constant and temporally varying environment has been reported in spatially coupled Rosenzweig-MacArthur models Arumugam et al. (2021).

The above two coupled Rosenzweig-MacArthur model with $\alpha$ and $\beta$ limiting predator and prey dispersal rates, respectively, can also be written as a rescaled model with standard diffusive coupling represented as

	$\displaystyle\frac{dV_{i}}{dt}$	$\displaystyle=r_{e}V_{i}\left(1-\frac{V_{i}}{k_{e}}\right)-\frac{aV_{i}}{V_{i}+b}H_{i}+\epsilon_{v}(V_{j}-V_{i})\;,$		(2a)
	$\displaystyle\frac{dH_{i}}{dt}$	$\displaystyle=\frac{caV_{i}}{V_{i}+b}H_{i}-m_{e}H_{i}+\epsilon_{h}(H_{j}-H_{i})\;,$		(2b)

where, $r_{e}=r-d_{v}(1-\beta),k_{e}=\frac{r_{e}}{r}k=1-\frac{d_{v}}{r}(1-\beta),\epsilon_{v}=d_{v}\beta,m_{e}=m+d_{h}(1-\alpha)$ and $\epsilon_{h}=d_{h}\alpha$. Here, $r_{e}$ is the decreased growth rate, $m_{e}$ is the increased mortality, $k_{e}$ is the effective carrying capacity, and $\epsilon_{v}$ and $\epsilon_{h}$ are the rescaled coupling coefficients. Since the range of $\beta$ and $\alpha$ lies within $[1,0]$, $r_{e}\leq r$ and $m_{e}\geq m$ resulting in decreased prey growth rate and increased predator death rate in the rescaled version of coupled Rosenzweig-MacArthur model. Note that the results of both the rescaled model with standard diffusive coupling and the original Rosenweig-MacAurther model with limiting predator and prey dispersal rates in the coupling will be exactly similar to each other. In other words, the exact mapping between these two models can also be interpreted as that the observed results may not a pure manifestation of the topological features. However, we prefer to proceed our analysis without rescaling the original Rosenweig-MacAurther model but only by engineering the coupling strategies even though limiting the predator and prey dispersal rates effectively results in an increase in the predator’s death rate $m_{e}$ and decrease in the prey growth rate $r_{e}$, respectively.

In this section, we will unravel the effect of limiting the predator dispersal on the dynamical transitions of the metacommunity model as a function of spatial parameters (predator and prey dispersal rates). We fix $\beta=1$, so that the preys are allowed to disperse completely, whereas only the predator dispersal is controlled by decreasing the value of $\alpha$. Bifurcation diagrams, plotted using XPPAUT  Ermentrout (2002), illustrating the dynamical transitions as a function of the predator dispersal rate are depicted in Fig. 2 for the carrying capacity $k=0.5$. Extrema of stable homogeneous oscillatory state are plotted as filled circles, while that of unstable oscillatory states are represented by unfilled circles. Stable and unstable steady states are indicated by solid and dashed curves/lines, respectively. For $\alpha=1$, the dispersal among the metacommunity is the usual diffusive coupling between the patches. Only stable homogeneous (synchronized) oscillation among the populations is observed all along the low values of the predatory dispersal for $\alpha=1$, where the homogeneous steady state is unstable (see Fig. 2(a)). The former prevails in the entire explored range of the predator dispersal rate. Subcritical pitchfork bifurcation onsets at $d_{h}=0.48$ resulting in the coexistence of unstable homogeneous and inhomogeneous steady states along with the homogeneous oscillatory state in the range of $d_{h}\in[0.48,0.597)$. A Hopf bifurcation at $d_{h}=0.597$ leads to the onset of stable inhomogeneous steady (asynchronous) states, which coexist along with unstable inhomogeneous oscillatory states. The onset of the former gives rise to the emergence of alternative stable states of the populations thereby leading to an increase in the degree of persistence of the metacommunity in the range of $d_{h}\in[0.597,1.0]$.

Now, we control the dispersal rate of the predators using the limiting factor $\alpha$ and investigate the effect of limited dispersal on the metacommunity persistence. The dynamical transitions of the prey density are depicted in Fig. 2(b) for $\alpha=0.985$. It is evident from the figure that even a feeble decrease in the limiting factor, that is a negligible fraction of decrease in dispersal, results in drastic changes in the metacommunity dynamics. The Hopf bifurcation, resulted in the onset of alternative stable states, now emerges at a higher dispersal rate at $d_{h}=0.7745$ reducing the tolerance of the metacommunity persistence to a narrow range of predator dispersal rate. This effect of limiting dispersal increases to a larger degree for further negligible decrease in the fraction of dispersal. As a consequence, the alternative stable states lose its stability in the explored range of $d_{h}$ even for $\alpha=0.98$, denoting a very small decrease in the degree of the dispersal, endangering the persistence of the metacommunity due to the presence of the homogeneous (synchronous) oscillatory state as the only dynamical state of the metacommunity, which is highly prone to extinction. In the following, we will explore the effect of controlling the predator dispersal on the metacommunity dynamics as a function of the dispersal rate of both prey and predator populations.

We have depicted the spread of the inhomogeneous steady (asynchronous, alternative) states in the two-parameter bifurcation diagram in Fig. 3(a) for different degree of the dispersal by reducing the limiting factor $\alpha$. Homogeneous oscillatory state is stable in the entire parameter space, while the inhomogeneous steady states are stable in the shaded region, where both the former and the latter coexists. The transition from homogeneous to inhomogeneous steady state is via the Hopf bifurcation as depicted in one parameter bifurcation diagrams in Fig. 2. Note that the Hopf bifurcation curves are obtained from XPPAUT. The spread of the inhomogeneous steady states is maximum for $\alpha=1$, which monotonically decreases for further decrease in the degree of the predator dispersal. Eventually, the inhomogeneous steady states destabilized completely in the same parameter space, where it was stable earlier, at a critical value of $\alpha$ resulting in the homogeneous (synchronous) oscillatory state as the only dynamical state of the metacommunity. The normalized area of the spread of the stable inhomogeneous steady states defined as $R=S(\alpha)/S(\alpha=1)$, where $S(\alpha)$ denotes the spread of the stable inhomogeneous steady states for $\alpha$, is depicted in Fig. 3(b) as a function of $\alpha$. It is evident from the figure that the spread of the stable inhomogeneous steady states decreases monotonically and eventually vanishes at the critical value of $\alpha=\alpha_{c}\approx 0.9818$ corroborating the dynamical transitions observed in the one- and two-parameter bifurcation diagrams.

It is extremely difficult to extract the inhomogeneous steady states and deduce its stability condition analytically for finite values of prey dispersal rate $d_{v}$. However, for $d_{v}=0$, the inhomogeneous steady states $(V_{1}^{\ast},H_{1}^{\ast},V_{2}^{\ast},H_{2}^{\ast})$ can be deduced as

	$\displaystyle V_{1}^{\ast}$	$\displaystyle=\frac{q_{1}+\sqrt{q_{3}-q_{1}^{2}}}{q_{2}}\;,$
	$\displaystyle V_{2}^{\ast}$	$\displaystyle=\frac{-\sqrt{q_{4}^{2}-4q_{5}\left(\alpha bkd_{h}-bkd_{h}-\alpha bV_{1}^{\ast}d_{h}-bkm+\alpha kV_{1}^{\ast}d_{h}-\alpha V_{1}^{{}^{\ast}2}d_{h}\right)}+q_{4}}{2q_{5}}\;,$
	$\displaystyle H_{1}^{\ast}$	$\displaystyle=-\frac{\left(b+V_{1}^{\ast}\right)\left(kV_{1}^{\ast}d_{v}-kV_{2}^{\ast}d_{v}-krV_{1}^{\ast}+rV_{1}^{\ast 2}\right)}{akV_{1}^{\ast}}\;,$
	$\displaystyle H_{2}^{\ast}$	$\displaystyle=-\frac{\left(b+V_{2}^{\ast}\right)\left(kV_{2}^{\ast}d_{v}-kV_{1}^{\ast}d_{v}-krV_{2}^{\ast}+rV_{2}^{\ast 2}\right)}{akV_{2}^{\ast}}\;,$

where,

	$\displaystyle q_{1}$	$\displaystyle=\left(d_{h}(2m-2ac)+(m-ac)^{2}-\left(\alpha^{2}-1\right)d_{h}^{2}\right)\left(ack+(\alpha+1)(b-k)d_{h}+m(b-k)\right),$
	$\displaystyle q_{2}$	$\displaystyle=2\left(ac+(\alpha-1)d_{h}-m\right)\left(ac-(\alpha+1)d_{h}-m\right){}^{2},$
	$\displaystyle q_{3}$	$\displaystyle=q_{2}\left(ac-(\alpha+1)d_{h}-m\right)\left(a^{2}c^{2}k^{2}+2d_{h}\left(ack(\alpha b-k)+m\left(b^{2}+k^{2}\right)\right)-2ack^{2}m-\left(\alpha^{2}-1\right)\left(b^{2}+k^{2}\right)d_{h}^{2}+m^{2}\left(b^{2}+k^{2}\right)\right),$
	$\displaystyle q_{4}$	$\displaystyle=-ack-bd_{h}-bm+kd_{h}+km,\qquad\text{and}$
	$\displaystyle q_{5}$	$\displaystyle=\left(-ac+d_{h}+m\right).$

The corresponding Jacobian matrix can be given as

$A=\left(\begin{array}[]{cccc}-d_{v}+r\left(1-\frac{2V_{1}^{\ast}}{k}\right)+s_{1}-s_{2}&-s_{5}&d_{v}&0\\ cs_{2}-cs_{1}&cs_{5}-d_{h}-m&0&\alpha d_{h}\\ d_{v}&0&-d_{v}+r\left(1-\frac{2V_{2}^{\ast}}{k}\right)+s_{3}-s_{4}&-s_{6}\\ 0&\alpha d_{h}&cs_{4}-cs_{3}&cs_{6}-d_{h}-m\end{array}\right),$

where, $s_{1}=\frac{ah_{1}V_{1}^{\ast}}{(b+V_{1}^{\ast})^{2}}$, $s_{2}=\frac{ah_{1}}{b+V_{1}^{\ast}}$, $s_{3}=\frac{ah_{2}V_{2}^{\ast}}{\left(b+V_{2}^{\ast}\right){}^{2}}$, $s_{4}=\frac{ah_{2}}{b+V_{2}^{\ast}}$, $s_{5}=\frac{aV_{1}^{\ast}}{b+V_{1}^{\ast}}$, and $s_{6}=\frac{aV_{2}^{\ast}}{b+V_{2}^{\ast}}$. One can deduce the characteristic equation from the Jacobian as

$\displaystyle\lambda^{4}-Tr(A)\lambda^{3}+a_{2}\lambda^{2}-a_{3}\lambda+Det(A)=0.$

(3)

The condition for the stability of the inhomogeneous steady (alternative) states can be obtained in terms of the system parameters as

	$\displaystyle F_{\alpha}(a,b,c,k,m,r,\alpha,d_{h})$	$\displaystyle\,=a_{3}(Tr(A)a_{2}-a_{3})$
		$\displaystyle\,-Tr(A)^{2}Det(A)>0.$		(4)

The inhomogeneous steady state is stable for $F_{\alpha}>0$ and unstable for $F_{\alpha}<0$. $\hat{F}=F_{\alpha}/F_{1}$ is depicted in Fig. 4 as a function of the limiting factor for three predator dispersal rates $d_{h}=0.8,0.9$ and $1$. The value of the limiting factor $\alpha$ at which $\hat{F}=0$ corresponds to the critical value $\alpha=\alpha_{c}$, where there is a switch in the stability of the inhomogeneous steady state. The Hopf bifurcation curves plotted using XPPAUT in Fig. 3(a) for the critical values of the limiting factor $\alpha_{c}=0.984$ and $0.982$ coincide with the predator dispersal rate ($x$-axis where $d_{v}=0$) nearly at $d_{h}=0.8$ and $0.9$, respectively, corroborating the critical value of the limiting factor, at which there is change in the stability of the inhomogeneous steady state, obtained using the analytical stability curve in Fig. 4.

In this section, we investigate the influence of the limited predator dispersal and a local parameter (carrying capacity) on the dynamics of the metacommunity for $\beta=1$. The predator density is depicted in the bifurcation diagram as a function of the carrying capacity $k$ (see Fig. 5). The prey and predator dispersal rates are fixed as $d_{v}=0.001$ and $d_{h}=1$, respectively. The dynamical states are the same as in Fig. 2. The trivial steady state lose its stability via a transcritical bifurcation at $k=k_{TB}$ leading to a nontrivial homogeneous steady state, which undergoes a supercritical Hopf bifurcation (HB1) resulting in stable homogeneous oscillation and an unstable steady state at $k=0.4$ for $\alpha=1$ (see Fig. 5(a)). Using the standard linear stability analysis around the homogeneous steady state $\left(\frac{kr+kd_{v}(\beta-1)}{r},0,\frac{kr+kd_{v}(\beta-1)}{r},0\right)$, the transcritical bifurcation at $k=k_{TB}$ can be deduced as

$$k_{TB}=\frac{b(d_{h}+m-\alpha d_{h})}{ac+d_{h}(\alpha-1)-m}.$$

Similarly, one can obtain the Hopf bifurcation

$$k_{HB}=\frac{b\left[ac+m+(1-\alpha)d_{h}\right]}{ac+(\alpha-1)d_{h}-m},$$

by performing a linear stability analysis around the homogeneous steady state $(V^{*},H^{*},V^{*},H^{*})$, where $V^{*}=\frac{bm+bd_{h}(1-\alpha)}{ac-m-d_{h}(1-\alpha)}$, $H^{*}=\frac{bc\left[ack(r+d_{v}(\beta-1))-\left(br+kr+kd_{v}(\beta-1)\right)\left(m+d_{h}(1-\alpha)\right)\right]}{k(ac-m+d_{h}(\alpha-1))^{2}}$. Note that for $\beta=\alpha=1$, the Hopf and pitchfork bifurcation points/curves can be reduced to that in Ref. Arumugam et al. (2021).

The stable homogeneous oscillatory state prevails in the entire explored range of $k\in[0.4,2]$. Subcritical Hopf bifurcation (HB2) emerges at $k=0.458$ resulting in unstable inhomogeneous oscillations and stable inhomogeneous steady states in the range of $k\in[0.458,1.043)$ leading to the existence of alternative stable states of the metacommunity. The latter lose its stability via a supercritical Hopf bifurcation (HB3) at $k=1.043$ resulting in stable inhomogeneous oscillations in the range $k\in[1.043,1.226)$ and unstable inhomogeneous steady states in the range $k\in[1.043,2.0]$, leaving behind the synchronized oscillatory states as the only stable state of the metacommunity A saddle-node bifurcation emerges at $k=1.226$. In the following, we will unravel the effect of limited dispersal as a function of the carrying capacity.

The bifurcation diagram elucidating the dynamical transitions of the metacommunity for $\alpha=0.95$ is depicted in Fig. 5(b). The dynamical states and their transitions are almost similar to that discussed in Fig. 5(a) for $\alpha=1.0$. Nevertheless, it is evident from the figure that the spread of the stable inhomogeneous steady (asynchronous alternative states) states increased to a larger extent, thereby increasing the degree of persistence of the metacommunity. However, this tendency, enhancing the spread of stable inhomogeneous steady states, of the degree of the dispersal persists up to a critical value and then the steady states lose its stability eventually resulting in the homogeneous (synchronous) oscillatory state as the only dynamical state as illustrated in Fig. 5(c) for $\alpha=0.9$. The spread of the stable inhomogeneous steady states for different values of $\alpha$ is depicted in Fig. 6 as a function of the carrying capacity $k$. It is also evident from the figure that the degree of the spread of inhomogeneous steady states increases initially for a small decrease in the dispersal and eventually decreases for further decrease in the dispersal beyond a critical value of $\alpha$.

The spread of the homogeneous and inhomogeneous states are shown in Fig. 7 as two parameter bifurcation diagrams in $(k,d_{h})$ parameter space for different degree of the predator dispersal. The homogeneous state prevails in the entire parameter space, where as the stable inhomogeneous steady states prevail only in the shaded regions leading to multistability. It is also evident from the figures that the spread of the stable inhomogeneous steady states increases for a small decrease in the degree of the predator dispersal up to a critical value of $\alpha$ and then it starts decreasing as a function of $\alpha$. The stable inhomogeneous steady states eventually destabilized from the entire parameter space, where it was stable, resulting in the homogeneous states as the only dynamical states of the metacommunity. Thus, the degree of limiting the predator dispersal results in increasing the persistence of the metacommuntiy until a critical value, above which it is endangered by the presence of homogeneous (synchronous) states as the only dynamical states of the metacommunity.

Now, we will investigate the effect of limiting the prey dispersal on the metacommunity dynamics as a function of the spatial parameter. For complete dispersal of the predators, that is for $\alpha=1$, the effect of limiting the prey dispersal on the metacommunity dynamics can be seen only in a narrow range of the parameter space. Hence, we have fixed $\alpha=0.985$, and depicted the two-parameter bifurcation diagram in the $(d_{h},d_{v})$ parameter space, where the effect of limiting the prey dispersal is well pronounced as evident from Fig. 8. Homogeneous oscillatory state is stable in the entire parameter space, while the inhomogeneous steady states are stable only in the shaded region. The transition from homogeneous oscillator state to inhomogeneous steady states is via the transcritical Hopf bifurcation (HB1) curve, which is obtained from XPPAUT. It is to be noted that for $\beta=1$, that is for uncontrolled prey dispersal, the spread of the inhomogeneous steady states is rather limited to smaller values of $d_{v}$. However, decreasing $\beta$, that is limiting the prey dispersal, manifests the stable inhomogeneous steady states in a larger region of the $(d_{h},d_{v})$ parameter space. For instance, the spread of the stable inhomogeneous steady states are depicted in Fig. 8 for $\beta=0.7$ and $0.5$. Thus, limiting the prey dispersal results in a large degree of the metacommunity persistence in contrast to the effect of limiting the predator dispersal. It is also to be noted that for different degrees of local parameter (carrying capacity), limiting the prey dispersal results in an increase in the metacommunity persistence.

To elucidate the generic nature of our results, we consider an arbitrary network consisting of $N=20$ patches, whose governing equation is represented as

	$\displaystyle\frac{dV_{i}}{dt}$	$\displaystyle=rV_{i}\left(1-\frac{V_{i}}{k}\right)-\frac{aV_{i}}{V_{i}+b}H_{i}+\frac{d_{v}}{d_{j}}\sum_{j=1}^{N}g_{ij}\left(\beta V_{j}-V_{i}\right)\;,$		(5a)
	$\displaystyle\frac{dH_{i}}{dt}$	$\displaystyle=\frac{caV_{i}}{V_{i}+b}H_{i}-mH_{i}+\frac{d_{h}}{d_{j}}\sum_{j=1}^{N}g_{ij}\left(\alpha H_{j}-H_{i}\right)\;,$		(5b)

where $i=1,\ldots,N$. $g_{ij}$ encodes the topology of the underlying network. $g_{ij}=g_{ji}=1$, if $i$th and $j$th oscillators are connected, otherwise $g_{ij}=g_{ji}=0$. $d_{j}=\sum_{j=1}^{N}g_{ij}$ corresponds to the degree of the $j$th oscillator. The other parameters are the same as in Eq. (1). We have fixed the values of the parameters as $\beta=1.0,r=0.5$, $a=1$, $c=0.5$, $b=0.16$, $m=0.2,d_{h}=1$ $k=0.5$ and $d_{v}=0.001$. The dynamics of a random (Erdos-Renyi) network is depicted as spatiotemporal and time series plots in Fig. 9 for different degree of the predator dispersal. The network exhibits a ten cluster state, corresponding to ten stable steady states, for $\alpha=1$ as evident from the spatiotemporal and time series plots in Fig. 9(a). Nevertheless, decreasing the degree of dispersal to $\alpha=0.99$ results in increase in the number of clusters (see Fig. 9(b)). Note that the multi-cluster state corresponds to the alternative states of the metacommunity elucidating a high degree of their persistence. However, further decrease in the degree of the dispersal decreases the number of multi-cluster and eventually manifests as a single cluster state, as depicted in Fig. 9(c) for $\alpha=0.9$, at a critical value of $\alpha$ signaling a low degree of the metacommunity persistence.

The number of clusters is depicted as a function of the limiting factor $\alpha$ in Fig. 10(a) for the same values of the parameters as in Fig. 9. The size of the multi-cluster states initially varies for small decrease in the limiting factor and eventually decreases resulting in a single cluster state (homogeneous state) at a critical value of $\alpha=\alpha_{c}=0.907$. Thus, it is evident that the observed results of the two patch metacommunity hold for a random network of metacommunity thereby corroborating the generic nature of the results.

Denoting $\rho_{j}$’s ($j=1,2,\ldots,N$) as the eigenvalues of the connectivity matrix $g_{ij}$, characterizing an arbitrary network, which can be ordered as $1.0=\rho_{1}\geq\rho_{2}\geq\ldots\geq-1/(N-1)\geq\rho_{N}\geq-1.0$ Atay (2006). The smallest bounding eigenvalue of the connectivity matrix $g_{ij}$ corresponding to the random network employed in Figs. 9 and 10(a) is found to be $\rho_{N}=-0.987$. It is known that the role of coupling topology on the stability of the homogeneous steady state is determined by $\rho_{N}$ Zou et al. (2013, 2021), which completely characterizes the effect of the connection topology corresponding to an arbitrary network. We have depicted the critical value of the limiting factor $\alpha_{c}$ as a function of $\rho_{N}$ met in Fig. 10(b), which elucidates that our results are generic to any arbitrary network. Note that the entire range of $\rho_{N}\in[-1,0]$ captures all possible coupling topologies. It is to be noted that limiting the prey dispersal rate always results in multi-cluster states (inhomogeneous steady state) promoting metacommunity persistence and hence the above analysis of calculating $\alpha_{c}$ cannot be applied to this case.

Further, different networks can be constructed for the same probability of rewiring ‘$p$’ from complex networks perspective. Noted that the range of $p\in[0,1]$ covers the entire class of network topologies from regular, small world to random (scale-free) networks Watts and Strogatz (1998). The critical value of $\alpha$ is depicted in Fig. 11 as a function of $p$ with error bars characterizing the effect of various network sizes $N$. Specifically, we have chosen $N=20$ to $500$ in steps of $10$. The mean of all the critical values of the limiting factor corresponding to each $N$ is indicated by the filled circle for each $p$, while their variance is denoted by the error bars. Inset in Fig. 11 clearly depicts the error bars. The critical value of the limiting factor in the entire range of $p$ again corroborates the generic nature of our results.