Open Access
Open Peer Review

This article has Open Peer Review reports available.

How does Open Peer Review work?

Spatial model for risk prediction and sub-national prioritization to aid poliovirus eradication in Pakistan

  • Laina D. Mercer1Email authorView ORCID ID profile,
  • Rana M. Safdar2,
  • Jamal Ahmed2, 3,
  • Abdirahman Mahamud2, 3,
  • M. Muzaffar Khan2, 3,
  • Sue Gerber4,
  • Aiden O’Leary2, 5,
  • Mike Ryan2,
  • Frank Salet2,
  • Steve J. Kroiss1,
  • Hil Lyons1,
  • Alexander Upfill-Brown6 and
  • Guillaume Chabot-Couture1
BMC Medicine201715:180

https://doi.org/10.1186/s12916-017-0941-2

Received: 7 April 2017

Accepted: 6 September 2017

Published: 11 October 2017

Abstract

Background

Pakistan is one of only three countries where poliovirus circulation remains endemic. For the Pakistan Polio Eradication Program, identifying high risk districts is essential to target interventions and allocate limited resources.

Methods

Using a hierarchical Bayesian framework we developed a spatial Poisson hurdle model to jointly model the probability of one or more paralytic polio cases, and the number of cases that would be detected in the event of an outbreak. Rates of underimmunization, routine immunization, and population immunity, as well as seasonality and a history of cases were used to project future risk of cases.

Results

The expected number of cases in each district in a 6-month period was predicted using indicators from the previous 6-months and the estimated coefficients from the model. The model achieves an average of 90% predictive accuracy as measured by area under the receiver operating characteristic (ROC) curve, for the past 3 years of cases.

Conclusions

The risk of poliovirus has decreased dramatically in many of the key reservoir areas in Pakistan. The results of this model have been used to prioritize sub-national areas in Pakistan to receive additional immunization activities, additional monitoring, or other special interventions.

Keywords

Disease mappingPolio eradicationRisk mappingSpatial epidemiologyHurdle modelsPakistanRisk prioritizationVaccination campaignsSupplementary immunization activities

Background

The Global Polio Eradication Initiative (GPEI) has seen great success since its launch in 1988. At the time of this writing, only Pakistan, Afghanistan, and Nigeria remain endemic for polio with only 37 cases of wild poliovirus serotype 1 (WPV1) recorded in 2016. The 20 cases reported in Pakistan in 2016 represent a historically low case count for a calendar year and a 63% reduction in cases compared to 2015. However, transmission is still occurring on a considerable geographic scale, with four of eight provinces reporting WPV1 cases in 2016. Pakistan has approximately 25 million children under the age of 5 years [1], which presents challenges for allocating limited resources. As the program approaches the goal of zero cases, identifying the districts that are most likely to be infected and prioritizing those districts for interventions is a priority of the program and should accelerate the path towards eradication.

In the Eastern Mediterranean Region of the World Health Organization (WHO) key strategies of the of GPEI included (1) achieving high coverage of at least three doses of oral polio vaccine (OPV), (2) implementation of supplementary immunization activities (SIAs), and (3) the development of sensitive epidemiological and laboratory surveillance using standard WHO definitions [2]. Eradication efforts began in Pakistan in 1994, when the first SIA was conducted, and in 1995, when acute flaccid paralysis (AFP) surveillance commenced [3]. The case burden in Pakistan has dramatically decreased since the 1990s, but sustained implementation of the first two GPEI strategies has been challenging due to poor rates of routine immunization (RI) and security issues. As of 2012, only 53% of children in Pakistan were receiving all basic vaccines, including Bacillus Calmette–Guerin, measles and three doses of polio and diphtheria, pertussis, and tetanus, with provincial rates as low as 16% and 29% in Balochistan and Sindh provinces, respectively [4]. Additionally, vaccination bans and security limitations in the Federally Administered Tribal Areas (FATA) and Khyber Pakhtunkhwa (KP) and violence against vaccinators in FATA, Balochistan, and Sindh provinces have periodically limited the program’s efforts to consistently implement SIAs with high population coverage in conflict-affected areas since 2008 [5]. These programmatic challenges have resulted in pockets of underimmunized (fewer than three OPV doses) children and have allowed transmission to persist.

