RESULT
Genetic differences, population structure, ancestry, and infection complexity
All five microsatellite loci were highly polymorphic, with the number of alleles ranging from 8‐12 and an average of 9.8 alleles per locus. Pm_09 was the most polymorphic microsatellite, while Pm_11 had the lowest variability (Table 4). The expected genetic diversity (Nei’s genetic divergence) averaged across all loci was 0.72, ranging from 0.61 (Pm_11) to 0.85 (Pm_02). This parameter estimates the probability that two randomly selected genotypes are different, with a score ranging from 0 (no genotypes are different) to 1 (all genotypes are different). The average genotype evenness was 0.65; This evenness index measures the frequency distribution of genotypes, where a population dominated by a single genotype has a value close to 0 and a population with genotypes at equal frequencies has a value close to 1.18.
Table 4. Summary statistics of microsatellite loci showing richness, diversity and genotype evenness at each locus.
| Locus | Allele | 1‐D | Hexp | Evenly |
| Pm_09 | 12 | 0,75 | 0,76 | 0,63 |
| Pm_34 | 9 | 0,67 | 0,68 | 0,63 |
| Pm_11 | 8 | 0,6 | 0,61 | 0,52 |
| Pm_02 | 10 | 0,84 | 0,85 | 0,84 |
| Pm_47 | 10 | 0,7 | 0,71 | 0,63 |
| Medium | 9,8 | 0,71 | 0,72 | 0,65 |
1‐D: Diversity Index Simpson, Hexp: Nei’s unbiased genetic diversity.
At the population level within each country, only two pairs of multi‐locus genotypes (MLGs) were found in isolates from Burkina Faso and Nigeria. The diversity and richness of the genotypes were high, as demonstrated by the high values of the Lambda index or Simpson index and Nei’s expected heterozygosity‐unbiased genetic diversity (Table 5).
The lowest heterozygosity was observed in isolates from Cameroon (0.47), while the highest was observed in isolates from Ghana (0.80), both comprising only 3 isolates. In total, 20 distinct alleles were recorded across all populations, of which 50% were found in Tanzania. Genotypic uniformity within populations was found to be high (overall E.5 = 0.972).
The degree of linkage disequilibrium between loci in highly variable populations, as determined by the Index of Association (IA), ranged from 0 in Ghana to 1 in isolates from Cameroon.
Table 5. Summary statistics of microsatellite loci showing richness and diversity
and the level of genotype uniformity in each population
| Population | N | MLG | eMLG | SE | PA | H | G | lambda | E.5 | Hexp | IA |
| Burkina Faso | 17 | 16 | 9,67 | 4,71e−01 | 3 | 2,75 | 15,2 | 0,934 | 0,969 | 0,746 | 1,28e−01 |
| Cameroon | 3 | 3 | 3,00 | 0,00e + 00 | 0 | 1,10 | 3,0 | 0,667 | 1,000 | 0,467 | 1,00e + 00 |
| Ghana | 3 | 3 | 3,00 | 0,00e + 00 | 2 | 1,10 | 3,0 | 0,667 | 1,000 | 0,800 | −1,11e−16 |
| Guinea | 4 | 4 | 4,00 | 0,00e + 00 | 1 | 1,39 | 4,0 | 0,750 | 1,000 | 0,667 | − 5,56e−01 |
| Mali | 9 | 9 | 9,00 | 0,00e + 00 | 2 | 2,20 | 9,0 | 0,889 | 1,000 | 0,733 | 5,16e−01 |
| Nigeria | 18 | 17 | 9,71 | 4,56e−01 | 2 | 2,81 | 16,2 | 0,938 | 0,970 | 0,703 | 2,83e−01 |
| Southeast Asia | 1 | 1 | 1,00 | 0,00e + 00 | 0 | 0,00 | 1,0 | 0,000 | NaN | NaN | NaN |
| Tanzania | 19 | 19 | 10,00 | 2,51e−07 | 10 | 2,94 | 19,0 | 0,947 | 1,000 | 0,722 | − 1,11e−01 |
| Total | 74 | 70 | 9,93 | 2,54e−01 | 20 | 4,23 | 66,8 | 0,985 | 0,972 | 0,722 | 1,01e−01 |
Number of observed individuals (N), number of observed multilocus genotypes (MLG);
Expected number of MLGs with minimum sample size ≥10 based on rarefaction method (eMLG);
Standard error based on eMLG (SE); Number of distinct alleles (PA); Shannon‐Wiener index of MLG diversity (H);
Stoddart and Taylor index of MLG diversity (G); Simpson index (lambda), and E.5 evenness index (E5E5).
Hierarchical heatmap based on genetic distances from microsatellite and SNP data showing relationships between samples P. malariae regardless of geographic origin (Fig. 2), with samples from West, Central, and East Africa generally clustering into the same few branches. Three genetic clusters were observed in the SNP data, indicating lower overall genetic distances due to low frequencies of minor alleles, which was more evident with the filtered and denoised datasets (Fig. 2b,c). The distribution of genetic distances between the SNP and microsatellite data differed, with wider distances between a small number of samples observed from the microsatellite data (Supplementary Fig. 1).
This is understandable given that microsatellites have many alleles, are likely to be neutral, and are more easily differentiated between pairs of samples. Wright fixation index (FST) estimates of differentiation between populations were affected by small sample sizes but were generally low, suggesting poor differentiation between national populations due to intrapopulation genetic variation. However, SNP data showed relatively high genetic distance between Cameroon (which had the smallest sample size) and Burkina Faso (FST = 0.437) (Supplementary Figure 2).

