Comparative Analysis of Plasmodium falciparum Genotyping via SNP Detection, Microsatellite Profiling, and Whole-Genome Sequencing.

Research efforts to combat antimalarial drug resistance rely on quick, robust, and sensitive methods to genetically characterize Plasmodium falciparum parasites. We developed a single-nucleotide polymorphism (SNP)-based genotyping method that can assess 33 drug resistance-conferring SNPs in dhfr, dhps, pfmdr1, pfcrt, and k13 in nine PCRs, performed directly from P. falciparum cultures or infected blood. We also optimized multiplexed fragment analysis and gel electrophoresis-based microsatellite typing methods using a set of five markers that can distinguish 12 laboratory strains of diverse geographical and temporal origin. We demonstrate how these methods can be applied to screen for the multidrug-resistant KEL1/PLA1/PfPailin (KelPP) lineage that has been sweeping across the Greater Mekong Subregion, verify parasite in vitro SNP-editing, identify novel recombinant genetic cross progeny, or cluster strains to infer their geographical origins. Results were compared with Illumina-based whole-genome sequence analysis that provides the most detailed sequence information but is cost-prohibitive. These adaptable, simple, and inexpensive methods can be easily implemented into routine genotyping of P. falciparum parasites in both laboratory and field settings.

INTRODUCTION

In 2019, there were an estimated 229 million malaria cases and 409,000 deaths, spanning 87 malaria-endemic countries, with most deaths resulting from Plasmodium falciparum infections in African children (1). The emergence and spread of multidrug-resistant P. falciparum strains complicate efforts to reduce malaria morbidity and mortality. Crucial components of malaria research include tracking the emergence of drug-resistant parasites, identifying and investigating genetic markers that correlate with treatment failure, and discovering effective antimalarials with novel modes of action.

These research initiatives have incentivized a long-standing interest in developing robust, sensitive, and cost-effective genotyping methods. Accuracy is important to ensure correct strain identification, and sensitivity is essential for assessing clonality by detecting lower frequency alleles such as those found in mixed infections or cross-contaminated lab cultures. Even a 10% contamination of a drug-sensitive strain with a drug-resistant strain can impact drug susceptibilities (2, 3). Cost-effectiveness is also necessary to ensure routine use. These methods should also allow the study of parasitized blood dried on filter paper, a routine means of storing samples in the field.

Historically, P. falciparum parasites were genotyped by assessing length- and sequence-polymorphic molecular markers such as the merozoite surface proteins msp1 and msp2 and the glutamine rich protein glurp (4–9). These markers can distinguish new infections from recrudescence, and help assess treatment efficacy (10). These three genes, however, are under host immune selective pressure and their analysis has been shown to underestimate treatment failure rates and insufficiently capture parasite genetic diversity (11). Thus, microsatellites and single-nucleotide polymorphisms (SNPs) are commonly used for genotyping parasites, and these analyses can be complemented by whole-genome sequencing (WGS) (12).

Microsatellites form when the DNA replication complex slips off the template strand of one repeat unit, and re-attaches at a different unit, creating “errors” or new strands with longer or shorter sequences (13). In P. falciparum, microsatellites are mostly (TA)n or (T or A)n, ranging from 10 to 30 repeats (14). SNPs most commonly form due to errors during DNA replication and after environmental exposures to DNA-damaging agents. Nonsynonymous SNPs or copy number variations (CNVs) can mediate resistance through amino acid substitutions and/or expression level changes in resistance determinants (15). Approaches to identify resistance determinants include molecular genotyping, genome-wide association studies, and positional mapping using genetic cross progeny (16–18).

Resistance to quinine, a drug used to treat severe malaria, has been partially associated with SNPs in the P. falciparum chloroquine resistance transporter (pfcrt) and P. falciparum multidrug resistance-1 (pfmdr1) genes (19–21). Resistance to chloroquine (CQ), the former front-line antimalarial, is mediated primarily by polymorphisms in PfCRT, of which K76T is critical, and can be modulated by polymorphisms in PfMDR1 (22–24). Novel SNPs encoding for mutations in PfCRT (including T93S, H97Y, F145I, I218F, M343L, and G353V) that emerged on the CQ-resistant Dd2 PfCRT variant, as well as amplifications in plasmepsins (pm) 2 and 3 genes, were recently shown to be the primary mediators of resistance to piperaquine (PPQ), a partner drug used in artemisinin (ART)-combination therapy (ACT) (25–27). PfMDR1 N86Y and D1246Y, as well as select PfCRT variants, have also been associated with in vitro and in vivo resistance to amodiaquine, another ACT partner drug (24, 28–31). Increased pfmdr1 copy number in vitro and in vivo has been associated with low-grade resistance to mefloquine, which is partnered with artesunate (1, 32, 33). SNPs encoding for mutations in the folate biosynthesis enzymes dihydrofolate reductase (DHFR) (N51I, C59R, S108N, I164L) and dihydropteroate synthase (DHPS) (A437G, K540E) confer resistance to pyrimethamine and sulfadoxine, respectively (15, 34). DHPS A581G associates with enhanced in vitro resistance to both drugs (35).

