Learn to Vaccinate: Combining Structure Learning and Effective Vaccination for Epidemic and Outbreak Control

We thank the reviewer for their time and thoughtful comments. We address each of the points raised below and will clarify them in the revised manuscript.

Stationarity

We refer the reviewer to our response to Q6 of Reviewer Gdj4 for a more detailed response. In short, either the process dies out quickly, in which case vaccination is unnecessary, or it enters a meta-stable regime where it behaves like a Markov chain with a stationary distribution, justifying Assumption 3.1.

Comparison with the SIR model and sample complexity

We thank the reviewer for this observation. Prior work on SIR-based structure learning relies on observing multiple independent cascades, such as Gomez-Rodriguez et al. 2012, since in SIR models each node can be infected at most once. As a result, interactions between nodes are limited within each cascade, and restarting the process is necessary to collect sufficient signal.

In contrast, the SIS model permits reinfection, enabling a single persistent cascade to generate rich temporal correlations over time. Our learning algorithm leverages this property and thus does not require multiple independent cascades. Consequently, our sample complexity guarantees given in Theorem 3.1 are expressed in terms of the number of times specific infection patterns (e.g., $I(Y_j = 0, Y_i = 1, Y_S = y_S)$) appear in the data, rather than the number of cascades.

Empirically, as shown in Figure 3 (Appendix A.2), SISLearn achieves an F1 score of 0.80 with just 400 rounds of observation, demonstrating strong learning performance.

Observation of edges vs infection states

We agree that combining partial edge and infection information is a realistic and important direction. We note that the setting of unknown edges but observed infection states is standard in prior work on SIR-based structure learning, like Gomez-Rodriguez et al., 2012 or Netrapalli & Sanghavi, 2012, which motivates exploring the analogous setup in the SIS setting.

That said, our approach is highly flexible. Since SISLearn learns the neighborhood of each vertex independently, it can easily incorporate partial edge knowledge. Similarly, it can tolerate missing infection data: as long as some infection history is observed for every vertex and pairs of consecutive observations (e.g., at times $t$ and $t{+}1$) are available, the algorithm remains effective, though naturally requiring more samples.

Specific parameters

We emphasize that our parameter choices were guided by the goal of evaluating non-trivial regimes of the SIS model—where vaccinations are both necessary and impactful. If the infection-to-recovery ratio is too low, the disease naturally dies out without intervention, making vaccination redundant. If it is too high, the process becomes effectively unstoppable without unrealistic levels of vaccination. The intermediate regime is where extinction dynamics are sensitive to targeted interventions—precisely the setting where intelligent strategies matter.

To demonstrate robustness beyond this regime, we conducted additional experiments on both low (Figure 5) and high (Figure 6) infection-to-recovery ratios; results are provided here. As can be observed, our combined learning-and-vaccination approach still outperforms the other baselines. If the reviewer has a specific parameter regime in mind, we would be glad to run additional experiments.

Tree width of inferred networks

The China flu graphs (Sec. 5) and USA flu graphs (App. A) have treewidths of 8 and 58, respectively, computed exactly using the state-of-the-art treewidth solver from the PACE 2017 competition (Tamaki, 2019). For graphs with large treewidth (tw > 20), our proposed Greedy algorithm is preferable due to its polynomial-in-n complexity, making it significantly more efficient. Please also see our answer to the computational feasibility comment of Reviewer WHtL, where we discuss new experiments (new Figures 1 & 2).

Usage of data & observation horizon

Unlike SIR-based inference methods that rely on multiple independent cascades, SISLearn uses infection data from a single, ongoing epidemic, needing only the infection states of the vertices over time—even this can be relaxed, as discussed above. While our theoretical guarantees require that the process has reached meta-stability, we show in Appendix A.2 that SISLearn achieves an F1 score of 0.80 with just 400 rounds of data, without needing to wait for meta-stability to be reached.

References

Gomez-Rodriguez M, Leskovec J, Krause A. Inferring Networks of Diffusion and Influence. ACM Trans Knowl Discov Data. 2012.

Netrapalli P, Sanghavi S. Learning the graph of epidemic cascades. SIGMETRICS Perform Eval Rev. 2012.

Tamaki H. Positive-instance driven dynamic programming for treewidth. J Comb Optim. 2019.

We thank the reviewer for their positive assessment and their questions. We answer them below.

Adaptivity

