DISCUSSION
Parasite genetic diversity has been proposed as a key indicator of malaria transmission dynamics to track control and elimination progress (Stuckey et al., 2018). This demands an evaluation of whether genotyping provides insights about changes in transmission following intensified malaria control efforts. Here we have monitored the population genetics of P. falciparum and P. vivax over an eight-year period of LLIN-induced transmission decline in PNG, an area of high year round transmission. Despite large reductions in parasite prevalence and multiple infections from very high to low/moderate levels (Kattenberg et al., 2020; Koepfli et al., 2017; Koepfli et al., 2015), we show that population genetic changes were minimal with populations remaining diverse and unstructured. P. falciparum diversity decreased somewhat, though remained high relative to other malaria endemic areas outside Africa (Anderson et al., 2000; Branch et al., 2011; Chenet et al., 2012; Chenet et al., 2015; dalla Martha, Tada, Ferreira, da Silva, & Wunderlich, 2007; Noviyanti et al., 2015; Orjuela-Sanchez et al., 2009; Orjuela-Sanchez et al., 2013; Pava et al., 2017; Susomboon et al., 2008), whereas P. vivax diversity was unchanged throughout the study period. Surprisingly, P. falciparum populations that were structured pre-LLIN (2005, 2006) (Jennison et al., 2015; Schultz et al., 2010), were unstructured post-LLIN (2010, 2014), although a reduction in multiple infections and an increase in multilocus LD due to clonal haplotypes were detected. For P. vivax , there was also no evidence of population structure after LLIN, however increasing pairwise genetic differentiation within and between provinces was observed and clonal transmission and inbreeding had emerged in at least one village. These results demonstrate that declining transmission does not result in an immediate decrease in overall population diversity nor an increase in population structure. Sustained low transmission may be needed to observe these changes using these small panels of microsatellite markers. However, changes in multiple infection prevalence and multilocus LD indicate increasing heterogeneity in transmission within populations. These results have implications with respect to monitoring changing transmission for any pathogen using population genetic approaches.
Nkhoma and colleagues previously reported a limited impact on P. falciparum diversity following a decrease from moderate to low transmission on the Thai-Myanmar border over 10 years, however as we have observed in PNG, there was an increase in multilocus LD which and decreasing proportions of multiclonal infections (Nkhoma et al., 2013). Daniels et al . have also reported decreasing multiclonal infections, increasing proportions of repeated genotypes and multilocus LD, and long-term persistence of particular clones in Senegal (Daniels et al., 2013). These studies utilised 96 and 24 biallelic SNP markers respectively, compared to our study using a small number of microsatellite markers - similar panels of SNPs may reveal additional insights in the PNG context. As for P. vivax , Gunawardena and colleagues also reported sustained high P. vivax microsatellite diversity and multiclonal infections during a five-year period as the country progressed towards malaria elimination (Gunawardena et al., 2014). Population genetic signals of declining transmission might take longer to emerge for P. vivax due to new blood-stage infections from the hypnozoite reservoir and could explain why we only observed mLD in one village of East Sepik.
In PNG, the limited decline in P. falciparum diversity and loss of parasite population structure after LLIN distribution indicates that there may be increased gene flow between the sampled parasite populations, which was unexpected. Population structure prior to intensified control was thought to be the result of limited historical human migration due to the rugged terrain and lack of direct transport connections (Mueller, Bockarie, Alpers, & Smith, 2003; Riley, 1983) or population bottlenecks due to past control efforts or emergence of drug resistance (Anderson et al., 2000; Jennison et al., 2015; Schultz et al., 2010; Tessema et al., 2015). Changes in first-line treatment policies, for example the introduction of sulphadoxine/pyrimethamine (SP) in the early 2000’s and artemether-lumefantrine in 2008, might have played a role in shaping parasite population structure (Mu et al., 2005). Chloroquine (CQ) and/or SP resistance (near fixation of resistantpfmdr1 and pfdhps resistant alleles were observed in the same areas (Barnadas et al., 2015; Koleala et al., 2015; Mita et al., 2018)) and may have contributed to the observed bottleneck and mLD in pre-LLIN P. falciparum populations, with consequent reductions in effective population size, while outcrossing due to high transmission preserved within-population genetic diversity as the resistance mutation spread (Mita et al., 2018). From 2000 to 2011 the PNG population increased by over two million people (National Census Report, 2011), and local observations suggest that large numbers of migrants from East Sepik have moved into Madang in the last decade seeking employment opportunities. As a result, post-LLIN, greater connectivity between human populations may have broken down P. falciparum population structure and maintained high gene flow between P. vivaxpopulations (Fola et al., 2018).
For the post-LLIN East Sepik P. vivax population where prevalence dropped to below 5%, significant mLD was observed resulting from very closely related parasite strains circulating in a residual pocket of relatively high transmission within a single village. This suggests considerable inbreeding of parasites in that village, in a background of high genetic diversity resulting in a lack of evidence of a bottleneck at the population level. This paradoxical signature of significant mLD with high diversity and a considerable proportion of multiple clone infections of P. vivax appears to be a hallmark of lower transmission areas (Barry et al., 2015; Batista et al., 2015; Chenet et al., 2012; Delgado-Ratto et al., 2016; Delgado-Ratto et al., 2014; Ferreira et al., 2007; Hong et al., 2016; Noviyanti et al., 2015; Orjuela-Sanchez et al., 2013). Similar to P. falciparumpopulations though, there was a correlation between mLD and prevalence of infection for P. vivax . This shows that reductions in multiclonal infections and mLD is an earlier indicator of reduced transmission than genetic diversity and population structure (for both species).
Multilocus LD in post-LLIN P. vivax populations was explained by both focal clonal transmission and inbreeding, as similarly observed in other studies in Peru, Vietnam, and Papua Indonesia (Delgado-Ratto et al., 2014; Hong et al., 2016; Noviyanti et al., 2015). The explanation for this observation will most likely lie in unique P. vivaxcharacteristics related to hypnozoites, relapses and transmission dynamics (Abdullah et al., 2013; Bright et al., 2014; Chen, Auliff, Rieckmann, Gatton, & Cheng, 2007; Delgado-Ratto et al., 2014; Ferreira et al., 2007; Fola et al., 2018; Iwagami et al., 2012; White, 2011). At high transmission (e.g. pre-LLIN) with high prevalence and high multiplicity of infection, these clusters of similar haplotypes might also occur, but could be obscured due to sampling limitations and the large number of different haplogroups circulating at the same time. As transmission declines, infections have fewer clones and the diversity of the hypnozoite reservoir decreases, resulting in fewer circulating haplogroups, lower levels of recombination between distinct genomes and more frequent clonal transmission and inbreeding, resulting in increased mLD as in this and other studies (Barry et al., 2015; Batista et al., 2015; Chenet et al., 2012; Delgado-Ratto et al., 2016; Delgado-Ratto et al., 2014; Ferreira et al., 2007; Hong et al., 2016; Noviyanti et al., 2015; Orjuela-Sanchez et al., 2013). These signatures are more likely to be seen in a population with more sustained low transmission such as was the case for the East Sepik Province. In this region, malaria transmission is heterogeneous between villages. Besides the national malaria control program, other initiatives were also distributing LLINs in East Sepik Province prior to the nationwide distribution (Hetzel et al., 2014; Hetzel et al., 2012; Hetzel et al., 2016) suggesting that longer term sustained control efforts have been effective.
Considerable variance in the impact of the LLIN program was observed in the two provinces. In Madang, whilst P. falciparum rates steadily declined over the entire study period, there was a resurgence of submicroscopic P. vivax infections in 2014 (Koepfli et al., 2017). Although multiplicity of infection remained low, the lack of mLD and population structure suggests that this is not due to an outbreak, but more likely the residual P. vivax population was able to gain a foothold once again despite the intensive control efforts. In addition, an increase in the prevalence of pvdhfr triple and quadruple mutants, related with SP resistance, were observed in Madang province in 2010 (Barnadas et al., 2015), and a high proportion of resistant parasites could be a possible explanation for the higher infection prevalence by 2014. Different studies have shown that selective pressure of drugs such as CQ and/or SP have had an impact on shaping worldwide P. vivax populations in recent history (Hupalo et al., 2016; Pearson et al., 2016). However, a population bottleneck as seen in P. falciparum populations (Mita et al., 2018) did not occur in P. vivax populations of PNG. Malaria control had a greater impact on P. vivax prevalence in East Sepik and population structure was observed in one village post-LLIN.
This study has some limitations related to sampling and the genotyping approach used. The samples were collected in serial cross-sectional surveys over a period of malaria control initiated at different times in the two provinces. Fewer years of sustained control pressure compared to East Sepik might explain why, despite substantial prevalence decline in Madang Province by 2010, we did not observe any signature of changing population structure. Whilst, in East Sepik 2012, a cluster of closely related parasites was observed in one village suggesting more focal transmission than previous years. The microsatellite panels were selected as these have been the gold standard genotyping tool for large-scale malaria population genetic studies for many years (Anderson et al., 1999; Imwong et al., 2006; Karunaweera et al., 2006). However, they are limited in number (less than one per chromosome), rapidly evolving and prone to error. Although these markers are extremely useful for measuring parasite population structure on local and global scales (Auburn & Barry, 2017; Barry et al., 2015; Koepfli & Mueller, 2017; Sutton, 2013), they are not ideal for cross-study comparisons due to the difficultly in calibrating alleles from different data sources. The lack of raw data from the previously published P. falciparum study (Schultz et al., 2010), prevented the direct comparison of haplotypes and thus the analysis of population structure between the earlierP. falciparum time points for the East Sepik II (Wosera) and Madang populations. Furthermore, dominant haplotypes derived from multiple clone infections can be reconstructed incorrectly, thus inflating diversity values (Messerli, Hofmann, Beck, & Felger, 2017). However, only haplotype-based analyses such as multilocus LD and phylogenetic analysis were vulnerable to these possible artefacts, and would result in an overestimate of diversity and underestimate of LD. All other analyses were conducted using mean values across markers or allele frequencies and thus limit the impact of such errors. Moreover, the dataset included a large proportion of single infection haplotypes in all populations, and the detected clones included dominant haplotypes suggesting that allele calling was reliable. Finally, the highly polymorphic and rapidly evolving nature of microsatellite markers (Ellegren, 2004) may limit the resolution of the population genetic parameters such as population level diversity and population structure in high transmission areas (Branch et al., 2011). This may both lead to false assignment of unrelated parasites (e.g. from different provinces) as related and reduce the detection of truly related parasites (identical by descent), both of which would result in a lack of population structure. Other more stable markers, such as SNPs (Baniecki et al., 2015; R. Daniels et al., 2008) or whole genome sequencing (Hupalo et al., 2016; Miotto et al., 2013; Mu et al., 2005; Pearson et al., 2016; Volkman et al., 2007) may be more sensitive and specific to detect changes in parasite population structure.
According to the data presented, in two high transmission provinces of PNG, recent reductions in malaria transmission result in limited changes in parasite population diversity and structure. The high parasite diversity and lack of population structure are consistent with both species maintaining a large and evolutionarily fit population immediately after prevalence decline suggesting the gains made are fragile. However, fewer multiclonal infections, and the emergence of significant mLD for both species indicates there is focally interrupted transmission and suggests that these parameters may be useful markers for measuring control impact and early changes in parasite population structure with intervention. The results were in contrast to our expectations of decreasing diversity and increasing population structure (Jennison et al., 2015) and show that long term sustained control efforts may need to be maintained to observe significant change in population structure, at least as measured by the microsatellite panels used in this study. The PNG national malaria control program has made an impact on the malaria burden through the substantial reductions in infections circulating in the community (Hetzel et al., 2016; Kattenberg et al., 2020; Koepfli et al., 2017), however it appears that this has not been sustained long enough for the underlying parasite population to fragment. Finally, this study demonstrates how parasite population genetics can inform whether malaria intervention has reduced and fragmented the parasite population or if significantly more control effort will be required to do so.