Analysing Disease Spread on Complex Networks Using Forman-Ricci Curvature

1. Introduction

Network theory has become a pivotal framework in epidemiology, offering profound insight into the patterns and dynamics of disease spread in human populations [1,3,6,11]. The essence of applying network theory in this field lies in its ability to model complex interactions between individuals as nodes connected by edges [12,17], representing potential paths for pathogen transmission [15,16,18]. This approach contrasts markedly with traditional compartmental models such as SIR that often assume homogeneous mixing within populations [2,10,13], thus missing critical structural nuances that can significantly influence infection outcomes.

The integration of network theory into epidemiological modeling can provide a more granular perspective by considering variations in individual contact patterns, which can vary widely on social, spatial, and temporal scales [19]. For example, scale-free networks, characterised by their heavy-tailed degree distributions, highlight the presence of super-spreaders, that is, individuals linked to an exceptionally high number of contacts, thus disproportionately affecting the spread of diseases [5]. In contrast, small-world networks, known for their high clustering coefficient and short average path lengths [24], exemplify how tightly knit community structures can rapidly spread an infection locally while maintaining the potential for rapid global transmission through occasional long-range connections.

Erdős-Rényi networks, which can be generated by connecting nodes randomly with a fixed probability [22], offer a baseline for studying epidemic dynamics in homogeneous populations. These types of networks lack the hierarchical structure of the hubs, resulting in relatively uniform connectivity. As such, they are valuable for exploring epidemic spread in systems where all individuals are equally likely to interact.

On the other hand, power-law cluster networks combine the scale-free property of heavy-tailed degree distributions with high local clustering, simulating the structure of tightly knit communities with influential hubs [23]. This dual property of power-law cluster networks makes them particularly relevant for understanding disease dynamics in close-contact settings, where local clusters can amplify infections, while hubs facilitate rapid spread across broader networks.

By modelling these complex interaction patterns, network theory will allow epidemiologists to predict more accurately which nodes (i.e. individuals or groups) are at the highest risk of infection, identify potential outbreak sources, and evaluate the efficacy of intervention strategies designed to disrupt pathogenic pathways [14].

The Forman-Ricci curvature, adapted from its original mathematical context by Robin Forman, is a topological measure used to assess the "robustness" or "fragility" of interactions within a network [4]. In the setting of network theory, curvature is used to quantify how well connected a network is around a particular node or along an edge. Networks with positive curvature are typically dense and well-connected [21], indicating that nodes mutually reinforce each other’s connectivity, which can lead to rapid spread of disease within these clusters. Conversely, areas of negative curvature suggest bottlenecks or weak links in the network, which can act as potential barriers to disease transmission.

Incorporating curvature into epidemiological models provides a nuanced view of disease dynamics, highlighting potential areas where interventions could be most effective. For instance, targeting nodes or edges with high positive curvature for vaccination or quarantine measures could theoretically disrupt the spread more efficiently than strategies that do not consider network topology.

Despite the promising applications of curvature in network analysis [20], a significant gap exists in the literature regarding its use in epidemiology-specific network models. Both undirected and directed networks are vital for modelling infectious diseases. Undirected networks often represent symmetric interactions, such as household or community contact networks, while directed networks capture inherently directional interactions, such as sexual contacts, vector-borne transmission, or transmission in transportation networks.

However, studies investigating the relationship between curvature metrics, such as Forman-Ricci curvature, and disease dynamics in these networks have not been carried out to the best of our knowledge. In particular, directed networks have not received sufficient attention in epidemiology despite their importance in understanding asymmetric interactions. Existing studies on epidemics in networks largely focus on undirected models and overlook how the directionality of edges might influence the implications of curvature on disease spread [7,8,9].

Thus, empirical studies correlating the Forman-Ricci curvature with real-world disease outcomes have not been conducted. This gap is especially evident in complex networks, where the nature of edge connections can significantly alter the network’s topology and, therefore, its epidemiological behaviour.

In this research, we investigate how the Forman-Ricci curvature influences disease dynamics in networks, which could improve our ability to predict outbreak patterns and design more effective public health interventions.

The rest of the paper is structured as follows: Materials and methods are described and formulated in Section 2. In Section 3, validation and simulation studies were presented, and in Section 4, we conclude with a summary and discussion, highlighting key findings and potential directions for future research.

2. Materials and Methods

2.1. Networks Generation and Modelling

The development of network models for this study involves the construction of networks that capture distinct topological properties, including scale-free, small world, random graph, and power-law cluster networks. Each of these models is used to provide unique insights into the dynamics of infectious disease transmission within networks and their implications for public health interventions.

2.1.1. Scale-Free Networks

Scale-free networks exhibit a power-law degree distribution, where a few highly connected nodes, or hubs, dominate. These networks, generated by the Barabási-Albert (BA) model, use preferential attachment, where new nodes are more likely to connect to existing high-degree nodes. The degree distribution follows:

$$P\left(k\right)\propto k^{-\gamma }$$

(1)

where $\gamma$ typically ranges from 2 to 3. Scale-free networks are robust to random failures but highly vulnerable to targeted attacks on hubs, making them crucial for studying super-spreader events in epidemics.

2.1.2. Small-World Networks

Small-world networks combine high clustering coefficients with short path lengths. The Watts-Strogatz (WS) model generates these networks by rewiring edges in a regular lattice with probability p, introducing shortcuts that significantly reduce path lengths, and is given by:

$$L\left(p\right)\approx {\displaystyle \frac{logN}{log(\langle k\rangle p)}}$$

(2)

where $L\left(p\right)$ is the average path length, N the number of nodes, and $\langle k\rangle$ the average degree. These properties facilitate rapid transmission of diseases, providing insight into epidemic dynamics.

2.1.3. Erdős-Rényi Networks (Random Graphs)

The Erdős-Rényi (ER) model creates random graphs by connecting n nodes with probability p. The degree distribution follows a Poisson distribution:

$$P\left(k\right)={\displaystyle \frac{{\lambda }^k{e}^{-\lambda }}{k!}}$$

(3)

where $\lambda =\langle k\rangle$ is the mean degree. ER networks lack hubs, resulting in homogeneous connectivity. Their simplicity makes them ideal for modelling baseline epidemic dynamics in non-hierarchical networks.

2.1.4. Power-Law Cluster Networks

The Power-Law Cluster (PLC) model builds on the Barabási-Albert model by incorporating a clustering probability p, encouraging triangular connections. This creates networks with a power-law degree distribution and high clustering coefficients, mirroring real-world network structures. The clustering probability that governs triangle formation is given by:

$$P\left(\mathrm{triangle}\right)\propto p$$

(4)

where p balances preferential attachment with local clustering. PLC networks are particularly effective for studying disease spread in close-knit communities.

2.2. Curvature Calculation

2.2.1. Forman-Ricci Curvature for Undirected and Directed Networks

Forman-Ricci curvature extends the concept of Ricci curvature to both undirected and directed networks, mathematically represented as $G=(V,E)$, where V is the set of nodes (vertices) and E is the set of edges connecting the nodes. The curvature provides a measure of robustness and structural significance, which varies depending on the local topology of the graph and the weights assigned to its edges.

For an undirected network, the Forman-Ricci curvature for an edge e connecting nodes $v_1$ and $v_2$ is given by the equation:

$$F\left(e\right)=w_e\left({\displaystyle \frac{w_{v1}}{w_e}}+{\displaystyle \frac{w_{v2}}{w_e}}-\sum _{e_{v1}\sim e}{\displaystyle \frac{w_{v1}}{\sqrt{w_e{w}_{e_{v1}}}}}-\sum _{e_{v2}\sim e}{\displaystyle \frac{w_{v2}}{\sqrt{w_e{w}_{e_{v2}}}}}\right)$$

(5)

where $w_e$ denotes the weight of the edge, and $w_{v1}$ and $w_{v2}$ are the weights of the nodes connected by edge e. This formulation provides insights into network robustness and identifies structurally significant regions [20]. This equation takes into account the neighboring edge weights and node weights to provide a comprehensive evaluation of local connectivity and resilience.

For a directed network, the Forman-Ricci curvature for an edge e from node $v_1$ to node $v_2$ is computed using the equation:

$$F\left(e\right)=w_e\left({\displaystyle \frac{w_{v1}}{w_e}}-\sum _{e_{v1}\sim e}{\displaystyle \frac{w_{v1}}{\sqrt{w_e{w}_{e_{v1}}}}}\right)+w_e\left({\displaystyle \frac{w_{v2}}{w_e}}-\sum _{e_{v2}\sim e}{\displaystyle \frac{w_{v2}}{\sqrt{w_e{w}_{e_{v2}}}}}\right)$$

