Echinococcus multilocularis and Echinococcus shiquicus in a small mammal community on the eastern Tibetan Plateau: host species composition, molecular prevalence, and epidemiological implications

Study area

Field studies were conducted at Yunbo Gou (33°11'N, 97°39'E) in northwestern Shiqu County (Ganze Tibetan Autonomous Prefecture, Sichuan Province, China) with an elevation of 4200–4700 m above sea level. Habitat vegetation is primary Kobresia meadow, with shrubs, mainly Potentilla fruticosa and Salix cupularis, from the middle to the summit of some hills. Wang et al. [22] classified the vegetation in this region into four categories: grassland; grassland and shrub; shrub; and disturbed areas. Grassland was the main vegetation type, covering > 90% of the study area [23]. The warm season extends from late June to mid-August, and is the suitable time for small mammal capturing.

Sampling of small mammals

Small mammals were collected between July and August 2014, during the annual wildlife plague (Yersinia pestis) surveillance, conducted by the Shiqu County Center for Disease Control (Shiqu CDC). Four 50 × 50 m trapping plots were randomly set on grassland at Yunbo Gou. Small mammals in the plot were trapped using break-back traps (size: 12 × 6.5 cm) set at the entrance of their dens. In total, 400 traps were set in each trapping plot. The trapping period of each plot was set for 24 h (10:00 h to 10:00 h the next day).

Each small mammal was sexed and its body measured; the head was stored in a 50 ml capped tube with 95% ethanol for further species identification in the laboratory. The bodies were then dissected, and any lesions of Echinococcus spp. in the thoracic and peritoneal cavities and organs were carefully checked. When a lesion was detected, a small portion was cut and checked under a microscope for presence of protoscoleces. To a typical Echinococcus lesion, protoscoleces can be checked out, while to those atypical lesions, ones that were either too small (e.g. with a diameter < 1 mm) or calcified, protoscoleces could not be observed. Therefore, to further confirm the existence of Echinococcus spp., samples of each typical and atypical lesion and from livers of small mammals without visible lesions were stored separately in 2 ml storage tubes in 95% ethanol and stored at -20 °C for further PCR analysis.

Small mammal species identification

Small mammal species identification was based on pelt color patterns, body measurement data, and skull-mandible morphological characteristics (e.g. of the molars) according to Luo et al. [24] and Smith & Xie [25]. To further confirm our identification, the specimens were compared with small mammal species specimens of the museum of the Northwest Institute of Plateau Biology, Xining.

DNA extraction and PCR

Other samples for molecular analyses

DNA of Echinococcus spp. from other host species involved in local transmission cycles was also used in this study. These samples included: three Echinococcus-positive Tibetan fox fecal samples previously used by Jiang et al. [15]; four positive domestic dog fecal samples infected by E. multilocularis previously used by Boufana et al. [13]; and AE lesion samples from six human patients living in Shiqu County between 2002 and 2007 [28]. The pretreatment and copro-DNA extraction protocols of Tibetan fox and dog fecal samples followed Jiang et al. [15] and Boufana et al. [13], respectively. Treatment of the human AE samples followed the small mammal sample pretreatment and DNA extraction protocol as described above.

PCR product cloning and sequencing

