We propose a two-scale epidemic model that distinguishes intra-community infection dynamics from intercommunity transmission. At the micro-layer, an epidemic outbreak within a community is simulated with the Gillespie SIR engine with additional infections. At the macro scale each community collapses to a single node in a mobility-weighted meta-population graph; stochastic transmissions between nodes are timed with a hazard-integral implementation of Gillespie's SSA allows dynamically varying infectivity and susceptibility parameters to account for behavior changes, interventions, or emerging variants. We aim to reduce the complexity of visual analysis and simulation by introducing techniques to capture key aspects of two-scale network dynamics. For visualization, we propose Sankey diagrams to depict virus strain infection rates across communities. Our two-scale simulation approach considers cumulative infection exposure, allowing nodes to continue receiving infection pressure from external sources after initial outbreaks, capturing the efect of ongoing inter-community interactions. By abstracting macro-edges from point-to-point contacts a complexity reduction by a multiplicative factor proportional to the squared average size of a community can be achieved. Barabási-Albert model of complex scale-free network Continuous-Time Markov Chain Intensive-Care Unit (special Infectious state) Susceptible-Infectious-Removed epidemic model Stochastic Simulation Algorithm (Gillespie)
BA CTMC ICU SIR SSA
The aim of the paper is to develop a stochastic simulation and visual analysis technique for epidemic outbreaks on complex networks that would be possible to scale up to a country with millions of inhabitants. We condense a community into a single node (metanode) to model transmission between metanodes and use an agent-based simulation with infections from the outside to depict infection dynamics within communities.
Pandemic preparedness depends on our ability to anticipate how quickly and in what way an infection will spread, so that intensive care capacity, test-and-trace teams, and vaccination logistics can be mobilized in time and in the right places. In practice, the most critical policy indicators are how high the simultaneous case load climbs and how long the authorities have before that peak arrives. Classical compartmental models capture long-term averages, but they smooth out exactly those local surges that overwhelm hospitals. A finer-grained lens is, therefore, essential.
One principled way to incorporate such fine -grained detail is to represent every individual as a node
in a contact network and let infection propagate along its edges. Our previous work [
To overcome this limitation, we introduce a two-layer simulation model. By isolating the infection process within each community (micro layer) from the inter-community spread (macro layer), we efectively eliminate a large fraction of the network edges from the simulation space. This reduces overall complexity: the simulation time scales roughly linearly with the number of communities, in contrast to the quadratic growth in a full agent-based model. Moreover, representing each community as a node on the macro layer opens up new avenues for localized intervention. Policies such as mobility restrictions, targeted vaccination, or school closures can be applied to individual communities and their impact on national-level epidemic dynamics can be eficiently simulated and analyzed.
Out work is mainly inspired by Colizza et al. [
Multi-scale epidemic models extend classical compartmental models by embedding them in
metapopulation frameworks, where nodes (metanodes) represent spatially distinct populations (e.g., cities or
regions), and edges (metalinks) represent interactions such as travel or migration. The foundational
work by Colizza and Vespignani [
Recent developments focus on integrating high-resolution mobility data, heterogeneous contact patterns,
and agent-based models at the local level while maintaining coarse-grained population-level transmission
dynamics between communities. Notable approaches include the GLEaM (Global Epidemic and Mobility
model) [
Advances in computational tools have also enabled scalable simulations across multiple levels of
granularity. For instance, EpiGraph [
One more work that describes multiscale agent-based model is [
Ongoing challenges include aligning model granularity with data availability, handling nonMarkovian transitions, and calibrating heterogeneous models under uncertainty.
National-level epidemic simulations must reconcile two conflicting requirements: (i) microscopic realism inside each locality and (ii) computational tractability across millions of individuals. Our Two-Scale (a) Microscopic view. Agent-level contact graph (b) Two-layer abstraction. The inter-community comprising four densely connected communities links have been abstracted to appear as single links (colours) and a sparse set of inter-community links. (black edges) instead of point-to-point connections.
Epidemic Model achieves that balance by running two inter-connected layers whose time lines advance in perfect synchrony (Fig. 1b): Micro layer. Every community is treated as an independent population and simulated with a stochastic Gillespie-based engine. At the end of each micro step micro (typically a few hours) the model outputs two scalars: contagiousness () and susceptibility ().
Macro layer. Communities are the nodes of a mobility-weighted graph. A continuous-time Markov chain (CTMC)whose time-dependent rates () = () () capture cross-community mobility generates the moments at which one community seeds a new infection in another. Event timing is obtained with the hazard–integral method described in Section 3.2.
Layer synchronisation. The two layers are updated in strict lock-step: 1. During each micro interval micro every community integrates its own SIR dynamics and refreshes (, ); those numbers are immediately passed to the macro layer to update all (). 2. The macro CTMC accumulates hazard () = ∫︀ Θ( ) on the same grid = 0 + micro; when () reaches its random budget ∼ Exp(1) it triggers exactly one inter-community infection, seeds an infectious individual in the destination micro model, sets 0 ← ⋆, draws a new ′, and resumes.
This bidirectional handshake guarantees that the two layers remain numerically consistent at every discrete step.
Sankey Flow diagram. The diagram in Fig. 2 summarizes the simulation output of a simulated epidemic outbreak on the ’network of communities’ scale. Figure 2 arranges the communities in two vertical bands. The bigger a node that represents a community is, the more infections the community transmitted or received. Nodes on the (left) represent the source community of an infection event, while those on the central and right marks the recipient communities. The flow from left to central column counts all inter-community (macro-layer) infection events, while the frow from central to right column tallies the accumulated intra-community (micro-layer) infections. The thickness of each ribbon is scaled as
= √︀ , where is the total number of transmissions and a global scale factor chosen so that the thinnest visible edge remains ≈ 0.4 pt. This square-root mapping preserves small but epidemiologically relevant lfows without allowing the largest outbreaks to dominate the canvas.
Start by following the thickest ribbons leaving the macro-layer. The most prominent stream pinpoints the community that exports the bulk of infections. Next, compare the width of its return flow: a narrow link back to itself but broad outgoing links to others signals a net exporter of infection pressure. In contrast, symmetric back-and-forth ribbons between two nodes indicate reciprocal mixing. The microlayer further reveals where those imports amplify, because the right-hand node width is proportional to the cumulative intra-community burden. Thus, a modest import flow that swells into a thick micro-layer node flags a vulnerable community with intense local spread. quick visual estimate of the proportional contribution of each external source to a given outbreak.
Because the underlying contact graph (in which the simulation 2 was performed) is generated with the Barabási–Albert (BA) model, it contains high-degree hubs that contact many nodes in almost every community. As a result, the in the network are not strongly separated communities, therefore the totals of inter-community and intra-community infections are similar. Real-world contact networks are also approximately scale-free but typically have hubs of lower degree—households, workplaces, or schools rather than global “super-connectors.” We therefore expect real data to show a larger share of intra-community transmission than the BA benchmark depicted here.
Complexity advantage. Let be the total number of agents (nodes) split into communities. Denote by their local sizes, by the number of edges in the intra-community contact graph of community , and by inter the number of inter-community mobility edges in the macro graph. A naïve agent-level SSA stores and updates one channel per susceptible–infectious edge, i.e. (tot) memory and (tot ) run time, with tot = ∑︀ + inter, T - total number of events. In dense networks with nodes, up to ( 2) edges may occur which can significantly strain the memory and time resources, complicating simulations. Assuming the population can be partitioned into equally sized communities, and that all inter-community connections between the same pair of communities are aggregated into a single weighted macro-link, the space complexity of the model can be reduced to (︀ )︀ 2)︁ . Simplifying the expression, we obtain that the space complexity reduces by a multi︁( plicative factor ︁( 2 )︁
sizes and Gillespie techniques as discussed in [
. The time complexity can be precisely determined using specific community ⟨⟩ are the current number of infectious nodes and their average intra-degree; • each Gillespie micro engine costs () memory and (⟨⟩) time per infection, where and • the macro layer needs only (inter) channels and advances them in (inter) per macro event by the hazard–integral rule. inter︀) , yielding near-linear scaling in the number of communities.
Overall storage scales as (︀ ∑︀ +inter ≪ ︀)
(tot) and the run time per global step as (︀ ∑︀ ⟨⟩+ Targeted–intervention capability. Because each community evolves in its own micro simulator, the platform naturally supports spatially targeted measures—lock-downs, surge vaccination or travel restrictions—are applied by simply altering that community’s parameters. An intervention applied to community simply modifies its micro parameters (e.g. infection and recovery rates, edge weights or community structure), which in turn alters () and (). Those updated scalars automatically feed into the macro rates () and (), letting the model capture knock-on efects such as reduced export of cases from a locked-down community or heightened import risk for under-vaccinated comminuty.
The national (or regional) population is partitioned into discrete communities—cities, districts, university campuses, etc.— which become the nodes of a directed, weighted meta-population network. The weight on edge (, ) measures how easily infection can pass from community to community ; it can be estimated from commuting statistics, transportation schedules, mobile-phone mobility data, or other social-interaction proxies.
Every node carries two time-dependent scalars, obtained from a within-community • Susceptibility () – proportional to the remaining susceptible individuals in the community . Inter-community transmission rate. The instantaneous hazard (rate) that one infection is seeded from community to community is
() = () () , where is infection rate or the biological transmission rate of the pathogen and is a constant calibration factor that matches the macro layer infection rate to the micro layer one. Only one kind of stochastic event exists at the macro layer: “one additional infection appears in community ”. After such an event, both () and () are updated by the micro-layer models and all () are refreshed. (micro-layer) simulator:
community ; • Contagiousness () – proportional to the current number of infectious individuals in the
Gillespie algorithm samples the waiting time as ∆ = − ln Θ , ∼ (0, 1) .
This relies on the density () = Θ − Θ of the Exp(Θ) distribution. In our application the rate is not constant because () and () evolve continuously. The correct total rate is therefore a function of time: With a constant total rate Θ = ∑︀̸= , the Non-homogeneous Poisson formulation. For a time-dependent rate, the probability that the next macro event has not occurred by real time is Θ( ) = ∑︁ ().
̸= (︀ event > )︀ = exp(︁ − ∫︀0 Θ( ) ,
︁)
⏟
()
⏞
where 0 is the moment the previous macro event happened and () is the cumulative hazard.
This survival function generalises the familiar − Θ of the homogeneous case. For more information
see [
0 () =
Θ( ) , where Θ( ) = ∑︀̸= and convert it through the inverse cumulative distribution function (CDF) of the unit-rate exponential: () is the instantaneous network rate. Draw one uniform deviate ∼ (0, 1) = − 1( ) = − ln =⇒ ∼ Exp(1).
Thus is a sample obtained from the CDF of the Exp(1) distribution. The next inter-community infection occurs at the first calendar time ⋆ > 0 satisfying
(⋆) = .
Intuitively, acts as a random “risk budget”; the event fires when the accumulated hazard reaches this threshold. If Θ( ) is constant, () = Θ( formula ∆ = − ln /Θ .
− 0) and the rule collapses to the familiar Gillespie
Because Θ( ) is available only at discrete micro-layer update times, we approximate the integral in rectangles of width (e.g. = 0.1 day): +1 = + Θ (︀ + 2 )︀ , +1 = + . [, +1); assuming Θ is nearly constant there, we interpolate to obtain the exact calendar time We accumulate these rectangles until +1 ≥
. The event must occur inside the last interval ⋆ = +
Θ (︀ − +2 )︀ .
Selecting the transmitting edge. At ⋆ the probability that edge (, ) is the one that fires is
Pr{︀ (, ) fires at ⋆}︀ = Θ( (⋆⋆)) .
We seed one infection in community , update the micro models for and , set 0 ← ⋆, draw a new ′ ∼ Exp(1), and repeat the cycle.
Accuracy and step-size choice. The local error of the midpoint rectangle is ( 3); choosing one order of magnitude smaller than the reporting interval (e.g. 0.1 day when epidemiological data are daily) keeps timing errors below a few minutes—negligible relative to epidemiological uncertainty. Constant-rate sanity check. If Θ( ) ≡ is constant, the integral collapses to ( − 0) and the method reproduces the classical formula ∆ = − ln / . Thus the hazard-integral procedure is a true generalisation of Gillespie’s original algorithm.
This layer runs an exact continuous-time SIR process on the contact graph of a single community. It
is implemented with the incremental Gillespie algorithm described in our earlier work [
Reaction channels.
2 →− 3
∅→−− − − # Transition Propensity
1 + →−− − micro + 1 = micro
imp() 2 = 3 = imp()
Channel 1 captures within-community contagion, channel 2 recovery, and channel 3 the import process: every macro-layer jump seeds one new infectious individual, drawn uniformly at random from the current susceptible list (alternative heuristics—e.g. degree-weighted seeding— are possible). Synchronising with the macro layer. At the end of each micro step of length micro the engine exports () = (), () = (), which refresh the inter-community rates (). Conversely, when the macro CTMC fires an intercommunity infection at time ⋆, the micro model of the destination community receives one additional via channel 3, after which all propensities and the next-event clock are recomputed. Hence the two layers remain fully synchronised at every discrete tick.
Impact of imports on epidemic curves. Figure 3 illustrates how two external seeding events (red arrows) alter the virus dynamics within a single community: the first import at = 7 (orange curve) raises and advances the peak, while a second import at = 9 (black curve) produces an even higher and earlier epidemic crest. The plot confirms that the micro–macro simulator is supposed to accurately capture the compound efect of imported infections on local dynamics.
erty that underpins the classic Gillespie exponential draw:
• The internal hazard Θ int() = micro/ + evolves continuously. • The import hazard imp() changes asynchronously whenever the macro layer injects a case, producing instantaneous jumps in the total hazard.
Two technical issues follow: 1. Interrupted exponentials. If an import lands before the next internally scheduled event, the pending exponential draw becomes invalid and must be discarded; a new draw is generated from the updated state. 2. Cumulative-hazard evaluation. A draw-once strategy keeps a running integral () = ∫︀0 Θ int() and fires when () first exceeds ∼ Exp(1), but the integral must be reset each time imp() jumps.
These complications motivate a dedicated time-estimation scheme, presented in the following subsection, that merges internal CTMC dynamics with macro imports without biasing event times.
The intra–community clock must decide which happens first: 1. an internal infection or recovery generated by the SIR Gillespie channels (propensity Θ int() = micro/ + ), or 2. an import scheduled by the macro layer at a deterministic calendar time imp, ( = 1, 2, . . . ).
Because the SIR propensities are constant between events, Θ int does not vary on [︀ 0, min{imp,1, ⋆})︀ . This permits an exact, (1) waiting–time draw without numerical integration.
1. Initial draw. At time 0 sample one exponential budget
= − ln , ∼ (0, 1).
If no import were pending, the internal waiting time would be ∆ int = /Θ int. 2. Compare with next import. Let bar be the next deterministic import time (or +∞ if the queue is empty).
bar − 0 {︃< ∆ int import arrives first ,
≥ ∆ int internal event wins. 3. Case A: internal event wins. Fire the infection or recovery at ⋆ = 0 + ∆ int, update (, , ), set 0 ← ⋆, and return to step 1. 4. Case B: import arrives first. Advance the clock to bar = 0 + with < ∆ int; insert one infectious node (channel 3), which changes the rate to Θ i′nt.
Compute the residual hazard budget and resume the clock with res = − Θ int, ∆ after = res . Θ i′nt
Pop the next import from the queue (if any) and go back to step 2.
We introduced a two-scale modeling framework for large-scale epidemic simulation that decouples intra-community infection dynamics from inter-community transmission. The core idea is to model each community independently at the micro-layer — using stochastic agent-based algorithm — and aggregate inter-community interactions via a macro-level meta-population graph. This separation of scales enables significant computational savings while preserving the essential structure of disease spread across regions.
Our analysis shows that, assuming a balanced partition into equally sized communities, simulation and storage complexity can be reduced by a factor of , from ( 2) to (( )2 · ), where is the total population. This reduction is achieved both through an optimized two-level Gillespie algorithm and by learning intra-community dynamics via machine learning surrogates, enabling near-quadratic gains in some configurations. Similar principles apply to visualization: abstracting macro-transmissions and aggregating infection curves allows epidemic flow diagrams to remain interpretable even at national scale.
The important next step would be a thorough analysis of the approximation error introduced by our abstraction to macro-edges and the approximation of inter-infection error introduced by infection of a random node inside a community. We could use Gillespies model on the full network as a reference model, at least for moderate size graphs, and use empirical data or paralllel computing clusters for sampling micro-simulation reference data for large scale reference networks.
Looking ahead, several extensions are promising. First, improving the fidelity of intra-community dynamics through hybrid machine learning–mechanistic models could boost both accuracy and speed. Second, using data from real-world mobility and contact data for macro-layer weighting will enhance policy relevance. Finally, we envision this framework as a foundation for real-time decision support systems that simulate interventions interactively across spatial and social layers. During the preparation of this work, one of the authors (Y. Kuryliak) used ChatGPT-4o and ChatGPT-o3 in order to improve grammar and text readability. Further, the author used ChatGPT to generate ifgure 1(b). After using these tool, the author reviewed and edited the content as needed and take full responsibility for the publication’s content.