Delayed parasite clearance times after treatment with the current front-line antimalarial, ART, and its derivatives are driven by mutations in the Kelch13 (K13) propeller domain, of which C580Y dominates in the Greater Mekong Subregion (GMS) (36–38). Isolates harboring C580Y have also been observed in Guyana and Papua New Guinea (39, 40). ART resistance-conferring K13 mutations in the GMS (including R539T and C580Y) are not yet prevalent in Africa, although the emergence of the ART resistance-conferring R561H mutation in Rwanda and recent evidence that it can confer resistance in vitro raises important concerns (41–45). This worsening situation indicates an urgent need to monitor k13 as well as genetic markers of resistance to ACT partner drugs.

Microsatellites have also been useful for studying evolution, migration, and interrelatedness between strains, including analysis of linkage disequilibrium and selective sweeps in drug-resistant P. falciparum isolates (46). Microsatellite genotyping has revealed a common genetic signature surrounding the K13 C580Y variant (PfPailin lineage) in a majority of isolates across the GMS, indicating a transnational ART selective sweep that originated in Cambodia (47). This PfPailin lineage subsequently gave rise to the KEL1 (K13 C580Y)/PLA1 (multicopy pm2-3) co-lineage, which more recently includes new PfCRT variants and which is primarily responsible for DHA-PPQ treatment failure in the GMS (48–53). Sets of 12 microsatellite markers have also been used to genotype P. falciparum clinical isolates from finger-prick blood samples or to distinguish genetic cross progeny (54–56).

Here, we provide rapid and inexpensive tools for classifying strains, identifying drug resistance haplotypes, and assessing clonality. We describe a rapid and adaptable method to genotype SNPs from cultured parasites or infected blood, and compare this platform with optimized genotyping using either whole-genome sequencing or microsatellite markers that can distinguish geographically diverse laboratory strains.

RESULTS

PCR-based SNP genotyping from parasite cultures.

We selected a panel of 12 well-characterized laboratory-adapted P. falciparum strains representing diverse genomes, drug susceptibility profiles, and geographic origins, and used next-generation sequencing to analyze their genomes (Table 1; Fig. S4C). We then developed PCR and Sanger sequencing reaction protocols that can test directly from asexual blood stage parasites, without a genomic DNA (gDNA) extraction step, and that can be completed in less than a day (Fig. 1). Our nine PCR conditions and 16 Sanger sequencing reactions cover 33 amino acid changes in five genes associated with drug resistance: dhfr, dhps, pfmdr1, pfcrt, and k13 (Fig. 2A to C; Table S1). Primers and PCR conditions are listed in Table S1, and the protocol is listed in the supplemental methods. All 33 variant codons were assessed in six strains: NF54, Cam3.II, RF7, GB4, Dd2, and 7G8 (Fig. 2B and C). All of the amino acid sequences determined by PCR/Sanger sequencing correctly matched the whole-genome sequence data, validating our methods.

TABLE 1

All single-nucleotide polymorphisms (SNPs) and their corresponding amino acid mutations assessed by PCR/Sanger sequencing of parasite cultures for a panel of 15 P. falciparum strains

Six other SNPs can be assessed by our parasite culture PCR conditions but are not listed in this table, as all 15 strains have the wild-type allele: (pfmdr1 G293, and pfcrt T93, F145, I218, C350, G353).

Gene IDs are dhfr (PF3D7_0417200), dhps (PF3D7_0810800), pfmdr1 (PF3D7_0523000), pfcrt (PF3D7_0709000), and k13 (PF3D7_1343700).

The full names for the antimalarials are PYR (pyrimethamine), CYC (cycloguanil), SDX (sulfadoxine), MFQ (mefloquine), QN (quinine), ADQ (amodiaquine), CQ (chloroquine), PPQ (piperaquine), and ARTs (artemisinins).

3D7 is a clone of NF54; as such, these are considered to be the same strain.

FIG 1

Workflow of the single-nucleotide polymorphism (SNP) and microsatellite (MS) genotyping methods described herein. SNPs can be identified using PCR and Sanger sequencing of parasite cultures or whole-genome sequencing (WGS) of purified parasite DNA. PCR-based SNP genotyping utilizes asexual blood stage cultures whereas MS genotyping and WGS require genomic DNA. MS sizes can be identified visually by gel electrophoresis or quantitatively by fragment analysis (FA) or by WGS. FA can be conducted on a DNA analyzer such as the ABI3730xl or SeqStudio Genetic Analyzer (Applied Biosystems). WGS requires a quality control step using the Bioanalyzer or TapeStation systems (Agilent). Protocol durations are for groups of up to ∼24 samples and minimum amounts of DNA typically used per protocol are indicated. Image was created with BioRender.com (refer to Tables S1 and S2 and the supplemental methods).

FIG 2