Since 2011, the Pakistan program has been implementing and enhancing a National Emergency Action Plan (NEAP) for polio eradication to improve management and accountability strategies, highlight core reservoirs of transmission, and to ensure the program is creating and using high quality data [5, 6]. To aid in the prioritization of sub-national areas for programmatic interventions, we developed a spatial model to estimate the risk of future WPV1 cases for the 155 districts of Pakistan. Previous studies highlight the utility of spatial risk models for guiding programmatic interventions for polio, such as the 86% accuracy for predicting districts at risk for future WPV1 cases in Nigeria [7]. These models guided the prioritization of sub-national areas for immunization planning and allocation of technical and administrative field personnel. The use of a spatial risk model, which is statistically evaluated based on its accuracy for predicting locations of cases, represents a methodological departure from the common approach of compiling programmatic indicators of disease risk, assigning weights based on expert opinion, and linearly combining into a risk score [8, 9]. Unfortunately, a spatial risk model has not previously been applied to model the risk of WPV1 in Pakistan. In this paper, we will describe our efforts to model the risk of future WPV1 cases in districts of Pakistan and describe how these efforts have been incorporated by the National Emergency Operating Centre (N-EOC) in Islamabad to the 2016–2017 NEAP to prioritize the districts of Pakistan.

Methods

Description of data

Our modeling efforts rely on the Pakistan AFP surveillance data, which is managed by WHO and is the source for information on paralytic cases of polio. AFP cases are identified through the extensive surveillance network in Pakistan, which relies on both passive and active (sweeps and case search) surveillance with oversite at the district, provincial, and national level [6]. However, AFP can be caused by many viruses other than WPV [10]. For each childhood AFP case, the initial investigation includes two stool samples (to be tested for poliovirus), demographics, date of paralysis onset, history of OPV through RI as well as SIAs, and clinical symptoms. For the purposes of our analyses, AFP cases that are found to be positive for WPV1 are treated as the outcome of interest and AFP found not to be infected with polio, non-polio AFP (NPAFP), are treated as a random sample from the population [1113]. A surveillance system that achieves an annual rate of 2 or more NPAFP per 100,000 children under age 15 is considered sufficient to detect circulating WPV1 [14].

Statistical methods

Using the NPAFP dose histories we constructed estimates of district-level vaccination rates, such as zero-dose RI (received zero OPV doses from RI) and underimmunization (defined as children who have received three or fewer doses of OPV from birth until AFP investigation). However, the number of NPAFP observations in each district within a 6-month period can be quite small, ranging from 0 to 74 in our data, resulting in differences during a 6-month period that are implausible when considering children under 5 years of age. To alleviate this problem, we used a hierarchical Bayesian space–time model with a temporally structured space–time interaction [15] to generate smoothed estimates of district-level vaccination rates. Briefly, the space–time model borrows information over space and time to estimate the underlying rates from which the observed data was drawn. Additional details about the model specifications and interaction selection can be found in Section 1 of Additional file 1.

A dynamic immunity model, originally developed for Northern Nigeria [16], was implemented using Pakistan-specific vaccine efficacies [13] to estimate district-level population immunity. The immunity model uses a hierarchical Bayesian model to estimate district-level, annual, age-specific SIA coverage using the ages and SIA doses reported by NPAFP cases. A hypothetical cohort then progresses through the true SIA calendar experiencing the estimated yearly age-specific SIA coverage rates and vaccine type. The result is monthly district-level serotype-specific estimates of under-five population immunity which have been averaged over each 6-month period.

Polio incidence has been rare and spatially heterogeneous in Pakistan since 2003, with 91% (3804 out of 4185) of 6-month district-level observations reporting zero WPV1 cases. To account for this large number of zeros we implemented a spatial Poisson hurdle model [17] using a hierarchical Bayesian framework. The Poisson hurdle model explicitly models the excess observations of zero cases we would expect from areas without circulating virus with a Bernoulli component as well as the total number of cases given at least one case with a truncated Poisson. The probability of at least one WPV1 case (Bernoulli) and the total number of WPV1 cases given at least one case (truncated Poisson) in a district during a 6-month period were jointly modeled as a function of covariates from the previous period, namely a set of independent and spatially structured random effects, also known as the convolution model [18], and an observation-level random effect to account for overdispersion. Notably, a bivariate prior distribution was assigned to the independent district random effects to allow for a correlation between the district-level Bernoulli and truncated Poisson models. The expected number of WPV1 cases, which is defined as the product of the probability of at least one WPV1 case and the expected number of WPV1 cases given at least one case, was used as the measure of risk. A similar modeling approach has been used to predict WPV1 and WPV3 in Nigeria [7]. Full modeling details, including prior distribution specifications, can be found in Section 2 of Additional file 1.

