Sequence polymorphisms ofmsp1block 2 inP. falciparumisolates from STP
Twenty-two haplotypes of msp1 block 2 were identified in 118 samples collected between 2010 and 2016. The sequence accession nos. for the msp1 haplotypes were MW001371- MW001392. After translating nucleotides into protein sequences, nine haplotypes belonged to the K1 family, nine belonged to the MAD20 family, and four belonged to the RO33 family. The K1 and MAD20 allelic types contained tri-peptide repeat units, showing length polymorphisms from 312 to 366 bp (104–122 amino acids) in the K1 types, and 285–411 bp (95–137 amino acids) in the MAD20 types. In the K1 haplotypes, flanked by family-specific regions (E1 and E2 in Figure S2), the tri-peptide repeat region (R1) started from SAQ and terminated with SGT, and could be divided into two subgroups on the phylogenetic tree according to the central repeat region (Fig. 3). One group included KH2, 4, 5, 6, and 9, which contained SAQ and SGA-rich repeats, and was more similar to the P. falciparum NF7 clone isolated from West Africa. Another group included KH1, 3, 7, and 8 consisting of SGT and SGP-rich repeats, and was more similar to the 3D7 reference clone. Among the total 118 MSP1 sequences, K1 haplotypes were identified in 41 isolates (34.7%, 41/118). The prevalent haplotype in the K1 allelic family was KH1 (36.6%, 15/41), followed by KH3 (26.8%, 11/41).
In the MAD20 allelic types, the tri-peptide repeat region (R1) started with SKG, SSG, or SGG motifs, followed by the SVT motif (Figure S2). The central repeat region majorly contained SGG and SVA-rich repeats and terminated with SGG. The predominant type in the MAD20 family was MH3 (63.5%, 33/52), which was translated from two nucleotide sequences with a synonymous SNP, represented as MH3-1 and MH3-2 on the phylogenetic tree (Fig. 3). Except for MH3, each MSP1 haplotype was translated from a single nucleotide sequence. Two subgroups were found in the MAD20 family in the phylogenetic tree. One included MH1, 5, 6, 7, and 8, and was grouped together with the HB3 reference sequence, and sequences identified in Tanzania and Oceania islands. Another subgroup included MH2, 3, and 4, which was similar to sequences identified in Southeast Asia and Ghana.
The RO33 allelic type had the least isolates (21.2%, 25/118) and haplotype diversity in STP across the sampling period. The RO33 allele did not encode any tandem repeats. The RO33 haplotypes had identical lengths (285 bp) with four nonsynonymous SNPs (Figure S2). The predominant haplotype in the RO33 family was RH1 (84%, 21/25), isolated from samples collected after 2012 across all the districts in São Tomé. Phylogenetic analysis showed that RH1 and RH2 were grouped together with the RO33 clone isolated from Ghana, while RH3 and RH4 were grouped together with the sequence isolated from Malawi.
Sequence polymorphisms ofmsp2block 3 inP. falciparumisolates from STP
Twenty-one haplotypes were identified from 118 samples, of which 11 belonged to the 3D7/IC family (IH1-IH11, accession nos. MW001393-MW001403), and 10 belonged to the FC27 family (FH1-FH10, accession nos. MW001404-MW001413). Two tandem repeat regions of MSP2 block 3 (R1 and R2) were flanked by three family-specific regions (E1-E3, Figure S3). The R1 region in 3D7/IC haplotypes consisted of GSA-rich repeats differing in the lengths and numbers of repeat units, resulting in diverse length polymorphisms (from 430 to 600 bp). The R2 region of 3D7/IC was composed of 5–14 poly-threonine (T) residues, and only IH2 and IH9 presented additional copies of the TPA motif before the poly-threonine segments. Several nonsynonymous amino acid substitutions and deletions were found in the family-specific regions (E1-E3) of the 3D7/IC types. At the bottom of the phylogenetic tree of 3D7/IC (Fig. 4A), IH5 and IH11 were grouped together with reference sequences from Asian countries. The middle part of the tree showed that IH1, 2, 8, 9, and 10 were similar to sequences identified in Gambia and Nigeria. Finally, IH3, 4, 6, and 7 were grouped together with the 3D7 clone and sequences identified from Asian and Oceanian countries.
Among the MSP2 sequences, the majority (72.9%, 86/118) belonged to the FC27-allelic type. Ten haplotypes were identified in the FC27 family, and FH6 and FH1 were the most prevalent types, accounting for 30% each. The FC27 haplotypes were divided into two subgroups. One subgroup included FH1, 2, and 3, containing three copies of the family-specific repeats (32 amino acids) in the R1 region, without showing tandem repeats in the R2 region. In contrast, another subgroup included FH4 to FH10, showing tandem repeats of 12 amino acids (aa) in the R2 region, without presenting any repeats in the R1 region (Figure S3). At the top of the phylogenetic tree (Fig. 4B), the prevalent haplotype FH6 was most similar to the reference sequences from the neighboring country, Gabon. On the other hand, haplotypes FH1 to FH3 located at the bottom of the tree were grouped together with the sequences found in Asian countries, Myanmar and India.
Temporal changes ofmsp1andmsp2alleles before and after the recent peak incidence
Both MSP1 and MSP2 sequences were acquired from 118 pre-treatment samples collected from 2010 to 2016, which was a period of low transmission leading to the pre-elimination phase. During this period, a slight increase in annual malaria incidence was found in 2012 (incidence rate = 6.9%, Fig. 1). The incidence declined after 2012 and decreased to ~ 1% after 2014, reaching the pre-elimination criteria. This study tracked the temporal changes in MSP1 and MSP2 in infective parasites before and after the recent incidence peak. The results showed significant changes in both MSP1 and MSP2 allelic types before and after the peak incidence at 2012 (MSP1: χ² = 7.94, df = 2, p = 0.02; MSP2: χ² = 23.46, df = 1, p < 0.001), with the predominant alleles in MSP1 changing from K1 to MAD20, and those in MSP2 changing from 3D7/IC to the FC27 type (Fig. 5).
The major MSP1 allelic types before 2012 were found to be the haplotypes KH1 and KH8, which were not detected in the samples isolated after 2012. An increase in MAD20 and RO33 allelic types was found in parasites collected after 2012, particularly in the MH3 and RH1 haplotypes (Fig. 5). Haplotype MH3 accounted for 12% of the parasites isolated before 2012, but became the predominant type after 2012, accounting for 27%-35% of the prevalence. Haplotype RH1 was not identified in the isolates collected before 2012, but was found to constitute 20%-27% of the parasite population isolated after 2012.
Before 2012, 67% (18/27) of the parasite isolates were identified as the MSP2-3D7/IC allelic types, of which IH1 and IH4 were the prevalent types, and the minor types were IH3 and IH5-IH10 (Fig. 5). However, after 2012, only IH1, IH2, IH3, and a new IH11 were found in the MSP2-3D7/IC allelic family, and the prevalence of 3D7/IC alleles dropped by 50%. The predominant MSP2 types were replaced by FC27 alleles, since their prevalence increased by 50% after 2012. Haplotypes FH1 and FH6 were present throughout the sampling period in our sample set, and increased their frequencies after 2012. Haplotype FH5 was isolated only after 2014, and became the predominant type in 2016.
In our sample set, 6.8% (8/118) of the samples showed polyclonal infections by nested PCR and gel purification (Table S2). Two patients were infected with more than one haplotype in both MSP1 and MSP2, while the other six patients were infected with 2–3 haplotypes in MSP1, but a single haplotype in MSP2. The majority of the samples were identified as monoclonal infections, being infected by one haplotype in either MSP1 or MSP2. Overall, the prevalent circulating parasite strains before 2012 were MSP1-K1 and MSP2-3D7/IC, and those after 2012 were MSP1-MAD20 and MSP2-FC27. The changes occurred in the year of peak incidence, showing that the prevalent strains had changed from the low-transmission to the pre-elimination period in STP.
Prevalence of anti-malarial drug resistance mutations inP. falciparumisolates
Antimalarial drug-resistance markers were genotyped in 60 pre-treatment samples and 51 post-treatment samples collected from 2014 to 2016. The dominant alleles of pfmdr1 and pfcrt were pfmdr1 86Y (82.9%, 92/111), pfmdr1 184F (62.2%, 69/111), pfmdr1 D1246 (96.4%, 107/111), and pfcrt 76T (92.8%, 103/111) (Table 1). Thus, the dominant allelic types of pfmdr1 86 + 184 + 1246-pfcrt 76 were YFD-T (51.4%, 57/111), followed by YYD-T (27.0%, 30/111). A total of 93 pfK13 sequences (artemisinin resistance marker, 986 bp) were also amplified, but no mutations were found in this genetic marker.
The post-treatment samples were selected from 41 patients showing inadequate treatment outcomes, of which nine (22.0%, 9/41) of them had insufficient clearance on day 3, nine (22.0%, 9/41) showed early treatment failure on day 3, 13 showed late parasitological failure after day 7 (31.7%, 13/41), and ten (24.4%, 10/41) showed gametocytes during follow-up. Thirty-one (75.6%, 31/41) of these patients were treated with quinine and 10 (24.4%, 10/41) were treated with ASAQ on day 0. Five patients were originally treated with quinine and showed positive parasitemia for more than one follow-up day, and the alternative drugs were administered afterwards: ASAQ was administered to three patients, AL was administered to one patient, and both ASAQ and AL were administered to one patient after day 7. Fisher’s exact test showed no significant differences in the pfmdr1 and pfcrt polymorphisms after treatment with either quinine or ASAQ at baseline.
Matched-pair comparisons of pfmdr1 and pfcrt polymorphisms were performed between pre-treatment (day 0) and post-treatment samples (day R) by McNemar’s test, and no significant changes were found. More than 80% (41/51) of the parasite isolates on day R shared identical pfmdr1 and pfcrt polymorphisms with their initial infections on day 0. This was also confirmed by the genotyping results of msp1 and msp2 to discriminate recrudescence from new infections. The results showed that 83% (34/41) of the recurrent patients shared the same msp1, msp2, pfmdr1, and pfcrt alleles between pre- and post-treatment samples, indicating that the positive parasitemia at follow-up was mainly due to the incomplete clearance of the initial infections (recrudescence). However, seven patients (17%, 7/41) showed different pfmdr1 or pfcrt patterns between day R and day 0 (Table S3). Among these seven patients, five showed substitutions one week after treatment, while two showed early substitutions at day 3. The substitutions mostly changed from pfmdr1 86Y to N86, 184F to Y184, and pfcrt 76T to K76 after treatment with quinine (6 patients) and ASAQ (1 patient) at baseline. One patient who was treated with AL by day 7 showed substitutions from pfmdr1 86Y to N86 and pfcrt 76T to K76 on days 7, 9, and 14. Moreover, through confirmation by nested PCR in these seven patients, four were confirmed to have polyclonal infections at day 0, and new infections were also detected in four recurrent cases. Thus, substitutions of pfmdr1 or pfcrt types may be due to the polyclonal infections at the beginning or new infections during follow-ups.
The pfmdr1 N86Y and Y184F were correlated with the MSP1 allelic families using Fisher’s exact test (p < 0.001, Figure S4). The pfmdr1 N86 was only found in the MSP1-RO33 type in our samples, whereas 86Y was mainly present in MSP1-MAD20 (48.8%, 41/84) and MSP1-K1 (42.9%, 36/84) allelic types. The combination of pfmdr1 N86 + 184F + D1246 type (resistance to LF) was only found in the MSP1-RO33 and MSP2-FC27 allelic types. In contrast, pfmdr1 Y184 was mainly found in the MSP1-K1 family (73.3%, 22/30), while 184F was predominant in the MSP1-MAD20 type (56.1%, 37/66). The combination of pfmdr1 86Y + Y184 + pfcrt 76T, which was known to be associated with resistance to AQ, was present mainly in the MSP1-K1 type (81.5%, 22/27).
Finally, temporal changes were only found in pfmdr1 Y184F among the 60 pre-treatment samples from 2014 to 2016 (χ² = 8.25, df = 2, p = 0.02). The pfmdr1 184F type showed a significant increase of 30%-40% in 2016 (Figure S5). In conclusion, from 2014 to 2016, pfmdr1 86Y (75%-95%) and pfcrt 76T (95%-100%), which were associated with the reduced sensitivity of AQ, the first-line ACT partner drug used in STP, were constantly present in the local parasite population. The pfcrt 76T mutation was nearly fixed in the local parasite population, and the pfcrt 72–76 mutant type was identified as the CVIET haplotype from sequencing results. Among these markers, pfmdr1 Y184F showed the largest variation in our sample set from 2014 to 2016, while the other markers showed the same prevalent types across the period.
Characteristics of infective parasite density, treatment outcomes, and factors related with recurrence
To examine the host and temporal factors that may affect parasite density in initial infections, a linear regression model was used. The outcome variable was log10- transformed parasite density (parasites/µL) at the date of initial infection. The predictive host-related variables were age, gender, residential area, and time-related variables, including period and endemic season. Results showed that “age” is a significant factor (coefficient = -0.012, p < 0.001), which may decrease the infective parasite density with an increase in age (Table S4). The levels of parasitemia in infected hosts were significantly higher (coefficient = 0.40, p < 0.001) in the low-transmission period (2010–2013) than in the pre-elimination period (2014–2016).
A total of 7,482 malaria cases with complete records were enrolled in the analysis, and 86.7% (6,490/7,482) of the cases showed adequate clinical and parasitological responses (ACPR) after treatment. The ACPR rate in the ACT treatment group ranged from 88.2%-97.6%, which was slightly higher than that in the quinine treatment group (73.6%-94.0%). The proportion of cases with inadequate treatment outcomes (positive parasitemia at days 3–28 after treatment) after quinine treatment (17.3%, 847/4,906) was higher than that after ACT (5.6%, 145/2,576). In particular, patients treated with quinine were mostly hospitalized patients with higher parasitemia levels at enrollment, and their follow-up assessments showed higher probabilities of insufficient clearance on day 3 (0.7%-14.2% in patients receiving quinine and 0.3%-4.2% in those receiving ACT) and late parasitological failure after day 7 (3.6%-7.8% in patients receiving quinine and 0.8%-4.7% in those receiving ACT) (Fig. 6). In the low-transmission period, more than half of the patients (61.5%-81.5%) were treated with quinine, while a lower proportion of patients were treated with ACT (18.5%-38.5%). After the decrease in parasitemia levels in the pre-elimination period, cases treated with ACT were greater or nearly equal to the quinine-treated numbers, and the overall ACPR rate was higher in the pre-elimination phase (2014–2016) than in the low-transmission period (2010–2013).
Logistic regression models (Table 2) and mediation analysis (Figure S6) were performed to characterize the factors that were correlated with early and late parasitological failure. Results from logistic regression models showed that age was a significant protective factor against recurrence (OR = 0.97 for early recurrence, and OR = 0.98 for late recurrence, both p < 0.001). Thus, younger children had higher odds of recurrence than the elderly. In the low-transmission period, the odds of early and late recurrence were significantly higher than those in the pre-elimination period (OR = 1.63–5.82, p < 0.01). Both parasite density (log10 transformed, OR = 1.44, p < 0.001) and treatment selection (quinine vs. ACT, OR = 1.96, p < 0.001) were significant risk factors for early recurrence; however, only treatment selection (quinine vs. ACT, OR = 1.91, p < 0.001) was a significant risk factor for late recurrence. Since treatment selection would be partially determined by the infective parasite density at enrollment, the mediation pathway of “initial parasite density→ treatment selection→ recurrence” was further tested in the model. The direct effect of parasite density on early and late recurrence was 0.46 (p < 0.001), and 0.16 (p = 0.009), respectively. After partial mediation by treatment selection, the coefficient of parasite density declined, and its indirect effect on both early and late recurrence was 0.09. Moreover, log10 parasite density was no longer a significant factor for late parasitological failure after the mediation effect of treatment selection. In conclusion, quinine treatment was a significant risk factor correlated with early and late recurrence, and a partial mediator between infective parasite density and recurrence during follow-up.