Parasite culture-based PCR genotyping to assess 33 single-nucleotide polymorphisms (SNPs) in dhfr, dhps, pfmdr1, pfcrt, and k13 (see Table S1 and the supplemental methods). (A) Amino acid mutations in these drug resistance determinants can be assessed by each PCR (separated by vertical lines) and are listed with forward, reverse, and sequencing primers that are color-coded by the sequencing primer. (B) Sanger sequencing electropherograms of these amino acid positions for NF54 and Cam3.II parasites. (C) Sensitivity test results of the parasite culture PCR/Sanger sequencing method for the dhfr genotype in nine mixtures containing defined ratios of NF54 and Cam3.II cultures (see Fig. S1 for other positions). Representative electropherograms showing S108 (NF54) or S108N (Cam3.II) alleles.

To assess sensitivity, we assayed nine mixtures of synchronized blood stage cultures of NF54 and Cam3.II strains at different ratios (0:100; 5:95; 10:90; 30:70; 50:50; 70:30; 90:10; 95:5; 100:0) and visually examined the electropherograms (Fig. 2C; Fig. S1). PCR/Sanger sequencing was sensitive up to a 70:30 or 90:10 polyclonal mix, with sensitivity as high as 95:5 in one instance. This finding suggests that PCR/Sanger sequencing is useful for quickly validating strains and profiling drug resistance markers, but has limited ability to detect minor alleles in a mixed sample.

Evaluation of microsatellite-based gel electrophoresis and fragment analysis.

We next assessed microsatellite-based genotyping, and identified a set of five microsatellite markers (Table 2) that uniquely identify our panel of 12 P. falciparum strains (Fig. 1). Markers were selected from a published list (56) based on their distribution on distinct chromosomes, size difference relative to 3D7 (mostly ≥ 6 bp), and the minimum set of markers required to distinguish these strains. For this PCR-based approach, we compared agarose gel electrophoresis with capillary-based fragment analysis (FA). We also optimized microsatellite typing methods and assessed their resolution and sensitivity (see Table S2 and the supplemental methods). The KAPA HiFi DNA polymerase (Roche) was chosen due to its documented low error rate and efficient amplification across highly AT-rich genomes (57).

TABLE 2

Microsatellite sizes determined from whole-genome sequencing (WGS) of the 15 P. falciparum strains assayed herein

Indicates the five MS markers selected for genotyping laboratory strains.

For each strain, we observed expected size differences relative to the 3D7 reference genome (Fig. 3A; Fig. S2A). However, for the FA method, the absolute sizes of the microsatellite products varied occasionally between runs. Thus, we recommend that for both methods a reference strain (e.g., 3D7) should be included in every run. The correlation of the relative microsatellite size from FA versus WGS yielded a mean R2 = 0.98 (Fig. S2B) across these five microsatellite markers. Gel electrophoresis had a minimum resolution of 9 bp difference between strains, while FA was accurate to within a 2 bp difference.

FIG 3

Selection of five microsatellite (MS) markers to distinguish 12 P. falciparum strains by gel electrophoresis or fragment analysis (FA) (see Table 2). (A) Gel electrophoresis image capturing the MS PCR product sizes that vary between strains. (B) Sensitivity of gel electrophoresis MS genotyping for C4M30 tested in 11 mixtures containing different ratios of HB3:Dd2 genomic DNAs (see Fig. S2 for remaining markers). (C) Sensitivity of FA genotyping of C4M30 in these mixtures to determine the proportion of HB3 and Dd2 fragment sizes. Bars represent the average fragment size determined from the fluorescence signals in two experiments. Error bars represent SEM values. L, ladder. (D) FA electropherograms obtained from fluorescently labeled PCR products multiplexed in two groups for NF54 and Cam3.II parasites. FAM, fluorescein amidite; RFU, relative fluorescence unit.

We also assessed the sensitivity of gel electrophoresis-based and FA-based microsatellite genotyping. This was tested with 11 mixtures containing different ratios of HB3:Dd2 input gDNA (0:100; 2.5:97.5; 5:95; 10:90; 30:70; 50:50; 70:30; 90:10; 95:5; 97.5:2.5; 100:0). Mixtures were genotyped at the five microsatellite markers (Fig. 3B and C; Fig. S3A and B). FA-based detection of the minor population was sensitive to as low as 2.5% for all microsatellite markers, whereas gel electrophoresis could only detect down to 10% for C4M30, C3M67, and TA1; and 30% for TAA81 and TA127. To increase throughput, we also multiplexed FA-based microsatellite genotyping. Forward primers were fluorescently labeled with 6-FAM (blue) and ATTO550/NED (yellow) and used in two separate PCRs for C2M18 and BM5, respectively, for each strain (Fig. 3D; Multiplex Group A). We also labeled forward primers with 6-FAM, ATO550, and ATTO532/VIC (green) for C4M30, C13M87, and TA127, respectively (Multiplex Group B). PCR products were mixed at either a 2:1 ratio for Multiplex Group A or at a 1:1:1 ratio for Multiplex Group B. These multiplex runs showed 100% accuracy for the five markers.

Classification of genetic cross progeny by whole-genome sequencing versus fragment analysis-based microsatellite genotyping.