(6)

where $w_e$ denotes the weight of the edge, and $w_{v1}$ and $w_{v2}$ are the weights of the nodes. This formulation incorporates both node and edge weights and adapts to the asymmetry introduced by edge directions, making it highly relevant for analyzing directed networks [20].

The summation terms account for the edges adjacent to the nodes connected by the edge e. If we denote the set of incoming and outgoing edges for node v by $E_{I,v}$ and $E_{O,v}$, then the In-Forman curvature $F_I\left(v\right)$ and Out-Forman curvature $F_O\left(v\right)$ are given as follows:

$$F_I\left(v\right)=\sum _{e\in E_{I,v}}F\left(e\right)$$

(7)

$$F_O\left(v\right)=\sum _{e\in E_{O,v}}F\left(e\right)$$

(8)

Here, the summations are taken over only the incoming and outgoing edges, respectively. Hence, the total amount of flow through node v is given as:

$$F_{I/O}\left(v\right)=F_I\left(v\right)-F_O\left(v\right)$$

(9)

The formulation for directed networks accounts for the asymmetry introduced by edge directions, making it suitable for networks where directionality plays a significant role.

2.2.2. Algorithm Development

We implemented the above formulae of Forman-Ricci Curvature (FRC) in both undirected and directed networks using Python with the NetworkX library. The computation process involves iterating over all edges in the network, calculating curvature values based on node and edge weights, and accounting for directionality, in the case of directed networks. The algorithms are detailed below.

Algorithm (1) computes the Forman-Ricci curvature for undirected networks by iterating through each edge $e=(u,v)$ in the graph $G=(V,E)$. It calculates the curvature $F\left(e\right)$ based on the edge weights, node weights, and the contributions from adjacent edges. The resulting curvature provides insights into the structural robustness and connectivity properties of the network.

Algorithm (2) computes the Forman-Ricci curvature for directed networks by iterating through each directed edge $e=(u,v)$ in the graph $G=(V,E)$. The calculation incorporates the contributions from incoming edges to node u and outgoing edges from node v. By considering edge weights, node weights, and directionality, the algorithm quantifies the structural robustness and flow characteristics of the network.

Algorithm 1 Forman-Ricci Curvature for Undirected Networks

Require:  Undirected graph $G=(V,E)$ Ensure:  $F_e$ for all $e\in E$ and node curvatures $F_v$ 1: for each edge $(u,v)\in E$ do 2:     $w_e\Leftarrow 1$, $w_u\Leftarrow 1$, $w_v\Leftarrow 1$ 3:     ${\mathrm{sum}}_u\Leftarrow 0$, ${\mathrm{sum}}_v\Leftarrow 0$ 4:     for each neighbor $u^{\prime }$ of u, $u^{\prime }\ne v$ do 5:         $w_{u{u}^{\prime }}\Leftarrow 1$ 6:         ${\mathrm{sum}}_u\Leftarrow {\mathrm{sum}}_u+w_u/\sqrt{w_e·w_{u{u}^{\prime }}}$ 7:     for each neighbor $v^{\prime }$ of v, $v^{\prime }\ne u$ do 8:         $w_{v{v}^{\prime }}\Leftarrow 1$ 9:         ${\mathrm{sum}}_v\Leftarrow {\mathrm{sum}}_v+w_v/\sqrt{w_e·w_{v{v}^{\prime }}}$ 10:     $F_e\Leftarrow w_e·\left({\displaystyle \frac{w_u}{w_e}}+{\displaystyle \frac{w_v}{w_e}}-{\mathrm{sum}}_u-{\mathrm{sum}}_v\right)$ 11:     Assign $F_e$ to edge $(u,v)$ 12: for each node $v\in V$ do 13:     $F_v\Leftarrow$ average of $F_e$ over all edges incident on v

Algorithm 2 FRC for Directed Networks

Require:  Directed graph $G=(V,E)$ Ensure:  $F_e$ for all $e\in E$ and node curvatures $F_v$ 1: for each directed edge $(u,v)\in E$ do 2:     $w_e\Leftarrow 1$, $w_u\Leftarrow 1$, $w_v\Leftarrow 1$ 3:     ${\mathrm{sum}}_u^{\mathrm{in}}\Leftarrow 0$, ${\mathrm{sum}}_v^{\mathrm{out}}\Leftarrow 0$ 4:     for each incoming neighbor $u^{\prime }$ of u with edge $(u^{\prime },u)\in E$, $u^{\prime }\ne v$ do 5:         $w_{u^{\prime }u}\Leftarrow 1$ 6:         ${\mathrm{sum}}_u^{\mathrm{in}}\Leftarrow {\mathrm{sum}}_u^{\mathrm{in}}+{\displaystyle \frac{w_u}{\sqrt{w_e·w_{u^{\prime }u}}}}$ 7:     for each outgoing neighbor $v^{\prime }$ of v with edge $(v,v^{\prime })\in E$, $v^{\prime }\ne u$ do 8:         $w_{v{v}^{\prime }}\Leftarrow 1$ 9:         ${\mathrm{sum}}_v^{\mathrm{out}}\Leftarrow {\mathrm{sum}}_v^{\mathrm{out}}+{\displaystyle \frac{w_v}{\sqrt{w_e·w_{v{v}^{\prime }}}}}$ 10:     $F_e\Leftarrow w_e·\left({\displaystyle \frac{w_u}{w_e}}-{\mathrm{sum}}_u^{\mathrm{in}}\right)+w_e·\left({\displaystyle \frac{w_v}{w_e}}-{\mathrm{sum}}_v^{\mathrm{out}}\right)$ 11:     Assign $F_e$ to directed edge $(u,v)$ 12: for each node $v\in V$ do 13:     $F_v\Leftarrow$ average of $F_e$ over all directed edges incident on v (incoming and outgoing)

2.3. Integration of Disease Spread Models with FRC

Classical mathematical epidemiological models such as the SIR (Susceptible-Infectious-Recovered) model offer a foundational framework for studying disease dynamics under the assumption of a homogeneously mixed population. In this classical setting, every individual is assumed to interact equally with every other individual, leading to the following system of ordinary differential equations:

$$\left\{\begin{array}{c}{\displaystyle \frac{dS}{dt}}=-\beta SI,\hfill \\ {\displaystyle \frac{dI}{dt}}=\beta SI-\gamma I,\hfill \\ {\displaystyle \frac{dR}{dt}}=\gamma I\hfill \end{array}\right.$$

(10)

Here, $S\left(t\right),I\left(t\right),R\left(t\right)$ denote the proportions of susceptible, infected, and recovered individuals at time t; $\beta$ is the transmission rate; and $\gamma$ is the recovery rate.

However, in real-world scenarios, populations are not well-mixed. Individuals interact through structured, often heterogeneous social or physical contact patterns. This motivates the shift from classical compartmental models to network-based SIR models, where the population is modeled as a graph $G=(V,E)$. Each node $i\in V$ represents an individual, and each edge $(i,j)\in E$ represents a potential interaction or contact through which infection can spread. Each node maintains compartmental state variables $S_i\left(t\right),I_i\left(t\right),R_i\left(t\right)$, which represent the local probability or proportion of the node being in the susceptible, infected, or recovered state, respectively. The local conservation condition holds, that is, $S_i\left(t\right)+I_i\left(t\right)+R_i\left(t\right)=1\forall i\in V.$ To enhance realism in simulating transmission dynamics, we integrate Forman-Ricci Curvature (FRC) into the proposed model to capture features such as edge robustness, node connectivity, and network fragility. We define the pairwise transmission rate between nodes i and j as ${\beta }_{ij}:={\beta }_0·exp\left[F\left(e\right)\right]$, where ${\beta }_0$ is the baseline transmission rate, and $F\left(e\right)$ is the Forman-Ricci curvature of the edge $e=(i,j)$. This modification allows the transmission rate ${\beta }_{ij}$ to adapt based on the structural properties of the edge. Edges with higher curvature ($F\left(e\right)>0$) suggest more robust connections, increasing transmission likelihood. Conversely, edges with lower curvature ($F\left(e\right)<0$) indicate structural fragility or weaker connectivity, reducing transmission probability. The dynamics of the proposed system are governed by the following set of differential equations for each node $i\in V$:

$$\begin{array}{cc}\hfill {\displaystyle \frac{d{S}_i}{dt}}& =-\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j\hfill \end{array}$$

(11)

$$\begin{array}{cc}\hfill {\displaystyle \frac{d{I}_i}{dt}}& =\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i\hfill \end{array}$$

(12)

$$\begin{array}{cc}\hfill {\displaystyle \frac{d{R}_i}{dt}}& =\gamma I_i\hfill \end{array}$$

(13)

$\gamma$ is the recovery rate and was assumed uniform across all nodes, captures how infected individuals transition to the recovered state, and $N\left(i\right)$ is the neighbourhood of node i, defined as $N\left(i\right)=\{j\mid (i,j)\in E\}.$ In the case of directed networks, we refine the notion of neighbourhood into In-neighbourhood $\left(N_{\mathrm{in}}\left(i\right)=\{j\mid (j,i)\in E\}\right)$, and Out-neighbourhood $\left(N_{\mathrm{out}}\left(i\right)=\{j\mid (i,j)\in E\}\right)$. This distinction allows the model to capture directionality in contact, e.g., where the disease may be transmitted from j to i, but not vice versa. The neighbourhood $N\left(i\right)$ defines the local connectivity of each node i, determining potential sources and targets of infection. This modified SIR model, enriched with curvature-dependent dynamics, provides a refined tool for studying epidemic propagation in realistic and structured populations. This bridges geometric and epidemiological perspectives by linking transmission dynamics to the underlying network geometry, making it particularly useful for modeling diseases in complex systems with heterogeneous interaction patterns. We present the detailed mathematical properties and analysis of the proposed model in Appendix A.

3. Validation and Simulation Studies

Studies were carried out to validate the computed curvature values by comparing the curvature distributions with theoretical expectations in well-understood network configurations. Furthermore, the computed values were benchmarked against established network measures, ensuring consistency and robustness in capturing the structural properties of both undirected and directed networks. These studies provide a solid foundation for understanding how curvature metrics impact disease propagation and offer insights for optimising public health interventions.

In addition, we conducted computational simulations to observe how diseases spread through networks with varying curvature metrics. The analysis focused on how changes in curvature influence disease dynamics, such as the size and speed of outbreaks, in both undirected and directed networks. These results were then compared with those derived from the traditional SIR model to explore the added predictive power of incorporating curvature-based measures. With this approach, we highlight the potential of curvature metrics to understand and mitigate the spread of infectious diseases in networked systems.

3.1. FRC Distributions in Undirected and Directed Networks

Figure 1 presents the FRC distributions for four undirected network models: Erdős–Rényi (1000 nodes, edge probability $p=0.003$), Watts–Strogatz (1000 nodes, mean degree $k=10$, rewiring probability $\beta =0.1$), Barabási–Albert (1000 nodes, $m=3$ edges per new node), and Power-law Cluster (1000 nodes, $m=5$, and clustering probability $p=0.05$).

These models demonstrate distinct structural characteristics reflected in their FRC distributions. The Erdős–Rényi model exhibits a symmetric distribution, reflecting its randomness, while the Watts–Strogatz model shows tighter clustering indicative of small-world properties. With their scale-free and clustered structures, the Barabási–Albert and Power-law Cluster models display highly skewed FRC distributions dominated by negative curvature values, highlighting the prominence of hubs and clustering. These differences suggest that the network topology, captured through FRC, can profoundly influence disease dynamics, with implications for predicting epidemic patterns and tailoring public health interventions to effectively disrupt transmission.

Figure 2 shows the distribution of in-degree and out-degree FRC for directed networks generated using four network models: the Erdős–Rényi model (1000 nodes, edge creation probability $p=0.003$), the Scale-Free model (1000 nodes, power-law degree exponents ${\lambda }_{in}={\lambda }_{out}=2.1$), the Watts-Strogatz model (1000 nodes, each connected to $k=10$ nearest neighbors, rewiring probability $p=0.1$), and the Power-Law Cluster model (1000 nodes, attachment parameter $m=3$, triangle formation probability $p=0.05$).

The Erdős–Rényi model shows a relatively uniform FRC distribution with clustering near zero, reflecting its random and homogenous connectivity. The Scale-Free model displays an extreme skew in curvature values, with highly negative curvatures driven by a few dominant hubs and numerous sparsely connected nodes, characteristic of its power-law degree distribution. The Watts-Strogatz model demonstrates moderate variability in curvature values, highlighting a balance between local clustering and global randomness typical of its small world topology. The Power-Law Cluster model exhibits a distinctive curvature distribution with sharply negative values for both in-degree and out-degree, emphasising its hierarchical structure and strong transitivity.

3.2. Analysis of the Correlation Between FRC and Measures of Centrality

Figure 3 presents the relationship between FRC and three centrality measures, that is, node degree, clustering coefficient, and betweenness centrality, in undirected networks. Figure 3(a) reveals an inverse correlation, where nodes with higher curvature values tend to have lower degrees, indicating that less connected nodes often exhibit higher curvature. Figure 3(b) shows a weak positive relationship, suggesting that nodes involved in tightly knit communities exhibit slightly higher curvature. Lastly, Figure 3(c) highlights that nodes with higher betweenness centrality, critical for network flow, tend to have lower curvature values. These findings emphasise the nuanced role of FRC in capturing structural and functional properties of undirected networks, with potential implications for understanding network resilience and disease dynamics.

Figure 4 illustrates the relationship between FRC and node degree (in-degree and out-degree) in directed network models. The plots reveal consistent negative correlations between out-Forman and out-degree across all models, with the Scale-Free model showing the strongest correlation ($r=-0.90$), followed by Power-law Cluster ($r=-0.88$), Watts–Strogatz ($r=-0.55$) and Erdős–Rényi ($r=-0.52$). These findings indicate that nodes with high out-degrees, typically potential super-spreaders, tend to exhibit highly negative curvature, suggesting structural fragility or influence within the network. In contrast, correlations between in-Forman and in-degree are generally weak or negligible, reflecting directional asymmetry in network dynamics. This analysis reinforces the role that the FRC can play as a structural marker for identifying key nodes that facilitate transmission and thus are vital for targeted epidemiological interventions.

3.3. SIR Dynamics with Undirected Forman-Ricci Curvature

Figure 5,Figure 6,Figure 7 present SIR model simulations incorporating undirected FRC for four different undirected network models, with the same configurations as in Section 3.1, with initial conditions of $90\%$ susceptible, $10\%$ infected, and $0\%$ recovered, under varying transmission rates ($\beta =0.1,0.3,0.5,0.7$). For each model, two sets of subplots are shown: the first set illustrates the complete epidemic progression across the susceptible, infected, and recovered compartments, while the second highlights only the infected population over time to emphasize peak intensity and timing.

In the Erds-Rényi network, a randomly connected structure, we observe that as $\beta$ increases, infections peak earlier and affect more individuals. At low transmission ( $\beta =0.1$), the outbreak remains minimal and prolonged, while higher values such as $\beta =0.7$ induce rapid, acute epidemics. The inclusion of FRC emphasizes how weak or sparsely connected areas structurally accelerate disease spread as transmissibility increases. This showcases the influence of network homogeneity and curvature on the magnitude of the epidemic.

In Watts–Strogatz network, which exhibits small-world topology and high clustering, we see similar acceleration in epidemic spread with increasing $\beta$ as that of an Erdős-Rényi network, but maintains delayed peaks and smoother transitions. Curvature-based adjustments to transmission rates help capture the interplay between local clusters and occasional long-range connections. This highlights the importance of structural awareness when designing interventions in environments with both dense neighbourhoods and shortcut links.

Additionally, the Barabási–Albert (Scale-Free) network, characterised by a few highly connected hubs, the epidemic spreads rapidly even at lower $\beta$. The FRC-informed model captures these vulnerabilities, as peak infection occurs early and resolves quickly due to the influence of super-spreader nodes. Despite high connectivity, the infection stabilises fast, emphasising how negative curvature around hubs can serve as early diffusion accelerators while structural saturation quickly exhausts susceptible paths.

Finally, the Power-Law Cluster network demonstrates moderate epidemic escalation with increasing $\beta$. We observed that at lower transmission, infection is well contained due to strong clustering and local loops. However, at higher $\beta$, the presence of hubs combined with clustered neighbourhoods facilitates quicker but more contained outbreaks. FRC effectively captures these hybrid dynamics, showing fast initial spread through hub nodes followed by localised diffusion within tight communities, suggesting strategies that combine hub targeting and cluster-focused containment.

These simulations reveal that integrating FRC into epidemic models enhances network sensitivity to structural nuances. Negative curvature highlights transmission-prone zones, while clustered or bridge regions indicate strategic control points. Such curvature-informed insights support more precise and adaptive public health interventions tailored to network topology.

3.4. SIR Dynamics with Directed Forman-Ricci Curvature

Figure 8,Figure 9,Figure 10,Figure 11 present SIR dynamics under directed FRC for directed network models with the same configurations for directed networks covered in Section 3.1, with initial conditions as those of undirected networks. In the simulations, we explore how infection spreads under varying infection rates, highlighting the structural influence of directionality on epidemic trajectories.

In the simulation of the directed Erdős–Rényi Model, as captured in Figure 8, the infection curve rises gradually and decays steadily, with peak magnitude increasing proportionally with ${\beta }_0$. Directionality introduces slight asymmetries but does not significantly alter the spread compared to the undirected case. However, recovery occurs marginally faster, indicating that directional constraints reduce long transmission loops.

For directed Watts–Strogatz Model, in Figure 9, the infection peaks occur more rapidly and decay faster, especially at higher ${\beta }_0$. The directed small-world topology facilitates short-range, clustered transmission while limiting extensive spread due to one-way paths. Compared to the undirected model, the tail of the infection curve is steeper, reflecting quicker containment.

In the setting of the directed Barabási-Albert model (scale-free) (Figure 10), highly connected nodes (hubs) play a dominant role in transmission. Directionality affects whether hubs act as spreaders or sinks. Although infection increases sharply with increasing ${\beta }_0$, it decreases faster than in the undirected counterpart, likely due to the directional isolation of the peripheral nodes.

Finally, for the directed Power-law Cluster Model, captured in Figure 11, the model shows the most pronounced epidemic surges, particularly at ${\beta }_0\ge 0.5$, with sharp infection spikes followed by rapid resolution. The directionality intensifies structural bottlenecks, resulting in explosive yet short-lived outbreaks, reflecting the dynamics of clustered hubs infecting local neighborhoods efficiently.

Undirected networks, which allow for bidirectional transmission, tend to produce more prolonged epidemics with flatter infection curves and slower recoveries. Directed networks, on the other hand, yield sharper peaks and faster resolutions. Forman-Ricci curvature adjusts transmission likelihood in both types, but in directed settings it becomes sensitive to edge orientation, thus reshaping transmission pathways and influencing bottlenecks.

Incorporating FRC into both directed and undirected networks offers a nuanced tool for modeling real-world disease transmission. Directional networks represent asymmetric transmission scenarios (e.g., sexual, airborne, vector-borne diseases), where certain nodes may serve as unidirectional superspreaders. Strong negative curvature combined with directionality marks critical transmission corridors. This highlights the importance of targeted interventions that account for both topology and edge orientation, allowing more effective containment and prevention strategies.

3.5. Dynamics of Traditional SIR Model

Figure 12 illustrates the classical SIR model dynamics under varying transmission rates. The top figure shows the progression of each compartment over time, while the bottom figure isolates the infected class to highlight its peak and duration.

We observed that at lower transmission rates (${\beta }_0=0.1$), the spread of infection is slow and contained, with only a modest fraction of the population becoming infected before the disease fades. As ${\beta }_0$ increases, the infection propagates faster and affects a larger portion of the population. At ${\beta }_0=0.3$, a pronounced peak emerges, indicating a rapid surge in infections followed by a decline as recovery progresses. This trend becomes more dramatic at ${\beta }_0=0.5$ and ${\beta }_0=0.7$, where the epidemic spreads swiftly and infects a significant majority of the population before recovery dominates.

The infected curve in the second plot reveals that higher values of ${\beta }_0$ lead to earlier and higher infection peaks, suggesting greater epidemic severity. This behavior is consistent with classical epidemiological theory, which states that higher transmission rates result in faster and more extensive disease spread.

While the classical SIR model offers foundational insights, it lacks structural features like network topology and curvature that can affect transmission pathways, limiting its ability to capture heterogeneous interactions observed in real populations.

4. Conclusions

This study presented a novel approach to modelling disease spread by integrating Forman-Ricci Curvature (FRC) into network-based SIR dynamics, evaluated across both undirected and directed complex networks. Through extensive simulations, we demonstrated that FRC serves as a meaningful structural descriptor that captures both local connectivity and global transmission potential across diverse network topologies.

Our results show that networks with high clustering and hierarchical structures, such as Power-Law Cluster and Barabási–Albert models, display skewed FRC distributions dominated by highly negative curvature, especially around hub nodes. These negatively curved regions correlate strongly with high-degree nodes, revealing structurally vulnerable zones where disease can spread rapidly. SIR simulations further confirmed that incorporating FRC into transmission dynamics can lead to more accurate predictions of epidemic progression, including peak timing and magnitude. Notably, directed networks displayed sharper epidemic curves with faster recoveries, highlighting the role of edge directionality in modulating disease pathways.

Correlation analysis between FRC and classical centrality measures, such as degree, clustering coefficient, and betweenness centrality, revealed that FRC provides a complementary perspective. Particularly, strong negative correlations between FRC and out-degree centrality in directed networks emphasise FRC’s potential in identifying super-spreader nodes and critical transmission corridors.

The curvature-adjusted SIR model successfully captured heterogeneities in the network structure and improved the predictive capacity over the traditional SIR framework, especially in non-homogeneous populations. These findings support the use of curvature as a strategic variable in epidemiological modelling, offering actionable insights for designing targeted interventions, such as hub vaccination or localised containment.

Public health strategies should account for network topology, and especially curvature-informed metrics, when designing control interventions. This is particularly crucial for diseases transmitted in structured populations, such as schools, workplaces, or transportation systems.

Building on this foundation, future research will explore the integration of FRC in real-world epidemiological datasets. We will also extend the framework to include temporal networks, mobility-driven dynamics, and integration with Graph Neural Networks (GNNs) to enable data-driven learning of transmission structures. Validating the model using empirical outbreak data will further establish its applicability in real-time disease surveillance and intervention planning.

Appendix A.1. Assumptions of the Modified Network-Based SIR Model

The following are the assumptions for the modified network-based SIR model, similar to those in traditional compartmental models but incorporating network-specific dynamics:

• The total population at each node $N_i=S_i+I_i+R_i$ remains constant. The model assumes a closed population where no individuals enter or leave, and $N_i$ does not vary over time.
• In both undirected and directed networks, the interactions within a local neighbourhood $N\left(i\right)$ are influenced by the weights of nodes and edges. For undirected networks, this means that susceptible individuals within the local neighbourhood are more likely to interact with infectious individuals based on the weighted strength of their connections. For directed networks, interaction likelihood depends on the direction of edges and their weights. Therefore, not all susceptible individuals in the local neighbourhood have equal interaction probabilities with infectious individuals. Instead, interaction probability is proportional to the edge weight, ${\beta }_{ij}$ and the edge direction in the case of directed edges.
• The network structure remains fixed throughout the duration of the epidemic. Edges, which represent potential contacts between individuals at different nodes, do not change over time. Nodes and edges are weighted, but there are no dynamic rewiring processes or edge weight adjustments over time.
• Disease transmission occurs when a susceptible node $S_i$ is connected to an infectious node $I_j$ through an edge in the network. The transmission rate ${\beta }_{ij}$ between nodes depends on the weights of the edges and nodes, which is adjusted by the FRC to account for the local network structure and geometry.
• Infected individuals at each node recover at a constant rate $\gamma$. Once an individual recovers, they cannot become susceptible again, implying permanent immunity. The recovery rate $\gamma$ is not influenced by edge weights or network structure.
• The baseline transmission rate ${\beta }_0$ and the recovery rate $\gamma$ are assumed to be constant over time. While ${\beta }_{ij}$ between nodes is adjusted using FRC and can vary depending on local network properties, the baseline transmission rate and recovery rate themselves are fixed.
• The model assumes instantaneous transmission and recovery processes. When a susceptible individual comes into contact with an infectious individual, infection happens immediately at a given rate, and recovery from infection also occurs without delay at a fixed rate.
• The model does not consider natural births or deaths in the population. Individuals are only classified as susceptible, infected, or recovered; these states evolve without demographic changes.
• The contact rate between individuals is influenced by the network structure using FRC. High-curvature regions represent tightly-knit areas where individuals have a greater chance of interacting, leading to potentially higher transmission rates. Low-curvature or negative-curvature regions indicate bottlenecks or sparse connections, leading to lower transmission rates.
• Individuals are fixed at nodes within the network. There is no movement of individuals between nodes, meaning that transmission is restricted to direct network connections as defined by edges. The model does not include the spatial mobility of individuals across the network.
• The model does not account for external interventions such as vaccinations, quarantine, or social distancing measures. The spread of the disease is entirely driven by interactions inherent to the network’s structure, as well as the transmission and recovery processes.

Appendix A.2. Feasibility of the Model

The feasibility of the model ensures that the solutions $S_i\left(t\right)$, $I_i\left(t\right)$, and $R_i\left(t\right)$ remain biologically meaningful for all nodes i and time $t\ge 0$. This involves verifying the model’s non-negativity, conservation, and boundedness properties.

Appendix A.2.1. Non-Negativity of Solutions

Recall that the equations that governed the dynamics of the system are given by equations 11 - 13. To ensure $S_i\left(t\right),I_i\left(t\right),R_i\left(t\right)\ge 0$, we observe that; from equation 11, the change in $S_i$ depends on $-{\beta }_{ij}{S}_i{I}_j$, which is non-positive for all $S_i,I_j\ge 0$. Thus, $S_i\left(t\right)$ decreases or remains constant. Also, from equation 12, $I_i$ is increased by ${\sum }_{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j$, a non-negative term, and decreased by $\gamma I_i$, which is proportional to $I_i$. This guarantees $I_i\left(t\right)\ge 0$. And finally, from equation 13, $R_i\left(t\right)$ is increased by $\gamma I_i$, a non-negative term, ensuring $R_i\left(t\right)\ge 0$. Thus, the state variables $S_i\left(t\right),I_i\left(t\right),R_i\left(t\right)$ remain non-negative for all $t\ge 0$.

Appendix A.2.2. Conservation of Population

The conservation equation is given by:

$$S_i\left(t\right)+I_i\left(t\right)+R_i\left(t\right)=1,\forall i,t\ge 0,$$

ensures that the total proportion of the population at each node remains constant over time. This is derived as follows: Summing the governing equations for $S_i$, $I_i$, and $R_i$, we have that:

$${\displaystyle \frac{d}{dt}}(S_i+I_i+R_i)={\displaystyle \frac{d{S}_i}{dt}}+{\displaystyle \frac{d{I}_i}{dt}}+{\displaystyle \frac{d{R}_i}{dt}}.$$

Now, substituting the equations, we obtained:

$${\displaystyle \frac{d{S}_i}{dt}}+{\displaystyle \frac{d{I}_i}{dt}}+{\displaystyle \frac{d{R}_i}{dt}}=-\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j+\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i+\gamma I_i.$$

Simplifying, we get:

$${\displaystyle \frac{d}{dt}}(S_i+I_i+R_i)=0.$$

Thus, the total population proportion at each node remains conserved.

Appendix A.2.3. Boundedness of Solutions

For the boundedness of solutions to occur:

$$0\le S_i\left(t\right),I_i\left(t\right),R_i\left(t\right)\le 1,\forall i,t\ge 0.$$

This follows from the conservation of the total population at each node:

$$S_i\left(t\right)+I_i\left(t\right)+R_i\left(t\right)=1,\forall i,t\ge 0.$$

Now, the term $S_i\left(t\right)$ can only decrease due to infection ($-{\sum }_{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j$). Also, the term $I_i\left(t\right)$ increases due to infection and decreases due to recovery ($-\gamma I_i$). Additionally, the term $R_i\left(t\right)$ can only increase due to recovery ($\gamma I_i$).

Since $S_i\left(t\right)$, $I_i\left(t\right)$, and $R_i\left(t\right)$ are interdependent through the conservation equation, their sum is constant and bounded within $[0,1]$.

Appendix A.2.4. Feasible Domain

The solutions $S_i\left(t\right),I_i\left(t\right),R_i\left(t\right)$ are constrained to the feasible domain:

$$\mathcal{D}=\{(S_i,I_i,R_i)\mid S_i,I_i,R_i\ge 0,S_i+I_i+R_i=1\}.$$

The system dynamics ensure that the solutions remain within $\mathcal{D}$ for all $t\ge 0$, as non-negativity ensures $S_i,I_i,R_i\ge 0$, and conservation ensures $S_i+I_i+R_i=1$.

The boundedness, non-negativity, and conservation properties of the modified SIR model ensure that the solutions $S_i\left(t\right),I_i\left(t\right),R_i\left(t\right)$ are biologically meaningful and the solutions remain within the feasible domain $\mathcal{D}$ for all $t\ge 0$. Hence, the model is well-posed and provides realistic predictions of disease dynamics.

Appendix A.3. Existence and Uniqueness of Solutions

The existence and uniqueness of solutions for the modified SIR model was established using the Picard-Lindelöf theorem. This theorem provides conditions under which a system of ordinary differential equations (ODEs) admits a unique solution for given initial conditions. Specifically, the theorem requires the vector field defined by the ODEs to be both Lipschitz continuous and locally bounded [25].

Recall that the modified SIR model for each node i is represented by the system of equations 11, 12and 13. Hence, the vector field defining the system is given by:

$$\mathbf{F}=(f_{S_i},f_{I_i},f_{R_i}),$$

where $f_{S_i},f_{I_i},f_{R_i}$ are the right-hand sides of the equations. The initial conditions for the system are:

$$S_i\left(0\right)=S_{i,0},I_i\left(0\right)=I_{i,0},R_i\left(0\right)=R_{i,0},\forall i,$$

with $S_{i,0}+I_{i,0}+R_{i,0}=1$ and $S_{i,0},I_{i,0},R_{i,0}\ge 0$.

Appendix A.3.1. Lipschitz Continuity

To verify Lipschitz continuity, we consider each function $f_{S_i},f_{I_i},f_{R_i}$ separately, that is, $f_{S_i}$ is a linear combination of terms of the form $-{\beta }_{ij}{S}_i{I}_j$. $f_{I_i}$ is a combination of terms of the form ${\beta }_{ij}{S}_i{I}_j$ and $-\gamma I_i$ and $f_{R_i}$ is a linear term $\gamma I_i$. Each of these functions is a polynomial in $S_i,I_i,R_i$ with coefficients dependent on the bounded transmission rate ${\beta }_{ij}$ and recovery rate $\gamma$. Since polynomials are differentiable and their derivatives are bounded over the domain $[0,1]$, the vector field $\mathbf{F}$ satisfies the Lipschitz condition:

$$\parallel \mathbf{F}\left(\mathbf{X}\right)-\mathbf{F}\left(\mathbf{Y}\right)\parallel \le {L\parallel \mathbf{X}-\mathbf{Y}\parallel ,\forall \mathbf{X},\mathbf{Y}\in [0,1]}^n,$$

where L is a Lipschitz constant, and $\mathbf{X},\mathbf{Y}$ are state vectors for the system.

Appendix A.3.2. Lemma (Local Boundedness of the Vector Field)

The vector field $\mathbf{F}=(f_{S_i},f_{I_i},f_{R_i})$ governing the modified SIR model is locally bounded, i.e., for any bounded domain $[0,1]$, there exists a finite constant $M>0$ such that:

$$\parallel \mathbf{F}\left(\mathbf{X}\right)\parallel \le M,\forall \mathbf{X}\in {[0,1]}^n,$$

where $\mathbf{X}=(S_i,I_i,R_i)$ are the state variables at each node i.

Proof.

• Boundedness of Variables: The variables $S_i,I_i,R_i\in [0,1]$ are proportions of the population at node i, representing susceptible, infectious, and recovered states, respectively. By definition, these variables satisfy:

  $$S_i+I_i+R_i=1,\forall i\mathrm{and}t\ge 0,$$

  ensuring their boundedness within $[0,1]$.
• Boundedness of Transmission Rate: The transmission rate ${\beta }_{ij}$ is defined as:

  $${\beta }_{ij}={\beta }_0·exp\left[F\left(e\right)\right],$$

  where ${\beta }_0>0$ is a constant baseline transmission rate and $F\left(e\right)$ is the Forman-Ricci curvature of edge $e=(i,j)$, which is finite for any edge in the network.
  Since both ${\beta }_0$ and $exp\left[F\right(e\left)\right]$ are finite, hence ${\beta }_{ij}$ is bounded. The upper bound of ${\beta }_{ij}$ can be expressed as:

  $$0<{\beta }_{ij}\le {\beta }_{\max }:={\beta }_0·exp\left[F_{\max }\right],$$

  where ${\beta }_0>0$ is the baseline transmission rate, $F\left(e\right)$ is the Forman-Ricci curvature of edge $e=(i,j)$, and $exp\left[F\right(e\left)\right]$ adjusts the baseline transmission rate based on the local geometric structure of the network. $F_{\max }$ is the maximum curvature value over all edges e in the network, ${\beta }_{\max }$ is the upper limit of ${\beta }_{ij}$, which guarantees that the transmission rate ${\beta }_{ij}$ remains finite. This ensures that ${\beta }_{ij}$ is bounded for all edges in the network, making it suitable for mathematical and computational modeling.
• Boundedness of Recovery Rate: The recovery rate $\gamma >0$ is a constant and hence finite.
• Boundedness of $\mathbf{F}$: Each component of $\mathbf{F}$ is bounded:

  $$\begin{array}{cc}\hfill f_{S_i}& =-\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j\mathrm{is}\mathrm{bounded}\sin \mathrm{ce}{S}_i,I_j\in [0,1]\mathrm{and}{\beta }_{ij}\mathrm{is}\mathrm{finite},\hfill \\ \hfill f_{I_i}& =\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i\mathrm{is}\mathrm{bounded}\mathrm{as}\mathrm{both}\mathrm{terms}\mathrm{involve}\mathrm{bounded}\mathrm{variables},\hfill \\ \hfill f_{R_i}& =\gamma I_i\mathrm{is}\mathrm{bounded}\sin \mathrm{ce}\gamma \mathrm{and}{I}_i\mathrm{are}\mathrm{finite}.\hfill \end{array}$$

Since each component is bounded, the vector field $\mathbf{F}$ is locally bounded, the vector field $\mathbf{F}$ is locally bounded over the domain ${[0,1]}^n$. This ensures the modified SIR model satisfies the boundedness criterion required for well-posedness. □

Appendix A.3.3. Proof of Lipschitz Continuity

To verify that the vector field $\mathbf{F}=(f_{S_i},f_{I_i},f_{R_i})$ is Lipschitz continuous, we show that:

$$\parallel \mathbf{F}\left(\mathbf{X}\right)-\mathbf{F}\left(\mathbf{Y}\right)\parallel \le {L\parallel \mathbf{X}-\mathbf{Y}\parallel ,\forall \mathbf{X},\mathbf{Y}\in [0,1]}^n,$$

where $\mathbf{X}=(S_i,I_i,R_i)$ and $\mathbf{Y}=(S_i^{\prime },I_i^{\prime },R_i^{\prime })$ are state vectors of the system, and $L>0$ is a Lipschitz constant.

Step 1: Define the Functions: The vector field $\mathbf{F}=(f_{S_i},f_{I_i},f_{R_i})$ is defined as:

$$f_{S_i}=-\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j,$$

$$f_{I_i}=\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i,$$

$$f_{R_i}=\gamma I_i.$$

Here, $S_i,I_i,R_i\in [0,1]$, ${\beta }_{ij}:={\beta }_0·exp\left[F\left(e\right)\right]$, where ${\beta }_0>0$ and $F\left(e\right)$ is the bounded Forman-Ricci curvature.

Step 2: Difference of Functions: Consider two state vectors $\mathbf{X}=(S_i,I_i,R_i)$ and $\mathbf{Y}=(S_i^{\prime },I_i^{\prime },R_i^{\prime })$. For each component:

The difference in $f_{S_i}$ is:

$$f_{S_i}\left(\mathbf{X}\right)-f_{S_i}\left(\mathbf{Y}\right)=-\sum _{j\in N\left(i\right)}{\beta }_{ij}\left(S_i{I}_j-S_i^{\prime }{I}_j^{\prime }\right).$$

The difference in $f_{I_i}$ is:

$$f_{I_i}\left(\mathbf{X}\right)-f_{I_i}\left(\mathbf{Y}\right)=\sum _{j\in N\left(i\right)}{\beta }_{ij}\left(S_i{I}_j-S_i^{\prime }{I}_j^{\prime }\right)-\gamma \left(I_i-I_i^{\prime }\right).$$

The difference in $f_{R_i}$ is:

$$f_{R_i}\left(\mathbf{X}\right)-f_{R_i}\left(\mathbf{Y}\right)=\gamma \left(I_i-I_i^{\prime }\right).$$

Step 3: Boundedness of Differences: Now, we give the boundedness of these differences.

• Transmission Term: For ${\beta }_{ij}$, since ${\beta }_0$ and $F\left(e\right)$ are bounded, we have:

  $${\beta }_{ij}\le {\beta }_{\max }={\beta }_0·exp\left[F_{\max }\right].$$
  Thus, ${\beta }_{ij}$ is bounded for all edges e.
• State Variables: Since $S_i,I_i,R_i\in [0,1]$, their differences are bounded:

  $$|S_i-S_i^{\prime }|\le 1,|I_i-I_i^{\prime }|\le 1,|R_i-R_i^{\prime }|\le 1.$$
• Neighborhood Size: Let $|N(i)|$ denote the size of the neighborhood of node i. Since real-world networks typically have finite degrees, $|N(i)|$ is bounded.

Step 4: Lipschitz Condition for Each Component: For $f_{S_i}$:

$$|f_{S_i}\left(\mathbf{X}\right)-f_{S_i}\left(\mathbf{Y}\right)|\le \sum _{j\in N\left(i\right)}{\beta }_{ij}\left(|S_i||I_j-I_j^{\prime }|+|I_j^{\prime }||S_i-S_i^{\prime }|\right).$$

Using the boundedness of ${\beta }_{ij}$, $|S_i|$, and $|I_j|$:

$$|f_{S_i}\left(\mathbf{X}\right)-f_{S_i}\left(\mathbf{Y}\right)|\le {\beta }_{\max }|N\left(i\right)|\left(|S_i||I_j-I_j^{\prime }|+|I_j^{\prime }||S_i-S_i^{\prime }|\right).$$

For $f_{I_i}$:

$$|f_{I_i}\left(\mathbf{X}\right)-f_{I_i}\left(\mathbf{Y}\right)|\le \sum _{j\in N\left(i\right)}{\beta }_{ij}\left(|S_i||I_j-I_j^{\prime }|+|I_j^{\prime }||S_i-S_i^{\prime }|\right)+\gamma |I_i-I_i^{\prime }|.$$

Again, by bounding ${\beta }_{ij}$ and $|S_i|$, we find:

$$|f_{I_i}\left(\mathbf{X}\right)-f_{I_i}\left(\mathbf{Y}\right)|\le {\beta }_{\max }|N\left(i\right)|\left(|S_i||I_j-I_j^{\prime }|+|I_j^{\prime }||S_i-S_i^{\prime }|\right)+\gamma |I_i-I_i^{\prime }|.$$

For $f_{R_i}$:

$$|f_{R_i}\left(\mathbf{X}\right)-f_{R_i}\left(\mathbf{Y}\right)|=\gamma |I_i-I_i^{\prime }|.$$

Step 5: Combine Components: The total difference is:

$$\parallel \mathbf{F}\left(\mathbf{X}\right)-\mathbf{F}\left(\mathbf{Y}\right)\parallel \le L\parallel \mathbf{X}-\mathbf{Y}\parallel ,$$

where L depends on ${\beta }_{\max },\gamma$, and $|N(i)|$. Since these quantities are finite, L is a finite constant, proving that $\mathbf{F}$ is Lipschitz continuous.

Hence, the vector field $\mathbf{F}$ satisfies the Lipschitz condition, ensuring the existence and uniqueness of solutions for the modified SIR model by the Picard-Lindelöf theorem which states that a system of ODEs:

$${\displaystyle \frac{d\mathbf{X}}{dt}}=\mathbf{F}\left(\mathbf{X}\right),\mathbf{X}\left(0\right)={\mathbf{X}}_0,$$

has a unique solution if $\mathbf{F}\left(\mathbf{X}\right)$ is Lipschitz continuous and locally bounded. Since these conditions are satisfied by the modified SIR model, the system admits a unique solution for any given initial condition.

Thus, the modified SIR model is well-posed because it admits unique solutions for any admissible initial condition. The solutions are biologically feasible, that is, ($S_i,I_i,R_i\in [0,1]$) and conserve the total population. This guarantees the model’s reliability for studying epidemic dynamics on complex networks.

Appendix A.4. Stability Analysis

To study stability, we linearize the system around its equilibrium points and analyze the eigenvalues of the Jacobian matrix.

Appendix A.4.1. Disease-Free Equilibrium (DFE)

The Disease-Free Equilibrium (DFE) represents the state of the system where there is no infection in the population. At this equilibrium, all individuals are susceptible and there are no infectious or recovered individuals.

$$S_i^{*}=1,I_i^{*}=0,R_i^{*}=0,\forall i.$$

Appendix A.4.2. Linearization Around the DFE

To analyze the stability of the DFE, we linearize the system around this equilibrium. Substituting $S_i^{*}=1$, $I_i^{*}=0$, and $R_i^{*}=0$ into the modified equations 11, 12and 13. At the DFE:

$${\displaystyle \frac{d{S}_i}{dt}}=0,{\displaystyle \frac{d{R}_i}{dt}}=0,$$

and the equation for ${\displaystyle \frac{d{I}_i}{dt}}$ reduces to:

$${\displaystyle \frac{d{I}_i}{dt}}\approx \sum _{j\in N\left(i\right)}{\beta }_{ij}{I}_j-\gamma I_i.$$

Appendix A.4.3. Jacobian Matrix

The linearized dynamics can be represented in matrix form using a Jacobian matrix denoted as J. For the infectious states $I_i$, to get the Jacobian matrix J, we computed the partial derivatives of the infectious dynamics equation for $I_i$:

$${\displaystyle \frac{d{I}_i}{dt}}=\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i.$$

Now, the Jacobian matrix J is defined as:

$$J_{ij}={\displaystyle \frac{\partial \left(\frac{d{I}_i}{dt}\right)}{\partial I_j}}.$$

Case 1: When $i=j$ (Diagonal Elements of J)

For the diagonal elements, we compute:

$${\displaystyle \frac{\partial }{\partial I_i}}\left(\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i\right).$$

Breaking this into terms:

• From the first term ${\sum }_{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j$, only the $j=i$ component contributes, which makes:

  $${\displaystyle \frac{\partial }{\partial I_i}}\left({\beta }_{ii}{S}_i{I}_i\right)={\beta }_{ii}{S}_i.$$
• From the second term $-\gamma I_i$:

  $${\displaystyle \frac{\partial }{\partial I_i}}\left(-\gamma I_i\right)=-\gamma .$$

Thus, for $i=j$:

$$J_{ii}={\beta }_{ii}{S}_i-\gamma .$$

Case 2: When $i\ne j$ (Off-Diagonal Elements of J)

For the off-diagonal elements, we compute:

$${\displaystyle \frac{\partial }{\partial I_j}}\left(\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i\right).$$

Here, only the term ${\beta }_{ij}{S}_i{I}_j$ contributes, as all other terms are independent of $I_j$:

$${\displaystyle \frac{\partial }{\partial I_j}}\left({\beta }_{ij}{S}_i{I}_j\right)={\beta }_{ij}{S}_i.$$

Thus, for $i\ne j$:

$$J_{ij}={\beta }_{ij}{S}_i.$$

Combining the diagonal and off-diagonal elements, the Jacobian matrix then becomes:

$$J_{ij}=\left\{\begin{array}{cc}{\beta }_{ij}{S}_i,\hfill & i\ne j,\hfill \\ -\gamma +{\beta }_{ii}{S}_i,\hfill & i=j.\hfill \end{array}\right.$$

Here:

• ${\beta }_{ij}$ is the transmission rate between nodes i and j, which incorporates the Forman-Ricci Curvature $F\left(e\right)$ as:

  $${\beta }_{ij}:={\beta }_0·exp\left[F\left(e\right)\right].$$
• $\gamma$ represents the recovery rate, which is constant across all nodes.

In matrix notation, the Jacobian can be expressed as:

$$J=\mathcal{B}-\gamma I,$$

where:

• $\mathcal{B}$ is the transmission matrix with elements $\mathcal{B}={\beta }_{ij}$.
• I is the identity matrix of appropriate dimensions.

The transmission matrix $\mathcal{B}$ is an $n\times n$ matrix, where n is the number of nodes in the network. Each element ${\beta }_{ij}$ represents the transmission rate from node j to node i, given by:

$$\mathcal{B}=\left[\begin{array}{ccccc}{\beta }_{11}& {\beta }_{12}& {\beta }_{13}& \cdots & {\beta }_{1n}\\ {\beta }_{21}& {\beta }_{22}& {\beta }_{23}& \cdots & {\beta }_{2n}\\ {\beta }_{31}& {\beta }_{32}& {\beta }_{33}& \cdots & {\beta }_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ {\beta }_{n1}& {\beta }_{n2}& {\beta }_{n3}& \cdots & {\beta }_{nn}\end{array}\right].$$

Properties of $\mathcal{B}$:

• If nodes i and j are not directly connected, ${\beta }_{ij}=0$.
• Diagonal elements (${\beta }_{ii}$) are set to zero since self-transmission is not possible in standard SIR frameworks:

Hence,

$$\mathcal{B}=\left[\begin{array}{ccccc}0& {\beta }_{12}& {\beta }_{13}& \cdots & {\beta }_{1n}\\ {\beta }_{21}& 0& {\beta }_{23}& \cdots & {\beta }_{2n}\\ {\beta }_{31}& {\beta }_{32}& 0& \cdots & {\beta }_{3n}\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ {\beta }_{n1}& {\beta }_{n2}& {\beta }_{n3}& \cdots & 0\end{array}\right].$$

The stability of the DFE depends on the eigenvalues $\lambda$ of the Jacobian J. The eigenvalues satisfy the characteristic equation:

$$det(J-\lambda I)=det(\mathcal{B}-(\gamma +\lambda \left)I\right)=0.$$

This simplifies to:

$$\lambda =\rho \left(\mathcal{B}\right)-\gamma ,$$

where $\rho \left(\mathcal{B}\right)$ is the spectral radius (largest eigenvalue) of the transmission matrix $\mathcal{B}$.

The eigenvalues $\lambda$ of the Jacobian matrix J are computed to assess the stability of the equilibrium. The eigenvalues satisfy the characteristic equation:

$$det(J-\lambda I)=0.$$

The equilibrium is locally stable if all eigenvalues $\lambda$ satisfy $\Re \left(\lambda \right)<0$, where $\Re \left(\lambda \right)$ is the real part of $\lambda$.

Appendix A.4.4. Stability Condition for the DFE

The DFE occurs when:

$$S_i^{*}=1,I_i^{*}=0,R_i^{*}=0.$$

This represents a fully susceptible population with no infections. The DFE is stable if all eigenvalues $\lambda$ of J are negative. From the eigenvalue equation:

$$\lambda =\rho \left(\mathcal{B}\right)-\gamma ,$$

the stability condition becomes:

$$\lambda <0\Rightarrow \rho \left(\mathcal{B}\right)<\gamma .$$

Appendix A.4.5. Interpretation of the Stability Condition

• Basic Reproduction Number ($R_0$): The spectral radius $\rho \left(\mathcal{B}\right)$ reflects the effective reproduction number in the network. If $\rho \left(\mathcal{B}\right)>\gamma$, infections can spread through the network. Defining $R_0={\displaystyle \frac{\rho \left(\mathcal{B}\right)}{\gamma }}$, the DFE is stable if:

  $$R_0<1.$$
• Network Impact: The stability condition highlights the role of network structure, encoded in $\mathcal{B}$, in determining disease spread. Highly connected networks (with large $\rho \left(\mathcal{B}\right)$) are more likely to destabilize the DFE, leading to an epidemic outbreak.
• Role of Forman-Ricci Curvature: The curvature-based transmission rate ${\beta }_{ij}$ modifies $\mathcal{B}$, influencing $\rho \left(\mathcal{B}\right)$. Regions of high curvature (low ${\beta }_{ij}$) contribute to stabilizing the DFE, while regions of low curvature (high ${\beta }_{ij}$) increase the risk of destabilization.

Hence, the Disease-Free Equilibrium (DFE) is stable if:

$$\rho \left(\mathcal{B}\right)<\gamma \mathrm{or}\mathrm{equivalently}{R}_0<1.$$

This condition ensures that infections die over time, preventing an epidemic outbreak. The spectral properties of the transmission matrix $\mathcal{B}$, influenced by both network topology and Forman-Ricci Curvature, are critical in determining the stability of the DFE.

Appendix A.5. Endemic Equilibrium

The endemic equilibrium represents a state where the infection persists within the population, with all compartments ($S_i$, $I_i$, and $R_i$) having non-zero values.

$$S_i^{*}\ne 0,I_i^{*}\ne 0,R_i^{*}\ne 0.$$

At equilibrium, the rates of change of all compartments vanish, that is:

$${\displaystyle \frac{d{S}_i}{dt}}=0,{\displaystyle \frac{d{I}_i}{dt}}=0,{\displaystyle \frac{d{R}_i}{dt}}=0.$$

Substituting these conditions into the system of equations for the modified SIR model, we obtain:

• Susceptible Dynamics:

  $${\displaystyle \frac{d{S}_i}{dt}}=-\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i^{*}{I}_j^{*}=0.$$
  This implies:

  $$\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i^{*}{I}_j^{*}=0.$$
  Since $S_i^{*}>0$ and $I_j^{*}>0$ for the endemic equilibrium, this equation holds automatically when the infectious equilibrium condition is satisfied.
• Infectious Dynamics:

  $${\displaystyle \frac{d{I}_i}{dt}}=\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i^{*}{I}_j^{*}-\gamma I_i^{*}=0.$$
  Rearranging:

  $$\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i^{*}{I}_j^{*}=\gamma I_i^{*}.$$
  This balance equation ensures that the new infection rate equals the equilibrium recovery rate.
• Recovered Dynamics:

  $${\displaystyle \frac{d{R}_i}{dt}}=\gamma I_i^{*}=0.$$
  This is consistent with the equilibrium condition for $I_i^{*}$.
• Population Conservation: The total population at each node remains conserved:

  $$S_i^{*}+I_i^{*}+R_i^{*}=1.$$

Combining the above equilibrium conditions, the endemic equilibrium is defined by the following system of equations:

$$\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i^{*}{I}_j^{*}=\gamma I_i^{*},\forall i,(\mathrm{Infectious}\mathrm{equilibrium}\mathrm{condition}).$$

$$S_i^{*}+I_i^{*}+R_i^{*}=1,\forall i,(\mathrm{Population}\mathrm{conservation}).$$

These equations are nonlinear due to the product $S_i^{*}{I}_j^{*}$, and their solution depends on the network structure ($N\left(i\right)$), the transmission rates (${\beta }_{ij}$), and the recovery rate ($\gamma$).

Given the nonlinearity of these equations, analytical solutions are generally infeasible, necessitating numerical methods for their solution.

Appendix A.5.1. Existence of the Endemic Equilibrium

To establish the existence of an endemic equilibrium, these conditions must hold:

(1)

Non-Zero Infection Condition: For $I_i^{*}>0$ to hold, the reproduction number $R_0$ must satisfy:

$$R_0={\displaystyle \frac{\rho \left(\mathcal{B}\right)}{\gamma }}>1,$$

where $\rho \left(\mathcal{B}\right)$ is the spectral radius of the transmission matrix $\mathcal{B}$. This ensures that the infection can sustain itself within the population.

(2)

Positivity of Solutions: The feasible domain ($[0,1]$) for $S_i^{*}$, $I_i^{*}$, and $R_i^{*}$ guarantees biologically meaningful solutions for the equilibrium.

Jacobian Matrix Calculation for Endemic Equilibrium

The Jacobian matrix $J^{*}$ of endemic equilibrium captures the sensitivity of the system near the endemic equilibrium. Its entries are:

$${\displaystyle \frac{\partial F_{I_i}}{\partial S_j}}=\left\{\begin{array}{cc}\sum _{j\in N\left(i\right)}{\beta }_{ij}{I}_j^{*},\hfill & \mathrm{if}j\in N\left(i\right),\hfill \\ 0,\hfill & \mathrm{otherwise},\hfill \end{array}\right.$$

$${\displaystyle \frac{\partial F_{I_i}}{\partial I_j}}=\left\{\begin{array}{cc}\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i^{*}-\gamma ,\hfill & \mathrm{if}j=i,\hfill \\ {\beta }_{ij}{S}_i^{*},\hfill & \mathrm{if}j\in N\left(i\right),\hfill \\ 0,\hfill & \mathrm{otherwise}.\hfill \end{array}\right.$$

In conclusion, the endemic equilibrium:

$$S_i^{*}\ne 0,I_i^{*}\ne 0,R_i^{*}\ne 0,\forall i,$$

exists if the basic reproduction number $R_0>1$. Its stability depends on the spectral properties of the Jacobian matrix $J^{*}$, influenced by network topology and the curvature-based transmission rates (${\beta }_{ij}$). Numerical methods are essential for solving and validating the endemic equilibrium conditions in practical scenarios.

Appendix A.6. Global Stability of the Disease-Free Equilibrium (DFE)

We analyze the global stability of the Disease-Free Equilibrium (DFE) of the network-based SIR model given by equations (11)–(13) using the Lyapunov method.

Theorem A1

(Global Stability of the DFE). Consider the network-based SIR model defined by equations (11)–(13), where $S_i,I_i,R_i$ represent the proportions of the susceptible, infectious, and recovered individuals at node i, respectively. Let the transmission rate be ${\beta }_{ij}:={\beta }_0·exp\left[F\left(e\right)\right]$, and let $\gamma >0$ denote the recovery rate. Then, the Disease-Free Equilibrium (DFE),

$$S_i^{*}=1,I_i^{*}=0,R_i^{*}=0\forall i,$$

is globally asymptotically stable.

Proof. 

To prove the global stability of the DFE, we construct the following Lyapunov candidate function to capture the total "infectious energy" in the network:

$$\mathcal{V}=\sum _{i=1}^n{\displaystyle \frac{I_i^2}{2}},$$

where n is the number of nodes in the network. Clearly, $\mathcal{V}\ge 0$ for all i, and $\mathcal{V}=0$ if and only if $I_i=0$ for all i, which corresponds to the DFE.

Taking the time derivative of $\mathcal{V}$ along the trajectories of the system:

$${\displaystyle \frac{d\mathcal{V}}{dt}}=\sum _{i=1}^n{I}_i{\displaystyle \frac{d{I}_i}{dt}}.$$

From the infectious dynamics in the model, we substitute:

$${\displaystyle \frac{d{I}_i}{dt}}=\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i,$$

so that:

$${\displaystyle \frac{d\mathcal{V}}{dt}}=\sum _{i=1}^n{I}_i\left(\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_j-\gamma I_i\right).$$

Expanding the sum:

$${\displaystyle \frac{d\mathcal{V}}{dt}}=\sum _{i=1}^n\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_i{I}_j-\sum _{i=1}^n\gamma I_i^2.$$

We analyze the terms:

• Transmission Term: ${\sum }_{i=1}^n{\sum }_{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_i{I}_j$. At the DFE, $S_i^{*}=1$, so $S_i\le 1$ for all i. Since ${\beta }_{ij}\ge 0$ and $I_i,I_j\ge 0$, this term is non-negative, and vanishes only when $I_i=0$ for all i.
• Recovery Term: $-{\sum }_{i=1}^n\gamma I_i^2$. This term is non-positive since $\gamma >0$ and $I_i^2\ge 0$. It strictly decreases $\mathcal{V}$ unless all $I_i=0$.

Thus, the derivative of $\mathcal{V}$ satisfies:

$${\displaystyle \frac{d\mathcal{V}}{dt}}\le \sum _{i=1}^n\sum _{j\in N\left(i\right)}{\beta }_{ij}{S}_i{I}_i{I}_j-\sum _{i=1}^n\gamma I_i^2.$$

At the DFE, $I_i=0$ for all i, so both terms vanish and ${\displaystyle \frac{d\mathcal{V}}{dt}}=0$.

To apply LaSalle’s Invariance Principle, we observe:

• $\mathcal{V}$ is positive definite.
• ${\displaystyle \frac{d\mathcal{V}}{dt}}\le 0$, i.e., negative semi-definite.
• The largest invariant set in $\left\{{\displaystyle \frac{d\mathcal{V}}{dt}}=0\right\}$ is the DFE, since the only way for the derivative to vanish is for $I_i=0\forall i$.

Hence, according to LaSalle’s invariance principle [26], all trajectories converge to the largest invariant set where $I_i=0$ for all i. Therefore, DFE is globally asymptotically stable. □