Model selection was carried out in two stages. In the first stage, we fit models based on the 64 combinations of seasonality (low season is January through June), vaccine-derived immunity for type 1 poliovirus, underimmunized fraction, zero dose RI fraction, recent neighboring cases (defined as paralytic WPV1 cases in neighboring areas in the previous 6 months and square root transformed), and recent cases (defined as paralytic WPV1 cases and square root transformed) for the Bernoulli portion, with only random effects in the truncated Poisson model, and selected the model that minimized the deviance information criteria (DIC) [19]. In the second stage, we included the selected covariates for the Bernoulli portion and considered all combinations of covariates for the truncated Poisson portion of the model and then again selected the covariates which minimized the DIC.

The predictive accuracy of our model was assessed by comparing our predictions to held out data. The selected model was first fit using data from 2003 until the first 6 months of 2009 and then used to predict cases in the second half of 2009. This procedure was repeated, using all historical data, for each subsequent 6-month period with the final predictions for the second half of 2016 relying on data from 2003 through the first 6 months of 2016. For each set of predictions, the area under the curve (AUC) for the receiver operator characteristic (ROC) curve was calculated. The AUC for the ROC curve is a single value that summarizes the shape and position of the ROC curve and has the useful interpretation as the probability that a randomly selected district with a case will have a higher risk score than a randomly selected district without a case [20, 21]. AUCs between 0.7 and 0.9 suggest a moderate predictive power and AUCs above 0.9 suggest a strong predictive power [22].

All modeling was completed in R 3.2.1 [23]. The SIA coverage models were fitted via Markov chain Monte Carlo using RStan [24, 25]. The space–time models used for smoothing vaccination rates and the spatial Poisson hurdle model were fit using the Integrated Nested Laplace Approximation (INLA) [26, 27] as implemented in the INLA package [28, 29], which provides a fast and accurate alternative to Markov chain Monte Carlo methods for these type of space–time models [15, 3032].

Results

Estimated dose history and vaccine-derived immunity

AFP surveillance in Pakistan collected data on 43,301 NPAFP cases between January 2003 and June 2016, with an average annual rate increasing from 4.3 to 11.4 NPAFP per 100,000 children under the age of 5 years from 2003 to 2016. Space–time smoothing models fit to the NPAFP vaccination dose history data indicated that zero dose RI and underimmunized rates (fewer than three doses) are highly heterogeneous across Pakistan (Figs. 1 and 2). Both zero dose RI and underimmunization rates were high in most of Punjab, Sindh, and KP provinces, and lowest in the western provinces, Balochistan and FATA. Zero dose RI rates were generally higher than underimmunization rates, suggesting apparent gaps in RI. Our dynamic immunity model indicates that immunity is spatially and temporally variable across Pakistan (Fig. 3) due to differences in SIA coverage and the vaccine efficacies for each serotype, which vary among the vaccines [13].
Fig. 1

A map of smoothed estimates of zero routine immunization (RI) doses for January through June, 2016 (left) and an example of the smoothing models for observed zero RI rates for Khyber district in the Federally Administered Tribal Areas (FATA) of Pakistan from 2003 to 2016 (right)

Fig. 2

A map of smoothed estimates of underimmunized fraction for January through June, 2016 (left) and observed underimmunized fraction with smoothed estimates for Khyber district from 2003 to 2016 (right)

Fig. 3

A map of the vaccine-derived population immunity for type 1 poliovirus estimated by the dynamic immunity model as of June 30, 2016 (left), and the dynamic immunity traces and 95% credible interval for immunity in Khyber with supplementary immunization activity (SIA) calendar from 2003 to 2016 (right). The dashed marks along the horizontal axis show the timing of the SIAs and the color represents the vaccine used in the SIA, where mOPV1 represents monovalent oral polio vaccine (OPV) for serotype 1, bOPV represents bivalent OPV for serotypes 1 and 3, tOPV represents trivalent OPV (serotypes 1-3), and IPV is the inactivated poliovirus vaccine

Association between covariates and cases of WPV1

