Evolution and Transmission Dynamics of Wild Poliovirus in Pakistan and Afghanistan (2012-2023): implications for global polio eradication

doi:10.21203/rs.3.rs-5331855/v1

Download PDF

Article

Evolution and Transmission Dynamics of Wild Poliovirus in Pakistan and Afghanistan (2012-2023): implications for global polio eradication

https://doi.org/10.21203/rs.3.rs-5331855/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Despite concerted global vaccination efforts, wild poliovirus remains endemic in two countries in 2024, Pakistan and Afghanistan. This study uses phylogeographic analysis of poliovirus genetic and epidemiological data from clinical and wastewater surveillance to identify the causes of poliovirus persistence and routes of spread over the last decade (2012 to 2023). Poliovirus genetic diversity declined post-2020, with the die-out of one of two major genetic clusters and recent detections are now close genetic relatives. We find that high-risk and hard-to-access regions have been key to persistent polio transmission over the past decade, supporting transmission when it had stopped elsewhere. Karachi, one of the most densely populated cities globally, has acted as a major hub for the amplification and spread of poliovirus to other regions, many of which we show to be dead-end for onwards transmission despite frequent virus detection. This information has been used to guide polio vaccination and surveillance.

Biological sciences/Evolution/Phylogenetics

Health sciences/Medical research/Epidemiology

Since the launch of the Global Polio Eradication Initiative (GPEI) in 1988, the world has come close to eliminating wild poliovirus (WPV), with cases plummeting from hundreds of thousands annually to just a handful today. Two of the three WPV serotypes have been eradicated, and only serotype 1 (WPV1) remains in circulation, endemic to just two countries: Pakistan and Afghanistan.

Despite significant progress, these countries have never completely interrupted WPV1 transmission. The longest period without any reported cases began in early 2021, with Afghanistan going eight months and Pakistan fourteen months without new poliomyelitis cases.^1,2 During these periods, poliovirus continued to be detected in wastewater, indicating that the virus was still silently circulating in under-immunized populations, even in the absence of clinical cases. Reaching the goal of a polio-free world hinges on stopping WPV1 transmission in its final strongholds.

Between 2012 and 2023, there were two major outbreaks of WPV1 in Pakistan and Afghanistan, around 2013-14 and 2019–2020 (Fig. 1A). As of the end of September, there have been 58 reported paralytic cases in 2024 across the two countries, a significant expansion on the 12 reported cases in 2023 and the highest number in a single year since 2019. This highlights the continued threat of WPV1 transmission for the health of children in the region and globally. As recently as 2022, Mozambique and Malawi reported a 9-case WPV1 outbreak, linked to transmission in Southern Pakistan.³ The persistence of wild poliovirus in Pakistan and Afghanistan, alongside vaccine-derived cases in Africa and Asia, has resulted in polio being designated a Public Health Emergency of International Concern (PHEIC) by the World health Organization since 2014.⁴

Understanding WPV1 circulation in Pakistan and Afghanistan is complex due to several interacting factors which both help to support viral circulation and make eradication efforts difficult. Insecurity in the region, including armed conflict and targeted attacks on vaccinators, have made certain regions and communities difficult or impossible to access for health workers, allowing polio transmission to persist despite overall progress in vaccination coverage.^5,6 Additionally, the porous border between the two countries and movement of populations in the geographically isolated border regions (both nomadic and displaced populations) complicate efforts to contain the spread of WPV1, with some regions hard to reach for both vaccinators and healthcare services.⁷ Societal mistrust and misinformation around polio vaccination have also driven vaccine refusals, although overall reported coverage of vaccination campaigns remains high. Alongside issues in reaching rural communities, sanitation issues in dense urban populations as well as widespread flooding and monsoon weather exacerbate seasonal variation in polio transmission.

Understanding patterns of viral spread and persistence is further hampered by a vast majority of poliovirus infections being asymptomatic and by some areas having poor surveillance for paralytic cases. Environmental surveillance (ES) based on testing wastewater samples for the presence of poliovirus has been widely implemented in the region to help detect circulation of poliovirus in the absence of cases. ES data have been shown to provide early warning of circulation approximately 4 months prior to surveillance using AFP alone in Pakistan.⁸ In order to target vaccination responses to endemic circulation of WPV1 it is important to understand past transmission patterns in the region, identifying which areas support long-term persistence of viral circulation and act as sources for infection elsewhere. Here, we analyze routinely collected poliovirus genetic sequences from the VP1 region of the viral capsid collected over more than a decade (2012–2023) from both patients expressing symptoms of acute flaccid paralysis (AFP) (a rare and severe outcome of poliovirus infection) and ES.

As the VP1 region is sequenced routinely for genotyping and serotyping of detected viruses, a wealth of poliovirus VP1 sequences are produced globally.^8–15 So far, WPV1 sequences have primarily been used in the GPEI to routinely identify “genetic clusters” based on pairwise genetic relatedness of isolates.¹¹ Geographic analysis of the movement of these clusters has been published previously for short time periods, but an overall historic look and inference of virus movement over time and space is missing from the published literature.¹¹ Here we use a phylogeographic approach to investigate poliovirus spread in Pakistan and Afghanistan, identifying key transmission routes and describing trends in viral diversity and long term persistence over time.^16,17

Included genetic sequence data and association with reported cases