Our learning algorithm, SISLearn, relies on data drawn from the meta-stable distribution of the SIS process, where the probabilities of configurations do not change anymore. However, it can accommodate time-varying infection probabilities, provided the changes are infrequent enough for the system to reach a new meta-stable state between shifts. As shown in Figure 3 of our additional experiments here (also see our response to Q6 of Reviewer Gdj4), the meta-stable state is reached in as few as 50 time steps depending on the SIS parameters.

Importantly, since SISLearn sequentially learns each node’s neighborhood independently, even if the infection parameters shift mid-process, it suffices to re-estimate the new infection probabilities and allow the process to reach the new meta-stable state. Under these conditions, the theoretical guarantees would hold with minor changes. It should be noted that in practice, SISLearn performs well even without waiting for full stabilization. All experiments (main paper and appendix) use observations starting from the first round, and the algorithm still achieves high learning performance (see experiments in Appendix A.2).

Moreover, SISLearn is robust to parameter misspecification. It only requires an estimate of the infection probability, and even if this input is significantly incorrect, the algorithm performs well in practice. In new experiments where the true infection probability was 0.3 but the input was set to 0.7, SISLearn still achieved an F1 score of 0.80 (compared to 0.97 with the correct value), showing graceful degradation under substantial misspecification. These new experiments are shown in Figure 4 here.

Computational feasibility

We conducted two additional experiments to explicitly evaluate the computational feasibility of our DP algorithm. First, we measured runtime on graphs with fixed treewidth $\omega = 10$ and increasing number of vertices (up to $n=2000$); results provided in Figure 1 here. These confirm that DP scales cubically in $n$, consistent with our theoretical analysis (Section D.3.3), and remains practical even for graphs with as many as $2000$ vertices!

Second, we analyzed runtime on random Erdős–Rényi graphs with varying $n$ and naturally increasing treewidth; results provided in Figure 2 here. This experiment highlights that while DP runtime grows exponentially with increasing treewidth, our proposed heuristic (Algorithm 4) remains computationally efficient (scaling cubically as shown in Section 4.3) and achieves vaccination performance close to DP, as demonstrated in Section 5.

Thus, in practice, for graphs with known small-to-moderate treewidth ($\omega \lesssim 20$), the DP method is recommended even for large graphs. Otherwise, the greedy heuristic is an effective and scalable alternative, offering near-optimal performance with significantly lower runtime.

Regarding GNNs

While GNNs have shown success in epidemic modeling (Liu, 2024), to the best of our knowledge, no existing work directly employs GNNs for learning or vaccination strategies in SIS models. Developing such a GNN-based baseline would itself constitute a substantial research effort beyond a simple baseline comparison. However, if the reviewer is aware of specific GNN-based approaches applicable to our setting, we would gladly include them in our evaluation.

Link to new figures: https://drive.proton.me/urls/5N1PC0M6ZW#oCweKFODmBNn

References

Liu Z, Wan G, Prakash BA, Lau MSY, Jin W. A Review of Graph Neural Networks in Epidemic Modeling. Proceedings of the 30th ACM SIGKDD Conference on Knowledge Discovery and Data Mining. KDD ’24.

We sincerely thank the reviewer for their careful reading and attention to detail. The reviewer is absolutely right—Lemma D.1, as originally stated, was incorrect and developed too hastily. We have corrected the result below by applying Theorem 1.2 from Kontorovich & Ramanan (2008). We also address the reviewer's other concerns below.

Lemma D.1 (Replacement)

Given observations $\\{Y^{(t)}\\}_{t \in T_D}$ from the SIS process where $T_D$ are not necessarily consecutive time indices and $Y^{(t)} \in\\{0,1\\}^V$, under Assumption 3.1, for any subset $U \subseteq V$, any state $y_U \in \\{0,1\\}^{|U|}$, and any positive $\varepsilon$ and $\delta$, we have the following deviation bound on the unbiased estimator

$$\mathbb{P}\left(\left| \mathbb{P}(Y_U = y_U) - \hat{\mathbb{P}}(Y_U = y_U) \right| \geq \varepsilon \right) \leq\delta,$$

whenever $|T_D| \geq \frac{2\log(2 / \delta)}{\varepsilon^2 (1-\theta)^2}$, where $0 < \theta < 1$ is a constant depending on the graph structure and model parameters. Here, $\mathbb{P}(Y_U = y_U)$ is the marginal probability of the vertices in $U$ being described by state vector $y_U$ and $\hat{\mathbb{P}}(Y_U = y_U) = \varphi_U$ is the empirical estimate over $T = |T_D|$ samples.