The model that achieved the lowest DIC for the Bernoulli and truncated Poisson portions of the hurdle model are shown in Table 1. Associations are in the expected directions with high season, higher underimmunization rates, recent cases, and recent neighboring cases being positively associated with at least one case (odds ratios above 1) and higher rates of population immunity being associated with a lower probability of a case. The number of cases given at least one case is best described by high season, underimmunization rate, zero-dose RI rate, recent cases, and recent neighboring cases. As expected, season, underimmunization rates, zero-dose RI rate, and recent neighboring cases are positively associated with cases (relative rates above 1).
Table 1

Posterior medians and 95% credible intervals of covariates selected for final risk model

 

Bernoulli

Truncated Poisson

Indicator

Odds ratio

95% CI

Relative rate

95% CI

Season (reference Jan–June)

3.13

(2.44–4.03)

2.52

(1.69–3.82)

Type 1 immunity (10% difference)

0.89

(0.81–0.99)

  

Under immunization (10% difference)

1.12

(1.01–1.24)

1.28

(1.09–1.50)

Zero RI doses (10% difference)

  

1.24

(1.11–1.38)

Sqrt. recent cases

1.78

(1.44–2.19)

1.08

(0.90–1.30)

Sqrt. recent neighbor cases

1.18

(1.05–1.32)

1.31

(1.15–1.50)

Odds ratios are reported for the Bernoulli portion of the model and relative rates for the truncated Poisson

Predictive performance

The predictive accuracy of our selected model was high (Fig. 4). AUC values were consistently above 0.8 and often above 0.9. Over the past 3 years, the AUC ranged from 0.84 to 0.97, with a mean of 0.90. Furthermore, the sensitivity of a list made up of the districts with the 50 highest risk scores, that is the proportion of districts with a case that are ranked in the top 50 risk scores, ranged from 0.78 to 0.97 over the past 3 years. Both measures show good performance overall, aside from the second half of 2012, which had more districts in KP with cases than our model would have predicted.
Fig. 4

Area under the curve, the probability that a randomly selected district with a case will have a higher risk score than a randomly selected district without a case, for prediction of WPV1 cases by district as predicted based on model data from 2003 until 6 months prior to observed data (left) and sensitivity, or true positive rate, of a list containing the top 50 high risk districts for each time point (right)

Predicted risk for WPV1

The probability of at least one WPV1 case is highest along the regions of KP, FATA, and Balochistan, which border Afghanistan, in northern Sindh, and near Karachi (Fig. 5). The expected number of cases given at least one case, which is driven primarily by population size and population immunity, is highest in areas with low immunity and high population. The final risk score, defined as the expected number of cases and described in Section 2 of Additional file 1, indicates that the highest risk areas are found along the border with Afghanistan in the Quetta block (the districts of Quetta, Pishin, and K. Abdullah), in Northern Sindh, and near Karachi.
Fig. 5

The probability of at least one case (left), expected number of cases given at least one case (center), and the overall risk score (expected number of cases) (right), for July through December, 2016

Programmatic impact of the risk model

The N-EOC, in coordination with Provincial EOCs, have used the outcomes of the model to assess the overall risk profile of the program. To operationalize the risk modeling results, a combination of the modeling output described here, genetic sequencing results (including virus isolated from environmental and healthy children stool samples), and local knowledge of access and security status of districts and Union Councils, as well as sociocultural links within Pakistan and across the border with Afghanistan, were used to classify the 155 districts into four tiers (Fig. 6). Tier 1 represents the top 11 districts considered ‘core reservoirs’ of poliovirus in Pakistan. Through genetic sequencing data, there is strong evidence of persistent circulation of local lineages of WPV for at least 2 years in the core reservoirs [6]. Tier 2 consists of the next 33 districts, considered high risk districts, while tier 3 districts are considered vulnerable districts. Tier 4 includes all other districts. Table 2 shows the programmatic implications of the tier classification on the immunization and operational strategies deployed.
Fig. 6

The final risk tier classification to be included in the National Emergency Action Plan for polio eradication in Pakistan for 2016–2017 with Karachi in the inset. This final list incorporates modeling output, genetic sequencing results, and local knowledge of access and security

Table 2

Programmatic implications of risk classifications of districts for July 2016 to June 2017

Tier

Number of districts

Target population (%)

Goal

Strategy

1

11

4,042,214 (11%)

Interrupt endemic and/or persistent local transmission using multiple strategies

NID + SNID + CBV in selected UCs + Priority 1 for combined bOPV/IPV SIA + RI service delivery support and other auxiliary support

2

33

5,746,129 (16%)

Interrupt transmission if transmission is ongoing, decrease vulnerability