SNP and microsatellite genotyping methods from parasitized blood are less expensive, faster, and more accessible than WGS. In addition to characterizing lab strains, microsatellite genotyping can also be used to classify P. falciparum genetic cross progeny without existing WGS data (Fig. 4). This is of increasing relevance given the recent availability of humanized mouse models to generate genetic crosses (58), and the need to ensure clonal integrity during continuous cell culture. To assess feasibility, we used our five microsatellite markers to genotype parents of the HB3×Dd2 genetic cross and its representative progeny 1BB5, 3BA6, and GC03 (59). From the WGS data, heatmaps were generated to identify the location and length of microsatellite insertions and/or deletions, relative to the 3D7 reference (Fig. 4A and B). Microsatellite sizes relative to 3D7 were also obtained from the FA data (Fig. 4C). Results correctly assigned each microsatellite marker in the progeny to parental HB3 or Dd2 origins, consistent with WGS data for the progeny. Thus, in addition to ensuring clonal integrity of cross progeny clones, FA typing with microsatellite markers can also be used to screen for recombinant progeny as long as they are polymorphic between the genetic cross parents.

FIG 4

Comparison of microsatellite (MS) typing by fragment analysis (FA) and whole-genome sequencing (WGS) analysis to classify parents and progeny of the HB3×Dd2 genetic cross. (A) Heatmap of the WGS MS typing data indicating the fragment insertion or deletion size for each parasite relative to 3D7. The MS sequence for 3D7 (pink) and the position at which the insertion or deletion occurred within the MS fragment window are shown. The percentage of sequence reads covering the MS region that have base insertions or deletions are shown by the blue or red scale, respectively. The net MS sizes relative to the 3D7 reference are shown from (B) WGS-based or (C) FA-based genotyping methods.

Application of SNP and microsatellite genotyping approaches to clinical samples.

We tested whether our optimized PCR/Sanger sequencing methods using parasite cultures could be applied to simulated fresh and frozen whole-blood samples, as well as dried blood spots (DBSs), by using P. falciparum laboratory cultures mixed with whole-blood (final hematocrit 40–50%) (Table S3). While PCR/Sanger sequencing using simulated fresh and frozen whole-blood samples were able to detect dhfr, dhps, pfmdr1, pfcrt, and k13 down to 0.4% parasitemia, those using gDNA extracted from the DBSs were the most robust, with sensitivity down to 0.02% parasitemia.

We further investigated the applicability of these genotyping methods using gDNA extracted from DBSs collected from 16 symptomatic malaria patients in eastern Uganda in 2019 (Table S4). We were able to robustly determine the drug resistance alleles for all five genes. Almost all samples had multiple mutations in DHFR (N51I/C59R/S108N) and DHPS (K540E/A437G), 6/16 were pure PfMDR1 Y184F mutants (two more were mixed infections), and one sample was a PfCRT K76T mutant. No samples harbored any nonsynonymous mutation in the K13 beta-propeller domain. We also characterized these samples using microsatellite genotyping and fragment analysis. Results demonstrated greater genetic diversity among parasites using microsatellites compared to SNP genotyping. PCR/Sanger sequencing and FA analysis identified polyclonal infections in six and nine samples, respectively (Table S4).

SNP versus microsatellite-based clustering of geographically diverse parasites.

We also assessed whether 27 discriminatory SNPs in drug resistance genes or eight microsatellite markers could segregate our set of 12 laboratory strains according to their geographical origin (Fig. 5A and B). Both SNP and microsatellite genotyping approaches grouped the African strains. The 7G8 and HB3 strains from South and Central America, respectively, segregated correctly only with the SNP method. They did, however, consistently segregate closer to Southeast Asian lines than to African lines with both methods. Overall, these data suggest that the drug resistance-conferring SNPs and a small set of microsatellite markers can broadly segregate parasites by continent.

FIG 5

Hierarchical clustering of laboratory strains and Ugandan dried blood spot samples by their respective geographical origins using microsatellite (MS) or single-nucleotide polymorphism (SNP) markers. (A) Heatmap showing the clustering of 12 laboratory lines using 27 discriminatory SNPs in dhfr, dhps, pfmdr1, pfcrt, and k13. The other six SNPs were identical across all strains tested. (B) Heatmap depicting these laboratory strains clustered by eight MS markers. Parasite geographical origins are indicated by font color. Major versus minor alternate alleles, and the size of MS insertions or deletions, are indicated by the color scale. (C) Heatmap showing the clustering of 16 Ugandan samples alongside two control lines using 18 SNPs in dhfr, dhps, pfmdr1, pfcrt, and k13. (D) Heatmap of Ugandan samples and two control lines clustered by five MS markers (also see Table S4).

To identify clades within the Ugandan clinical samples, we also conducted SNP- and microsatellite-based clustering (Fig. 5C and D). Using SNPs, Ugandan parasites clustered into three major groups segregating by pfmdr1 and pfcrt status, whereas microsatellite clustering revealed five major groups. Interestingly, TO-04-22 appeared identical to U659, another Ugandan parasite isolated in 2007, by both SNP and microsatellite clustering, suggesting a shared ancestry.

