Conditional-path Monte Carlo for rare stochastic dynamics on networks: Details and derivations

Abstract: The simulation of rare macroscopic events in stochastic network dynamics, such as widespread epidemic outbreaks, cascading failures in communication networks, or the escape from metastable states in many-body systems, is severely hindered by methodological challenges like catastrophic rejection rates, weight degeneracy, genealogical correlations, and critical slowing down inherent to standard forward-time algorithms, splitting methods, and transition-path sampling. Conditional-path Monte Carlo (CPMC) overcomes these limitations by employing non-local Swendsen-Wang-like cluster updates that operate directly on full-system trajectories. Serving as the technical companion to [Sun, Moody, and Barthel, arXiv:2608.16171], this paper provides the rigorous mathematical foundations and algorithmic details underlying the CPMC framework. We formally define the joint path-graph probability weights and derive the transition and uniformization sum rules that guarantee detailed balance. Applying the framework to susceptible-infectious-susceptible (SIS) models, we systematically construct and optimize single-node and edge graph vertex sets specifically designed to prevent lock avalanches and maintain the structural mobility of the epidemic trunk. Furthermore, we detail a dynamic programming scheme to exactly implement complex boundary conditions - including patient-zero and macroscopic outbreak-size constraints - enabling the rejection-free generation of valid trajectories. Finally, we assess the computational complexity of the algorithm, describe parallelization strategies, and validate CPMC against exact solutions for dynamics on small networks.

• The paper establishes conditional-path Monte Carlo as a mathematically rigorous rare-event sampler whose joint path-graph weights satisfy marginalization and strict detailed balance.
• The method uses nonlocal cluster updates and a complete SIS edge-vertex construction, including E3, which supports arbitrary infection-rate asymmetries and one-directional edges.
• The paper introduces dynamic programming that enforces patient-zero and outbreak-size constraints exactly with zero rejection, while identifying autocorrelation optimization and broader model extensions as open problems.

This manuscript provides the rigorous mathematical foundations and algorithmic details for conditional-path Monte Carlo (CPMC), a rare-event sampling method for stochastic dynamics on networks, serving as the technical companion to the framework introduction by Sun, Moody, and Barthel (Sun et al., 17 Aug 2026). The central object of study is the sampling of full-system trajectories subject to macroscopic constraints—such as a single patient zero at $t=0$ combined with a large outbreak size at $t=T$ in susceptible-infectious-susceptible (SIS) epidemics—under which forward-time stochastic simulation algorithms (SSA), splitting methods, and transition path sampling (TPS) suffer catastrophic rejection rates, path degeneracy, genealogical correlations, or critical slowing down.

Motivation and relation to existing methods

The authors situate CPMC against four established families of techniques. The weighted SSA biases local propensities and corrects via likelihood ratios, but heuristic biasing choices inflate estimator variance, weight degeneracy destroys effective sample size, and the method remains blind to future topological bottlenecks such as low-rate bridging edges. Splitting methods (forward flux sampling, RESTART) and the weighted ensemble rely on order parameters that mask hidden barriers—for instance, ten infected nodes clustered in a dense core differ drastically from ten scattered peripherally in their outbreak potential—and suffer from critical slowing down and founder effects under cloning. TPS operates directly on trajectories but its local shooting and shifting moves cannot substantially mutate the topological core of infection trees without violating physical rules. CPMC departs from all of these by performing non-local, Swendsen-Wang-like cluster updates directly on complete trajectories, generating a Markov chain whose every element strictly satisfies the imposed constraint.

Path probabilities and the cluster-update construction

The formalism starts from the constrained path probability density $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$, where $C$ encodes trajectory constraints and $\Lambda^t$ is the total escape rate. The update decomposes as $\omega \to g \to \omega'$ through an intermediate graph $g$: first, graph vertices are assigned to every physical event of $\omega$ and inserted as background vertices via an inhomogeneous Poisson point process; second, free spacetime clusters are identified and their states randomly reassigned.

The key theoretical contribution is the joint path-graph weight $J(\omega,g)$ and two local sum rules: the transition sum rule, requiring each physical event rate to equal the summed rates of compatible graph vertices, and the uniformization sum rule, requiring state-independent uniformization rates to exactly fill the gap between local escape rates and background-vertex insertion. These guarantee the marginalization $P(\omega)=\sum'_{g\ni\omega}J(\omega,g)$, proved in an appendix, and imply strict detailed balance for the induced path transition density, since $t=T$0 reduces to a manifestly symmetric expression in $t=T$1. Graph construction follows from the Poisson-process structure of the background density, while cluster identification uses union-find with path compression and union-by-rank, yielding near-linear cost in the number of graph vertices.