NID + SNID + CBV in selected UCs + Priority 2 for bOPV/IPV SIA + RI service delivery support + other auxiliary support

3

24

7,246,474 (20%)

Decrease vulnerability

NID + SNID

4

87

19,638,741 (54%)

Maintain high population immunity

NID only

Tier classifications dictate inclusion in National Immunization Days (NIDs), Sub-National Immunization Days (SNIDs) and employment of community-based vaccination (CBV) strategies in selected Union Councils (UCs). Reproduction of Panel 1 in the National Emergency Action Plan ([6], p. 19)

IPV inactivated polio vaccine, OPV oral polio vaccine, SIA supplementary immunization activity, RI routine immunization

Discussion

Our results indicate that seasonality, immunity, underimmunization rate, recent cases, and recent cases in a neighboring district are most predictive of at least one WPV1 case. The number of cases given at least one case is similarly predicted by seasonality, underimmunization rate, zero RI dose rate, recent cases, and recent neighboring cases. The point estimates of all associations were in the expected direction of lower immunity and dose history as well as recent cases being associated with increased risk of cases and larger outbreaks.

Our modeling efforts suggest that the large outbreaks in 2014 and the recent improvements over the past 2 years can be described by population immunity driven primarily by SIAs in FATA and KP provinces. For example, the dramatic decrease in serotypes 1 and 2 immunities in KP beginning in 2010 were partially the result of declining SIA coverage rates (despite the high frequency of campaigns). Furthermore, type 2 immunity declined due to infrequent trivalent OPV campaigns, whereas type 1 immunity declined because SIAs were primarily using bivalent and trivalent OPV, which has a lower efficacy for type 1 than monovalent OPV type 1 vaccine. Finally, improvements beginning in 2015 were due to improved vaccinator access driven by military intervention in FATA [33].

Based on the recommendations in the 2016–2017 NEAP, finalized in May of 2016, Tier 1–3 districts participated in four bOPV SIAs in addition to the five national bOPV SIAs that covered all districts. Additionally, Tier 1 districts scaled up the community-based vaccination strategy, which employs local individuals, primarily women who are thought to have better access to children within homes, as permanent vaccinators within their communities. In the 6-month period between July and December, 2016, only two of the 44 Tier 1 and 2 districts experienced cases (one each), which reflects well on the efforts focused on those districts. Two of the three other districts, which reported cases during this time frame, were classified as Tier 4 districts (four of five cases in Tier 4), although they were ranked in the top 30 per the risk model; the remaining Tier 4 district, in northern KP province, would be considered relatively surprising from a modeling and programmatic perspective. We emphasize that even Tier 4 districts received considerable programmatic attention, with five planned SIAs across 2016.

Our approach does have several limitations. We have developed a model on a 6-month time scale, which is programmatically relevant but does not align with the approximately 1 month infectious period estimated for poliovirus [34]. Additionally, we have modeled observed WPV1 paralytic cases that only represent approximately 0.5% of WPV1 infections [35]. This absence of cases could be misleading if circulation is silent due to surveillance failure, waned mucosal immunity among older children or adults [36, 37] or, as observed in Israel, high rates of humoral immunity due to exclusive inactivated polio vaccine use leading to low mucosal immunity [38]. Finally, areas that are not explained well by the covariates will have large residual risk captured by the random effects and, as these are invariant in time, we will likely overestimate risk in areas with a long history of WPV1 cases, despite improvements in indicators.

In 2009, the Pakistan program initiated environmental surveillance (ES) to compliment AFP surveillance. Since 2009, the program has grown from 47 samples across 6 sites to 648 samples across 62 sites covering 33 districts in 2016 [14, 39]. Due to the selective and expanding deployment of ES over time and the unique interpretation of ES positives, which signify at least one infection, ES data it is not easily included as a predictor in our currently modeling framework. Similarly, as sampling locations have been selected based on a history of persistent infection [39], we would expect positive samples to primarily reinforce the high-risk status of areas with a history of cases if included as an outcome. It is a limitation of our case-based statistical model that there is not a straightforward way to incorporate the ES data as either a predictor or an outcome. Alternative approaches, such as a transmission model that includes genetic information or a statistical model that incorporates the relative sensitivity of ES and AFP, may be better suited to incorporate the ES data. However, although it is not explicitly included in the risk model for WPV1 cases, ES is used extensively within the program to identify infected areas when transmission is low and is essential for assessing progress towards eradication.

