Plasmodium knowlesi: Reservoir Hosts and Tracking the Emergence in Humans and Macaques

1Until recently, it was believed that malaria in humans was caused by only four species of parasite (Plasmodium falciparum, P. vivax, P. malariae and P. ovale). However, this perception changed when we discovered a large focus of human infections with P. knowlesi in the Kapit Division of Sarawak, Malaysian Borneo. These infections had predominantly been mistakenly identified as P. malariae by microscopy, since both species have similar morphological characteristics,. With subsequent reports of human infections in other parts of Malaysia,, and in Thailand,, Myanmar, Singapore,, the Philippines, Vietnam and Indonesia,, P. knowlesi is now recognized as the fifth species of Plasmodium responsible for human malaria. It causes a wide spectrum of disease and can lead to high parasite counts, severe complications and death,. In a recent study, we found that approximately 1 in 10 knowlesi malaria patients at Kapit Hospital developed potentially fatal complications, comparable to P. falciparum malaria, which is considered to be the most virulent type of malaria in humans.

2P. knowlesi is primarily a simian malaria parasite and was first isolated from a long- tailed macaque (Macaca fascicularis) imported to India from Singapore in 1931. Subsequently, P. knowlesi has been detected in wild long-tailed macaques of Peninsular Malaysia, and the Philippines, in pig-tailed macaques (M. nemestrina) of Peninsular Malaysia and in banded leaf monkeys (Presbytis melalophus) in Peninsular Malaysia. There has been no documented evidence of P. knowlesi or any other malaria parasites in monkeys in Malaysian Borneo, and although a monkey source for the hundreds of human P. knowlesi infections that have been described in the Kapit Division of Sarawak,, appeared likely, it remained to be proven.

3Prior to our report in 2004 of the large focus of human infections in Sarawak, Malaysian Borneo, when we utilized molecular methods for characterisation and PCR assays for detection of P. knowlesi, there had been only one confirmed case of a naturally-acquired P. knowlesi infection in a human. That person got infected with P. knowlesi while spending a few weeks in the forest of Pahang, Peninsular Malaysia in 1965. It is not known whether the large focus in Malaysian Borneo and subsequent recent reports of human knowlesi malaria in Southeast Asia represent a truly newly emergent malaria parasite in humans or whether human infections have been occurring for a relatively long period, but have gone undetected due to unavailability of molecular detection methods to distinguish between P. knowlesi and P. malariae. In order to identify the reservoir hosts and extend our understanding of the epidemiology and evolutionary history of P. knowlesi, we examined blood samples from wild macaques for malaria parasites, and analyzed the circumsporozoite protein (csp) gene and the mitochondrial (mt) genome of P. knowlesi isolates derived from humans and macaques.

