Genetic characterization of nodular worm infections in Asian Apes

Parasitic nematodes of Oesophagostomum spp., commonly known, as 'nodular worms' are emerging as the most widely distributed and prevalent zoonotic nematodes. Oesophagostomum infections are well documented in African non-human primates; however, the taxonomy, distribution and transmission of Oesophagostomum in Asian non-human primates are not adequately studied. To better understand which Oesophagostomum species infect Asian non-human primates and determine their phylogeny we analysed 55 faecal samples from 50 orangutan and 5 gibbon individuals from Borneo and Sumatra. Both microscopy and molecular results revealed that semi-wild animals had higher Oesophagostomum infection prevalence than free ranging animals. Based on sequence genotyping analysis targeting the Internal transcribed spacer 2 of rDNA, we report for the first time the presence of O. aculeatum in Sumatran apes. Population genetic analysis shows that there is significant genetic differentiation between Bornean and Sumatran O. aculeatum populations. Our results clearly reveal that O. aculeatum in free-ranging animals have a higher genetic variation than those in semi-wild animals, demonstrating that O. aculeatum is circulating naturally in wildlife and zoonotic transmission is possible. Further studies should be conducted to better understand the epidemiology and dynamics of Oesophagostomum transmission between humans, non-human primates and other wild species and livestock in Southeast Asia.