Despite the limitations, our modeling approach provides a principled framework for ranking districts for risk classification that performs well as measured by AUC and sensitivity. In practice, the risk analysis generally identifies the same Tier 1 districts as the aggregation of programmatic knowledge and scientific intuition of the N-EOC members and the greatest impact of the modeling approach is the promotion of districts that appear susceptible based on immunity profile, but have not yet had WPV1 cases, to a lower (higher risk) tier. These promotions impacted the allocation of resources by deploying the community-based vaccination strategy in additional districts and broadening the geographic scope of the four sub-national SIAs. Additionally, the method provides a metric for quantifying the absolute risk and changes in risk over time, which is not always captured well by cases or intuition exclusively.

Conclusion

This study serves as the first use of a spatial model for risk prediction and sub-national prioritization to aid in polio eradication in Pakistan. This risk modeling approach has been applied to the WPV1 case history and NPAFP dose histories to generate risk predictions and a ranked list of districts within Pakistan, which were subsequently used by the N-EOC to help assign districts to risk tiers as part of their NEAP for polio eradication. This approach will help maximize the impact of the resources available for polio eradication efforts in Pakistan.

Abbreviations

AFP: 

Acute flaccid paralysis

AUC: 

Area under the curve

DIC: 

Deviance information criteria

ES: 

Environmental surveillance

FATA: 

Federally Administered Tribal Areas

GPEI: 

Global Polio Eradication Initiative

INLA: 

Integrated Nested Laplace Approximation

KP: 

Khyber Pakhtunkhwa

NEAP: 

National Emergency Action Plan

N-EOC: 

National Emergency Operation Centre

NPAFP: 

Non-polio acute flaccid paralysis

OPV: 

Oral polio vaccine

RI: 

Routine immunization

ROC: 

Receiver operator characteristic

SIA: 

Supplementary immunization activity

WHO: 

World Health Organization

WPV1: 

Wild poliovirus type 1

Declarations

Acknowledgements

The authors would like to thank Bill and Melinda Gates for their sponsorship through the Global Good Fund. We would also like to thank all those in Pakistan involved in the AFP surveillance and laboratory testing.

Funding

LDM, SK, HL, and GCC are supported by Bill and Melinda Gates through the Global Good Fund. The funders had no role in study design, data collection, data analysis, the decision to publish, or preparation of the manuscript.

Availability of data and materials

The data that support the findings of this study are available from the National Emergency Operating Centre (N-EOC) in Islamabad, Pakistan, but restrictions apply to the availability of these data, which were used with permission for the current study, and are therefore not publicly available. However, data are available from the authors upon reasonable request and with permission of the N-EOC.

Authors’ contributions

All authors contributed to the conceptual design of the study. LDM completed all statistical analyses and drafted the manuscript. MRS, JA, AM, SG, AO, and MR contributed to the literature search, data interpretation, and writing. FS and MMK contributed to the data collection, data interpretation, and writing. SK, HL, AUB, and GCC contributed to data analysis and writing. All authors have read and approved the final manuscript.

Ethics approval and consent to participate

Permission to access the acute flaccid paralysis and vaccination campaign data was granted by the National Emergency Operation Centre in Islamabad, Pakistan.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Authors’ Affiliations

(1)
Institute for Disease Modeling
(2)
National Emergency Operations Centre for Polio Eradication
(3)
World Health Organization
(4)
Bill and Melinda Gates Foundation
(5)
United Nations Children’s Fund (UNICEF)
(6)
David Geffen School of Medicine, University of California