PCR products were subjected to agarose gel electrophoresis and stained with ethidium bromide (EB); positive screening results indicated that the target gene fragments had been amplified. Positive amplicons were excised carefully from the gel and purified with the TIAN gel Midi Purification Kit (TIANGEN, Beijing, China). Purified products were cloned into the T-Vector pMD 19 (TaKaRa) in strict accordance with product instructions and transformed into competent Escherichia coli cells. DNA sequencing was conducted by Shanghai Majorbio Bio-Pharm Technology Co. Ltd. The results were compared with the NCBI database. (http://www.ncbi.nlm.nih.gov/BLAST).

Statistics

Percentages of each trapped mammal species were calculated and a Chi-square test was used to test the sex bias among them. Plateau pikas and vole species were main species of trapped mammals. Therefore, the different distribution patterns between voles and pikas were also tested using a Chi-square test.

Given that PCRs with different primer pairs might have different results, the prevalence of the same Echinococcus spp. detected by different primers could be inconsistent. Therefore, all the visually identified (i.e. individuals with typical and atypical lesions) and Echinococcus DNA-detected individuals were analyzed using a Chi-square test for trends in proportions to determine if the detection of the same Echinococcus spp. by different primers (i.e. cox1 and nad1) was significantly different. Meanwhile, we defined the conservative prevalence of Echinococcus spp. by the percentage of positive samples detected by both cox1 and nad1 primers, and the maximum prevalence by the percentage of positive samples detected using at least one of the two genes.

To study how body condition might influence the detection of echinococcosis, all the visually identified (i.e. individuals with typical and atypical lesions) and Echinococcus DNA-detected individuals were analyzed using logistic regression models with four variables: (i) relative body weight (RW): the body weight of each individual divided by the heaviest one of the same species collected in this study; (ii) relative head-body length (RHBL): the head-body length of each individual divided by the longest one of the same species collected in this study; (iii) cross effect of relative weight and head-body length (CWL), expressed using the product of (ii) and (iii) of each individual; and (iv) lesions: three types, (typical, atypical, and no obvious lesion). The generalized R² [29] and AICc [30] of each model were calculated. The model with the lowest AICc was selected as the best model, while models with ΔAICc < 2 compared with the AICc of the best model were also selected.

The development of Echinococcus lesions is known to be positively related with age of the host [31–34]. Morris [35] recommended using the dry weight of eye lenses to evaluate age, as practiced by Burlet et al. [34, 36]. We could not use this method because no scrutinized age-eye lens weight analyses of our studied species have ever been reported in China. Thus, we used weight and body length to indicate the relative age of each trapped individual [35]. Logically, it should be easier to identify Echinococcus infection in larger, heavier and older individuals [35]. Given that body size can differ significantly among species, relative measurements (i.e. RW and RHBL) of each individual were calculated by dividing the measurement with the data of the largest or heaviest individual of its own species collected in this study.

All statistical analyses were conducted using R 3.4.0 (http://www.r-project.org).

Phylogenetic analyses

To analyze the phylogenetic relationships between E. multilocularis collected from small mammal hosts and from other host species of both local and wider geographical transmission cycles, maximum likelihood trees (ML trees) and Bayesian inference trees based on cox1 gene fragment haplotypes were constructed. Haplotypes of cox1 sequences were acquired from this study and 18 E. multilocularis cox1 sequences were selected from GenBank. One E. shiquicus (from one of the three Tibetan fox fecal samples, F12033) and one E. granulosus (s.s.) cox1 sequences (GenBank accession ID: KJ628374.1) were used as outgroups (Additional file 1: Table S1). When selecting sequences online, only data from indigenous host species were used. We did not build trees based on E. multilocularis nad1 gene sequences because the nad1 amplicons in E. multilocularis were too short (Table 1) and online data from different geographical regions were insufficient. Phylogenetic trees based on E. shiquicus sequences were not given in this study because only limited molecular data from this species from only a small area of the eastern Tibetan Plateau are available on GenBank, which would result in trees of E. shiquicus being less informative than trees of E. multilocularis.

Before construction of the phylogenetic trees, sequences were edited (Bioedit 7.0.9 [37]) and aligned (ClustalX2 [38]). Haplotypes were summarized, and Tajima’s D and Fu’s Fs tests were conducted by DnaSP v.5 [39]. Substitution saturation of the sequence matrix was tested by DAMBE 5 [40]. jModelTest v.2.1.4 [41] was used to test for the best-fit models of nucleotide substitution. Finally, ML trees were constructed using MEGA 7 by setting a ‘GRT+I’ substitution model with 1000 bootstrap replications. Bayesian trees were constructed using MrBayes 3.2.4 [42] by setting the ‘TIM3+I’ substitution model, using Markov Chain Monte Carlo (MCMC) posterior probability estimation for 2,000,000-generation with a 1000-generation sampling interval, and discarding the first 25% aging samples when summing up each tree. The best Bayesian tree was then compiled and processed by FigTree 1.4.3 (http://tree.bio.ed.ac.uk/software/figtree). Finally, a network diagram was drawn using Network 5.0 (http://www.fluxus-engineering.com).

Small mammal species composition

Prevalence of Echinococcus spp. in the small mammal community

There was a significant difference between the use of cox1 vs nad1 in the detection of Echinococcus infection. For E. shiquicus, significantly more infections were detected using nad1 primers than with cox1 (χ² = 10.480, df = 1, P = 0.001), whereas the use of cox1 primers detected more E. multilocularis infections than with nad1 (χ² = 7.415, df = 1, P = 0.006).

Regression model simulation revealed that Echinococcus infection was detected in individual hosts with significantly longer RHBL (log odds ratio, 19.134), while their RBWs were significantly lighter than those without infection (log odds, -16.010) (Table 4). Moreover, although typical lesions were useful signs of infection (with a log odds ratio of -0.556 relative to the molecular results), atypical lesions were misleading and caused significantly more misidentification of infections (log odds ratio of -21.982, Table 4).

Phylogenetic relationships between E. multilocularis cox1 gene haplotypes

In total, 45 haplotypes of a 768 bp long E. multilocularis cox1 gene fragment were acquired from 23 small mammals (76 sequences), six human AEs (6 sequences), four dog feces (4 sequences), and two Tibetan fox fecal samples (5 sequences) in this study, with an additional 17 sequences from online sources. Among the 76 sequences from Shiqu County, 33 haplotypes of E. multilocularis (i.e. Hap06-Hap38) were confirmed (Additional file 1: Table S1), of which Hap06 was the dominant haplotype, with 58 sequences covering all the human AE, dog, and Tibetan fox fecal samples, and 21 small mammal samples from L. fuscus, M. limnophilus and O. curzoniae (Fig. 1, Additional file 1: Table S1). The haplotype and nucleotide diversities for the E. multilocularis cox1 gene from the small mammal community in this study were 0.655 ± 0.064 and 0.0014 ± 0.0002, respectively. Significant negative values in both Tajima’s D (D = -2.83311, P < 0.001) and Fu’s Fs (Fs = -46.942, P < 0.001) tests indicated that the E. multilocularis population in the small mammal community of Shiqu County is currently expanding, and the dominant Hap06 haplotype may be under strong purifying selection.

Both ML and Bayesian trees were constructed based on 56 sequences, including 54 E. multilocularis sequences and two E. granulosus (s.s.) and E. shiquicus sequences as outgroups (Additional file 1: Table S1). Both trees demonstrated identical topological relationships between haplotypes, but only the Bayesian tree is included here (Fig. 2) because of its more concise structure. The Tibetan Plateau cluster was distinctive from other geographical regions of the world, and mainly comprised the Shiqu haplotypes Hap06-Hap38, the Qinghai Province haplotype Hap05 (also from the Tibetan Plateau), and two other sequences (Hap04 and Hap06-Vole-NO.RUS, Additional file 1: Table S1) from neighboring geographical regions (Fig. 2).

Small mammal host species community composition

Although all the six small mammal species reported in this study have been recorded previously on the Tibetan Plateau, information about small rodent species such as arvicolids and cricetids was largely lacking from Shiqu before the 21st century [45, 46]. He et al. [20] reported alveolar echinococcosis in the lagomorphs Ochotona curzoniae and Lepus oiostolus, and also in the commensal rodents N. irene, and Mus musculus in Shiqu. Raoul et al. [19] studied the habitat ecology of small mammal communities in Shiqu, reporting for the first time P. leucurus and C. kamensis in Sichuan Province. However, neither of these two studies reported L. fuscus in Shiqu. The distribution area of L. fuscus was judged to be restricted to southern Qinghai Province, and this was not considered to be a species distributed to Sichuan [24, 25, 43, 45]. However, in local plague surveillance studies [47, 48] and pest control schemes [49, 50], L. fuscus was judged to be a dominant rodent species in the high plateau pasture areas of northwest Sichuan. Lasiopodomys fuscus is morphologically similar to other vole species, such as P. leucurus, and thus could be included in reports in error if an accurate morphological identification is lacking [51]. Shiqu is located at the northwestern point of Sichuan Province, on the southern border of Qinghai; therefore, it is possible that what is thought to be the southern limit of the range of L. fuscus is incorrect [46]. The current study confirmed the presence of L. fuscus in Shiqu County, where it was the most abundant species of the five trapped rodent species and is likely to have a larger population than that of O. curzoniae in this region (Table 2).

Importance of rodent species in the transmission of echinococcosis

Although the potential importance of rodent species as intermediate hosts of E. multilocularis on the eastern Tibetan Plateau has been mentioned previously elsewhere [10, 14, 51, 52], O. curzoniae were the most frequently examined small mammal host species with large sample sizes [8, 17, 18, 20, 21, 52–54]. For example, when evaluating the prevalence of E. multilocularis, Zhao [21] reported a prevalence of 15.2% in O. curzoniae (34/224) and 20% in L. fuscus (1/5 individuals). Similarly, Zhang & Wang [18] reported a prevalence of 5.3% in O. curzoniae (62/1177), but only 0.4% in L. fuscus (1/269). Thus, much smaller sampling sizes might be an important reason for the reported low prevalence of E. multilocularis in rodent species. In terms of E. shiquicus, O. curzoniae was previously the only confirmed intermediate host species [12, 54]. By contrast, our data revealed that, among the three most-abundant small mammal species, both M. limnophilus and L. fuscus had a significantly higher molecular prevalence of E. multilocularis and E. shiquicus than did O. curzoniae (Table 3). Moreover, E. shiquicus DNA was detected not only in the three most-abundant small mammal species, but also in the other two rodent species sampled (P. leucurus and C. longicaudatus) (Table 3). Although trapping data might not reflect the exact abundance of each species, the higher molecular prevalence of both E. multilocularis and E. shiquicus (Table 2) in rodents than in pikas suggests that rodent species are probably more important intermediate host species than pikas during the transmission of both E. multilocularis and E. shiquicus.

Traditionally, epidemiological studies of echinococcosis in western China mainly focused on human and domestic animal populations because of their obvious direct interactions especially regarding transmission of human CE and a role for dogs in risk of human AE. Phylogenetic analyses revealed that all 33 cox1 gene haplotypes from Shiqu County were closely related (Fig. 1), and clustered with haplotypes from other studies to form the Tibetan Plateau group, which is distinct from haplotypes from other geographical regions (Fig. 2). Hap06 was the dominant haplotype of E. multilocularis discovered in humans, domestic dogs, and almost all the wildlife host species tested in this study, including the three most dominant small mammal species (i.e. O. curzoniae, M .limnophilus and L. fuscus) (Fig. 2, Additional file 1: Table S1). These results confirm that small mammal species and dogs together comprise a wildlife and peri-domestic ecosystem for transmission of E. multilocularis. Thus, understanding the epidemiology in wildlife host species is pivotal to understanding the life-cycles and transmission ecology of these parasites.

For a parasite species such as E. multilocularis, which can utilize multiple host species, understanding the interspecific dynamics between host species is essential for understanding its transmission [55]. Both Tajima’s D and Fu’s Fs tests revealed the population expansion of E. multilocularis in the small mammal community in the study region, as supported by Nakao et al. [28]. Could this expansion be the result of the frequently reported increases in small mammal populations and their interspecific dynamics in western China? Ochotona curzoniae and several vole species are the main prey of the Tibetan fox (V. ferrilata) on grasslands of the eastern Tibetan Plateau [16]. The significant differences in prevalence of E. multilocularis among voles and O. curzoniae (Table 3) suggest that interspecific dynamics among these small mammal species could be essential factors influencing the predator-prey food chain, which would affect the epidemiology of E. multilocularis across all pasture areas on the Tibetan Plateau. For example, Wang et al. [56] reported a high density of both O. curzoniae and voles in open grassland areas within 2 km of villages. Ochotona curzoniae is usually blamed for degrading the grassland ecosystem of the eastern Tibetan Plateau and, consequently, has been a main target for poisoning to protect the grasslands [49, 50, 57]. Areas around villages are usually where poisoning programs are located. The fact that there was a significantly higher prevalence of E. multilocularis in vole species compared with O. curzoniae (Table 3) suggests that pika may have a lesser role compared to voles in the transmission of E. multilocularis. Furthermore the deliberate poisoning of O. curzoniae could result in increased densities of vole species, accompanied by a higher prevalence of E. multilocularis in wildlife, dogs and humans in local areas.

Small mammal body condition and detection of Echinococcus DNA

The regression model revealed that individuals with larger RHBL were more likely to be detected with Echinococcus, based on the molecular tests (Table 4). Linear body dimensions are more or less correlated with age, especially in animals with short lifespans [35]. Thus, the results of the current study showed that older small mammals were more likely to have E. multilocularis or E. shiquicus infections. Similarly, Burlet et al. [34] reported that older Arvicola terrestris had a higher prevalence of E. multilocularis than younger animals. In small mammal hosts, E. multilocularis requires several months to grow from an oncosphere to a fertile metacestode [58]. Consequently, the development of metacestode lesions can be synchronized with the aging process of its hosts, such as rodents and pikas, further confirming the importance of understanding the host population structure and dynamics when evaluating the present infection burden and predicting the future trend of a specific parasitic species.

Although typical lesions were indicative of infection with specific Echinococcus species, as revealed by the regression model (Table 4), they were harder to find than were the frequently detected atypical lesions (Additional file 1: Table S2). However, low values of atypical lesions for determining Echinococcus infection (Table 4, Additional file 1: Table S2) required the use of molecular analyses. Nevertheless, such analyses provide the molecular-positive rate of the parasite, which is not equivalent to its true prevalence. For example, E. shiquicus DNA was detected in several rodent species; however, if no typical lesions with metacestodes were observed (Additional file 1: Table S2), then it was not possible to conclude that there was an established parasitic infection. Thus, the function of rodent species as natural intermediate hosts of E. shiquicus still requires further study. Moreover, many primers have been designed to test various nuclear [59, 60] and mitochondrial [15, 26, 27, 61] genes of Echinococcus species. In this study, the cox1 and nad1 genes were used. However, the significantly inconsistent results achieved with both genes (Additional file 1: Table S2) suggested that multiple parasite genes should be tested in the same epidemiological study. In addition, a prevalence interval that is based on the different levels of detection as determined by the available genes (e.g. the maximum and conservative prevalence defined in this study) would be more objective than using a single prevalence.