Figure 2. Clustering and heatmaps of pairwise genetic distances were calculated using the pHeatmap_1.0.12 package in R version 4.1.13, including (a) Bruvo distances from microsatellites (Msat), (b) Nei genetic distances between individual samples from SNPs in potential resistance loci (unfiltered data), and (c) Nei genetic distances calculated from SNPs after cleaning and filtering the read data. The origin of each sample by country or geographic region is represented by sidebars at the branch tips of the hierarchical trees.
Additionally, the microsatellite data identified five genetic clusters, but when visualized by DAPC scatter analysis, these clusters were not identified based on the geographic origin of the isolates. Similarly, SNP‐based clustering using the unfiltered dataset grouped the isolates into five clusters, each consisting of isolates from different geographic origins (Figure 3). Although the geographically independent clustering pattern was still observed in the filtered and cleaned datasets, the clusters became less clear as the number of isolates decreased, especially in the cleaned data (Figures 3c,d).
Using STRUCTURE software, the optimal number of clusters by the Evanno method was determined to be K = 3 (Figure 4a, Figure 4b), although additional peaks were also detected on the Evanno plot at K = 6 and K = 8. With a 70% probability threshold for assigning an individual to a particular cluster, the admixture model grouped the isolates into three ancestral clusters (Figure 4c), none of which were characteristic of isolates from any geographical area.

Figure 3. Discriminant analysis based on scatterplot principal component (DAPC) display of cluster populations of samples P. malariae using (a) microsatellites, (b) unfiltered SNPs, (c) SNPs filtered only for missing data, and (d) denoised and filtered SNPs.

Figure 4. Population and origin discrimination by STRUCTURE method showing the optimal K value (K = 3) according to Evanno’s method: (a) Mean value of K and (b) ΔK. (c) Bar chart showing individual Bayesian distribution probability of microsatellites for P. malariae from different countries in a mixed pattern.
Infections with mixed genotypes (FWS value < s 0.95) were identified in the SNP data from all populations, although most isolates had high FWS values from the denoised data, suggesting a single dominant genotype (Figure 5a,b). Isolates from Nigeria and Ghana had the lowest mean FWS values, thus indicating greater complexity.