Here, we analyze 3843 VP1 sequences from wild-type polioviruses detected in samples collected between January 2012 and June 2023. During this period, there were a far greater number of genetic sequences from ES detections than AFP cases (2727 ES to 1116 AFP) (Table S1), and the number of ES sites expanded significantly through the time-period covered (Fig. S6). Throughout this analysis, we divide Afghanistan and Pakistan into 16 regions, defined based on previous epidemiological analysis of polio, other diseases and climatic regions (see map in Fig. 1C and Methods). We provide alternative estimates of key reported parameters obtained with the full dataset (including AFP and ES sequences; presented in the main figures) and using only AFP data, considered a priori less biased than ES detections, which tend to be targeted to large sewage catchments and key high-risk areas for polio transmission.^18,19 Figures for the analysis of the AFP data alone are presented in the Supplementary Information. Nevertheless, there was a strong correlation between the number of sequences and the number of reported WPV1 AFP cases across regions (R = 0.75, p = 0.007) and across years (R = 0.82, p = 0.002). This suggests that ES data provide a representative sample across the study area and time period, although this is likely biased at the regional level by differences in ES sampling intensity.

Genetic diversity through time

A phylogeny of all sequences shows two clear lineages (here labelled A and B), which drove all WPV1 transmission during the most recent outbreak, in 2019–2020 (Fig. 2). The root of these two poliovirus lineages was found to be in mid-2009 (2009-08-22; 95% HPD: 2009-02-27 to 2010-02-02) when including all sequences, with no significant difference when using AFP data alone (2009-04-16; 95% HPD: 2008-02-23 to 2010-01-21). Lineage A died out in early 2021 and has not been detected since. Three distinct sub-lineages of the B lineage persisted into 2023 (Fig. 2).

Bayesian skyline analysis of AFP data inferred increases in viral diversity during the two major outbreaks, followed by a reduction between the outbreaks and a significant decline since early 2021 (Fig. 2C). In contrast, a lineages-through-time (LTT) analysis of both AFP and ES data demonstrated the early-warning capabilities of ES data, detecting an increase in viral lineages ahead of both major outbreaks, even when AFP data showed no such signal. Towards the end of the study period, the LTT analysis indicated expanding transmission outside the Pakistan Central Corridor and Afghanistan North Corridor, detectable through ES but not AFP. Both methods revealed that viral diversity in early 2021 dropped to very low levels, with only a few lineages present, supporting the three clusters detected in 2023.

Phylogeographic assessment of viral movement

Movement of poliovirus across the 16 geographic regions over the past decade was inferred using a discrete trait analysis (Methods), and showed repeated cyclical movement of poliovirus between the southern regions of both countries. This particularly affected the South Corridor regions and Karachi (Fig. 2, purple and dark green regions). The virus was also inferred to circulate and expand in the south linked to importation from the North Corridor regions (Fig. 2, deep blue regions). Interestingly, cyclical movement of the B lineage between Karachi and the South Corridor regions during the 2019–2020 outbreak was observed with the full dataset (Fig. 2), but not picked when restricting the analysis to AFP data alone (Fig. S3). Indeed, almost all AFP detections in the South Corridor Afghanistan region (72/73) and the Karachi region (6/7) in 2019 and 2020 were linked to the now extinct A lineage.

Overall, 895 (95% HPD 867–922) transitions were inferred between regions over the full time-period when including all data (Fig. 3), compared to only 268 (95% HPD 253–283) without the inclusion of ES data (Fig. S3). Splitting numbers of inferred transitions by region suggests the city of Karachi has been the primary region seeding transmission in other regions over the decade of the study (Fig. 3B). This association is found both with and without the inclusion of ES data, although a significantly greater number of exportations are estimated when both AFP and ES data are considered (all data: 240 exportations; 95% HPD: 212–266, AFP data: 63 exportations; 95% HPD 40–82, Fig. 3B, Fig. S3C). This analysis also highlights major historical importers of poliovirus, with the North and South Corridors of Pakistan being consistently high importers (Fig. 3A). Sindh and the East Pakistan region are found to be major importers only when ES data are included (Fig. 3A and S3B), suggesting that many of these importations, which are only detected through ES, do not circulate for a long enough time to result in paralysis cases.

Local transmission lineages

We define local transmission lineages (LTLs) as clusters of phylogenetically linked detections in the same geographic region (Fig. 4A shows descriptive diagrams of the splitting process and classification into distinct categories). Using these LTLs, the estimated mean time to detection of a virus imported into a region was 93 days (CI 85–101) overall across the whole period. Looking at the regions supporting persistent transmission in 2023 (the Central Corridor Pakistan and North Corridor Afghanistan regions), this was 112 days (CI 89–134), compared to 92 days (CI 84–101) in the other regions. Without the early warning given by ES, we estimate an overall mean time to detection of an imported cluster of 166 days (CI: 144–188). This corresponds to an estimated average improvement of 73 days (equivalent to around 2 and a half months), in the time to detection of an emerging lineage with the ES system, over AFP surveillance alone.