Oesophagostomum species are stout, white roundworms and have a direct life cycle. Eggs passed in the faeces, hatch and rapidly develop into L1 rhabditiform larvae. After 24 h of hatching, L2 become infective L3 within 3-4 days. L3 keep the protective cuticle of L2 and are capable of surviving long periods of adverse environmental conditions, e.g., hot-dry conditions of the dry season. Infection occurs by ingestion of filariform L3 larvae via vegetation eaten by the host. After ingestion, L3 pass to the cecum, where they exsheath within ~ 3 days of ingestion then invade the tunica mucosa, stimulating the formation of separate cysts around individual larvae in the gut wall. The larvae develop there to the L4 stage then once in the lumen, the larvae moult and reach the adult stage 26 .
For most nematodes parasites, species identification based on egg morphology is difficult, therefore recently molecular techniques are widely used to identify species. Accordingly, recent molecular genetic studies have demonstrated that the ITS nucleotide sequence of rDNA allows an unequivocal identification and discrimination of a range of strongylid nematode species, irrespective of the developmental stage of the parasites 27 .
Based on molecular findings, Oesophagostomum spp. are known to frequently infect domestic and wild pigs, ruminants and primates [28][29][30][31] . Eight species of Oesophagostomum have been observed to occur in African NHPs and three of them (O. bifurcum, O. stephanostomum and O. aculeatum) are also reported in humans [32][33][34] . Human cases have been attributed to an animal origin and NHPs have been proposed as a potential reservoir 33 . Infections in wild primates appear to be asymptomatic, although clinical signs and mortality due to Oesophagostomum have been recorded in captive settings 35 . Oesphagostomum spp. have been reported from orangutans based on morphological and coprological analysis 36 . Arizono et al. 37 molecularly identified O. aculeatum infections for the first time in Japanese macaques. Recently, molecular markers showed a widespread distribution of O. aculeatum in the Bornean primate community, including orangutans 38 , although there is a paucity of investigations in gibbons. Currently, there is no study done on O. aculeatum distribution in the Sumatran NHPs community. The aim of this study was therefore to genetically investigate Oesphagostomum parasitising in wild and semi-wild Asian apes (two orangutan species; P. pygmaeus and P. abelii and one gibbon species; Hylobates albibarbis) to confirm for the first time if Sumatran orangutans and wild Asian gibbons are being infected and to investigate the distribution and population structure of nodular worm infections in these species.
We obtained a 758-bp-length fragment of ITS1-5.8S-ITS2 ribosomal DNA (here after, rDNA) from a total of 27 different individuals, 3 from gibbons and 24 from orangutans. By comparison with published sequences using NCBI BLAST, our samples showed a 98% identity for O. stephanostomum (AB821014.1), 95% O. detantum (AJ619979.1) and 94% Oesophagostomum sp. MOS-2014 (AB908964.1). The resulting long DNA sequence overlapped with published sequences and contained no insertions or deletions, making alignment unequivocal. After analysis of the 27 Oesophagostomum sequences, we identified seven unique haplotypes. H1, H4, H2 and H3 were the most frequent haplotype detected in 48%, 14%, 11% and 11% of the infections, respectively. There is no shared haplotype between the two wild orangutan species P. pygmaeus; Bornean Sebangau and P. abelii; Sumatran Suaq (hereafter referred to as Seb and Suaq respectively). Only two haplotypes (H1 and H2) are shared between the Sumatran populations of Suaq and the semi-wild population in Bukit Lawang (hereafter referred to as the Buk). H3 is the only haplotype shared between wild P. pygmaeus and H. albibarbis in Sumatra and H5, H6 and H7 are unique haplotypes represented only in Seb (Fig. 2a). The nucleotide position of polymorphic sites are presented in Fig. 2b.
The best sequence evolution of Oesophagostomum sequences obtained in this study with published other Oesophagostomum species sequences including outgroup, were HKY86 + G6 by JModelTest 40 . Published reference sequences were included to identify putative species and only 258 bp sequence were available. Based on phylogenetic analyses of the Oesophagostomum ITS2 rDNA (258 bp) sequences, the phylogenetic tree resolved these sequences into 8 clades according to host species (Fig. 3 Tree topologies showed a strong pattern of geographical structure across the geographic range of Oesophagostomum spp. between Southeast Asian and African primates. Our results show that Southeast Asian NHPs are infected by O. aculeatum only (Fig. 3). Within the clade 6, we observed two groups. The first group contained two sequences from gibbons and one sequence from Sumatran orangutan and one sequence from a Bornean orangutan. The second group contained all of the sequences from semi-captive Sumatran orangutans. Six additional sequences; one sequence from a gibbon, two sequences from two Bornean orangutans, two sequences from two semi-wild Sumatran orangutans and one sequence from one free ranging Sumatran orangutan sorted into group two which were 99% similar to the same O. aculeatum reference sequences.
The total number of sequences obtained and analysed for each population is given in Table 2. Overall, the Bornean (Seb) population exhibited the greatest haplotype diversity (0.867) and nucleotide diversity (0.0089). However, the Sumatran populations showed almost the same haplotype diversity in both (0.6 for Suaq and 0.59 for Buk) but the Suaq population displayed four times higher nucleotide diversity than the Buk population ( Table 2).
At host intra species level, nucleotide diversity estimated from the rDNA marker was 3.6-fold higher in Suaq (Sumatran wild) then Buk (Sumatran semi-wild). At host inter species level, Seb (Bornean wild) samples had the highest nucleotide diversity; 7.5 and 2 times higher than semi-wild and wild populations from Sumatra, respectively.
Genetic differentiation between Sumatran and Bornean populations was very strong as shown by high pairwise Fst estimates at rDNA (Table 3). Regarding genetic differentiation between populations, the lowest Fst values was observed between Buk and Suaq populations (Fst = 0.06882, P > 0.05) while the highest significant genetic differentiation was observed between Buk and Seb (Fst = 0.49204, P < 0.05) ( Table 3). The analysis of the    Figure S1).