Figure 5. Infection complexity in different countries represented by Wright’s inbreeding coefficient (Fws), using (a) unfiltered SNP data and (b) filtered and denoised SNP data.
SNPs at homologous loci associated with drug resistance
Sequence alignment of homologous resistance genes using BOWTIE identified 5 more mutations from unmerged reads than merged reads using both variant calling tools, while merged reads identified significantly more variants with BWA sequence alignment (Supplementary Tables 2 and 3). Combining all sequence alignment and variant calling algorithms, 20 consistent SNPs from the two methods (although present at low frequency) were retained. Of these, 7 SNPs were retained from the denoised and filtered dataset (Table 6). Eight of these SNPs encode nonsynonymous variants, resulting in amino acid changes; notably, 3 of the 4 consistent SNPs reported in the Pmaat1 gene were nonsynonymous variants.
Non‐synonymous SNPs observed in the sulfadoxine resistance gene of P. malariae (Pmdhps) at position S452C located near the S436A/F mutation in P. falciparum.Similarly, non‐synonymous SNPs were detected in multidrug resistance genes. (Pmmdr1) at position V100L is also located near the N86Y/F mutation in P. falciparum.
Table 6. Single nucleotide mutations detected in homologous drug resistance genes P. malariae.
| Numerical order | Gene | Location | Chain Reference | Replacement String | Mutation type | Change acid codon | Amino acid changes | Impact | Frequency (%) | Retained after noise filtering |
| 1 | Pmcrt | 348 | C | G | Transposition | ACC/ACG | Thr/Thr | Synonymous | 2,53 | No |
| 2 | Pmmdr1 | 1389 | C | T | Transform | AGC/AGT | Ser/Ser | Synonymous | 10,13 | Have |
| 3 | Pmmdr1 | 1743 | A | T | Transposition | GGA/GGT | Gly/Gly | Synonymous | 1,27 | No |
| 4 | Pmmdr1 | 1845 | G | A | Transform | TTG/TTA | Leu/Leu | Synonymous | 6,33 | Have |
| 5 | Pmmdr1 | 298 | G | T | Transposition | GTA/TTA | Val/Leu | Not synonymous | 1,27 | No |
| 6 | Pmaat1 | 1201 | G | C | Transposition | GTT/CTT | Val/Leu | Not synonymous | 2,53 | No |
| 7 | Pmaat1 | 183 | C | T | Transform | AGC/AGT | Ser/Ser | Synonymous | 3,79 | No |
| 8 | Pmaat1 | 434 | A | G | Transform | AAT/AGT | Asn/Ser | Not synonymous | 1,27 | No |
| 9 | Pmaat1 | 451 | T | A | Transposition | TTG/ATG | Leu/Met | Not synonymous | 1,27 | No |
| 10 | Pmatp4 | 3129 | A | G | Transform | TTA/TTG | Leu/Leu | Synonymous | 2,53 | No |
| 11 | Pmatp4 | 603 | G | C | Transposition | GGG/GGC | Gly/Gly | Synonymous | 2,53 | No |
| 12 | Pmnhe | 591 | A | G | Transform | TCA/TCG | Ser/Ser | Synonymous | 1,27 | No |
| 13 | Pmdhps | 1879 | A | T | Transposition | AGC/TGC | Ser/Cys | Not synonymous | 1,27 | No |
| 14 | Pmcytb | 412 | C | T | Transform | CTA/TTA | Leu/Leu | Synonymous | 1,27 | No |
| 15 | Pmcytb | 479 | G | A | Transform | GGT/GAT | Gly/Asp | Not synonymous | 1,27 | No |
| 16 | Pmcytb | 528 | T | C | Transform | TAT/TAC | Tyr/Tyr | Synonymous | 3,79 | Have |
| 17 | Pmcytb | 668 | A | G | Transform | AAT/AGT | Asn/Ser | Not synonymous | 1,27 | Have |
| 18 | Pmcytb | 690 | T | C,A | Transform / Transpose | TTT/TTC,TTA | Phe/Phe,Leu | Synonyms, Non‐synonyms | 7,59 | Have |
| 19 | Pmcytb | 708 | A | T | Transposition | GCA/GCT | Ala/Ala | Synonymous | 3,79 | Have |
| 20 | Pmcytb | 819 | T | A | Transposition | ATT/ATA | Ile/Ile | Synonymous | 2,53 | Have |
The linkage disequilibrium (LD) heatmap showed significant LD patterns among SNPs on homologous drug resistance genes. High r² values were observed between most SNPs on the mitochondrial gene Pmcytb and between Pmcytb SNPs and SNPs on other target genes such as Pmdhps, Pmdhfr, and Pmmdr1 from the unfiltered dataset (Figure 6a). The high linkage disequilibrium observed in mitochondrial SNPs was maintained in both the filtered and denoised datasets (Figure 6b,c).