Vertex sets for SIS models

The systematic design of vertex sets constitutes the most substantive technical content. For single-node vertices, a one-parameter solution family exists, indexed by $t=T$2; minimizing the uniformization rate gives $t=T$3 (only the $t=T$4 vertex), but the authors report that intermediate values such as $t=T$5 are empirically preferable due to autocorrelations—a trade-off between per-update cost and mixing that they acknowledge is resolved only by numerical experience rather than theory.

For edge vertices, three sets are developed. Set E1 admits a closed-form solution space but carries two defects: it requires infection-rate asymmetry bounded by a factor of two ($t=T$6), and it produces upstream lock avalanches—a single downstream fixing vertex locks the patient zero with probability decaying exponentially in the infection-tree length, immobilizing the epidemic trunk. Set E2 introduces the swap vertex $t=T$7, whose inputs are temporally decoupled from its outputs, breaking upstream lock propagation, together with vertices that prevent downstream avalanches; E2 works for arbitrary rate ratios but degenerates on extremely asymmetric edges, losing downstream-flexibility vertices. Set E3 adds blocker vertices $t=T$8 and $t=T$9, restoring full flexibility for any asymmetry, including one-directional edges. An appendix proves that the nine enumerated edge vertices exhaust all possibilities compatible with an infection event and at least one free variable.

Exact constraint enforcement via dynamic programming

The treatment of boundary conditions is a distinctive strength. Conditional locks are simple but break detailed balance at the boundary and inflate autocorrelations. Cluster-level rejection sampling preserves detailed balance and succeeds when the required excess infections $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$0 are small relative to the mean cluster contribution, by Hoeffding-type concentration—but stalls otherwise. The fully general solution is a Knapsack/subset-sum-style dynamic programming table $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$1 counting valid assignments of $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$2 free clusters against both boundary constraints simultaneously, built in $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$3 time and followed by a backward pass that samples each valid configuration with equal probability and zero rejections. This yields exact, rejection-free satisfaction of joint patient-zero and outbreak-threshold conditions—an operationally severe requirement that typically destroys acceptance rates in SSA-based schemes. An initialization procedure bootstraps the chain by incrementally raising the outbreak threshold.

Cost scaling, parallelization, and validation

Per-update complexity scales as $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$4, i.e., linear in total recovery and infection rate integrals, plus the dynamic programming overhead. Autocorrelation times, by contrast, depend delicately on network structure and vertex-rate choices; the paper is candid that predicting them remains open. Parallelization options include trivially independent chains (limited by 20–40% burn-in fractions), temporal decomposition of graph generation, concurrent lock-free union-find, and parallel tempering over model and algorithmic parameters.

Validation compares CPMC against exact master-equation solutions $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$5 on two four-node networks—a homogeneous ring (reduced to six macrostates by translation symmetry, with $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$6, $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$7) and a hub–leaf diamond (nine macrostates, $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$8, $$P(\omega)=C(\omega)\prod_{x\in\omega}\omega_x(s_x)\, e^{-\int_0^T dt\,\Lambda^t}$$9)—showing excellent agreement for unconstrained dynamics. Larger-scale validation under actual outbreak constraints is deferred to the companion paper (Sun et al., 17 Aug 2026), so this document alone does not demonstrate performance gains on realistically sized networks.

Limitations and open questions

The paper concedes several points plainly. Ergodicity of specific vertex-rate choices is not proven in general; the $C$0 single-node solution may fail ergodicity depending on the constraint implementation. Optimal vertex rates ($C$1, $C$2, $C$3 in E2/E3) that minimize autocorrelation times remain an open, model-dependent optimization problem. The conditional-lock scheme strictly violates detailed balance as a boundary effect. Extensions to SIR/SEIR models with absorbing states, time-dependent rates, and temporal networks require reworking the background-vertex insertion process. Finally, the combination with parallel tempering and adaptive tuning of vertex rates is described only prospectively, with supporting work noted as in preparation.

Conclusion

The paper establishes CPMC on firm mathematical footing: joint path-graph weights with provable marginalization and detailed balance, a complete taxonomy of SIS-compatible edge vertices with explicit rate solutions, and an exact dynamic-programming mechanism for enforcing macroscopic trajectory constraints without rejection. Its principal open problems concern autocorrelation optimization and generalization beyond SIS-type reversible dynamics.