Monitoring Plasmodium falciparum and Plasmodium vivax using microsatellite markers indicates limited changes in population structure after substantial transmission decline in Papua New Guinea

Monitoring the genetic structure of pathogen populations may be an economical and sensitive approach to quantify the impact of control on transmission dynamics, highlighting the need for a better understanding of changes in population genetic parameters as transmission declines. Here we describe the first population genetic analysis of the major human malaria parasites, Plasmodium falciparum (Pf) and Plasmodium vivax (Pv) populations following nationwide distribution of long-lasting insecticide treated nets (LLIN) in Papua New Guinea (PNG). Parasite isolates from pre- (2005-6) and post-LLIN (2010-2014) were genotyped using microsatellite markers. Despite parasite prevalence declining substantially (East Sepik: Pf=54.9-8.5%, Pv=35.7-5.6%, Madang: Pf=38.0-9.0%, Pv: 31.8-19.7%), genetically diverse and intermixing parasite populations remained. Pf diversity declined modestly post-LLIN relative to pre-LLIN (East Sepik: Rs = 7.1-6.4, He = 0.77-0.71; Madang: Rs= 8.2-6.1, He = 0.79-0.71). Unexpectedly, population structure present in pre-LLIN populations was lost post-LLIN, suggesting that more frequent human movement between provinces may have contributed to higher gene flow. Pv prevalence initially declined but increased again in one province, yet diversity remained high throughout the study period (East Sepik: Rs=11.4-9.3, He=0.83-0.80; Madang: Rs=12.2-14.5, He=0.85-0.88). Although genetic differentiation values increased between provinces over time, no significant population structure was observed at any time point. For both species, a decline in multiple infections and increasing clonal transmission and significant multilocus linkage disequilibrium (mLD) post-LLIN was a positive indicator of impact on the parasite population using microsatellite markers. These parameters may be useful adjuncts to traditional epidemiological tools in the early stages of transmission reduction.
Methods
Whole blood samples were collected from participants in cross-sectional studies conducted between 2005 and 2014 along the North Coast of PNG (Mueller et al., 2009; Schultz et al., 2010). In Madang Province (MAD), the same three catchment areas were studied in 2006 (Schultz et al., 2010; Senn et al., 2012), 2010 (Koepfli et al., 2015) and 2014 (Koepfli et al., 2017).

For both species, a panel of 9-10 neutral microsatellite markers were amplified in the selected samples (Table S1) using a multiplex primary PCR followed by individual nested PCRs as previously described (Anderson et al., 2000; Jennison et al., 2015; Koepfli et al., 2013; Schultz et al., 2010). For P. falciparum, samples were genotyped at nine previously validated and commonly used, putatively neutral, microsatellite loci including TA1, TAA60, Polya, ARA2, Pfg377, TAA87, PfPK2, TAA81 and 2490 (Anderson et al., 2000; Schultz et al., 2010). For P. vivax, 10 putatively neutral microsatellites were genotyped as previously described: MS1, MS2, MS5, MS6, MS7, MS9, MS10, MS12, MS15, and MS20 (Jennison et al., 2015; Koepfli et al., 2013). All PCR products were sent to a commercial facility for fragment analysis on an ABI3730xl platform (Applied Biosystems) using the size standard LIZ500. Primers used were the same for all datasets (Jennison et al., 2015; Schultz et al., 2010). The electropherograms were analysed with Genemapper V4.0 (Applied Biosystems) with the same peak calling strategy as described previously (Jennison et al., 2015; Schultz et al., 2010). To avoid artefacts, precautions were taken to ensure allele calling was consistent (Jennison et al., 2015), and carefully reconstructing dominant haplotypes as per previously described methods (Anderson, Su, Bockarie, Lagog, & Day, 1999; Jennison et al., 2015; Schultz et al., 2010). For the Madang 2005 and Wosera 2006 P. falciparum data, previously published cleaned and rounded microsatellite allele repeat numbers for P. falciparum single clone infections (Schultz et al., 2010) were converted back to allele sizes using the known number of nucleotides/repeat, whereas for P. vivax the raw data (allele calls) was available (Jennison et al., 2015). These data were combined with the newly generated MS data from the other studies before binning the alleles using the TANDEM software (Matschiner & Salzburger, 2009). Allele frequencies of the entire dataset (incl. previously genotyped datasets) were investigated and outlying alleles (most likely caused by PCR artefacts) were removed. Samples with missing data at six (60%) or more MS loci were excluded from further analysis.