Discussion and conclusions
In this study, we extracted DNA directly from nematodes eggs, without coprocultures and genetic characterization of L3 larvae, to identify nodular worm infection in Asian NHPs. Based on the internal transcribed spacer 2 (ITS 2) region of ribosomal RNA gene (rDNA), we confirmed that O. aculeatum is infecting Asian apes in Borneo. We confirm for the first time that Sumatran NHPs are also infected with this parasitic nematode. We confirm by molecular techniques that gibbons can be infected and therefore should be included in future studies of nodular worms, such as investigating how their unique ecology effects their risk of infection. Mitochondrial and nuclear genes are used in phylogenetic studies on strongylid species 41 . Based on the mitochondrial cox1 gene, Frias et al. 38 showed that Bornean primates are subject to O. aculeatum infection. However, in our study we used ITS sequences to discriminate nematode species as the ITS regions are among the most variable nuclear loci, with a sufficient number of comparative sequences in databases. Our result demonstrated the same phylogenetic similarity observed from Frias et al. 38 . Accordingly, the rapid mode of ITS evolution including frequent indels makes alignment of ITS haplotypes more challenging compared to protein-coding markers, which can adversely affect the reliability of phylogenetic reconstruction 42 .
The evolution and dispersion of strongylid parasites and their relationships with other STHs species found in NHPs especially in hominids have been a point of attention 28,36,43 . Several studies have shown that Oesophagostomum infections of African primates are mainly O. stephanostomum, O. bifurcum and those belonging to the unclassified Oesophagostomum sp. and presented in three genetic lineages 27,28,44,45 . Later, Ota et al. 46 , showed the existence of a fourth clade, O. aculeatum clade from Asian primates, distinct from the three clades of African species. However, the molecular phylogenies have been limited by a deficiency of representative samples from Asian apes. In our study, based on ITS and 5.8S nucleotide sequence information, we identified Oesophagostomum sp. as one of the prevalent helminth parasites of Asian non-human primates. The Oesophagostomum sequences recovered from Asian apes clustered within the O. aculeatum clade. Interestingly, we did not detect any other nodular worm infections in our dataset, as it is known that the primer set used in our study are able to detect other Oesophagostomum species 46 , suggesting that O. aculeatum is the only or at least the predominant Oesophagostomum species that parasitizes orangutans and gibbons in our study area. Recently, Frias et al. 38 showed that O. aculeatum is widely distributed in Borneo and that it infected all studied members of host species, such as long-tailed macaques, silvered langurs, proboscis monkeys and Bornean orangutans, with little clear genetic differentiation among them.
Previous studies have indicated that captive and/or semi-wild orangutans are subject to a higher prevalence of infection than wild ones 12,47,48 . Studies of terrestrial NHPs found a high STH prevalence of 81.2% in captive orangutans with ground dwelling behaviour 8,47,49,50 . In addition to these studies, our results significantly demonstrate that more semi-wild orangutans were positive for O. aculeatum than wild orangutans. However, it should be noted that this may be because of a sampling bias or over-representation of semi-wild orangutan samples in the dataset. As the presence of nodular worms have been confirmed in Asian primates in this and further studies with larger sample sizes should be carried out to confirm these findings.
However, several behavioural and ecological factors may contribute to higher infection with nodular worm species among semi-captive primates, such as restricted ranges making it difficult to avoid contaminated areas, changes in social structure, increased population density, dietary changes and stress and more contact with other animal host species. Another significant factor could be decreased arboreality and increased terrestrial locomotion in semi-captive animals, which allows more contact with the soil and increases exposure to STHs infective larvae. Regarding wild NHPs in Asia, previous studies have reported that the percentage of wild orangutans infected by STHs were 33% for P. pygmaeus 8 and 47% for P. abelii 47 . In this study, we report that 30% of free ranging animals were infected with O. aculeatum alone. Apart from decreased arboreal life style and increased contact with soil, increasingly overlapping home ranges with other wildlife species and livestock may result in increased levels of infection and parasite occurrence 8,47,51 . Additionally, the habitat along the swamp forest seems to be an ideal environment for the development of soil-transmitted species as it provides the required humidity 7 .
In addition, Ghai et al. 28 found that, in the case of Oesophagostomum, smaller primate groups with large daily travel distances had a higher prevalence. The free ranging orangutans have large home ranges, for example wild orangutans home ranges are estimated to be 423 ha/month for males and 131 ha/month for females 52 . However, in our study the free ranging animals with larger habitat ranges were not observed to have higher prevalence's of infection than semi-wild animals. The highest prevalence of O. aculeatum in orangutans was found in Bukit Lawang (75%) despite it having the lowest orangutan density (1.8 individuals/km 2 ) 48 , suggesting orangutan to orangutan transmission may not be a common source of infection directly or indirectly. Although Frias et al. 38 suggested that orangutans may also acquire Oesophagostomum from other NHPs, in our study, this could be related to increased human activity in that area 48 and increased approaching of humans, their dwellings and domestic animals or by semi-wild orangutans who have been previously exposed to human contact during rehabilitation. Ghai et al. 28 observed that, some primates, particularly those traversing large distances in small groups, were most susceptible to nodule worm infection. If data were available to investigate the relationship between range size and prevalence in wild orangutans and gibbons, a similar pattern observed by Ghai et al. 28 would be interesting to test.
The observed low infection prevalence from free ranging NHPs may also point towards the hypothesis that great apes are able to control parasite infection and avoid the progression of the disease by self-medication, i.e. leaf swallowing 35,53,54 and wild animals have more access to these plants. Recently  www.nature.com/scientificreports/ swallowing in an Asian species, the white-handed gibbon, suggested a similar self-medicative function against nematode infection 55 .
In this study, the patterns of genetic variability observed among populations differ. A high level of genetic diversity was found across all three O. aculeatum populations. A significant difference in the number of polymorphic sites, nucleotide diversity and average number of nucleotide differences was evident between wild and semi-wild populations. Free ranging populations had higher genetic diversity than semi-wild ones. In wildlife, genetic diversity is probably greater than in anthropized environments, possibly due to a larger range of hosts 56,57 . Within free ranging populations, the Seb population had a higher genetic diversity then Suaq. This higher genetic diversity was due to observed unique haplotypes, geographical isolation, population size and different host species. Moreover, O. aculeatum, is a specialised organism that infects a specific host and the genetics and behaviour of a host can play an important role on the parasites genomic variation within the different host species (P. pygmaeus, P. abelii and H. albibarbis). Our assessment of population genetic differentiation calculating Fst index revealed non-significant differentiation among Sumatran populations and indicated low levels of genetic differentiation within the same host population. This pattern was also observed between the two free ranging populations, Seb and Suaq, from different host species. Interestingly our result showed significant genetic diversity between Buk and Seb populations. The isolation by distance analysis indicated that more genetically similar relationships were found for more distant populations then nearby populations (between Suaq and Seb then between two neighbour semi wild and wild Sumatran, Suaq and Buk population). Thus, O. aculeatum isolated from different geographical regions may differ because of a heterogeneous species, less population diffusion, complex genetic structure and high differentiation.
In our study, the molecular analysis was limited to a single gene. While several studies demonstrated a relatively high genetic diversity in the Oesophagostomum genus 27,28,46 , the sequencing of additional genes or other molecular markers should be applied to confirm species level differentiation and the population genetic structure of O. aculeatum in Southeast Asia.
The origin of O. aculeatum is unclear but it is clear that O. aculeatum is confined to Asia and infects more wild populations such as other primates and elephants 37,38,58 . Faecal examination of long-tailed macaques in Asia has uncovered Oesophagostomum spp. infections 59 38 .
Our results confirm that O. aculeatum should be considered a pathogen for non-human primates in Southeast Asia. Oesophagostomum in semi-wild animals were found at a higher prevalence but showed a lower genetic variability compared with the ones in fully wild animals, further demonstrating the influence of human activities on the parasite dynamics in wild primate species. Further studies should be conducted to better understand the epidemiology and dynamic of Oesophagostomum transmission between human, non-human primates, other wild species and livestock, in Southeast Asia, following a one-health approach. Of particular interest are changes in transmission and prevalence in the context of the increasing human perturbation of natural forested habitats in Southeast Asia.