Figure 6. Linkage disequilibrium between SNPs of homologous drug resistance genes calculated from (a) unfiltered data, (b) missingness‐filtered data, and (c) filtered and denoised data.
DISCUSSION AND CONCLUSION
Malaria elimination programs and tools are currently focused on eliminating the two major and dominant malaria parasite species, P. falciparum and P. vivax. However, other species such as P. malariae are also circulating in malaria‐endemic areas and require attention to achieve elimination goals. As knowledge of the diversity and impact of elimination tools on these small parasites has received little attention from the academic or public health community, this study utilized the Pathogen Diversity Network Africa (PDNA) to collect a small but extensive sample set of P. malariae samples from seven African and one Asian countries.
The study described the population structure, high genetic diversity, and genotypic richness and uniformity of this malaria parasite. High levels of genetic diversity are essential for long‐term survival of a population, and the level of variation determines the ability of a species to adapt to environmental challenges posed by nature or control interventions. The high levels of diversity in P. malariae found in this study are similar to those previously reported in Kenya and Malawi, although the number of samples analyzed from some countries was small and uneven. The different malaria transmission intensities and intervention histories among the countries included in the study may have influenced the results obtained.
High transmission intensity leads to frequent heterozygous recombination of malaria parasites in the mosquito vector, disrupting linkage disequilibrium between variable loci and increasing genetic diversity within populations. Overall, the probability of random mating within malaria parasite populations varied along the transmission intensity from high to low, in the direction from West to Central and East Africa, or sub‐Saharan Africa. However, heterozygosity remained high and genetic distances between samples were relatively low, regardless of differences in malaria transmission intensity, with the exception of three samples from Cameroon. The low genetic distances observed in the study need to be confirmed further as there are no studies on P. malariae for direct comparison.
Nevertheless, the results are consistent with a recent study of P. falciparum in Nigeria, although in marked contrast to an older study in Senegal, which found relatively high levels of free mating in the presently sampled populations, despite different transmission patterns. Indeed, multilocus genotypes were rare across all populations, an indicator of recombination and the absence of clonal expansion, which is often observed in some P. falciparum populations with low or seasonal transmission. Recurrent gene flow between malaria parasite populations across countries, through human or vector migration, may have led to a lack of geographical differentiation. This is reflected in low microsatellite differential indices between countries, although the small number of isolates within each national population limits the precision of the inferred indices. Furthermore, it is also possible that gene flow is not the only factor explaining the high genetic diversity or lack of geographic differentiation observed in P. malariae; other factors, such as the absence of bottleneck events or interventions to reduce local diversity, should also be considered.
Population structure analysis using neutral microsatellite loci (i.e. SSRs) identified five clusters, each comprising isolates from different countries. This further demonstrates the high intrapopulation diversity at these markers and the absence of population‐specific selection at the loci, which could have driven population differentiation. This structure is inconsistent with the model of isolation by distance observed in P. falciparum, in which genetic clusters can be distributed among geographical populations in West, Central and East Africa. Since P. malariae occurs mainly in co‐infections with P. falciparum, the factors driving such independent structure may be different, or it is possible that this structure was established by an earlier event that preceded any demographic or isolation‐driven population drift.
Using the mixture model in STRUCTURE software, three optimal ancestral clusters were identified and the results also showed that all isolates had members from each ancestral cluster, regardless of national origin. STRUCTURE implements a Bayesian algorithm to identify groups of individuals subject to Hardy‐Weinberg equilibrium and linkage equilibrium. However, the robustness of this method has been shown to be affected by small sample sizes and unevenness between subpopulations and/or hierarchical levels in the population structure.
SNP data from P. malariae homologs and P. falciparum resistance genes clustered the isolates into five subgroups, although only three less distinct clusters remained after rigorous screening procedures, and membership of these groups did not completely match those identified by microsatellites. The distribution of distances between SNPs and microsatellites was different, with larger distances between fewer isolates when using SSR data. This was not unexpected, as SSRs are multi‐allelic, potentially neutral, and often more divergent between pairs of isolates. Therefore, further investigation of the factors that may drive population divergence in this malaria parasite would improve our understanding of its complexity, particularly with regard to control and malaria eradication strategies. Although it appears that the majority of infections had mixed genomes (polygenomic‐polygenomic), as indicated by Fws, this result was affected by the cleaning and filtering procedures used in the analysis. Therefore, further research with adequate sample sizes is needed to clarify whether co‐transmission of different lineages and increased recombination capacity exist for this Plasmodium species. Of note is the high complexity of infections, as this is one of the indicators for monitoring the effectiveness of interventions. Unlike P. falciparum, this complexity does not appear to be higher in areas with relatively high malaria transmission rates and this may be part of the unique biology of this species, which requires further investigation. As drugs and other interventions reduce populations, selection and changes in complexity should be monitored for this species.
Antimalarial drugs have served as first‐line drugs for P. falciparum, with resistance associated with mutations in multiple genes and evidence of positive selection across the genome. Scientists identified 20 mutations in P. malariae at homologous genes responsible for drug resistance, using a combination of different sequence alignment and variant detection algorithms. These potential variants have not been described in previous studies of targeted or whole‐genome gene scans, possibly due to differences in the isolates used or the methods applied. In this study, only high‐quality variants were retained, supported by a combination of two mapping algorithms and two SNP detection algorithms.
The majority of potential variants were synonymous mutations, but a number of nonsynonymous mutations were found in seven genes, notably in the Pmcytb gene, which is homologous to the P. falciparum gene responsible for atovaquone resistance. Atovaquone belongs to the quinoline class, and resistance in P. falciparum has been associated with mutations in the multidrug resistance gene (Pfmdr1), the chloroquine resistance transporter (Pfcrt), and the amino acid transporter (Pfaat1). Most of the genetic linkage (LD) observed in the unfiltered dataset was not replicated in the processed and filtered dataset, with the exception of LD in mitochondrial gene SNPs. This could be the result of common ancestry or selection of dominant genotypes by the drug or other factors.
In addition, the scientists also identified potential variants in Pmdhfr and Pmdhps, which are homologous to folate resistance genes in P. falciparum populations. Antifolate drugs, such as sulfadoxine‐pyrimethamine, are still widely used as malaria chemoprophylaxis in pregnant women and in combination with amodiaquine for seasonal malaria prophylaxis in West Africa. These selective pressures may be driving the emergence of the identified variants. Although the nonsynonymous SNPs reported in this study occur at low frequencies, identification, characterization, and association of these SNPs with resistance traits will require extensive genome‐wide surveillance and phenotype‐association studies, based on in vivo and ex vivo treatment efficacy testing.
A limitation of the study, which calls for caution in interpreting the results, is the small number of samples analyzed from different countries, together with the lack of biological and clinical data on the samples. Larger population studies of P. malariae with appropriate epidemiological or clinical data are needed to validate the findings from the small‐scale studies reported here. Another limitation is the use of standard bioinformatic analysis protocols designed for P. falciparum but not yet available for P. malariae. Although these protocols are acceptable for preliminary analysis, custom protocols, considering the possibility of amplification and sequencing errors, may be more appropriate for P. malariae, especially given the scarcity of population data with high‐quality confirmed variants from this parasite. The different types of samples analyzed may also be a limitation of the study. Dried blood samples often yield poor quality results, likely due to the low prevalence and low parasite density of non‐P. falciparum species. Therefore, venous blood sampling and the application of more robust molecular techniques that take these factors into account will benefit future molecular surveillance studies of P. malariae.
Currently, advancing malaria elimination requires innovative strategies that target all species of malaria parasites. One approach could be to integrate genomic surveillance of all Plasmodium spp. into malaria control and elimination programs in sub‐Saharan Africa, drawing on lessons learned from the COVID‐19 response to refine approaches as new variants are identified and monitored. This study demonstrates the importance of this approach for P. malariae.
CN. Nguyen Thai Hoang & TS.BS. Huynh Hong Quang
IMPE‐QN