4Nested PCR examination of blood samples from 108 wild macaques (82 long-tailed, 26 pig-tailed), sampled from 17 different locations in the Kapit Division of Sarawak, showed that 101 (94%) of the macaques were infected with malaria parasites. Long-tailed macaques had a higher prevalence of infection (98%) than pig-tailed macaques (81%) (Fisher's Exact P = 0.009) (Table 1). By nested PCR assays, we detected 5 species of Plasmodium, with P. inui being the most common (prevalence of 82%), followed by P. knowlesi (78%), P. coatneyi (66%), P. cynomolgi (56%), and P. fieldi (4%). Multiple species infections were very common, with 91 of the 108 (84%) macaques being infected by two or more species of Plasmodium each. There was a higher prevalence of P. knowlesi among long-tailed macaques (87%) than pig-tailed macaques (50%) (P = 0.006).

5To compare the molecular identity of the parasites in macaques and humans, we first sequenced the P. knowlesi csp gene in blood samples from 31 patients admitted to Kapit Hospital and 16 wild macaques. Most macaques (10 of 16), but only a minority of humans (3 of 31) contained 2 or more csp alleles (Fig. 1A). Overall, we derived 48 csp allele sequences of P. knowlesi from the macaques and 34 from the human samples, with 61 different alleles observed in total. Three of these csp alleles were shared between human and macaques, three were shared by macaques, and the remaining alleles were detected in only macaques or humans (Fig. S1 and Fig. 1B). We found that the central region of the P. knowlesi csp was composed of highly polymorphic repeat sequences (Table S1). Analysis of the aligned non-repeat regions of csp showed 19 polymorphic sites, of which 14 were shared polymorphisms in samples from both host populations (Fig. S2). The nucleotide diversity of csp was similar in both hosts (π = 2.2×10−2 in humans and 2.4×10−2 in macaques), although the haplotype diversity was marginally higher in macaques (H = 0.82, SD = 0.03) than in humans (H = 0.73, SD = 0.06). There was no clustering of csp allele sequence type associated with either host (Fig. 1B).

6We also sequenced the ∼6-kilobase mtDNA genome of P. knowlesi parasites isolated from 25 malaria patients and 11 macaques. Each human sample had a single mtDNA haplotype, while all except one macaque sample contained multiple (2 to 6) haplotypes (Fig. 2A). In total, we generated 54 complete mtDNA genome sequences, representing 37 different mtDNA haplotypes, with a higher number of haplotypes in the macaques (23 haplotypes from 11 samples) than in the humans (17 haplotypes from 25 samples). Six of the haplotypes were found in more than 1 sample, and 3 of these were shared between the human and macaque hosts (Fig. 2B). Forty-five single nucleotide polymorphisms (SNPs) and a 4-base insertion/deletion within the P. knowlesi mtDNA genome were identified (Fig. S3), and the level of nucleotide diversity (π) of mtDNA was estimated as 7.5±0.7×10−4.

7The Bayesian coalescent approach was used to estimate the time to the most recent common ancestor (TMRCA) for P. knowlesi. A nucleotide substitution rate for the mtDNA genome of 3.13×10−9 (95% HPD, 1.94–4.45×10−9) substitutions per site per year was estimated by comparing mtDNA sequences of P. knowlesi, P. fragile, P. cynomolgi, P. simiovale (parasites of Asian macaques) with P. gonderi (a parasite of African mangabeys) and Plasmodium sp. (Mandrill), assuming parasite lineages separated when Asian Old World monkeys and African Old World monkeys diverged 10 million years ago (MYA). This derived rate yielded an estimate of 257,000 (95% HPD: 98,000–478,000) years before present as the TMRCA of P. knowlesi (Fig. 3).

8There was no evidence of recombination in the mtDNA of P. knowlesi (Table S2), and reconstruction of the haplotype genealogical network demonstrated that no distinct lineages of mtDNA were exclusively associated with either human or macaque hosts (Fig. 2B). Our phylogeny-trait association tests based on association index (AI) and parsimony score (PS) of host-parasite phylogenetic substructure do not reject the null hypothesis of no association between parasite and host (Table S3). Hence, this further indicates the absence of a distinct lineage of P. knowlesi parasites associated with either human or macaque hosts.

9We observed an excess of unique mtDNA haplotypes, which appear at the edges of a star-like structure of the haplotype genealogical network (Fig. 2B), and this is indicative of an evolutionarily recent population expansion of P. knowlesi. The signature of population expansion is also evident from the unimodal shape of the pairwise mismatch distribution (Fig. 4A), and this is supported by a low Harpending's raggedness index (r = 0.009, P = 0.87). In addition we used the Tajima's D, Fu and Li's D and F (with out-group) and Fay and Wu's H statistics to detect deviation from the model of neutral evolution considering that the deviation from neutrality could be due to demographic processes such as population expansion, population bottleneck or mutation rate heterogeneity. We obtained significant negative values for all these statistics (Tajima's D = −1.88, P = 0.001; Fu and Li's D = −2.60, P = 0.02; Fu and Li's F = −2.80, P = 0.009; Fay and Wu's H = −10.33, P = 0.044), thereby providing further evidence for an expansion of the P. knowlesi parasite population.

10To further investigate the demographic history of P. knowlesi, we estimated the changes in effective population size of the parasite through time, using a coalescent approach called the Bayesian skyline plot. The plot indicates that P. knowlesi underwent a rapid population growth between approximately 30,000 and 40,000 years before present (Fig. 4B). In addition, we performed independent analyses for P. knowlesi mtDNA sequences derived from humans and macaques. Similar trends were reflected in the Bayesian skyline plots for each host (Fig. S4), showing that there are no differences between the demographic history of P. knowlesi for either host.

11We analysed mtDNA cytochrome b sequence data for long-tailed macaques in Southeast Asia using a similar approach, but did not find any evidence for changes in population size between 100,000 and 10,000 years before present (Fig. S5).

12Our study shows that wild macaques in the Kapit Division of Sarawak, Malaysian Borneo are infected with the same 5 species of Plasmodium found in macaques of Peninsular Malaysia,, and that these macaques have a very high prevalence of P. knowlesi and P. inui. In previous studies, we found that P. knowlesi is the most common cause of hospital admission for malaria in the Kapit Division and there are approximately 90 knowlesi malaria admissions, predominantly adults, at Kapit Hospital per year,,. The actual annual incidence of knowlesi malaria for the Kapit Division is probably higher, because not all persons with P. knowlesi infections may have sought treatment in hospital and there may be asymptomatic infections and misdiagnoses. Nevertheless, the restricted number of knowlesi malaria cases in the human population of 109,000, contrasts with the extremely high prevalence of P. knowlesi we detected in the wild macaques of the Kapit Division. These findings contrast with the absence of P. knowlesi infections in a survey of 99 long-tailed macaques in one region in Thailand. In that study, the majority of macaques were trapped near a temple in a region where very few human knowlesi malaria cases have been reported, and the absence of detectable P. knowlesi there could be due to the low abundance of mosquitoes of the Anopheles leucosphyrus group, which have been shown to be the most competent vectors of knowlesi malaria. We previously identified one member of this group, Anopheles latens, as the vector for P. knowlesi in the Kapit Division. This mosquito feeds outdoors after dusk and is attracted to humans and macaques at ground level, but prefers to feed on macaques at a higher elevation. Our findings here, of the higher number of P. knowlesi csp alleles and mtDNA genome haplotypes detected per infection in macaques compared with humans, and the very high prevalence of P. knowlesi in macaques, suggest that presently there is a greater intensity of transmission of P. knowlesi by the vectors among wild macaques, than from macaques to humans. These results, including our observation that certain alleles of the P. knowlesi csp gene and mtDNA genome haplotypes are shared between macaque and human hosts, taken together with previous epidemiological,, and entomological data,, strongly indicate that knowlesi malaria is a zoonosis in the Kapit Division and that wild macaques are the reservoir hosts.

13Our estimated TMRCA for P. knowlesi (98,000–478,000 years ago) indicates that P. knowlesi is derived from an ancestral parasite population that predates human settlement in Southeast Asia,. Therefore macaques, which colonized Asia more than 5 million years ago, were the most likely hosts during the initial emergence of P. knowlesi in this region. Our estimate also indicates that P. knowlesi is as old as, or older than the 2 most common human malaria parasites, P. falciparum and P. vivax, for which the TMRCA has been estimated to be 50,000–330,000, years and 53,000–265,000 years,, respectively.

14Our analyses of the mtDNA data indicate that that P. knowlesi underwent a period of population expansion, estimated at 30,000–40,000 years ago, which coincides with a time when Borneo was part of mainland Southeast Asia and the possibility of increased parasite admixture between macaque troops. This period is concordant with a time of exceptional human population growth in Southeast Asia, based on mtDNA sequence analysis. We did not detect a similar population expansion of macaques, but this analysis was based on the cytochrome b gene alone. It would be preferable to analyze mtDNA sequences of macaques sampled in Borneo to determine whether they underwent a parallel historical population expansion. It is possible that the population expansion of P. knowlesi was not directly linked to expansion in any primate host, but was rather due to the expansion or adaptation of the mosquito vectors.

15In conclusion, our results indicate that P. knowlesi in Sarawak is zoonotic, with humans sharing parasites with the original and preferred hosts, the macaques, most likely since they first came into close contact in the forests of Southeast Asia. A multi-gene family (KIR) in P. knowlesi encodes proteins with sequence motifs mimicking host cell receptor CD99 in macaques, and the observation that the KIR motifs are less perfectly matched to the human CD99 sequence also supports the hypothesis that the parasite is particularly adapted to macaque hosts. Humans acquire knowlesi malaria on occasions when they enter the habitats shared by macaques and mosquitoes of the Anopheles leucosphyrus group,,,, which are forest-dwelling mosquitoes that feed outdoors after dusk,. There is no evidence yet to suggest a host-switch by P. knowlesi, unlike other human malaria parasites such as P. vivax and P. falciparum that might have been part of ancient zoonoses,, but have since adapted to humans. However, it is possible that the current destruction of the natural forest ecosystem, with associated increase of the human population, may alter the parasite, macaque host and mosquito population dynamics and lead to an adaptive host-switch of P. knowlesi to humans.

16Currently, Malaysia has no legislation governing the use of animals in research. Nevertheless, this study was carried out in strict accordance with the recommendations by the Sarawak Forestry Department for the capture, use and release of wild macaques. A veterinarian took blood samples from macaques following anesthesia by intramuscular injection of tiletamine and zolazepam. All efforts were made to minimize suffering by collecting blood from macaques at the trap sites and releasing the animals immediately after the blood samples had been obtained. The Sarawak Forestry Department approved the study protocol for capture, collection of blood samples and release of wild macaques (Permits Numbers: NPW.907.4.2-32, NPW.907.4.2-97, NPW.907.4.2-98, 57/2006 and 70/2007). A permit to access and collect macaque blood samples for the purpose of research was also obtained from the Sarawak Biodiversity Centre (Permit Number: SBC-RP-0081-BS). Human blood samples were taken after written informed consent had been obtained from patients admitted to Kapit Hospital. This study was approved by the Medical Research and Ethics Committee of the Malaysian Ministry of Health (Reference number: KKM/JEPP/02 Jld.2 [133]), which operates in accordance to the International Conference of Harmonization Good Clinical Practice Guidelines.

17A total of 108 macaques were sampled from 2004 to 2008. Ninety were from 5 major sites and the remainder from 12 different locations in the Kapit Division of Sarawak. All locations were within 2 km from longhouse communities where human knowlesi cases had previously been reported. After blood was obtained from anaesthetised animals, they were tagged with a microchip (to prevent re-sampling) and released. Human blood samples were obtained from 31 patients with knowlesi malaria admitted to Kapit Hospital between 2000 and 2006.

18DNA was extracted from macaque and human blood samples as described previously. DNA samples from macaques were examined using nested PCR assays with genus and species-specific primers based on the small subunit ribosomal RNA genes. PCR primer sequences (for P. knowlesi, P. coatneyi, P. cynomolgi, P. fieldi and P. inui) and annealing temperatures are provided in Table S4.

19We amplified and sequenced the complete P. knowlesi csp gene and performed sequence analysis as described previously. At least two clones from each of the two PCR amplifications per sample were sequenced and for macaque samples, at times 5 PCR amplifications and cloning procedures were necessary before P. knowlesi sequences could be obtained. For amplification of the P. knowlesi mtDNA, two back-to-back primers were designed (Pkmt-F1, 5′- GGACTTCCTGACGTTTAATAACGAT-3′ and Pkmt-R1, 5′-TGGACGTTGAATCCAATAGCGTA-3′) by using previously described mtDNA sequence of P. knowlesi. PCR amplification was performed separately for each sample to prevent cross-contamination of DNA using the Elongase Amplification System (Invitrogen). The PCR product for each isolate was gel purified, cloned into pCR-XL-TOPO vector (Invitrogen) and sequenced using BigDye Terminator Cycle Sequencing kit (Applied Biosystems) with 28 internal primers (sequences in Table S5) that enabled sequencing of both DNA strands.

20The diversity of P. knowlesi mtDNA from 6 human samples was initially characterized. These were selected based on the criteria that each patient originated from a different geographical area of Kapit Division and patients had no records of recent travel history. Only females were chosen assuming that females travel less than the males. Sequencing data was obtained from at least 2 clones originating from separate PCR amplifications. Both DNA strands were sequenced from each clone and any nucleotide conflicts found were resolved following a third PCR amplification, cloning and sequencing.

21The remaining 19 human samples were randomly chosen from patients admitted between 2000 and 2006. Following PCR amplification and cloning, these samples were haplotyped by sequencing single DNA strand of the mt genome and at least 2 plasmid clones were sequenced for each sample. Any single nucleotide polymorphisms (SNPs) or singleton polymorphisms detected were verified by sequencing the polymorphic regions in at least 2 plasmid clones originating from separate PCR amplifications, and both DNA strands were sequenced.

22The mt genome was selected to examine the evolutionary history of P. knowlesi, just as the mt genomes of P. vivax and P. falciparum were previously found suitable; it does not undergo recombination so intraspecific phylogenetic analysis can be performed and it shows no evidence of non-neutral polymorphism.

23DNA sequence data were aligned using the Lasergene package (DNASTAR). Measures of genetic diversity were conducted using DnaSP v5.10.00 software. A minimum spanning network connecting the mtDNA haplotypes of P. knowlesi based on statistical parsimony method was constructed using the TCS 1.21 software.

24Host-parasite association was assessed based on the association index (AI) and parsimony score (PS) statistics, which account for phylogenetic uncertainty in analysis of phylogeny-trait correlations. The values of AI and PS statistics were calculated based on the posterior samples of trees produced by BEAST using the BaTS program. The null distribution for each statistic was estimated with 1,000 replicates of state randomization.

25The demographic expansion of P. knowlesi was examined based on pairwise mismatch distribution using Arlequin v3.1 software. Observed mismatch distribution was compared with that estimated under the sudden demographic expansion model using a generalized least-square approach. The deviations from the population expansion model were tested using the Harpending's raggedness index with a parametric bootstrap of 1000 replicates. Tajima's D, Fu and Li's D, Fu and Li's F and Fay and Wu's H statistics were performed using the software DnaSP v5.10.00. These statistics were calculated using the mitochondrial genome of P. coatneyi (AB354575) as out-group.

26The evolutionary rate, time to the most recent common ancestor (TMRCA) and the past population dynamics of P. knowlesi were inferred using the Bayesian Markov Chain Monte Carlo (MCMC) method implemented in the BEAST package v1.5.4. The mean substitution rate of mtDNA and TMRCA of P. knowlesi were estimated based on a time-calibrated Bayesian phylogenetic analysis of non-human primate malarias (P. gonderi, Plasmodium sp. (Mandrill), P. simiovale, P. fragile, P. cynomolgi and P. knowlesi) and human malarias (P. falciparum, P. vivax, P. malariae and P. ovale) (Table S6) (Figure 3), assuming co-divergence of the parasites with their host lineages, Asian Old World monkeys - African Old World monkeys at 10 MYA. The accession numbers of sequences derived from GenBank database are as follows; P. falciparum (M99416), P. malariae (AB354570), P. vivax (NC007243), P. ovale (AB354571), P. gonderi (AB434918), Plasmodium sp. (mandrill) (AY800112), P. simiovale (AB434920), P. inui (AB354572), P. hylobati (AB354573), P. cynomolgi (AB434919), P. simium (AY800110), P. fragile (AY722799) and P. coatneyi (AB354575). A General Time Reversible (GTR) substitution model with gamma distribution of rate variation among sites and a proportion of invariable sites as determined using Modeltest v3.7, an uncorrelated log-normal relaxed molecular clock model and a Bayesian skyline coalescent model (10 coalescent-interval groups) were used for this analysis. One hundred million generations of the MCMC chains were run with sampling every 10,000 generations and the first 10 million generations were discarded as burn-in. The BEAST output was analyzed using the Tracer v1.5 program (available at http://tree.bio.ed.ac.uk/software/tracer/) and uncertainty in parameter estimates was expressed as values of the 95% highest probability density (HPD). The trees produced by BEAST were annotated using TreeAnnotator, and maximum clade credibility tree was visualized using the FigTree v1.3.1 program (available at http://tree.bio.ed.ac.uk/software/figtree/).

27Past population dynamics of P. knowlesi parasites in terms of the change in effective population size (Ne) through time were independently analyzed using the P. knowlesi mtDNA datasets for humans and macaque, and also by combining both human and macaque P. knowlesi mtDNA datasets. Using the estimated mean substitution rate and Bayesian skyline coalescent model, the MCMC chains were run for 100 million generations with sampling every 10,000 generations and the first 10 million generations were discarded as burn-in.

28The changes in effective population size (Ne) through time were also drawn for M. fascicularis and M. nemestrina based on the cytochrome b sequences obtained from GenBank (Table S7). A BEAST analysis to determine the mean rate substitution of the cytochrome b (cytb) gene of macaques was performed using cytb sequences of M. fascicularis, M. nemestrina and Papio anubis (GenBank accession EU885461), and assuming baboons and macaques diverged 6.6 MYA. An estimated mean substitution rate of 4.56×10−8 substitutions per site per year was used to infer the Bayesian skyline plot for M. fascicularis and M. nemestrina. For each species, 100 millions generations were performed, with sampling every 10,000 generations and 10 percent of the sampling were discard as burn-in.

29For all analyses implemented in BEAST, at least 2 independent runs were performed and convergence of all parameters was determined based on Effective Sample Size (ESS) values of >200.

30The sequences generated during this study have been deposited in GenBank: P. knowlesi mitochondrial genome sequence data under the accession numbers EU880446–EU880499 and P. knowlesi csp gene sequences under the accession numbers AY327558–AY327572, DQ350272–DQ350306, DQ641526–DQ641528 and GU002471–GU002533.