Integrated analysis of KEL1/PLA1/PfPailin (KelPP) status.

We also applied our k13 SNP genotyping method, and previously published k13-flanking microsatellite genotyping and pm2 copy number quantification methods (26, 47), to test whether parasites had a KelPP co-lineage haplotype associated with an increased risk of DHA-PPQ treatment failure (Fig. 6, Table S5, and the supplemental methods). As expected, the k13-flanking microsatellite sizes for African and South American strains differed substantially from those of the Cambodian strains. The latter all shared the same microsatellite sizes on the four upstream (–50 kb to −0.15 kb) markers, but differed at the +8.6 kb marker. Almost all sampled Cambodia isolates from 2010 onwards had mutant K13 C580Y, consistent with the emergence in ∼2008 of the multidrug-resistant PfPailin lineage (47). RF7 was the only Cambodian strain to have both multicopy pm2 and mutant K13 and corresponds herein to the KelPP haplotype. This finding is also consistent with the rapid replacement of the KEL1 lineage with the KEL1/PLA1 co-lineage in Cambodia (38, 51, 52, 60).

FIG 6

Testing for the KEL1/PLA1/PfPailin (KelPP) co-lineage in 10 geographically diverse P. falciparum strains, based on kelch13 (k13)-flanking microsatellite sizes relative to 3D7, K13 mutations, and plasmepsin 2 (pm2) copy number (see Table S5 and the supplemental methods). The sizes of the k13-flanking microsatellite insertion or deletion are indicated by the color scale. RF7 was the sole KelPP strain in this set.

DISCUSSION

Herein, we describe several methods to genotype P. falciparum parasites for a range of experimental needs that include validating parental or genetically modified lab strains, profiling new field isolates, or identifying genetic cross recombinant progeny. These methods include PCR/Sanger sequencing from parasite cultures, and optimization of a set of five polymorphic markers for multiplexed microsatellite genotyping. The PCR/Sanger sequencing-based SNP genotyping method, applied directly to parasite cultures, is quick and accurate, making it ideal to routinely validate lab strains and confirm gene editing in transfected parasites (Table 3). This approach can help profile drug resistance markers and identify known and novel mutations with culture-adapted isolates. Our genotyping was also optimized to identify the KEL1/PLA1/PfPailin (KelPP) co-lineage that has spread rapidly throughout the GMS. This PCR-based method of detecting K13 C580Y directly from P. falciparum-infected red blood cells (RBCs) can expedite the identification of ART-resistant clinical isolates, and help prioritize further screening for the KelPP co-lineage. However, our panel of SNP-based genotyping markers will not detect emerging and/or novel drug resistance mutations outside the amplified gene fragments. In contrast, WGS analysis provides comprehensive SNP genotyping that can capture novel mutations, elucidate genetic diversity, and provide insights into the selective pressures impacting parasite genomes (61, 62). However, WGS is cost-prohibitive, time-consuming, and requires specialized instrumentation and bioinformatic tools, limiting its adoption for routine genotyping.

TABLE 3

Comparison of the genotyping methods used in this study

	Typing method	Resolution	Sensitivity	Time	Cost	Pros	Cons	Ideal application
SNP	Parasite culture PCR + Sanger sequencing	1 bp	≥10%	∼16 h	Moderately expensiveb	Quick, accurate	Requires specialized equipment, genetic info only around SNP	Validate strains, drug resistance profiling, confirm gene-editing
MS	Gel electrophoresis	9 bp	≥10%	∼7 h	Inexpensive	Quick, no specialized equipment required	Lowest resolution, visual determination, genetic info only on MS	Validate strains, identify novel genetic cross progeny
MS	Fragment analysis	2 bp	≥2.5%	>2 days	Moderately expensive	Quantitative, moderate sensitivity	Requires specialized equipment, genetic info only on MS	Validate strains, detect cross-contamination, identify novel genetic cross progeny, infer origin
SNP + MS	WGS	1 bp	≥4%c	>5 days	Expensive	High sensitivity and accuracy, provides unbiased genotyping data	Expensive, time-consuming, requires specialized equipment and bioinformatics skills	Drug resistance profiling, detect contamination, identify novel genetic cross progeny, infer origin

SNP, single-nucleotide polymorphism; MS, microsatellite; WGS, whole-genome sequencing.

The cost per SNP using parasite culture PCR genotyping depends on the PCR reaction used.

Assuming a mean fold coverage of 50× or greater.