Proof: Let $\mathcal{S}$ be the state space of our Markov chain, i.e., $\mathcal{S} = \\{0,1\\}^V$, $T = |T_D|$ the number of observations, and $\\mathcal{S}^T$ be the product state space of $T$ observations. By definition, the process $\\{Y^{(t)}\\}\_{t \in \mathbb{N}}$ has the Markov property in time, hence so does a sub-sequence $\\{Y^{(t)}\\}\_{t \in T_D}$. We define $\varphi_U:\mathcal{S}^T \rightarrow \mathbb{R}$ as:

$$\varphi_U(Y^1, \dots, Y^T) = \frac{1}{T} \sum_{t=1}^{T} \mathbb{1} \\{Y_U^{(t)} = y_U\\},$$

where $\mathbb{1} \\{\cdot\\}$ is the indicator function.

Notice that this function is $1/T$-Lipschitz with respect to the Hamming metric $d(X,Y) = \sum_{t=1}^T \mathbb{1} \\{X^{(t)} \neq Y^{(t)}\\}$ on $\mathcal{S}^T$.

We apply Thm. 1.2 by Kontorovich & Ramanan (2008) for Markov chains, stating that for a $c$-Lipschitz function $\varphi$ on $\mathcal{S}^n$,

$$\mathbb{P}\\{|\varphi-\mathbb{E} \varphi| \geq t\\} \leq 2 \exp \left(-\frac{t^2}{2 n c^2 M_n^2}\right).$$

In our case, the sequence length is $T$ (replacing $n$ in the theorem), $c=1/T$, $M_n$ becomes $M_T$, and $t$ becomes $\varepsilon$. The bound then becomes

$$\mathbb{P}\\{|\varphi_U-\mathbb{E} \varphi_U| \geq \varepsilon\\} \leq 2 \exp \left(-\frac{\varepsilon^2}{2 T (1/T)^2 M_T^2}\right) = 2 \exp \left(-\frac{T \varepsilon^2}{2 M_T^2}\right),$$

where $M_T = (1-\theta^T)/(1- \theta)$ and $\theta$ is the Markov contraction coefficient given by

\theta = \sup\_{\\{Y^{\prime}, Y^{\prime \prime} \in \mathcal{S} \\}} \left\\|p\left(\cdot \mid Y^{\prime}\right)-p\left(\cdot \mid Y^{\prime \prime}\right)\right\\|_{\mathrm{TV}}.

We have $\theta <1$ since every configuration $Y$ can transfer into the all-zero configuration $\underline{0}$. Hence, for any two states $Y^{\prime}, Y^{\prime \prime} \in \mathcal{S}$ both $p(\underline{0}|Y^{\prime})>0$ and $p(\underline{0}|Y^{\prime \prime})>0$ and the support of $p(\cdot|Y^{\prime})$ and $p(\cdot|Y^{\prime \prime})$ is not disjoint which implies $\theta <1$.

Therefore, $M_T <1/(1-\theta)$ and the probability bound becomes:

$$\mathbb{P}\\{|\varphi_U-\mathbb{E} \varphi_U| \geq \varepsilon\\} < 2 \exp\left(-\frac{T \varepsilon^2 (1-\theta)^2}{2}\right).$$

Finally, setting the RHS to be $\leq \delta$ and solving for $T$ we get

$$T \geq \frac{2\log(2 / \delta)}{\varepsilon^2 (1-\theta)^2}. \quad \square$$

We emphasize that Lemma D.1 is not a central result of our paper. It is used solely to derive concentration bounds for estimating the direct and conditional influence quantities, thereby enabling our sample complexity analysis. Crucially, the correctness of the SISLearn algorithm itself is unaffected by this lemma. With the corrected bound, all downstream results remain valid (modulo minor constant adjustments), and in fact, Lemmas D.2–D.4 become cleaner, as we can eliminate the $1/m$ term in favor of a dependence on the (constant) Markov contraction coefficient $\theta$, which depends on the graph structure and the infection parameters.

Regarding rounds

The real-world interpretation of a "round" depends on the application. It may correspond to seconds in financial or communication networks, or days in epidemiological settings. The time to reach stationarity similarly depends on the timescale of the underlying process.

Regarding graph recovery

We provide experimental results on the learning performance of SISLearn in App. A.2. We also performed new experiments on the robustness of SISLearn, given in new Fig. 4 here (see Q4 of Reviewer Gdj4).

References

Kontorovich L., Ramanan K. Concentration inequalities for dependent random variables via the martingale method. The Annals of Probability. 2008.

