Another look at silent circulation of poliovirus in small populations

Background Silent circulation of polioviruses complicates the polio endgame and motivates analyses that explore the probability of undetected circulation for different scenarios. A recent analysis suggested a relatively high probability of unusually long silent circulation of polioviruses in small populations (defined as 10,000 people or smaller). Methods We independently replicated the simple, hypothetical model by Vallejo et al. (2017) and repeated their analyses to explore the model behavior, interpretation of the results, and implications of simplifying assumptions. Results We found a similar trend of increasing times between detected cases with increasing basic reproduction number (R0) and population size. However, we found substantially lower estimates of the probability of at least 3 years between successive polio cases than they reported, which appear more consistent with the prior literature. While small and isolated populations may sustain prolonged silent circulation, our reanalysis suggests that the existing rule of thumb of less than a 5% chance of 3 or more years of undetected circulation with perfect surveillance holds for most conditions of the model used by Vallejo et al. and most realistic conditions. Conclusions Avoiding gaps in surveillance remains critical to declaring wild poliovirus elimination with high confidence as soon as possible after the last detected poliovirus, but concern about transmission in small populations with adequate surveillance should not significantly change the criteria for the certification of wild polioviruses.


Introduction
The detection of a wild poliovirus (WPV) in Borno state in Northeast Nigeria in 2016 (Nnadi et al., 2017), at a time when many hoped WPV circulation had stopped in Africa, raised questions about the possibility of prolonged undetected poliovirus circulation in small populations. The serotype 1 WPV detected in 2016 in Borno occurred almost two years after the prior most Abbreviations: AFP, acute flaccid paralysis; CFP, case-free period; CNC, confidence about no circulation; CNCx%, time when the confidence about no circulation exceeds x%; DEFP, detected-event-free period; OPV, oral poliovirus vaccine; POE, Probability of eradication; TBC, time between detected cases; TUC, time of undetected circulation after the last detected-event; TUCx%, xth percentile of the TUC; WPV, wild poliovirus. recent case in Nigeria in 2014, although it likely originated from a distinct lineage last detected in Borno in 2013 (Nnadi et al., 2017). Due to civil unrest, polio eradication efforts could not reach parts of Borno with vaccination or surveillance for several years, and some areas remain inaccessible (Nnadi et al., 2017). Similar resurfacing of a WPV after multiple years without detection previously occurred in the Sudan and Chad (World Health Organization, 2005a, 2005b, and surveillance continues to occasionally detect WPVs with no closely linked ancestors in Pakistan and Afghanistan (i.e., orphan viruses) (Cowger et al., 2017). Poliovirus surveillance primarily relies on the detection of symptomatic cases of acute flaccid paralysis (AFP) and tests of their stool for the presence of poliovirus. However, only less than 1% of WPV infections in susceptible individuals result in paralytic symptoms (i.e., approximately 1/200, 1/2,000, and 1/1,000 for WPV serotypes 1, 2, and 3, respectively) (Nathanson & Kew, 2010). This implies that hundreds to thousands of infections can occur between successive cases, leading to a situation of silent circulation. While all documented instances of prolonged undetected circulation occurred in the context of likely gaps in surveillance, the possibility of silent poliovirus circulation contributes to the ability of transmission to go undetected.
Motivated by the re-emergence in Borno, a recent study used a stochastic model to examine whether simple, small, hypothetical populations not reached at all by vaccination can perpetuate WPV transmission while experiencing very long intervals between polio cases . The potential for silent circulation influences the confidence about the interruption of poliovirus transmission as a function of time after the last detected poliovirus, which informs the decision to certify wild poliovirus eradication in a country, region, and globally. Global certification of a wild poliovirus serotype affects when the world can coordinate the cessation of the use of homotypic oral poliovirus vaccine (OPV), which represents a necessary step to end all paralytic poliomyelitis disease (polio) caused by polioviruses, because OPV itself can in rare instances cause polio (Duintjer Tebbens et al., 2006;World Health Organization, 2005a, 2005b. Delaying OPV cessation leads to substantial additional costs (Duintjer Tebbens et al., 2011), but premature OPV cessation implies a risk of a WPV reemergence that could become very difficult to control. As of early 2018, of the six World Health Organizations (WHO) regions, only the African and the Eastern Mediterranean regions remain not certified as polio-free due to the possibility of continued WPV transmission in Borno (African region) and continuing polio cases in Pakistan and Afghanistan (Eastern Mediterranean region) (World Health Organization, 2018). Moreover, of the three WPV serotypes, the world certified serotype 2 WPV eradication in 2015 (Global Polio Eradication Initiative, 2015) and did not report any polio cases due to serotype 3 WPV since 2012 (Kew et al., 2014), with only serotype 1 WPV continuing to cause reported polio cases.
Based on epidemiological experience and modeling, certifying a region or the world as polio-free requires no observed polio cases for at least three years (Debanne & Rowland, 1998;Eichner & Dietz, 1996;Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;World Health Organization, 2017). Prior to , numerous modeling studies addressed the possibility of prolonged poliovirus circulation without any detected polio cases in different populations using different metrics, and spanning different model structures and assumptions about surveillance, vaccination, transmissibility, and demographics (Berchenko et al., 2017;Eichner & Dietz, 1996;Famulare, 2016;Houy, 2015;Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;Kalkowska, Duintjer Tebbens, & Thompson, 2012). These models generally agree that in most realistic situations, 3 years without any detected cases implies at least 95% confidence about no circulation, although our studies show that results depend on surveillance quality, serotype, seasonality, and vaccination strategies (Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;Kalkowska et al., 2012). In contrast to the prior literature,  prominently reported (in their abstract) a 22% chance of at least 3 years between successive polio cases in a population of 10,000 people despite a 100% case detection rate. In a correction, Vallejo, Keesling, Koopman, and Singer (2018) report an error in the calculation of the average age of infection that resulted in the use of unrealistically high values for the basic reproduction number (R 0 ), which led them to note the limited applicability of their findings to real-world situations. Nonetheless, the potential use of their findings in deliberations about polio certification motivated us to replicate the  analyses to take a closer look at the model behavior, the interpretation of the results (particularly in the context of the relevant prior literature that they did not consider), and the implications of some of their simplifying assumptions. Although they issued a correction that highlighted one unrealistic model input assumption (Vallejo et al., 2018), we highlight multiple other attributes that they did not consider that further limit the utility of their approach and model for supporting policy, including the unrealistic assumptions about the existence of completely isolated small populations that remain at the endemic equilibrium with no vaccination and yet benefit from perfect surveillance.

Material and methods
We used the model structure, notation, input values, and microsimulation algorithm descriptions provided by  to reconstruct their model in Java using the Eclipse open development platform (Eclipse Foundation, Inc., Ottawa, Ontario, Canada). The first two authors independently reconstructed the model. The model assumes a constant population (i.e., birth rate ¼ death rate) distributed over five compartments: susceptible (i.e., never infected), first infection (susceptible to paralysis), temporary full immunity to infection after recovery from infection, partial susceptibility to infection but not paralysis after waning, and reinfection of partially susceptible individuals. The paralysis-to-infection ratio determines the probability of a first infection leading to a polio case, and the detection rate determines the probability that surveillance would detect a polio case. The model ignores the 1e3-day latent period between virus exposure and becoming infectious to others (Duintjer Tebbens, Pallansch, Kalkowska, et al., 2013), and characterizes the infectious period using a single stage, which implies an exponential recovery process with an unrealistically long tail (Duintjer Tebbens, Pallansch, Kalkowska, et al., 2013;Lloyd, 2001). The model uses the Gillespie algorithm (Gillespie, 1976), which randomly generates when transition events between the compartments occur. Consequently, each stochastic iteration of the model results in different event times. This leads to difficulties in interpretation of the results in Figures 5 and 6 from  that they presented as a function of "event steps,"  and therefore we provided the corresponding results also as a function of time, which provides a more natural interpretation.
Briefly,  explored the model behavior in terms of time between detected cases (TBC), which represents a different metric than those developed in the prior literature (Berchenko et al., 2017;Eichner & Dietz, 1996;Famulare, 2016;Houy, 2015;Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;Kalkowska et al., 2012).  focused on the proportion of stochastic iterations of their model in which any TBC exceeded 3 years. They varied the detection rate, the population size, the basic reproduction number (R 0 ), and the pattern of waning immunity (i.e., "fast-shallow," "intermediate," or "slow-deep" , which considered the same R 0 values, and which Vallejo et al. (2018) identified as unrealistic).  also kept the paralysis-to-infection ratio constant at 0.005, which corresponds to the consensus estimate for serotype 1 WPV (Duintjer Tebbens, Pallansch, Kalkowska, et al., 2013;Nathanson & Kew, 2010). They initialized their simulations at the endemic equilibrium of the related differential equation based (DEB) model and reported results based on a single simulation run of 1,000 stochastic iterations. We emphasize that for polio, the assumption of a small, isolated population with no immunization that starts at the endemic equilibrium with perfect surveillance represents an unrealistic situation (i.e., in reality, such a population cannot indigenously sustain persistent transmission), but for purposes of this analysis we replicated the  model as published.
For each scenario (i.e., defined as a set of assumed values for population size, R 0 , waning immunity pattern , and detection rate), we ran the corresponding DEB model for a period of 100 years to reach its endemic equilibrium. Using this equilibrium, rounded to natural numbers, as the initial condition for the microsimulation model, we ran each stochastic iteration for 15 years or until infection died out, whichever came first. Instead of just running a single simulation and assuming that it represented the truth, we recognized the stochastic nature of the model and the potential for variability in the results as a function of the number of iterations. To characterize the stochastic variability as a function of number of iterations in the simulation, we explored the stochastic variability in the results by performing multiple simulation runs of different numbers of iterations (n i ). We recorded the variation around the mean results of n i ¼ 1,000, 10,000, or 100,000 iterations by running n s ¼ 10 simulations for each n i to demonstrate the range of results possible with a single simulation run of n i iterations. We considered the scenarios of a population (N) of 3,500 or 10,000, R 0 equals 10, a 100% detection rate, and all three patterns of waning immunity. Due to the relative instability we found in the results with relatively small numbers of iterations (i.e., n i ¼ 1,000), we reproduced the  analyses with 100,000 iterations.
We also explored different ways to interpret the results. Notably, the TBC metric informs the probability that consecutive detected cases occur at least 3 years apart in a stochastic iteration, given the detection of cases. Put differently, it answers the question: "given that we know some transmission is occurring in a small population, what is the probability of observing at least 3 years between successive (detected) cases?" However, a more relevant policy question is: "for a small population that we know supported transmission initially, what is the probability that we would observe a gap of at least 3 years between cases?" Answering the latter question requires dividing the number of iterations with an occurrence of a TBC of at least 3 years by the total number of stochastic iterations performed (instead of only those with at least two detected cases). Neither metric directly informs the question of interest for decision makers: "given that we haven't observed any cases for X amount of time, what is the probability that circulation continues?" The prior literature focused on answering this question. Characterizing the probability of circulation still occurring given 3 years without a detected case requires direct comparison of the frequency of 3 years without a detected case occurring despite continued transmission with the frequency of 3 years without a detected case when transmission truly stopped (as recognized by the seminal work by Eichner & Dietz, 1996). Prior studies used the following more appropriate metrics for measuring the confidence about no circulation of poliovirus transmission (Eichner & Dietz, 1996;Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;Kalkowska et al., 2012): POE e "the probability of eradication defined as the fraction of stochastic iterations in which die-out occurs (i.e., prevalence below the transmission threshold)" DEFP e "the detected-event-free period defined as the time in months since the last detected case" CNC e "confidence about no circulation given the DEFP approximated as (1 -the number of DEFPs equal to t months with ongoing WPV circulation, divided by all DEFPs of t months)" CNCx% e "the time when the confidence about no circulation exceeds x% (i.e., CNC95%, CNC99%)" TUC e "the time of undetected circulation after the last detected-event (for those iterations in which extinction occurs)" TUCx% e "the xth percentile of the TUC (i.e., TUC95%, TUC99%)" (Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015).
We note that models that only considered the occurrence of cases regardless of poliovirus detection (Eichner & Dietz, 1996;Kalkowska et al., 2012) used a case-free period (CFP) instead of the DEFP. Subsequent modeling extended the concept of the CFP to the DEFP to allow for the detection of transmission in environmental surveillance, because acute flaccid paralysis surveillance may not detect all cases and environmental surveillance can detect the presence of transmission without cases (Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;Kalkowska, Duintjer Tebbens, Pallansch, & Thompson, 2018). The estimation process evaluates the DEFP and CNC discretely using a time step of months.
Finally, while the assumed mean infectious period of 28 days appears realistic based on the evidence and expert opinion (Duintjer Tebbens et al., 2013;Duintjer Tebbens et al., 2013), the single-stage process implies unrealistic proportions of infected individuals of 34%, 12%, and 4% remaining infectious to others at 1, 2, or 3 months, respectively. Recognizing the issues associated with the use of a single stage infection process (Duintjer Tebbens, Pallansch, Kalkowska, et al., 2013;Lloyd, 2001), we explored the impact of breaking the infection process into either 2 or 4 even stages, which effectively narrows the distribution of the infectious period while maintaining the same mean infectious period.

Results
Fig . 1 shows the total number of infected individuals as a function of time for 10 different stochastic iterations of one model scenario (i.e., R 0 ¼ 10, N ¼ 3,500, and fast-shallow waning). We show 10 iterations to demonstrate the substantial variation between iterations in both the intensity of transmission and the time until elimination of the virus. With all stochastic models, the number of iterations performed will determine the stability of the results, with larger numbers of iterations required to capture the model behavior by averaging over all iterations. Table 1 reports our attempt to replicate the summary statistics for TBC distributions for different scenarios presented in Figs. 2 and 3 of the  paper using a single simulation of 100,000 iterations. Our results show systematically lower results (shorter TBCs) for the interquartile ranges. Table 2 summarizes the results of our analysis of the main results highlighted by . We note that they report (in their Tables 3 and 4) one metric as a proportion and one metric as a percentage, which we replicate while noting that this requires the reader to multiply or divide by 100 to get the metrics on the same scale. To facilitate direct comparison, the first two columns of our Table 2 show the results as reported in Table 3 of  of (i) the proportions of all iterations with more than 1 detected polio case (">1 case/n i ") and (ii) the percent of iterations with more than 1 detected polio case for which any TBC exceeded 3 years (">3 years/>1 case"). As mentioned in the methods section we recognized that the second of these metrics provided inflated estimates of the total probability of observing TBCs of greater than 3 years in the entire simulation given the conditioning on observing more than one case. Specifically, the ">3 years/>1 case" metric alone does not account for the possibility that polio cases may not occur in the first place, and it only provides information about the probability that we may detect silent circulation after more than 3 years in the case that we know that silent circulation is ongoing, which we would not know (and which makes this not a metric of interest to policy makers). Consequently, we added a more informative metric (shown in bold in Tables 2e4) that reports the percent of all iterations for which any TBC exceeded 3 years (">3 years/n i "). We inferred the value of this metric for the  results from the reported results.
Our results shown in Table 2 report the mean values and standard deviations (SDs) that we obtained by running 10 simulations of 1,000, 10,000, or 100,000 iterations. We find similar estimates for the ">1 case/n i " results reported by , but significantly lower values for the ">3 years/1 case" and ">3 years/n i " results. As noted in the methods, substantial variability can arise from stochastic iterations (e.g., Fig. 1), and this can lead to unstable estimates of relatively rare outcomes from the simulations (e.g., a TBC>3 years).  did not consider stochastic variability or explore the importance of running a sufficient number of iterations to obtain stable results (i.e., they used n s ¼ 1, n i ¼ 1,000). Comparison of our results in the right columns of Table 2 to the  results in the left columns show substantial stochastic variation around the means based on 1,000 stochastic iterations, especially for the ">3 years/>1 case" and the ">3 years/n i " metrics, which often involve very small numbers of iterations with a TBC exceeding 3 years. For example, the typical variation around the mean of the ">3 years/>1 case" metric for a population size of 3,500 exceeds the value of the mean, suggesting a long right tail of the distribution of this metric due to a bound of zero at the low end. Consistent with the central limit theorem (Durrett, 2010), increasing the number of iterations by a factor of 100 to 100,000 lowers the standard deviations roughly by a Fig. 1. The variation in the total infected population over time based on 10 iterations for R 0 ¼ 10, N ¼ 3,500, and fast-shallow waning. factor of 10. Given the relatively rare nature of the outcomes of interest for the TBC and other metrics associated with undetected circulation, we suggest that the results from Table 2 indicate the need for more than 1,000 iterations. Table 3 recreates the results from Tables 3 and 4 reported by . We found no significant difference between the proportions of all iterations with more than 1 detected polio case, as expected based on Tables 1 and 2. However, we found between 2 and 16 times lower percentages of iterations with a TBC exceeding 3 years given two or more detected polio cases, compared to , depending on the scenario. The maximum probability of 22% reported by  of at least 3 years between successive polio cases given >1 case in a population of 10,000 people despite a 100% case detection rate decreased to 12% in our independent replication of the model (see Table 3). We emphasize that this change is independent of the correction made by Vallejo et al. (2018), which highlighted the unrealistic nature of an R 0 of 20. Moreover, we emphasize that this metric only considers those iterations in which two or more polio cases occur, which represent a subset of all the possible realizations of endemic transmission that the model produces. If we instead focus on the more appropriate percentage of all iterations for which any TBC exceeds 3 years (representing all possible realizations of what could happen in the model), then our finding suggests a maximum probability of only 9% that more than 3 years go by between successive cases (Table 3, N ¼ 10,000, R 0 ¼ 20, fast-shallow waning, 100% detection, last column). We also believe that  omit the silent circulation from the start of simulation until the first case, which matters particularly when only 1 case occurs and is relevant to the question of silent circulation. Notably, if we count those iterations in which more than 3 years went by until the first case then this leads to a maximum probability of 11.5% (N ¼ 10,000, R 0 ¼ 20, fastshallow waning, not shown). Table 4 highlights the ambiguity about the appropriate probabilistic interpretation of the TBC-based metrics and the omission of information about how many iterations exhibit persistent transmission for 15 years by providing comparisons to Table 1 The distribution of time (years) between detected polio cases for different population sizes, R 0 values, and waning immunity patterns  (based on n i ¼ 100,000). Compare to Figures 1 and 2  a more comprehensive set of metrics, based on prior research (Eichner & Dietz, 1996;Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015;Kalkowska et al., 2012). The metrics from the prior literature represent a better and more appropriate characterization of the outcome of interest to decision makers because they reflect the actual situation by looking at the probability of transmission continuing given the observation of no cases as a function of time. While it may seem subtle, this differs significantly from the TBC. Table 4 shows the results of the POE, CNC95%, CNC99%, TUC95%, and TUC99% for comparison to the  metrics and estimates for N ¼ 3,500 or 10,000 with realistic R 0 ¼ 10. For all patterns of waning immunity, the model predicts POE ¼ 100% for N ¼ 3,500 and suggests time periods of 2.25 (3.08) years or less without polio cases required to achieve 95% (99%) confidence in the interruption of transmission. Also, the time of undetected circulation from the Table 2 The influence of the sample size on of the proportions of all iterations with more than 1 detected case (">1 case/n i "), the percent of iterations with more than 1 detected case and time between detected cases greater than 3 years (">3 years/>1 case"), and the percent of all iterations with time between detected cases greater than 3 years (">3 years/n i ") for different population sizes (i.e., N ¼ 3,500 and 10,000), R 0 ¼ 10, all three waning immunity patterns , and 100% detection rate. n s , n i n s ¼ 1, n i ¼ 1,000 n s ¼ 10, n i ¼ 1,000 n s ¼ 10, n i ¼ 10,000 n s ¼ 10, n i ¼ 100,000   mean (SD) mean (SD) mean (SD) mean (SD) mean (SD) mean (SD) mean (SD) mean (SD) mean (SD) Abbreviations: n s , number of simulations; n i , number of iterations; SD, standard deviation; yrs, years.

Table 3
The proportion of iterations with more than 1 detected case (">1 case/n i "), the percent of iterations with more than 1 detected case and time between detected cases greater than 3 years (">3 years/>1 case"), and the percent of all iterations with time between detected cases greater than 3 years (">3 years/ n i ") for different population sizes, R 0 values, and patterns of waning immunity patterns  with 100% and 50% detection rates (based on n i ¼ 100,000). last detected case until die-out does not exceed 2.07 (2.99) years with 95% (99%) certainty. Increasing the population size (N ¼ 10,000) decreases the POE by less than 1%, while the CNCs and TUCs rise to 3.42 and 2.57 years or less at the 95% confidence level. Table 4 also shows how the results change with the use of a 2-or 4-stage infection process. Applying a more realistic distribution for the infectious period increases the POE to 100% for the larger population (i.e., N ¼ 10,000), while reducing the time to reach 95% (99%) confidence about no circulation by up to 10 months (13 months).  2017), we do not observe drastic fluctuations in the average number of infectious people, which represents a different behavior of our implementation of their model than they reported. Recognizing that the time between events differs in each iteration, which complicates interpretation of the x-axis in Fig. 2aec, Fig. 2def shows the same results plotted on the more intuitive scale of time. Fig. 2def further show that circulation persists for the entire 15 years only in a very small fraction of the iterations, which we emphasize do not include any vaccination. Figs. 3 and 4 show the impact of varying the birth rate on the 100 event step and natural time scales, respectively. Given that the birth rate may influence the length and strength of virus circulation, we suggest that the initial conditions for Fig. 3ced and 4c-d should appropriately change with the introduction of new birth rates because the equilibrium changes (we note that  did not change them), consequently leading to different starting points and shapes of the curves (see long-dashed lines). Thus, we suggest that analyses assumed to start at the endemic equilibrium should use the correct endemic equilibrium for the inputs so as to avoid the effects of the initial conditions representing a different point than the endemic equilibrium.

Discussion
Despite our efforts to carefully follow the model and algorithm in , we could not reproduce their results of a relatively high probability of unusually long silent circulation of polioviruses in small populations. Although the two authors who independently reproduced the model obtained similar results to each other (within the stochastic variation), our results differ from the results reported by  in a way that suggests some difference between the approaches that we could not identify based on the information provided in the published article or in follow up with the authors. Our use of more simulation runs and more iterations does not explain this difference. We found some similar general trends in the length of times between polio cases as a function of R 0 or population size, but our analysis suggests much shorter times between polio cases and much lower probabilities of prolonged undetected circulation with perfect or sub-optimal casebased surveillance.
We also identified several issues with the presentation of the results by . First, focusing in the abstract on the percentage of runs with at least three years between successive polio cases among iterations in which two or more polio cases occurred ignores the inflation implied by conditioning on the probability that circulation continues for long Table 4 The proportion of iterations with more than 1 detected case (">1 cases/n i "), the percent of iterations with more than 1 detected case and time between detected cases greater than 3 years (">3 years/>1 cases"), the percent of all iterations with time between detected cases greater than 3 years (">3 years/n i "), POE, CNC95%, CNC99%, TUC95%, and TUC99% for different population sizes, R 0 ¼ 10, varying waning immunity patterns , and 100% detection rate (based on n i ¼ 100,000). enough in a small population to generate at least two polio cases. While Tables 3 and 4 in  provide the information required to get to the percentage of all runs with at least three years between successive polio cases, the reader must multiply the two numbers to obtain this result. Second, even with that adjustment, as noted above, the TBC metric does not represent the probability of circulation as a function of time without detections, which requires conditioning on the absence of detections. Eichner and Dietz (1996) developed a metric that truly informs the probability of circulation as a function of time without detections, which various other groups subsequently adopted to further address this question (Houy, 2015;Kalkowska, Duintjer Tebbens, Pallansch, et al., 2015). Third, the use of only 1,000 iterations remains insufficient to generate statistically robust findings for the relatively rare events that reflect the tails of highly skewed distributions. Fourth, the average prevalence of infections by 100 event steps remains challenging to interpret, because the time between events varies between individual stochastic iterations and depends strongly on the realized prevalence of infections. Results shown as a function of time ( Fig. 2def and 4) provide a more direct means for poliovirus epidemiologists to verify the plausibility of the findings, as does showing the behavior of individual realizations (Fig. 1).
Besides the questions about the numerical results and their presentation, the findings from  remain limited by the simplifying assumptions of their hypothetical model. Despite the stated focus on investigating the silent circulation of poliovirus in small populations unreached by vaccination, the model ignores the heterogeneity and conditions that affect the ability of polioviruses to circulate in sub-populations and among different age groups. While low-level heterogeneity remains challenging to characterize despite its importance for die-out behavior, consideration of a more realistic aging process and different levels of isolation of an unvaccinated subpopulation within the generally well-vaccinated broader population may offer different insights. Furthermore, given that the results clearly depend on the assumed R 0 , we note the importance of the Vallejo et al. (2018) correction that emphasized that the R 0 upper-bound of 20 does not produce results Fig. 2. The count of the total infected population for N ¼ 10,000 with fast-shallow waning immunity dynamics, and: (a) R 0 ¼ 10, (b) R 0 ¼ 15, (c) R 0 ¼ 20 averaged over n i ¼ 1,000 at every 100 event steps and (d) R 0 ¼ 10, (e) R 0 ¼ 15, (f) R 0 ¼ 20 averaged over n i ¼ 100,000 over time.
consistent with epidemiological observations in any realistic population. In our experience, the maximum R 0 in the context of approximate fast-shallow waning equals 13 among a broad set of populations considered, including northern India (Duintjer Tebbens, Pallansch, Kalkowska, et al., 2013;Kalkowska, Duintjer Tebbens, Grotto, et al., 2015;, although we note the dependence of R 0 values on  the model formulation. For an isolated, rural population in northern Nigeria, we suspect much lower R 0 values (e.g., maximum of 10), which explains why we focused on presenting our results for the R 0 value of 10. Thus, while we appreciate the importance of research that explores the robustness of assumptions, policy recommendations should draw from realistic models that represent the last WPV reservoirs. Moreover, poliovirus transmissibility varies seasonally, which implies greater fluctuations in the prevalence and leads to an increased probability of die-out during each low season compared to the  model, which does not include seasonality. The assumption of small populations that remain untouched by vaccination for up to 15 years also remains highly questionable in the context of ongoing global polio eradication efforts and ignores the possible secondary OPV benefits derived from massive OPV use in the surrounding population. With almost all (>99.7%, see POE in Table 4) simulations resulting in eradication, the assumption of an endemic equilibrium state for perfectly isolated small populations does not appear realistic. This suggests the need to allow at least a small chance of transmission between sub-populations. In a situation of many small sub-populations but with some interaction among them, the chance that any individual sub-population sustains prolonged transmission without experiencing cases increases, which would imply a higher overall probability of silent circulation. However, such prolonged circulation would also imply an increased probability that the virus spreads to other subpopulations where it can cause cases and get detected, which decreases the overall probability of silent circulation. Examining the net effect of these counteracting trends remains a topic of further research. Lastly, the assumed single-stage process implies unrealistic proportions of infected individuals. The precise shape of the duration distribution remains uncertain and does not matter much for questions about vaccination threshold, equilibrium behavior, and cumulative incidence, but it substantially affects poliovirus persistence in small populations ( Table 4).
Instead of concentrating solely on variation of the inputs using a model that simplifies much of the complexity related to poliovirus transmission, we suggest that future research of silent circulation in small populations should focus on realistic model structures and model inputs representative of the last WPV reservoirs.

Declarations of interest
None.