Alternative TaqMan- or GoldenGate-based methods described by other groups utilize a custom panel of informative SNPs to create parasite “barcodes,” identify lineages, and explore external selection pressures and parasite transmission (9, 63–66). These customized methods, which require gDNA extraction, were designed to discriminate between two separate DNA templates differing at a predefined SNP position. TaqMan-based allelic discrimination was also applied to genotype pfcrt K76T directly from whole blood (67). Another approach employs a high-resolution melting analysis, in which <250 bp amplicons harboring SNP positions of interest are melted post-PCR, and variants are detected if the shape of the melt curves and T values differ from those of the wild-type control (68, 69). While this method can detect the presence of more than one variant in each SNP position and several SNPs within a single amplicon, it requires Sanger sequencing to elucidate variants. Other genetic epidemiology tools include amplicon sequencing (70, 71), multiplexed mass spectrometry-based genotyping (71), high-throughput targeted sequencing using molecular inversion probes (72), and ligase detection reaction fluorescent microsphere assays (43). These methods, however, require gDNA extraction and specialized equipment. Amplicon sequencing and molecular inversion probes must be combined with next-generation sequencing while ligase detection reaction fluorescent microsphere assays can only detect specific alleles on a bead array laser scanner. Our method of SNP/Sanger sequencing directly from parasite cultures does not require gDNA extraction, can detect multiple alleles in the amplified fragment, is cost-effective, and is readily customizable. These attributes will be valuable in detecting novel K13 resistance-conferring mutations, especially in high-transmission regions in Africa (41). We show herein that this method can be successfully applied to define drug resistance genotypes directly from fresh or frozen parasite samples as well as gDNA extracted from DBSs from P. falciparum-infected patients.

Gel electrophoresis-based microsatellite genotyping can be readily implemented, is quick, inexpensive, and has sufficient accuracy and resolution for validating strains in the laboratory setting and/or in the field. This method is also amenable to identifying recombinant genetic cross progeny. However, given the low resolution and qualitative assessment, some strains will be difficult to genotype if the size difference between strains is <9 bp or if a minor population is <10%. We found the fluorescence-based capillary electrophoresis FA genotyping method to have a higher resolution (2 bp) and greater sensitivity (minor parasite population >2.5% of total). This method is also quantitative, making it suitable for detecting cross-contamination in a laboratory culture or polyclonal infections in field samples, as we showed with isolates from Uganda. Another advantage is that FA can be multiplexed either pre- or post-PCR to increase throughput. We recommend the latter as melting temperatures (T) may vary between primers. Also, microsatellites are often laden with AT-repeats, making it difficult to design compatible primer pairs that lack complementarity. While microsatellites have been previously multiplexed by hemi-nested PCR to assess field isolates (55), these markers may not be as relevant for multiplexed genotyping of common lab strains. Microsatellite genotyping can also be conducted directly from asexual blood stage P. falciparum cultures, as shown previously (73).

Although less expensive and quicker than WGS, microsatellite genotyping is inferior to WGS in a few aspects. WGS genotyping does not need a positive control strain (e.g., 3D7) for each run, can identify the genomic location of microsatellite insertions and deletions, and can assess a much larger number of microsatellites (74). However, errors can be introduced in this method because of the similarity in length between Illumina reads and microsatellite arrays that can complicate sequence alignment, and homopolymers that can cause sequencing errors (75, 76).

On average, microsatellites account for 3% of the human genome and ∼11% of the P. falciparum genome (76, 77). P. falciparum is thought to have an especially high number of microsatellites in part due to the AT-rich genome (81% overall; 90% in noncoding regions) that may increase the chances of polymerase slippage and AT-repeat microsatellite formation (12). Across the P. falciparum genome, an average of 80% of insertions and deletions (indels) are short tandem repeats that include microsatellites, and 83% of indels are located in noncoding regions (78). The indel:SNP ratio is ∼3:1 and ∼1:2 in noncoding and coding regions, respectively, reflecting greater constraints on nucleotide diversity in genes. Our observation that drug resistance-conferring SNPs are subjected to stronger selective pressures than microsatellites was evident in the geographical segregation of parasite lines through our phylogenetic analyses using SNPs or microsatellite markers (Fig. 5).

Despite their seeming potential for increased mutability, microsatellites are relatively stable in P. falciparum. In an in vitro evolution experiment, the mutation rate for 3D7 P. falciparum microsatellites was on average 3.1 × 10−7 to 2.2 × 10−8 per cell division, with longer microsatellite array lengths and motifs correlating with higher mutation rates (76). This mutation rate is low when compared with those of other organisms including humans (∼10−3 per locus per gamete), Escherichia coli (∼10−2), Saccharomyces (∼10−5), and Drosophila (∼10−6 to 10−4) (79, 80), possibly because Plasmodium spp. lack the nonhomologous end-joining pathway during DNA repair, which may reduce replication errors during DNA damage repair (76, 81, 82). Nonetheless, we cannot disregard the possibility that microsatellite sizes may change over long-term in vitro culture (83).

The five microsatellite markers that we selected span 3′ untranslated regions (UTRs), or exon-UTR or exon-intron boundaries of essential genes (Table S2) (84). The repeat motifs in these microsatellite markers were mostly di- to tetra-nucleotide repeats (Table S2). As larger motif sizes are associated with higher microsatellite mutation rates in vitro (76), the small repeat motif sizes in our microsatellites suggest a low propensity to mutate in culture. Thus, our minimal set of five microsatellite markers should be robust in genotyping strains of diverse geographical and temporal origin. This approach can also be applied to multiplexed microsatellite genotyping of genetic cross progeny.