We thank the reviewer for their detailed review and their insightful questions. We respond to each of them in detail below.

Q1

We propose three vaccination strategies, all of which solve the Spectral Radius Minimization (SRM) problem:

(1) An optimal polynomial-time algorithm for trees (Appendix B, Algorithm 5).

(2) An optimal (DP) algorithm that solves SRM on arbitrary graphs (Section 4.2, Algorithm 2), with polynomial runtime on graphs of bounded treewidth.

(3) A greedy polynomial-time heuristic (Section 4.3, Algorithm 4).

The abstract refers to the second method. By “optimal,” we mean that it exactly solves the SRM problem: given a budget $K$, it finds a subset of vertices of size $\leq K$ whose removal minimizes the spectral radius. This is stated and proven in Theorem D.1.

Q2

The VUG problem is formulated as an online problem: the agent observes the infection states over time and must decide when and whom to vaccinate, subject to a global budget $K$ (Section 2). In particular, the agent is not required to vaccinate all $K$ nodes at once.

That said, our proposed method instantiates this framework in an offline fashion: we collect $T_D$ rounds of observations, learn the graph using SISLearn, and then vaccinate all $K$ nodes at time $t = T_D + 1$. This choice reflects that, assuming the graph is learned perfectly, vaccinating all at once is optimal for minimizing extinction time.

Exploring adaptive or staggered vaccination—where the agent begins intervening while still learning—is a compelling direction for future work.

Q3

This assumption is not essential for our approach or theoretical results. All our lemmas, theorems, and algorithmic results hold under different initial infection states, whether they involve a single infected vertex, an arbitrary deterministic or probabilistic subset, or any other initial configuration. We adopted the standard uniform-seed initialization solely because it is common in the SIS literature.

Q4

The only disease parameter SISLearn requires is the infection probability, $p_\text{inf}$, which can be estimated from observational data, like done in Kirkeby et al. (2017) via a mean-field approximation. Importantly, SISLearn is robust to parameter misspecification, as evident by new experiments where the true infection probability was 0.3 but the input was set to 0.7, and SISLearn still achieved an F1 score of 0.80 (compared to 0.97 with the correct value). These new experiments are shown in new Figure 4 (see linked file below).

Q5

While real-world social networks may not exhibit small treewidth, our work addresses this directly: our DP algorithm is ideal for graphs with small to moderate treewidth (e.g., up to ~10-20; see new Figure 1, and our response to Reviewer WHtL), while our fast and scalable greedy heuristic (Algorithm 4) performs well on general graphs—including those with high density (see new Figure 2).

Q6

Our assumption of stationarity arises naturally from well-established theoretical results for SIS-type processes. It has been shown that if the infection parameters are above the epidemic threshold $\rho(\mathcal{G}) \geq p_{\text{rec}}/p_{\text{inf}}$, the process may reach a meta-stable distribution, resembling a stationary distribution (see, among others, Schonmann 1985, Liggett 1999, and Mountford et al. 2013). Conversely, if the threshold condition is not satisfied, the infection dies out rapidly, trivially resolving the vaccination problem.

More formally, following the coupling argument of Cator and Van Mieghem (2013), the SIS process can be related to a modified Markov chain that excludes the absorbing all-zero state. This modified chain is ergodic and therefore has a proper stationary distribution. Thus, until extinction occurs, the original SIS process behaves like an ergodic Markov chain with a stationary distribution.

In other words, either the SIS process (quickly) reaches a meta-stable state where our assumption holds, or it dies out. The former is backed up by new experiments, given in new Figure 3, where the SIS process reached the meta-stable state in an average of 100 rounds.

Thus, the stationarity assumption, while seemingly strong, is both theoretically well-founded and empirically justified.

New figures: https://drive.proton.me/urls/5N1PC0M6ZW#oCweKFODmBNn

References

Kirkeby C, et al. Methods for estimating disease transmission rates: Evaluating the precision of Poisson regression and two novel methods. Sci Rep. 2017.

Cator E, et al. Susceptible-infected-susceptible epidemics on the complete graph and the star graph: Exact analysis. Phys Rev E. 2013.

Liggett TM. Stochastic Interacting Systems: Contact, Voter and Exclusion Processes. Springer. 1999.

Mountford T, et al. Metastable densities for the contact process on power law random graphs. Electronic Journal of Probability. 2013.

Schonmann RH. Metastability for the contact process. J Stat Phys. 1985.