Materials and methods
Study site, sample collection and origin of faecal samples. Faecal samples were collected from three sites in Indonesia, one site in Borneo and two sites in Sumatra between 2004 and 2011. In Borneo, 11 faecal samples from orangutans (P. pygmaeus) and 5 samples from gibbons (H. albibarbis) were collected in Sebangau. Samples were collected within the LAHG (Laboratorium Alam Hutan Gambut: The Natural Laboratory for the Study of Peat Swamp Forest), an area of 500 km 2 in the north-east of the Sebangau Forest, which was designated for the purpose of scientific research in 1997, and is managed by CIMTROP (Centre for International Cooperation in Management of Tropical Peatland) at the University of Palangka Raya. Sebangau is located in Central Kalimantan and consists of peat swamp forest. 39 faecal samples from Sumatran orangutans (P. abelii) were collected from two localities in the Gunung Leuser National Park. Suaq Balimbing is a coastal swamp situated on the western border of Gunung Leuser National Park. Bukit Lawang is located in Northern Sumatra, on the eastern border of the Gunung Leuser National Park, in hill dipterocarp forest ( Fig. 1 and Table 1).
Except for orangutans in Bukit Lawang, all animals in this study can be considered wild. They have no physical contact with humans. In Bukit Lawang orangutans are considered semi-wild because they were released after a reintroduction process and therefore are less hesitant to have contact with humans 12 .
Wild and semi-wild animals in the field were followed from nest to nest, for several hours and/or days. All faecal samples were collected immediately after defecation from identified known individuals in conformity with the ethical treatment of non-human primates. Collected samples were preserved in ethanol (96%). Given codes, sex, animal category, population and location were recorded for each animal (Supp Table 1). www.nature.com/scientificreports/ Collection of eggs and microscopy. For each sample; two grams of faecal material was homogenised with 15 ml water and filtered thought a 100 µm sieve ( ~) for microscopic examination. The filtered material was then examined using Sheather's flotation solution 63 and the sedimentation method, recommended in the approved guidelines of the Clinical and Laboratory Standards Institute for identification of intestinal tract parasites 64 . The whole sediment sample examination was repeated several times after the same sedimentation. All observed eggs were examined under an Olympus BX50 light microscope under the magnification X100. Egg size, shape, colour and internal structures were recorded. Eggs tentatively identified as Oesophagostomum spp., based on morphological characteristics, were transferred using a micropipette to eppendorf tubes contain 200 μl PBS for molecular confirmation.
DNA extraction, PCR conditions and sequencing. The isolated eggs from one individual were combined into one tube and were washed and centrifuged twice with 200 μl PBS to remove any residual ethanol. DNA was extracted using the QIAamp DNA mini kit. DNA was purified from tissues (Qiagen, France) according to the manufacturer's recommendation with the following modifications: samples were subsequently mixed with 0.2 g sterile mixture of 0.1/0.5 mm glass beads (Bertin Technologies, www. prece llys. com) and 180 μl of ATL buffer (Qiagen). Samples were homogenised using a TissueLyser (Qiagen Retsch GmbH, Hannover, Germany) for 5 min at maximum speed. Hereafter, the suspensions were incubated over night at 56 °C. DNA was eluted in 70 μl of distillate water. The entire ITS1-5.8S-ITS2 region of the Oesophagostomum genus was amplified by a two-step semi nested PCR using forward and reverse primers (Supp Table 2) flanking the 3′ terminus of 18S rDNA and the 5′ terminus of 28S r DNA, respectively. PCR amplification was performed in 25 μl containing 12.5 μl AmpliTaq Gold Master Mix, 0.5 μl of 10 μM each primer and 5 μl of template. The second of the two nested PCR steps was done following the same protocol using 3 μl of the first PCR product as template. For each run, as a negative control, nuclease-free water was added to the PCR mix instead of the DNA sample.
The PCR I cycles were: 95 °C (10 min The amplified products were subjected to 2% agarose gel electrophoresis in TAE buffer. PCR products were commercially purified and sequenced in both directions using nested PCR primers by Macrogen (Macrogen, Amsterdam, Netherlands).

Statistical analysis.
Only PCR positive samples were taken for analysis. Statistical analysis was carried out using R version 3.5.1 65 . Generalized linear models are recommended for the analyses of parasite data but as full factorial or simpler models could not be fitted due to lack of convergence, χ2 tests were used. Tests carried out were two-tailed with Yates correction. Animals were deemed positive if any faecal sample collected from them during the study tested positive by PCR. A p ≤ 0.05 was considered statistically significant. For more details, see Stuart et al 12 .
Sequence alignment, sequences and phylogenetic analyses. All sequences were manually checked, assembled and concatenated using CodonCode Aligner software (www. codon code. com). Sequences were aligned using Clustal W with Bioedit software version 7.2.5.
For each population, aligned FASTA files were collapsed into variable sites and haplotypes for parsimony network were reconstructed using DNASP v.5.10.01 66 . A statistical parsimony network to infer relationships among sequences was created with Network 5.0 67 .
To examine the relationship of the all known Oesophagostomum species to date, a phylogenetic tree was constructed using a set of reference sequences belonging to different species from Genbank (Supplementary Table S4). For a more descriptive tree, only two of the multiple representative sequences of a haplotype were used in analysis. Phylogenetic analyses were then performed after multiple alignments of the obtained partial ITS1-5.8S-ITS2 sequences (258 nucleotides). We used both Bayesian (Fig. 3) and maximum likelihood (ML) (Supplementary Figure S2) methods for phylogenetic tree construction. The best-fitting ML model under the Akaike information criterion was HKY86 + G6 for nucleotides as identified by JModelTest v.0.1.1 40 .
For phylogenetic Bayesian inferences, the software MrBayes 3.2.6_1 68 was performed using the web service of Phylogeny.fr (https:// ngphy logeny. fr/ tools/ tool/ 281/ form). Markov Chain Monte Carlo (MCMC) parameters were set to sample a tree every 10 or 100 generations of 100.000 generations and the burn-in was set at 250 trees sampled. Resulted tree was accepted when the average standard deviation of split frequencies is lower than 0.05 and the average Potential Scale Reduction Factor (PSRF) for parameter values is around 1.0 (0.9-1.1), simultaneously.
The highest-likelihood DNA trees and corresponding bootstrap support values were obtained by PhyML 69 (freely available at the ATGC bioinformatics platform http:// www. atgc-montp ellier. fr/) using nearest neighbour interchange plus subtree pruning recrafting (NNI + SPR) branch swapping and 100 bootstrap replicates. Necator americanus (LC036563) was used to root the tree. The resulting tree was drawn using iTOL v.5.4 70 .
Genetic diversity and population differentiation. The genetic diversity estimates (N, Number of sequence obtained; S, Number of polymorphic sites; K, Number of haplotypes; Hd, Haplotype diversity; π, Nucleotide diversity; D, Average number of nucleotide difference between populations) were computed using DNASP v.5.10.01 66 . Finally, to evaluate whether the genetic differentiation between populations was associated with the geographical isolation of islands, analysis of molecular variance AMOVA was performed in Arlequin www.nature.com/scientificreports/ ver 3.1 71 . Mantel test for matrix correlation between genetic distance and geographic distance was performed by using IBD 72 with 1000 permutations.
Ethics approval. All the research reported in this manuscript adhered to the legal requirements of the country in which the work took place. Since the collection of faecal samples from orangutans and gibbons was non-invasive and did not involve interaction with or distress to the animals, the study was not reviewed by an animal ethics committee.

Data avaliability
The newly generated sequences were deposited in the GenBank database under the accession numbers: MW756966-MW756992.