Taking the above assessments into account, we propose the following (Table 3): Parasite strain identity can be validated using PCR/Sanger sequencing-based SNP genotyping of parasite cultures, combined with microsatellite genotyping. These two methods could potentially be used to characterize parasite strains in large-scale surveys. Laboratory strain cross-contamination can be readily detected using FA microsatellite genotyping. This is especially useful when phenotyping P. falciparum genetic cross progeny, as it is critical to verify the identity and clonal purity of individual progeny to accurately map quantitative traits. To determine the multiplicity of infection in a field sample and identify new resistance determinants across the genome, more extensive methods such as amplicon sequencing and WGS are necessary (70, 85–89). Given the greater accessibility, ease, and throughput nature of our combined genotyping methods, we envision that these will be useful in the lab and field for routine use. Our adaptable genotyping methods can also help monitor the spread of novel drug-resistant traits once their causal determinants have been identified.

MATERIALS AND METHODS

Parasite cultures and strains.

Asexual blood stage P. falciparum parasites were cultured in human RBCs obtained from the Interstate Blood Bank (Memphis, TN) at 3% hematocrit, using RPMI 1640 supplemented with 25 mM HEPES, 2.1 g/liter sodium bicarbonate, 10 μg/ml gentamicin, 50 μM hypoxanthine, and 0.5% (wt/vol) AlbuMAX II (Thermo Fisher Scientific). Parasites were maintained at 37°C under 5% O2/5% CO2/90% N2 gas conditions. The GB4 parasite was obtained through BEI Resources, NIAID, NIH (MRA-925, contributed by Dr. Thomas Wellems). Other P. falciparum lines were kindly provided by Drs. Philip Rosenthal, Didier Ménard, Rick Fairhurst, Sarah Volkman, and Lise Musset (Table S5). Isolate 3D7 is a clone of NF54, originally isolated in 1979 in the Netherlands, and likely of African origin (41, 90).

Simulated clinical samples.

NF54- or Cam3.II-infected packed RBCs were mixed with whole-blood obtained from the New York Blood Center and serially diluted 5-fold to generate four samples each with final parasitemias of 2.0% to 0.016% and hematocrits of 40%–50% (N = 2). Samples of 200 μl were used to make simulated dried blood spots (DBSs) on Whatman GB003 filter paper and processed 5 days later, while 100 μl each were used directly for simulated fresh whole-blood PCR, or frozen at−20°C and thawed 5 days later for simulated frozen whole-blood PCR.

Ugandan DBSs.

We studied 16 isolates from patients who presented with symptomatic P. falciparum infections in eastern Uganda from April 2019 to June 2019, as described previously (43). In brief, subjects diagnosed with P. falciparum malaria and who consented to this study provided blood samples prior to malaria treatment. Blood was dried on Whatman 3MM filter paper, and samples were stored in zip-lock bags with desiccant at room temperature before extraction of DNA and molecular analysis. This study was approved by the Makerere University Research and Ethics Committee, the Uganda National Council of Science and Technology, and the University of California, San Francisco Committee on Human Research.

SNP genotyping using PCR/Sanger sequencing.

The 40 μl PCRs consisted of 3 μl of 1%–5% parasitemia packed RBCs or simulated whole-blood samples, 0.125 μM final concentration of each forward and reverse primer, and 2× MyTaq Blood-PCR kit (Bioline) (supplemental methods). The following PCR conditions were used: 95°C for 3 min, 35 cycles: 95°C for 30 s, the appropriate annealing temperature and duration, and extension at 62°C for durations listed in Table S1. PCR products were submitted for Sanger sequencing (Genewiz) and analyzed with SeqMan Ultra (Version 17, DNASTAR). To generate NF54 and Cam3.II parasite mixtures, these lines were synchronized using 5% D-sorbitol (Sigma) and their parasitemias measured by flow cytometry (91) (N = 3). When starting with genomic DNA samples, the same annealing and extension PCR conditions were used, but the MyTaq polymerase was replaced with the 2× KAPA HiFi PCR Mix (Roche) in 15 μl total reaction volumes, and the denaturation step was performed at 98°C for 20 s. For each PCR reaction, we used 4 μl of gDNA from either Ugandan or simulated DBSs.

Genomic DNA (gDNA) extraction.

gDNA was obtained from 3%–7% parasitemia cultures by lysing with 0.1% saponin in 1× PBS, and extracting with the QIAamp DNA Blood Midi Kit (Qiagen), with a combined RNase and Proteinase K treatment. DNA concentrations were determined using the Qubit dsDNA BR assay kit (Thermo Fisher Scientific). gDNA was extracted from both clinical (3 mm punches) and simulated DBSs (6 mm punches) using the QIAamp DNA Investigator kit (Qiagen) according to the manufacturer’s protocol, and eluted in 40 μl of AE buffer.

Microsatellite PCR gel electrophoresis and FA genotyping.