We classified the reconstructed LTLs into four categories based on the number of detections, the duration of detection, and the number of exportation events associated with each lineage. These categories are ‘Dead end’ (circulating for less than 6 months with less than 5 associated sequences), ‘Persistent’ (circulating for greater than 1 year with greater than 10 associated sequences), ‘Export’ (Any lineage with onwards transmissions not falling into the previous categories and ‘Other’ (Any remaining lineages not classified). The median number of lineages falling into each category over time is shown in Fig. 4B, with the 95% HPD across 300 sampled phylogenies shown in the Supplementary Information (Fig. S8). Median numbers of lineages obtained with only AFP data are also provided in the Supplementary Information (Fig. S4). Over the full period, an average of 72% (median: 588/816) of all detected LTLs across reconstructed phylogenies were dead ends, not leading to onwards transmission in other regions or persistent circulation, suggesting good population immunity across Pakistan and Afghanistan. Only the Central Corridor Pakistan (median: 8; CrI: 7–9), Karachi (median: 8, CrI: 6–11) and South Corridor Pakistan (median: 8, CrI: 6–10) regions were found to have significantly greater than 5 persistent LTLs, highlighting their historic importance as long term harbors of poliovirus transmission, and key regions to target for poliovirus eradication. Long-term persistence has also been supported in the two North Corridor regions on either side of the border, although with a smaller number of distinct persistent lineages.

Grouping the LTLs into categories also shows that persistent transmission in most regions was interrupted following the 2019–2020 outbreak (Fig. 4). Each of the North Corridor of Afghanistan and the Central Corridor of Pakistan regions, the current focal regions for polio eradication efforts in the region, are found to have harbored a single persistent LTL each since early 2021. These regions supported the majority of the continued circulation of polio during the low transmission period (clusters 2 and 3 on Fig. 2A). Another persistent lineage was hosted in the city of Karachi, driven by a single detection in wastewater in mid-2023 (cluster 1 on Fig. 2A). This detection was highly divergent to other viruses in the region (5.4% divergent from any previously sequenced virus), with our analysis suggesting this isolate had a most recent common ancestor in May 2020, (2020-05-20, CrI 2020-01-25 to 2020-08-22), around 3 years before its detection. The reasons for this genetically distinct detection and any potential reservoir which supported its long-term persistence have yet to be identified.

All recent circulation (apart from the single detection in Karachi in mid-2023) has been linked to the two distinct persistent chains of transmission in the North Corridor Afghanistan and Central Corridor Pakistan regions (clusters 2 and 3 on Fig. 2A), with spill-over events into other regions seen towards the end of the included data. Full plots showing which regions poliovirus was imported from and exported to for each LTL are provided in the Supplementary Information (Fig. S2).

Recurring patterns of spread

The 2013-14 and 2019-20 outbreaks occurred primarily in the same set of focal regions, with both outbreaks emerging from exportations from the North Corridor (Fig. 2A, Fig. 3C and D). This exportation is a pattern that repeats towards the end of the included data, in early 2023. This is less clear without ES data, where the main stem (the most common ancestral line) of the phylogenetic tree is inferred to move between regions more regularly (Fig. S3A). Although large numbers of AFP cases are detected in the Central Corridor of Pakistan during both outbreaks, very little exportation is observed from this region into other areas (Fig. 3C and D). The role of the city of Karachi to amplify and repeatedly export poliovirus during outbreak periods, with cyclical movement between this region and the Southern Corridor on both sides of the border is also visible on the phylogenetic tree (Fig. 2A and Fig. 3C and D). This cyclical movement is again less clear without the inclusion of ES data, although we do infer exportations from both Karachi and the surrounding Sindh region into the South Corridor with AFP data alone.

Sensitivity to sampling

AFP data are considered here a proxy unbiased sample for comparison with the mixed AFP and ES analysis, and the results using these data are presented in the Supplementary Information. AFP surveillance is often considered a gold-standard method for poliovirus detection due to the high reporting likelihood of such a severe outcome. Despite this, the non-Polio AFP reporting rate, which is independent of poliovirus transmission in each of the regions, increased over time in all regions and varied widely between regions (Fig. S5). A method to incorporate sampling bias into estimates of movements between regions on the AFP tree was applied to polio data and presented in the Supplementary Information.²⁰ This method accounts for sampling rate differences between regions at the sampled tips, in a similar way to previous published work.²¹ The results with this correction support the findings presented using the AFP data only, with all estimated rates falling within the bounds of the corrected analysis.

Routinely collected polio genetic sequence data has the potential to provide crucial insights to the Global Polio Eradication Initiative. Genetic data, compared to case data, contains information on the relatedness of viral isolates, which can be used to reconstruct viral movement.²² In this study, we have unraveled historic patterns of polio spread and identified regions which act as importers and exporters of wild poliovirus in the last areas of endemic transmission over an 11-year period.

The approach used here was popularized by Dudas and colleagues in 2017, revealing the importance of human mobility and cross-border movement in the 2013-16 West African Ebola outbreak.¹⁷ This approach is interesting for polio as previous analyses of polio genetic data have constructed clusters of related viruses without a geographic component.^12,23,24 Including geographic data provides insights into the direction of transmission, movement of lineages across regions, and patterns of spread, critical for informing public health interventions. Our results underscore the critical role of known reservoirs in the North and Central Corridor regions in maintaining poliovirus persistence between major outbreaks. The data suggest that sustained circulation within the North Corridor specifically, poses a significant risk of transmission expansion to other regions. These corridor regions are typically characterized by under-immunized populations, hard-to-access areas and challenges to healthcare delivery driven by insecurity, geography and socioeconomic factors; all of which contribute to their importance to polio transmission.^7,25 Following major vaccination efforts which brought the last major outbreak under control these regions were key to sustaining endemic polio circulation, as is clearly illustrated by the LTLs presented in this work. ²⁶

The Karachi and South Corridor regions should also remain key targets for vaccination and surveillance, even when not supporting ongoing transmission, due to their influential role in amplifying and supporting the past two major outbreaks. Karachi, one of the world’s most densely populated urban areas with poor sanitation in many parts of the city, supports large mobile and migrant populations, under-immunized communities and high risk areas for vaccination teams, making it a key melting pot for viral spread.^27–29 The importance of these regions is widely supported by previous spatial and genetic analyses of polio transmission in the region.^9,25

The analysis highlights the importance of the East Pakistan region, which includes five major population centers, as both an importer and exporter of poliovirus, often acting as a link between circulation in the north and south of Pakistan. This role is not seen when looking at AFP data alone, suggesting that the burden is not high enough to result in significant numbers of paralytic polio cases in the region. This also underscores the utility of ES as a surveillance mechanism to reveal viral circulation even in the absence of clinical paralysis cases, which are a very rare outcome for polio-infected individuals.³⁰ Most detections in the East Pakistan region do not lead to sustained transmission, with viral lineages in the region circulating for an average of 143 days (mean; CI: 114–172). It is possible that a large amount of shedding in wastewater is due to individuals in transit through these major cities, which would be of interest as an early marker of viral dispersion from endemic areas.

Comparing the full dataset to data only from AFP cases allows us to show the early warning predictive ability of ES for poliovirus transmission, as previously investigated by Cowger et al. (2017).⁸ Our study supports their finding of an improvement in the detection of circulation with the inclusion of wastewater surveillance sites for polio monitoring over AFP reporting alone, although to a lesser extent. Here we find a 2.5-month improvement in time to detection of imported poliovirus compared to a 4-month improvement in the previous study. This is likely driven by the differences in how transmission chains are classified between the two studies, with the previous study looking at the detection time of a genetically similar virus before and after each polio case with and without the inclusion of ES data, where here we generate an inferred importation time of virus into a region.⁹ Using the date of the ancestral node as the date of importation, the most conservative estimate for time to detection in our model, a mixed ES and AFP sample provides a more comparable 3.7 month improvement over AFP sampling alone.

Incorporating ES data in a phylogenetic analysis for polio poses an issue for interpretation and lineage classification. Wastewater sites are in fixed locations and sampled at regular intervals, in contrast to stool samples from individuals with AFP, who are sampled directly from the human population. Differences in sampling intensity and targeting of ES sites to densely populated or high-risk areas need to be considered when interpreting the results of this analysis, as well as incomplete reporting of the viruses present in wastewater, particularly from sites with large and diverse catchment populations. Sampling-independent phylogeographic methods such as the structured coalescent could be used to perform this type of analysis, although the number of regions which can be considered and computational time and power required preclude their use for fast reporting, as required in an eradication context.^31–33

Incorporating data from wastewater samples is an increasingly important question for phylogenetics and phylodynamics, with widescale expansion in wastewater sampling following the global SARS-CoV-2 pandemic. The error estimation in this analysis is an improvement over current regularly reported poliovirus genetic analyses, but further work is needed to ensure robust results are provided to policymakers. Spatially biased data in simple continuous-time Markov chain methods such as the one reported here have previously been shown to lead to estimates with overly high confidence in the reported values, as well as decreasing accuracy with the degree to which data are biased.³⁴ A proposed adaptation to account for differential sampling over time is briefly presented in the Supplementary Information, allowing a correction for sampling differences in a simple phylogeographic model. Work to improve these methods and raise awareness of the assumptions made when treating location as a discrete trait is increasingly important as phylogeography becomes increasingly popular as an epidemiological tool.^35–37

With the expansion of wild polio transmission in late 2023 and 2024, the number of cases has reached its highest level in five years. This resurgence poses significant challenges for polio eradication efforts in the last endemic countries. To achieve the goal of a polio-free world, sustained commitment, community engagement and healthcare security are essential. This analysis underscores that the regions with known accessibility issues, sanitation deficiencies and healthcare shortcomings are key to the persistence of wild poliovirus. Continued international cooperation is vital to overcoming barriers to vaccination and surveillance efforts in the region.

In this analysis we infer phylogenies of wild serotype 1 poliovirus with genetic sequence data and use a Continuous Time Markov Chain (CTMC) method to estimate the likely locations of internal nodes in these phylogenies. ³⁸ To achieve this, we first reconstruct a set of empirical trees with BEAST 1.10.5 before implementing the Markov Jump method of Minin and Suchard to estimate ancestral locations and transitions between regions.^38–40 The R pipeline created for this analysis recreates aspects of the baltic python pipeline (https://github.com/evogytis/baltic) with R code to generate migration maps, phylogenetic trees with overlaid case data, and local transmission lineages from a set of posterior phylogenies output from BEAST. The code is built on several R packages, primarily the phylogenetic package ape and the ggplot2 and ggtree plotting packages and is archived on github (https://github.com/JorgensenD/WPV1_genetics_2023).^41–43

Data

We analyzed VP1 gene sequences of serotype 1 wild poliovirus detected across Pakistan and Afghanistan through stool samples and ES surveillance. These sequences were generated by the Regional Reference Laboratory for Poliomyelitis at the National Institute of Health, Islamabad, Pakistan as part of routine poliovirus surveillance. This laboratory performs virus culture and sequencing for both Pakistan and Afghanistan. The culture and sequencing of these viruses followed the WHO gold-standard protocol. Viruses which have a cytopathic effect are characterized by an intratypic differentiation reverse transcription quantitative PCR (ITD RT-qPCR) and poliovirus detections are referred for VP1 sequencing.¹⁵

The sequences are 906 bases long with mostly synonymous amino acid changes which occur according to a molecular clock at a high rate compared with most viral pathogens.^44–46 The full dataset used in this analysis includes 3,843 sequences which could be matched to location and date metadata and cover the period from 1st January 2012 to 28th August 2023 (Table S1).

The two data sources, AFP and ES, have different collection methodologies, with nubers of AFP cases generally matching well with intensity of transmission. AFP data are considered a-priori less biased than ES sampling, which is restricted to specific sampling sites.^18,30 AFP data themselves can be biased due to differences in reporting rates by region, with some regions better served by healthcare services than others. We show non-polio AFP rates as a measure of surveillance intensity independent of polio transmission in Fig. S5. ES samples can often return a large number of closely related or identical sequences due to the use of multiple cell culture flasks for each sample and the presence of viral mixtures, meaning there are duplicate sequences in the numbers included in Table S1. Here we retain duplicate sequences as they do not affect inferred transitions although they could be removed to reduce computational cost. Where we estimate possible numbers of transitions in each outbreak, we exclude these zero length branches between pairs of identical sequences. The number of ES sites sampled has expanded over time, from 23 in 2012 (all in Pakistan) to 221 in 2022 (189 in Pakistan and 32 in Afghanistan) (Fig. S6). These sites are targeted to known reservoirs of WPV1 transmission, transportation hubs and areas with large mobile populations.⁴⁷ To assess the potential impact on the results of using sequences from ES, we performed phylogeographic analysis on two datasets, the full dataset including all sequences presented in the main text and an AFP-only dataset presented in the supplementary material.

Outlying sequences

Prior to reconstruction of ancestral locations with discrete trait analysis, the sequence data were aligned with MAFFT and checked for outliers.⁴⁸ A maximum likelihood phylogeny was generated with IQ-Tree (version 1.6.12) and re-rooted with the R package ape (version 5.3) to give a phylogeny most compatible with a strict molecular clock.^41,49 Root-to-tip divergence was calculated in R with EpiGenR (version 1.1.0) (Fig. S1).⁵⁰

Discrete locations

In order to give a clear picture of geographic movement of poliovirus within the endemic countries, we divided them into 16 discrete geographic regions (Fig. 1B). These regions were informed by previous epidemiological analyses which highlighted the “corridor” areas between the two countries and the port city of Karachi in the south of Pakistan as key regions for poliovirus transmission and persistence.^{10,12,13,51,52} This corridor designation is routinely used by the GPEI and partners. The corridors all straddle the border between the two countries, with our analysis further separating these regions by country; they are: the South Corridor, centered around the major border crossing between Quetta and Kandahar; the North Corridor, centered around a major border crossing between Peshawar and Jalalabad linking the capitals of the two countries; and the Central Corridor, incorporating a more porous border with various smaller crossing points. We further split Afghanistan according to previously used ‘agro-climatic zones’, merging northern regions which report few poliovirus cases (Fig. 1A).^53,54 Outside of the corridor regions, Pakistan is subdivided primarily by province, with North Punjab and surrounding regions in Azad Jammu and Kashmir and Islamabad grouped together and separated from the less densely populated south Punjab region.⁵⁵

Bayesian phylogeny generation

Prior to Bayesian tree generation, a single maximum likelihood tree was inferred from the available sequence data with IQ-TREE.⁴⁹ From this initial tree, branches were collapsed and polytomies resolved randomly to give eight unique starting trees for Bayesian phylogeny generation. These starting trees were dated using the R package treedater with the molecular clock rate estimated using initial bounds of 0.005 and 0.015 substitutions per site per year and a starting value of 0.01, informed by previous analyses of the VP1 region.^23,56

Eight Markov chains, initiated with each of the unique starting trees were run in parallel for 15 million iterations with a phylogeny sampled every 1000 iterations. A strict clock with an informative normal prior centered on a clock rate of 0.01 was used, again informed by previous analyses and the inferred root to tip regression. Multiple chains were used to speed computation and to allow convergence to be assessed by ensuring effective sample sizes for each parameter were over 200 and 95% credible intervals were not significantly different across chains. After the iteration limit was reached, the initial 30% of each chain was removed as burn-in and the tree files combined with the BEAST utility logcombiner.⁴⁰ For the complete dataset of ES and AFP data, two of the eight chains failed to converge and were dropped from further analysis. For the AFP dataset alone, all chains were retained. The trees produced were down-sampled to give a set of 1000 trees to be used as an empirical tree set for DTA.

Markov Jumps and Rewards (MJ)

Discrete trait transitions over the phylogeny were estimated using the Markov Jumps and Rewards (MJ) method of Minin and Suchard, included in BEAST 1.10.5.³⁹ This method was chosen over simpler discrete trait CTMC as it can be used to estimate numbers of transitions between locations directly in a BEAST model. Despite the large number of parameters to estimate, the method can be sped up significantly by supplying a pre-estimated empirical tree set to estimate over as we did in this analysis. As the numbers of transitions are sampled directly as parameters in BEAST, statistical support and bounds for numbers of transitions can be trivially calculated in post-analysis of the Markov Chain Monte Carlo (MCMC) trace files. MJ was run for 10 million iterations with one of the 1000 empirical trees randomly sampled at each iteration. A 10% burn-in was removed after the iteration limit was reached.

Spread3 v0.9.7 was used to assign Bayes Factor support to the inferred number of transitions between pairs of discrete locations.⁵⁷ Transitions with Bayes Factor support above 10 are presented in the results section.

Measures of genetic diversity through time

We reconstructed viral diversity over time using two different approaches to account for the mixed sample of environmental surveillance (ES) and acute flaccid paralysis (AFP) data. For AFP data, which are sampled exclusively from the human population, we employed a Bayesian skyline model to estimate past population dynamics, reconstructing the effective population size over time (Ne(t)). This is produced using a stepwise-constant coalescent Bayesian skyline method on the available AFP sequences with 100 million iterations and 50 segments in BEAST 1.10.5. Each segment contains the same number of coalescent events, meaning that segments represent shorter time periods when there is more branching activity occurring in the tree. The mean and 95% CrI over time are calculated from the posterior after removing a 10% burn-in in the BEAST utility Tracer. However, for the full dataset, which includes wastewater samples, the Bayesian skyline approach was not appropriate. Instead, we used a lineages through time (LTT) analysis, based on the number of branching events from the posterior of the empirical tree generation procedure outlined previously. The LTT plot was censored six months before the end of the data as we cannot differentiate between right censoring and lineage extinction.

Consensus tree generation

To allow for the visualization of geographic clusters of transmission, the set of phylogenetic trees output from BEAST were combined to produce a representative consensus tree. A maximum clade credibility tree with common ancestor heights was generated with the BEAST utility TreeAnnotator, including location assignments at each internal node from MJ. Common ancestor heights were used to ensure only positive branch lengths were reconstructed in the final maximum clade credibility phylogeny. This consensus tree was used to produce the example representative phylogeny presented in the results section and example phylogenies presented in the supplementary Information.

Local transmission lineages

When inferring local transmission lineages (LTLs) we use a sample of 300 phylogenies from the posterior of the MJ analysis. These phylogenies are split wherever a movement between regions is inferred (two connected nodes of the phylogeny are inferred to be in different regions). Here we use the midpoint of the branch between connected nodes with different inferred locations as the estimated time of this movement (Fig. 4). This technique was based on visualizations presented by Dudas et al. for the 2013–2016 West African Ebola outbreak.¹⁷

LTLs across the sample of phylogenies are summarized by classifying them into three distinct types. The first, `dead-end' lineages, are defined as circulating for less than 6 months with fewer than 5 associated sequences and without onwards transmission to other regions. The second, `persistent' lineages are detected for greater than one year with greater than 10 associated sequences. The third are `export' lineages which include any remaining lineages not included in the previous two groups which seed onwards transmission in another region. This is visualized in Fig. 4.

We calculate a mean time to detection and 95% CI for each region using the inferred LTLs. Time to detection is calculated as the time from the midpoint of the branch connecting an LTL to a node in a different region to the first detection in the LTL by any surveillance method. As there will be edge effects due to censoring of the data at the end of the period covered, with any extant lineages which are yet to be detected at the endpoint not included, we do not assign any lineages extant in the last 6-months of the data to the dead-end group.

Mbaeyi C, Baig S, Safdar MR et al (2022) Progress Toward Poliomyelitis Eradication — Pakistan, January 2021–July 2022. MMWR Morb Mortal Wkly Rep 71:1313–1318
Mohamed A, Akbar IE, Chaudhury S et al (2022) Progress Toward Poliomyelitis Eradication ― Afghanistan, January 2021–September 2022. Morb Mortal Wkly Rep 71:1541
Malawi Ministry of Health, Global Polio Eradication Initiative, Davlantes E (2022) Notes from the Field: Initial Outbreak Response Activity Following Wild Poliovirus Type 1 Detection - Malawi, February 2022. MMWR Morb Mortal Wkly Rep 71:776–777
World Health Organization. WHO statement on the second meeting of the International Health Regulations Emergency Committee concerning the international spread of wild poliovirus. World Health Organization (2014) https://www.who.int/mediacentre/news/statements/2014/polio-20140505/en/ (accessed Nov 18, 2019)
Tharwani ZH, Shaeen SK, Arshad MS et al (2022) Polio amid a humanitarian crisis in Afghanistan: challenges and recommendations. Lancet Infect Dis 22:168–169
Verma AA, Jimenez MP, Tangermann RH, Subramanian SV, Razak F (2018) Insecurity, polio vaccination rates, and polio incidence in northwest Pakistan. Proc Natl Acad Sci U S A 115:1593–1598
Molodecky NA, Usman A, Javaid A et al (2021) Quantifying movement patterns and vaccination status of high risk mobile populations in Pakistan and Afghanistan to inform poliovirus risk and vaccination strategy. Vaccine 39:2124–2132
Cowger TL, Burns CC, Sharif S et al (2017) The role of supplementary environmental surveillance to complement acute flaccid paralysis surveillance for wild poliovirus in Pakistan – 2011–2013. PLoS ONE 12. 10.1371/journal.pone.0180608
Angez M, Shaukat S, Alam MM, Sharif S, Khurshid A, Zahoor Zaidi S (2012) Genetic relationships and epidemiological links between wild type 1 poliovirus isolates in Pakistan and Afghanistan. Virol J 9:51
Akhtar R, Mahmood N, Alam MM et al (2019) Genetic Epidemiology reveals three chronic reservoir areas with recurrent population mobility challenging poliovirus eradication in Pakistan. Clinical Infectious Diseases. ; published online Oct 26. 10.1093/cid/ciz1037
Mendes A, Whiteman A, Bullard K et al (2022) Spatial analysis of genetic clusters and epidemiologic factors related to wild poliovirus type 1 persistence in Afghanistan and Pakistan. PLOS Global Public Health 2:e0000251
Shaukat S, Angez M, Alam MM et al (2014) Molecular Characterization and Phylogenetic Relationship of Wild Type 1 Poliovirus Strains Circulating across Pakistan and Afghanistan Bordering Areas during 2010–2012. PLoS ONE 9:e107697
Alam MM, Sharif S, Shaukat S et al (2016) Genomic Surveillance Elucidates Persistent Wild Poliovirus Transmission During 2013–2015 in Major Reservoir Areas of Pakistan. Clin Infect Dis 62:190–198
Alam MM, Shaukat S, Sharif S et al (2014) Detection of multiple cocirculating wild poliovirus type 1 lineages through environmental surveillance: Impact and progress during 2011–2013 in Pakistan. J Infect Dis 210:S324–S332
World Health Organization. Polio laboratory manual. Geneva (2004) https://apps.who.int/iris/bitstream/handle/10665/68762/WHO_IVB_04.10.pdf?sequence=1&isAllowed=y (accessed June 19, 2019)
Grubaugh ND, Ladner JT, Kraemer MUG et al (2017) Genomic epidemiology reveals multiple introductions of Zika virus into the United States. Nature 2017 546:7658 ; 546: 401–5
Dudas G, Carvalho LM, Bedford T et al (2017) Virus genomes reveal factors that spread and sustained the Ebola epidemic. Nature 544:309–315
Asghar H, Diop OM, Weldegebriel G et al (2014) Environmental surveillance for polioviruses in the Global Polio Eradication Initiative. J Infect Dis ; 210 Suppl: S294–303
World Health Organization (2023) Field Guidance for implementation of environmental surveillance for poliovirus
Jorgensen D, Pons-Salort M, Grassly NC A simple correction to adjust for sampling biases in phylogeographic discrete trait analysis. bioRxiv 2023;: 2023.11.21.568020.
Faria NR, Rambaut A, Suchard MA et al (2014) The early spread and epidemic ignition of HIV-1 in human populations. Science 346:56
Faria NR, Suchard MA, Rambaut A, Lemey P (2011) Toward a quantitative understanding of viral phylogeography. Curr Opin Virol 1:423–429
Jorba J (2016) Phylogenetic Analysis of Poliovirus Sequences. Humana, New York, NY, pp 227–237
Gumede N, Jorba J, Deshpande J et al (2014) Phylogeny of imported and reestablished wild polioviruses in theDemocratic Republic of the Congo from 2006 to 2011. J Infect Dis ; 210 Suppl: S361–7
Mercer LD, Safdar RM, Ahmed J et al (2017) Spatial model for risk prediction and sub-national prioritization to aid poliovirus eradication in Pakistan. BMC Med 15:1–9
World Health Organization. Polio Information System (POLIS) (2014) extranet.who.int/polis/public/CaseCount.aspx
Kumar P, Arshad F, Shaheen SK, Nadeem A, Islam Z, Essar MY (2022) Water sanitation in Karachi and its impact on health. Annals Med Surg 77:103688
Abbasi FH, Shaikh AA, Mehraj J et al (2023) Vaccine Hesitancy and Perceptions of the Community about Polio in High-Risk Areas of Karachi, Sindh, Pakistan. Vaccines (Basel) 11. 10.3390/VACCINES11010070/S1
Wahid B, Kumari B, Saifullah KM, Idrees M (2023) The History and Current Killings of Polio Vaccinators in Pakistan: A Need for Targeted Surveillance Strategy. Asia Pac J Public Health 35:183
O’Reilly KM, Verity R, Durry E et al (2018) Field sensitivity of acute flaccid paralysis and environmental surveillance for poliovirus in Pakistan: an observational study. Lancet Infect Disease 18:176
Notohara M (1990) The coalescent and the genealogical process in geographically structured population. J Math Biol 29:59–75
Vaughan TG, Kühnert D, Popinga A, Welch D, Drummond AJ (2014) Efficient Bayesian inference under the structured coalescent. Bioinformatics 30:2272–2279
Müller NF, Rasmussen DA, Stadler T (2017) The Structured Coalescent and Its Approximations. Mol Biol Evol 34:2970–2981
Layan M, Müller NF, Dellicour S et al (2023) Impact and mitigation of sampling bias to determine viral spread: Evaluating discrete phylogeography through CTMC modeling and structured coalescent model approximations. Virus Evol 9. 10.1093/VE/VEAD010
Jorgensen D, Pons-Salort M, Grassly NC A simple correction to adjust for sampling biases in phylogeographic discrete trait analysis. bioRxiv 2023;: 2023.11.21.568020.
Inward RPD, Parag KV, Faria NR (2022) Using multiple sampling strategies to estimate SARS-CoV-2 epidemiological parameters from genomic sequencing data. Nat Commun 13. 10.1038/s41467-022-32812-0
De Maio N, Wu CH, O’Reilly KM, Wilson D (2015) New Routes to Phylogeography: A Bayesian Structured Coalescent Approximation. PLoS Genet 11:e1005421
Minin VN, Suchard MA (2008) Fast, accurate and simulation-free stochastic mapping. Philosophical Trans Royal Soc B: Biol Sci 363:3985
Minin VN, Suchard MA (2008) Counting labeled transitions in continuous-time Markov models of evolution. J Math Biol 56:391–412
Suchard MA, Lemey P, Baele G, Ayres DL, Drummond AJ, Rambaut A (2018) Bayesian phylogenetic and phylodynamic data integration using BEAST 1.10. Virus Evol 4. 10.1093/ve/vey016
Paradis E, Schliep K (2019) Ape 5.0: an environment for modern phylogenetics and evolutionary analyses in R. Bioinformatics 35:526–528
Yu G, Smith DK, Zhu H, Guan Y, Lam TT-Y (2017) ggtree: an r package for visualization and annotation of phylogenetic trees with their covariates and other associated data. Methods Ecol Evol 8:28–36
Wickham H (2009) ggplot2: Elegant Graphics for Data Analysis. Springer Science & Business Media
Koch F, Koch G (1985) History. The Molecular Biology of Poliovirus. Springer-, Vienna, pp 3–14
Minor PD, Ferguson M, Evans DMA, Almond JW, Icenogle JP (1986) Antigenic Structure of Polioviruses of Serotypes 1, 2 and 3. J Gen Virol 67:1283–1291
Ward CD, Stokes MA, Flanegan JB (1988) Direct measurement of the poliovirus RNA polymerase error frequency in vitro. J Virol 62:558
O’Reilly KM, Verity R, Durry E et al (2018) Population sensitivity of acute flaccid paralysis and environmental surveillance for serotype 1 poliovirus in Pakistan: An observational study. BMC Infect Dis 18. 10.1186/s12879-018-3070-4
Katoh K, Standley DM (2013) MAFFT multiple sequence alignment software version 7: Improvements in performance and usability. Mol Biol Evol 30. 10.1093/molbev/mst010
Nguyen L-T, Schmidt HA, von Haeseler A, Minh BQ (2015) IQ-TREE: A Fast and Effective Stochastic Algorithm for Estimating Maximum-Likelihood Phylogenies. Mol Biol Evol 32:268–274
Li LM, Grassly NC, Fraser C (2017) Quantifying Transmission Heterogeneity Using Both Pathogen Phylogenies and Incidence Time Series. Mol Biol Evol 34:2982–2995
Molodecky NA, Blake IM, O’Reilly KM et al (2017) Risk factors and short-term projections for serotype-1 poliomyelitis incidence in Pakistan: A spatiotemporal analysis. PLoS Med 14:e1002323
Polio Eradication (2021) Strategy 2022–2026: Delivering on a Promise. Geneva
Rehana S, Krishna Reddy P, Sai Bhaskar Reddy N et al (2022) Observed Spatio Temporal Trends Precipitation Temp Over Afghanistan. ;: 263–278
Osmani A, Robertson ID, Habib I, Aslami AA (2019) History and epidemiology of foot-and-mouth disease in Afghanistan: A retrospective study. BMC Vet Res 15. 10.1186/S12917-019-2119-Y
Pakistan Bureau of Statistics (2017) https://www.pbs.gov.pk/content/population-census (accessed Oct 5, 2022)
Volz EM (2012) Complex population dynamics and the coalescent under neutrality. Genetics 190:187–201
Bielejec F, Baele G, Vrancken B, Suchard MA, Rambaut A, Lemey P (2016) SpreaD3: Interactive Visualization of Spatiotemporal History and Trait Evolutionary Processes. Mol Biol Evol 33:2167–2169
Global Polio Eradication Initiative (2018) Best practice for monitoring the quality of polio eradication campaign performance. Geneva, http://polioeradication.org/wp-content/uploads/2018/02/best-practices-20180207-02.pdf
World Health Organization. Global Polio Surveillance Action Plan 2022–2024. Geneva (2022) https://apps.who.int/iris/bitstream/handle/10665/354479/9789240047310-eng.pdf?sequence=1&isAllowed=y (accessed Feb 7, 2023)
Amodan BO, Kisakye A, Okumu PT et al (2022) Trends of key surveillance performance indicators of acute flaccid paralysis: a descriptive analysis, Uganda, 2015–2020. BMC Public Health 22. 10.1186/S12889-022-14077-W
Soltani J, Esmailnasab N, Roshani D, Karimi M, Amjadi MJ (2014) Acute Flaccid Paralysis and Its Differential Diagnosis in in Kurdistan Province, Western Iran; an 11-Year Surveillance. Iran J Pediatr 24:131
Global Polio Eradication Initiative (2019) Polio endgame Strategy 2019–2023. Geneva
Faria NR, Suchard MA, Rambaut A, Streicker DG, Lemey P (2013) Simultaneously reconstructing viral crossspecies transmission history and identifying the underlying constraints. Philosophical Trans Royal Soc B: Biol Sci 368. 10.1098/rstb.2012.0196

There is NO Competing Interest.

supplementWPV11024.docx

Download PDF

Version 1

posted

You are reading this latest preprint version

Evolution and Transmission Dynamics of Wild Poliovirus in Pakistan and Afghanistan (2012-2023): implications for global polio eradication

Status:

Version 1

Abstract

Figures

Introduction

Results

Included genetic sequence data and association with reported cases

Genetic diversity through time

Phylogeographic assessment of viral movement

Local transmission lineages

Recurring patterns of spread

Sensitivity to sampling

Discussion

Methods

Data

Outlying sequences

Discrete locations

Bayesian phylogeny generation

Markov Jumps and Rewards (MJ)

Measures of genetic diversity through time

Consensus tree generation

Local transmission lineages

References

Additional Declarations

Supplementary Files

Status:

Version 1