References

  1. United Nations, Department of Economic and Social Affairs PD. World Population Prospects: The 2015 Revision, Key Findings and Advance Tables. Working Paper 2015. https://esa.un.org/unpd/wpp/publications/files/key_findings_wpp_2015.pdf. Accessed 1 Apr 2017.
  2. Centers for Disease Control and Prevention. Progress toward poliomyelitis eradication--eastern Mediterranean region, 1988-1994. MMWR Morb Mortal Wkly Rep. 1995;44:809–11. 817–8.Google Scholar
  3. Centers for Disease Control and Prevention. Progress toward poliomyelitis eradication--Pakistan, 1994-1998. MMWR Morb Mortal Wkly Rep. 1999;48:121–6.Google Scholar
  4. National Institute of Population Studies (NIPS) (Pakistan) and ICF International. Pakistan Demographic and Health Survey 2012-13. Islamabad, Pakistan, and Calverton, Maryland, USA: NIPS and ICF International; 2013.Google Scholar
  5. Alexander JP, Zubair M, Khan M, Abid N, Durry E. Progress and peril: poliomyelitis eradication efforts in Pakistan, 1994-2013. J Infect Dis. 2014;S152–61Google Scholar
  6. Government of Pakistan. National Emergency Action Plan for Polio Eradication 2016/2017. 2016. http://polioeradication.org/wp-content/uploads/2016/11/NEAP2016-2017_Pakistan.pdf. Accessed 1 Apr 2017.
  7. Upfill-Brown AM, Lyons HM, Pate MA, Shuaib F, Baig S, Hu H, et al. Predictive spatial risk model of poliovirus to aid prioritization and hasten eradication in Nigeria. BMC Med. 2014;12:92.View ArticlePubMedPubMed CentralGoogle Scholar
  8. Lowther SA, Roesel S, O’Connor P, Landaverde M, Oblapenko G, Deshevoi S, et al. World Health Organization regional assessments of the risks of poliovirus outbreaks. Risk Anal. 2013;33:664–79.View ArticlePubMedGoogle Scholar
  9. Lam E, Schluter WW, Masresha BG, Teleb N, Bravo-Alcántara P, Shefer A, et al. Development of a district-level programmatic assessment tool for risk of measles virus transmission. Risk Anal. 2017;37(6):1052–62.View ArticlePubMedGoogle Scholar
  10. Marx A, Glass JD, Sutter RW. Differential diagnosis of acute flaccid paralysis and its role in poliomyelitis surveillance definition of acute flaccid paralysis. Epidemiol Rev. 2000;22:298–316.View ArticlePubMedGoogle Scholar
  11. Jenkins HE, Aylward RB, Gasasira A, Donnelly CA, Abanida EA, Koleosho-Adelekan T, et al. Effectiveness of immunization against paralytic poliomyelitis in Nigeria. N Engl J Med. 2008;359:1666–74.View ArticlePubMedGoogle Scholar
  12. Jenkins HE, Aylward RB, Gasasira A, Donnelly CA, Mwanza M, Corander J, et al. Implications of a circulating vaccine-derived poliovirus in Nigeria. N Engl J Med. 2010;362:2360–9.View ArticlePubMedGoogle Scholar
  13. O’Reilly KM, Durry E, ul Islam O, Quddus A, Abid N, Mir TP, et al. The effect of mass immunisation campaigns and new oral poliovirus vaccines on the incidence of poliomyelitis in Pakistan and Afghanistan, 2001–11: a retrospective analysis. Lancet. 2012;380:491–8.View ArticlePubMedPubMed CentralGoogle Scholar
  14. Maes EF, Diop OM, Jorba J, Chavan S, Tangermann RH, Wassilak SGF. Surveillance systems to track progress toward polio eradication — Worldwide, 2015–2016. MMWR Morb Mortal Wkly Rep. 2017;66:359–65.View ArticlePubMedGoogle Scholar
  15. Schrödle B, Held L. Spatio-temporal disease mapping using INLA. Environmetrics. 2011;22:725–34.View ArticleGoogle Scholar
  16. Upfill-Brown AM, Voorman A, Chabot-Couture G, Shuaib F, Lyons HM. Analysis of vaccination campaign effectiveness and population immunity to support and sustain polio elimination in Nigeria. BMC Med. 2016;14:60.View ArticlePubMedPubMed CentralGoogle Scholar
  17. Neelon B, Ghosh P, Loebs PF. A spatial Poisson hurdle model for exploring geographic variation in emergency department visits. J R Stat Soc Ser A. 2013;176:389–413.View ArticleGoogle Scholar
  18. Besag J, York J, Mollié A. Bayesian image restoration, with two applications in spatial statistics. Ann Inst Stat Math. 1991;43:1–59.View ArticleGoogle Scholar
  19. Spiegelhalter DJ, Best NG, Carlin BP, der Linde A. Bayesian measures of model complexity and fit (with discussion). J R Stat Soc Ser B. 2002;64:583–616.View ArticleGoogle Scholar
  20. Hanley JA, McNeil BJ. The meaning and use of the area under a receiver operating characteristic (ROC) curve. Radiology. 1982;143:29–36.View ArticlePubMedGoogle Scholar
  21. Hajian-Tilaki K. Receiver operating characteristic (ROC) curve analysis for medical diagnostic test evaluation. Casp J Intern Med. 2013;4:627–35.Google Scholar
  22. Greiner M, Pfeiffer D, Smith R. Principles and practical application of the receiver-operating characteristic analysis for diagnostic tests. Prev Vet Med. 2000;45:23–41.View ArticlePubMedGoogle Scholar
  23. R Core Development Team. R: a language and environment for statistical computing, 3.2.1. Vienna, Austria; 2016. https://www.r-project.org/. Accessed 1 Apr 2017.
  24. Carpenter B, Gelman A, Hoffman M, Lee D, Goodrich B, Betancourt M, et al. Stan: a probabilistic programming language. J Stat Softw. 2017;76:Issue 1. doi:10.18637/jss.v076.i01.View ArticleGoogle Scholar
  25. Stan Development Team. Stan Modeling Language User’s Guide and Reference Manual, Version 2.8.0. 2015. http://mc-stan.org/. Accessed 1 Apr 2017.
  26. Rue H, Martino S, Chopin N. Approximate Bayesian inference for latent Gaussian models by using integrated nested Laplace approximations. J R Stat Soc Ser B. 2009;71:319–92.View ArticleGoogle Scholar
  27. Lindgren F, Rue H, Linström J. An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic differential equation approach (with discussion). J R Stat Soc Ser B. 2011;73:423–98.View ArticleGoogle Scholar
  28. Lindgren F, Rue H. Bayesian spatial modelling with R-INLA. J Stat Softw. 2015;63:Issue 19. doi:10.18637/jss.v063.i19.View ArticleGoogle Scholar
  29. Arab A. Spatial and spatio-temporal models for modeling epidemiological data with excess zeros. Int J Environ Res Public Health. 2015;12:10536–48.View ArticlePubMedPubMed CentralGoogle Scholar
  30. Ferkingstad E, Rue H. Improving the INLA approach for approximate Bayesian inference for latent Gaussian models. Electron J Stat. 2015;9:2706–31.View ArticleGoogle Scholar
  31. Fong Y, Rue H, Wakefield JC. Bayesian inference for generalized linear mixed models. Biostatistics. 2010;11:397–412.View ArticlePubMedPubMed CentralGoogle Scholar
  32. Held L, Schrödle B, Rue H. Posterior and cross-validatory predictive checks: A comparison of MCMC and INLA. In: Kneib T, Tutz G, editors. Statistical Modelling and Regression Structures. Berlin: Springer-Verlag; 2010. p. 91–110.View ArticleGoogle Scholar
  33. Ahmad SO, Bux AS, Yousuf F. Polio in Pakistan’s North Waziristan. Lancet Glob Heal. 2015;3(1), e15.View ArticleGoogle Scholar
  34. Alexander JP, Gary HE, Pallansch MA. Duration of poliovirus excretion and its implications for acute flaccid paralysis surveillance: a review of the literature. J Infect Dis. 1997;175:S176–82.View ArticlePubMedGoogle Scholar
  35. Nathanson N, Kew OM. From emergence to eradication: the epidemiology of poliomyelitis deconstructed. Am J Epidemiol. 2010;172:1213–29.View ArticlePubMedPubMed CentralGoogle Scholar
  36. Grassly NC, Jafari H, Bahl S, Durrani S, Wenger J, Sutter RW, et al. Asymptomatic wild‐type poliovirus infection in India among children with previous oral poliovirus vaccination. J Infect Dis. 2010;201:1535–43.View ArticlePubMedGoogle Scholar
  37. Koopman JS, Henry CJ, Park JH, Eisenberg MC, Ionides EL, Eisenberg JN. Dynamics affecting the risk of silent circulation when oral polio vaccination is stopped. Epidemics. 2017. Ahead of print. doi:10.1016/j.epidem.2017.02.013.
  38. Shulman L, Gavrilin E, Jorba J, Martin J, Burns C, Manor Y, et al. Molecular epidemiology of silent introduction and sustained transmission of wild poliovirus type 1, Israel, 2013. Eurosurveillance European Centre for Disease Prevention and Control (ECDC). 2014;19:20709.Google Scholar
  39. Alam MM, Shaukat S, Sharif S, Angez M, Khurshid A, Malik F, et al. Detection of multiple cocirculating wild poliovirus type 1 lineages through environmental surveillance: impact and progress during 2011-2013 in Pakistan. J Infect Dis. 2014;210:S324–32.View ArticlePubMedGoogle Scholar

Copyright

© The Author(s). 2017

Advertisement