PCRs with 15 μl volumes consisted of 12 ng of gDNA, 0.3 μM final concentration of each forward and reverse primer (with the forward primer fluorescently labeled for FA), and 2× KAPA HiFi PCR Mix (Roche) (supplemental methods). PCRs were run with the following conditions: 95°C for 3 min, 30 cycles: 98°C for 20 s, annealing temperature and duration as indicated in Table S2, and extension at 65°C for 15 s. For the gel electrophoresis method, 5 μl of each PCR product was run on a 2% agarose gel at 100V for 1.5 h–2 h with an Ultra Low Range DNA Ladder (Thermo Fisher Scientific) on both ends of the gel. The bands were visualized using a Gel Doc Imager (Bio-Rad) to identify the microsatellite sizes. For the FA method, fluorescently-labeled PCR products were treated with ExoSAP-IT Express Reagent (Thermo Fisher Scientific) and sent for capillary electrophoresis FA with the Liz500 ladder (Genewiz). Data were analyzed using Peak Scanner 2 (Thermo Fisher Scientific) to determine the microsatellite size of the fragment, which corresponds to the peak with the largest area in the electropherogram. The following parameters were used for analysis: peak sizes between 100 bp and 300 bp, peak height threshold of >500 RFU for the target channels (e.g., blue for 6-FAM) to eliminate PCR artifacts, and absolute and relative microsatellite sizes determined relative to the 3D7 reference strain. To multiplex FA, PCRs were set up individually for C2M18, BM5, C4M30, TA127, and C13M87, using forward primers fluorescently labeled with 6-FAM (blue), ATTO550/NED (yellow), or ATTO532/VIC (green) as listed in Table S2. PCR products from C2M18 and BM5 were mixed at a 2:1 ratio for each strain, and C4M30, TA127, and C13M87 were mixed at a 1:1:1 ratio, cleaned up with ExoSAP-IT, and sent for FA. Ugandan samples were run as 15 μl PCRs with 4 μl of gDNA for 35 cycles, cleaned up with EXO-SAP-IT, diluted 1:5, and sent for FA as singleplex reactions.

Clustering of parasite strains.

The 12 laboratory strains and 16 Ugandan isolates (alongside two laboratory strain controls) were hierarchically clustered by average linkage using Gene Cluster 3.0 (92) and the resulting dendrograms were visualized using Java TreeView (93). For SNP clustering of laboratory strains, the alleles were input as follows: wild-type reference = 0; most common “major” alternate = 2; less common “minor” alternate = 1. For Ugandan samples, the alleles were input as: wild-type homozygous reference = 0; homozygous alternate = 2; and heterozygous = 1. For microsatellite clustering of both laboratory strains and Ugandan samples, the microsatellite sizes relative to 3D7 were used and only the major peaks were used for clustering the Ugandan isolates (Fig. 5; Tables 1, 2, and S4).

Quantitative PCR (qPCR) assessment of plasmepsin 2 (pm2) copy number.

Multiplexed qPCR of pm2 and single-copy β-tubulin (as an internal control) used TaqMan assays on a QuantStudio 3 real-time PCR system (Thermo Fisher Scientific), as described previously (26). Five standards of gene fragments, mixed at 1:1, 2:1, 3:1, 4:1, and 5:1 molar ratios of pm2:β-tubulin, were included as copy number controls (49). The reaction efficiency was calculated from previous runs, and the final copy number was calculated using the standard curve method for relative quantification (n = 3). Refer to the supplemental methods for the full protocol.

Whole-genome sequencing (WGS) and microsatellite calling.

Parasite samples were subjected to WGS using the Illumina Nextera DNA Flex library preparation protocol. Briefly, 150 ng of gDNA was fragmented and tagmented using bead-linked transposons and subsequently amplified by five cycles of PCR to add dual-index adapter sequences to the DNA fragments. The libraries were quantified, pooled, and sequenced on Illumina MiSeq or NextSeq platforms. 3D7, NF54, GB4, 7G8, FCB, Cam3.II, RF12, RF7, Dd2, HB3, GC03, 1BB5, and 3BA6 were sequenced on a MiSeq instrument to generate paired-end 300 bp reads. F32, U659, and U815 were sequenced on a NextSeq 550 instrument to obtain 150 bp paired-end reads (Fig. S4C).

Sequence data were aligned to the P. falciparum 3D7 genome (PlasmoDB version 36.0) using Burrow-Wheeler Alignment (BWA). Reads that did not map to the reference genome and PCR duplicates were removed using SAMtools and Picard. The reads were realigned around indels using Genome Analyses Tool Kit (GATK) RealignerTargetCreator and base quality scores were recalibrated using GATK Table-Recalibration.

SAMTools mpileup was used to find SNPs in specified genomic regions (see Table 1) and were filtered based on quality scores (minimum base quality score ≥20, mapping quality >30, read depth ≥5) to obtain high-quality SNPs that were annotated with snpEFF. To determine the confidence of a microsatellite based on the number of reads, a custom Python script utilizing the pysam (94) module was written to call reads that harbored insertions or deletions at specified genomic loci within the respective windows (see Table 2). Integrative Genomics Viewer was used to verify the data.