Research | Open | Published:
Absence of genetic structure in Baylisascaris schroederi populations, a giant panda parasite, determined by mitochondrial sequencing
Parasites & Vectorsvolume 7, Article number: 1591 (2014)
Infection with the parasitic nematode, Baylisascaris schroederi (Ascaridida: Nematoda), is one of the most important causes of death in giant pandas, and was responsible for half of deaths between 2001 and 2005. Mitochondrial (mt) DNA sequences of parasites can unveil their genetic diversity and depict their likely dynamic evolution and therefore may provide insights into parasite survival and responses to host changes, as well as parasite control.
Based on previous studies, the present study further annotated the genetic variability and structure of B. schroederi populations by combining two different mtDNA markers, ATPase subunit 6 (atp6) and cytochrome c oxidase subunit I (cox1). Both sequences were completely amplified and genetically analyzed among 57 B. schroederi isolates, which were individually collected from ten geographical regions located in three important giant panda habitats in China (Minshan, Qionglai and Qinling mountain ranges).
For the DNA dataset, we identified 20 haplotypes of atp6, 24 haplotypes of cox1, and 39 haplotypes of atp6 + cox1. Further haplotype network and phylogenetic analyses demonstrated that B. schroederi populations were predominantly driven by three common haplotypes, atp6 A1, cox1 C10, and atp6 + cox1 H11. However, due to low rates of gene differentiation between the three populations, both the atp6 and cox1 genes appeared not to be significantly associated with geographical divisions. In addition, high gene flow was detected among the B. schroederi populations, consistent with previous studies, suggesting that this parasite may be essentially homogenous across endemic areas. Finally, neutrality tests and mismatch analysis indicated that B. schroederi had undergone earlier demographic expansion.
These results confirmed that B. schroederi populations do not follow a pattern of isolation by distance, further revealing the possible existence of physical connections before geographic separation. This study should also contribute to an improved understanding of the population genetics and evolutionary biology of B. schroederi and assist in the control of baylisascariasis in giant pandas.
Genetic markers are used increasingly to reconstitute the epidemiological history of human and veterinary diseases and to identify environmental factors underlying the spread of pathogens -. However, most of these studies focus on bacterial or viral diseases. For parasitic diseases, the complex life cycle of the causative organisms (such as helminthes), including different developmental stages and hosts, make their transmission heavily dependent on multiple environmental factors. Traditional approaches to study parasites, particularly parasitic nematodes, often concern their pathogenicity, prevalence and possible vaccines, but the complex evolutionary relationships of parasite populations are limitedly studied ,. Recently, DNA sequence data have been successfully used to characterize and understand the genetic diversity of parasites and to explain and predict their likely dynamic importance. Moreover, such genetic techniques have been also applied to identify recurrent co-evolutionary units within larger systems and to study how parasites survive short-term environmental perturbations and evolve in response to long-term environmental changes . Of course, a range of factors should be balanced when choosing suitable DNA markers, including cost and the purpose for which they are used. Mitochondrial DNA (mtDNA) sequences are powerful and reliable molecular markers for detecting population structures and inferring population differences as a result of its mutational rates which are proposed to be more rapid than nuclear DNA, and presumed lack of recombination and maternal inheritance. Additionally, mtDNA contains both rapidly diverging and highly conserved regions, thus making different mtDNA regions suitable for resolving distinct taxonomic issues in parasites -. For example, the genes for mitochondrial cytochrome c oxidase subunits I-III (cox1-3) and nicotinamide dehydrogenase subunit IV (nad4) have recently been used to investigate and compare the genetic variation and population structure of the parasitic nematodes Baylisascaris columnaris, Dictyocaulus viviparous, and Haemonchus contortus,. Moreover, Iñiguez et al. genetically characterized Ascaris spp. from humans and pigs using the mitochondrial cytochrome c oxidase subunit 1 (cox1) and NADH dehydrogenase subunit 1 (nad1) genes, and demonstrated that the haplotypes cox1 H14P3 and nad1 H12P17 are widely separated in Brazil .
The giant panda (Ailuropoda melanoleuca), a flagship species for wildlife conservation, is one of the world’s most recognized and threatened animals with a total wild population size estimated to be only 1,600; these individuals are restricted to six mountain ranges in China (Minshan, Qionglai, Qinling, Daxiangling, Xiaoxiangling and Liangshan), of which 44.36% are in Minshan, 27.38% in Qionglai and 17.23% in Qinling -. Infection with the parasitic nematode Baylisascaris schroederi (Ascaridida: Nematoda) is the most important cause of death in wild giant pandas and thus it poses a significant threat to these populations ,. Adult stages of this parasite usually inhabit the intestines of the giant panda, while migrating larval stages can disseminate into various body tissues. Damage to bodily tissues can include extensive inflammation and scarring of the intestinal wall and parenchyma of the liver and lungs (caused by larvae), as well as intestinal obstruction, inflammation and even death (mainly caused by adults) -. As with other ascaridoids, B. schroederi infection follows ingestion and the life cycle is complete without the need for an intermediate host. Infective second-stage (L2) larvae hatch in the small intestine, migrate into the liver and the lung, and finally arrive at the small intestine where they mature, mate and produce eggs . These eggs are highly resistant to decontamination and environmental degradation and can remain viable in moist soil for years, thus acting as a substantial reservoir for new infections. In nature, B. schroederi infection rate among wild pandas is often over 50% and can even be 100%, thus making it a leading cause of death in wild populations ,. Despite its detrimental health impact on wild pandas, this parasite was first described only in 1939 as Ascaris schroederi, before being renamed as B. schroederi in 1968 . Subsequent studies regarding B. schroederi focused mainly on morphological characteristics, fundamental biology, life cycle, pathogenicity, vaccine development and prevalence ,,. Recent advancements in PCR and sequencing technologies have led to the publication of several gene fragment sequences from B. schroederi (e.g., ITS-1, ITS-2 and 5.8S), as well as complete mtDNA and microRNA ,,,. These data have contributed to an improved understanding for the molecular biology and genetics of this parasite. However, studies to date that are involved to the genetic diversity of B. schroederi remain scarce, with only one report publically available in which a single mt cytb was utilized for this purpose .
Considering the limitation of using just a single gene to define the genetic structure of a parasite -, the aim of this present study was to combine two different mtDNA markers, cox1 and atp6, which possess a symmetrical evolutionary rate (a higher genetic variation of atp6 than cox1), to expand and improve our knowledge of the genetic diversity of B. schroederi. Our data analysis included plotting transition and transversion frequencies of these genes and determining the genetic distances between populations. These data provide more detailed information into the genetic diversity of B. schroederi populations in three endemic areas of China (Minshan, Qionglai and Qinling mountain ranges), and enable a reevaluation of the population structure. In addition, we compared the genetic diversity of the population and inferred the phylogenetic relationships between different isolates according to their geographical distribution.
Between two and three parasite specimens were collected per giant panda (10-15 pandas per mountain range), except for Qinling (see Figure 1 and Additional file 1: Table S1 for geographical sampling details). Fifty-seven worms originating from 28 dead, injured or rescued giant pandas from three mountain ranges (Minshan, Qionglai and Qinling) were obtained during 1994-2014 and all were included in this present study. All worms were washed in physiological saline buffer, identified according to morphology , and stored at -20°C for preservation.
Genomic DNA was extracted from 0.5 g of tissue from each individual B. schroederi by standard sodium dodecyl sulfate/proteinase K treatment and phenol/chloroform extraction, followed by column-purification using the Wizard® SV Genomic DNA Purification System (Promega, USA), eluted with 30 μL of elution buffer (pH 8.0) according to the manufacturers’ recommendations, and then stored at -20°C for PCR amplification. Two pairs of PCR primers were designed to amplify the complete mt atp6 and cox1 genes using Primer Premier 5.0 software (Premier Biosoft International, Palo Alto, CA, USA) and the published mtDNA sequence of B. schroederi (GenBank Accession No.: HQ671081). The primer pairs were as follows: atp6, forward: 5′-CGCGGATCCTTCGATATTCGTGGCCT-3′, and reverse: 5′-CGCAAGCTTCTAATATGGTGTCTTCGG-3′; cox1, forward: 5′-TTTAGAGGTTGGAATGTAGGGT-3′, and reverse: 5′- CCATCCCCTTAATCTGCAAT-3′. All PCR reactions contained ~20 ng of genomic DNA and were carried out in 50-μL reaction volumes containing 25 μL 2× Phusion High-Fidelity PCR Master Mix (Finnzymes OY, Espoo, Finland), 2 μL of each primer, 2 μL DNA and 19 μL of ddH2O under the following cycling conditions: an initial denaturation at 95°C for 5 min; then for cox1, 35 cycles at 95°C for 1 min, 50°C for 1 min, and 72°C for 1 min; but for atp6, 35 cycles of 95°C for 45 s, 50°C for 45 s, and 72°C for 45 s; followed by a final step at 72°C for 10 min. Negative controls were included in each PCR to test for possible contamination. All PCR products were examined on agarose gels (1%) to verify that they represented the target bands. The corrected gel-isolated fragments were subsequently column-purified and sequenced. Sequencing was performed by Invitrogen (Invitrogen Biotechnology Co. Ltd., Shanghai, China) using terminator-based cycle sequencing with BigDye chemistry (Applied Biosystems, Foster City, CA, USA) on an ABI 3730 DNA sequencer (Applied Biosystems). A “primer-walking” strategy (in both directions) was used and, to ensure maximum accuracy, each amplicon was sequenced twice independently. A third PCR product was sequenced in case of any discrepancies. Consensus sequences were deposited in the GenBank database.
Phylogenetic and population genetic analyses
Sequences with consensus lengths (600 bp for atp6 and 1578 bp for cox1) were aligned using the ClustalW program in MEGA v.5.0 . Nucleotide alignments were generated based on protein alignments according to codon alignment. Phylogenetic relationships and nucleotide/haplotype diversity were initially determined for each mitochondrial gene. Given the lack of variable and parsimony-information sites for each gene, the sequence data from atp6 and cox1 were concatenated into a single 2,178-bp fragment for each specimen for realignment according to the aforementioned method. During this procedure, any ambiguous regions within these alignments were filtered with Gblocks 0.91 b . Phylogenetic trees inferred from the atp6, cox1 and pooled datasets were constructed using the maximum parsimony (MP) method, with each node estimated using 1,000 bootstrap replicates. Data were plotted with MEGA v.5.0 . To further define the current population structure of B. schroederi, a Bayesian statistical method , coupled with the best-fit evolutionary model GTR + I + G determined by jModeltest  using the Akaike Information Criterion (AIC), was also used to analyze the genetic variation among populations based on the two mitochondrial genes separately and in combination (the model GTR + I + G was determined for both atp6 and cox1, and was directly applied in the following concatenated data). A total of 10 million generations were run and the trees were sampled every 10,000 generations. When the mean standard deviation of the split frequencies had reduced to less than 0.01, the first 250 trees were discarded as “burn-in” and the remaining trees were used to generate a 50% majority rule consensus tree for each data set. Homologous sequences obtained from Baylisascaris transfuga (Accession No. NC-015924.1 for both mitochondrial atp6 and cox1 genes) were used as outgroups in all phylogenetic analyses due to the close relationship between B. schroederi and B. transfuga. In addition, the Data Analysis in Molecular Biology and Evolution (DAMBE) program  was used to establish plots of transition (s) and transversion (v) frequencies against the K80 distance, and then the amino acid substitution propensity of these two genes was estimated.
Aligned FASTA files were collapsed into variable sites and haplotypes for parsimony network reconstruction using Dna SP 5.0 . A statistical parsimony network to infer relationships among sequences that had recently diverged was created with Network 4.6 . Haplotype (Hd) and nucleotide diversity (π) were calculated by Dna SP 5.0 . Furthermore, a Mantel test was conducted to calculate the correlation between pairwise geographies using Arlequin version 3.0 . Based on the same program, the coefficient of gene differentiation (G st) and index of genetic differentiation (F st) were also estimated, and the corresponding gene flows (N m) between populations were indirectly calculated as follows: N m = 0.5*(1-G st)/G st  or N m = (1-F st)/4F st . Of these parameters, the G st statistic has been used widely to evaluate population subdivision on the basis of the divergence of allele frequencies from random mating . Finally, to evaluate whether the genetic differentiation between populations was associated with the geographical isolation of the mountains, analysis of molecular variance (AMOVA)  was performed in Arlequin. The resulting data was used to examine population genetic structure and variation within and between the three sampling locations.
Historical demography of B. schroederi populations
Different approaches were used to infer historical demographic expansions, including those of Fu’s Fs, Tajima’s D and the mismatch-distribution tests -. To assess whether population expansion occurred in B. schroederi, 1,000 simulations in Arlequin  were performed to obtain the distribution of Fu’s Fs and Tajima’s D values under neutrality and to obtain P values for each population. Moreover, to measure frequency distributions of pairwise differences between sample sequences, mismatch-distribution analysis was also performed using Dna SP 5.0  for 1,000 simulations.
This study was reviewed and approved by the Research Ethics Committee in College of Veterinary Medicine, Sichuan Agricultural University (Ya’an, China; Approval No. 2011-028). All parasite samples were collected from giant pandas after the permissions of the Forestry Department of Sichuan Province and Shaanxi Rare Wildlife Rescue Breeding Research Center, with no specific permits being required by the authority for the parasite collection.
Haplotype sequence analysis
The 600-bp atp6 and 1,578-bp cox1 sequences from 57 B. schroederi isolates were deposited in GenBank under accession numbers KJ587749-KJ587805 and KJ587806-KJ587862, respectively. No insertion/deletions (indels) were detected. In total, 20 haplotypes (A1-20) were identified in atp6, 24 haplotypes (C1-24) in cox1 and 39 haplotypes (H1-39) in atp6 + cox1. Haplotypes A1 and A8 (inferred from atp6), C10 (inferred from cox1) and H11 (inferred from concatenated sequences) were shared by all samples collected from the three mountain ranges (Figure 2). Furthermore, A1, C10 and H11 haplotypes were used as alignment references to discover further variants. Interestingly, each mountain range population appeared to be associated with its own set of haplotypes. For example, the haplotypes C1, C2, C5, C6, C8, C9 and C12 deduced from cox1 were the typical haplotypes of the B. schroederi population in Qionglai. However, of the 39 haplotypes inferred from the concatenated sequences, just six haplotypes were shared by two or three of the mountain range populations (Figure 2C). Additionally, the genetic diversity of the 57 B. schroederi isolates from the three different regions was calculated (see Tables 1, 2, 3 and 4). A high level of haplotype diversity was maintained in the B. schroederi populations, but their nucleotide diversity was relatively low due to the richness of single-nucleotide substitutions. Of these populations, the Minshan mountain range population exhibited greatest haplotype diversity (0.84795 for atp6 and 0.92982 for cox1) and nucleotide diversity (0.00359 for atp6 and 0.00204 for cox1). However, the concatenated sequences seemed to have greater degree of haplotype diversity (all values >0.9) when compared with atp6 and cox1 (Table 4). In addition, nucleotide and deduced amino acid sequences were aligned. The maximum Kimura distance between the three B. schroederi populations was 0.003 for both atp6 and cox1. However, the estimated transition/transversion bias (R) differed between atp6 and cox1, and the R-value between the three populations was lower for cox1 (22.10) than for of atp6 (50.48) (data not shown).
To study the phylogenetic relationships between the B. schroederi isolates collected from different geographic regions, the atp6 and cox1 genes were analyzed, separately or in combination, to infer the phylogeny, and the corresponding constructed phylogenetic trees are shown in Figure 3. The samples in the three trees clustered into a mixed group with low posterior probability values, indicating that they would be difficult to distinguish from each other by phylogenetic analyses due to the small differences existing between individuals isolated from the different ranges. Therefore, it was impossible to elucidate evolutionary events based only on such undefined phylogenetic trees, and this meant that further network analysis was required.
In addition to phylogenetic trees, an entire network for all haplotypes was produced (Figure 2). Based on the 24 cox1 haplotypes, the network map revealed a star-like pattern around haplotype C5, which contained seven individuals, with three in the Qionglai and four in the Minshan populations. A total of 12 other haplotypes were found in the Qionglai population (C1–12), while there were 5 haplotypes present in the Qinling population (C3, C7, C10, C13 and 14). Interestingly, the haplotype C10 was shared by samples derived from all three mountain ranges, while haplotypes C4, C5 and C11 were shared by only the Qionglai and Minshan populations and C14 by only the Minshan and Qinling populations (Figure 2B, Additional file 2: Table S2). The remaining haplotypes appeared to be highly dispersed, and no obvious correlations were observed between sample clusters. For atp6, although the network analysis showed some patterns consistent with those revealed by cox1, lower haplotype resolution was observed when compared with cox1 (Figure 2A). Of these, haplotype A2, which originated from Qionglai and Minshan, was surrounded by the other haplotypes, but haplotype A1 was found in more individuals (19/57) than A2 (10/57). Within the concatenated sequence data, a more complex network was observed when compared with that produced for atp6 or cox1 genes alone. Further, 37 of the 39 haplotypes appeared to be around H11 and H6 haplotypes and formed two separate star-like patterns (Figure 2C). Of note, the haplotype H11 was surrounded by most of the isolates derived from the Qionglai population. To investigate the information of these sequences further, the number of haplotypes in each population was also calculated (see Additional file 2: Table S2).
Population genetic structure
Based on population genetic structure analysis, the genetic diversity of cox1 between the different B. schroederi isolates was highly consistent with that of atp6, but slightly lower than that of atp6 and cox1 in combination (Table 2, 3, 4 and data not shown). Interestingly, the corresponding nucleotide diversity (π) based on either atp6, cox1 or atp6 + cox1 exhibited significant stability between B. schroederi populations. In addition, for atp6, the genetic differentiation within populations was significant (P < 0.05) for all pairwise comparisons, which was in agreement with results from the phylogeny and network analyses. A similar result was also found for the concatenated sequence data (Table 4). In the cox1 data, 97.81% of the genetic variation was attributed to differences between individuals within populations, while the remaining 2.19% was due to differences between B. schroederi populations. These results confirmed that there was lower genetic differentiation between than within populations across the three B. schroederi populations examined here. This conclusion was validated further by the high rates of gene flow determined between different B. schroederi populations (Table 1). To test whether new B. schroederi isolates would constitute separate populations, G st and F st statistics were calculated to evaluate the variation in allele frequencies between B. schroederi populations. G st values varied from 0 to 1 with values > 0.25 and F st values ranged from 0 to 1 with values >0.27. As shown in Table 1, the highest G st and F st values for cox1 were observed between Minshan and Qinling populations (G st = 0.02541, P > 0.05; F st = 0.02669, P > 0.05), while the lowest G st and F st values were found between Qinling and Qionglai populations (G st = 0.01869, P > 0.05; F st = -0.02913, P > 0.05). Interestingly, a similar trend was also observed for F st values of atp6, but this was contrasted with G st values of atp6, where the highest G st value was detected between the Qinling and Qionglai populations (G st = 0.02126, P >0.05) (Table 1).
Based on cox1, atp6 and concatenated sequence datasets, Fu’s Fs and Tajima’s D statistics as well as the mismatch-distribution test were used separately to explore the demographic histories of the three B. schroederi populations. Of these, the results of the mismatch analysis inferred from either atp6 or concatenated data provided consistent insights into the demographic history of B. schroederi: all three populations (Minshan, Qionglai and Qinling) were unimodal and had undergone at least one expansion event (Figure 4A and 4C). Although the mismatch-distribution of the cox1 gene presented a multimodal outcome (Figure 4B), the statistical parameters RG and SSDs (not shown) meant that the hypothesis of population expansion could not be rejected. Additionally, the results of Fu’s Fs and Tajima’s D tests under neutrality generated by atp6, cox1 and concatenated sequences gave negative values (-13.735 for atp6, -11.065 for cox1 and -29.912 for atp6 + cox1 in Fu’s Fs, respectively; -2.11978 for atp6, -1.20886 for cox1 and -1.73588 for atp6 + cox1 in Tajima’s D, respectively) (Table 2, 3 and 4). Thus, these findings confirmed the existence of population expansion in the demographic history of all three B. schroederi populations and further supported the conclusions drawn from the mismatch analysis.
Parasites are a crucial component of natural communities, since they can directly or indirectly alter the structure of the communities by impacting the number of free-living species or their relative abundance . With advances in molecular genetics, rapidly increased knowledge concerning the genetic diversity of nematodes is becoming a prerequisite to elucidate basic biological and population characteristics of these parasites. More importantly, a better understanding of parasitic nematode population dynamics is fundamental to design new strategies to monitor and control these problematic organisms. Excitingly, many recent studies have used genetic markers (e.g., mitochondrial and nuclear DNA) to depict geographical movements of parasitic nematodes -. MtDNA markers are receiving increased attention for this purpose due to higher F st values than observed for nuclear DNA counterparts -,,. However, given the different evolutionary rates occurring in different regions of nematode mtDNAs  and the high degree of conservation of the cytb gene , in this present study two mitochondrial markers (atp6 and cox1) were combined to further explore the genetic variation and population structure of B. schroederi.
In this study, low G st and F st values indicated small genetic differentiation in the B. schroederi population (Table 1), and this was further supported by AMOVA analysis. Likewise, other factors also confirmed a non-significant level of genetic differentiation of the B. schroederi populations. For example, a high rate of gene flow (N m > 9) between the three populations was always observed for both atp6 and cox1, and it could be speculated that the effective number of B. schroederi individuals that exchanged genetic material was sufficient to overcome the effects of genetic drift. Thus, B. schroederi isolated from different geographical regions may be a homogenous species with strong population diffusion, little obstruction for gene flow, simple genetic structure and small differentiation. Moreover, the parasite is a specialized organism that parasitizes a specific host, and the genetics and behavior of the host can also have an important impact on parasite genomic variation. The population history of the giant panda, as characterized by Zhao et al. (2012), showed that pandas in the Qinling mountain range were distinct from other populations (including Minshan and Qionglai-Liangshan-Daxiangling-Xiaoxiangling) due to geographical isolation, human activities and global changes in climate . However, contrary to our expectations, B. schroederi showed only weak population subdivisions (Figure 2), particularly the Qinling population that shared frequent gene flow with the Minshan population. To a certain extent, these results implied that the evolutionary rate of B. schroederi may be disharmonious with its host. Similar conjectures have been proposed to interpret evolutionary relationships between the parasitic nematodes Syphacia obvelata and Trichuris muris and their host Mus musculus. Additionally, previous studies have indicated that a clear genetic structure would arise when inbreeding occurred, whereas outcrossing species tended to show little genetic structure ,. Therefore, we speculated that B. schroederi could be an ‘outcrossing’ species. This hypothesis is supported by knowledge of the genetic and reproductive behavior of the parasite and its host (for details, please see refs. ,,): (i) B. schroederi is a dioecious parasitic nematode with a direct developmental lifecycle; (ii) B. schroederi can be transmitted from one geographical origin to another because of the cross-regional reproductive behavior of pandas; (iii) thus, mating of B. schroederi could occur between individuals from different geographical areas. Alternatively, there could be another reason for possible mixed infections in the hosts, as B. schroederi specimens from the same individuals had low genetic differences and different haplotypes. Perhaps future complete mitochondrial genome-based population analysis involving a larger B. schroederi population would provide novel and/or further insights into this issue.
The other important finding of this present study was the low level of genetic diversity found across all three B. schroederi populations, which was similar to a previous report . Low nucleotide diversity was observed in B. schroederi and this revealed a relative lack of genetic variation across the ascaridoid, regardless of geographical origin and population size. Interestingly, similar findings were also described for other Ascaris species, including Ascaris galli, Ascaris suum, and Parascaris equorum. Together, these results implied that low level of diversity may be a common feature of ascaridoid nematodes. In addition, the leptokurtic distribution of the haplotypes showed several frequent haplotypes (C5, C10, C14; A1, A2, A8; H6, H11) as well as numerous rare haplotypes (Figure 2). These results indicated that the patterns of genetic differentiation were similar between the two mitochondrial genes (atp6 and cox1), which was probably due to the small sample size; however, the genetic diversity in the Qinling population was less clear. The analysis of atp6 and cox1 genes showed low geographic separation between the three mountain range populations, which was consistent with the findings of a previous report . Furthermore, network maps showed that the three populations only shared haplotypes A1 and A8 for atp6, C10 for cox1 and H11 for atp6 + cox1. The most frequent haplotypes were A1 of atp6, C10 of cox1 and H11 of atp6 + cox1 in the examined populations (Additional file 2: Table S2). This finding indicated that A1, C10 and H11 might be the most ancient haplotypes, as these often displayed a high frequency and tended towards a widespread geographic distribution . Of note, the haplotype A2 of atp6 was surrounded by nine clades while haplotype A1 was surrounded by just four, and this same phenomenon was observed for haplotypes C5 of cox1 and H6 of atp6 + cox1 (C5 and H6 had more clades than C10 and H11). Interestingly, haplotypes A2, C5 and H6 were shared only by Minshan and Qionglai mountain range populations. Thus, it could be speculated that these three haplotypes differentiated after from the ancient haplotype and would display a trend of widespread geographic distribution and may evolve to form a single haplotype across all three mountain ranges in the future. Taken together, our data indicated that the B. schroederi species was not genetically differentiated.
In general, certain aspects of parasites often closely track counterparts of their host’s biology due to long-term interplay between them ,. In the B. schroederi-panda co-evolutionary system, the host has recently experienced two population expansions, two bottlenecks and two divergences . Thus, it is reasonable to assume that the demographic history of the panda may influence that of its specialized parasite B. schroederi. Encouragingly, negative Fu’s Fs (high) and Tajima’s D (low) values (Table 2, 3 and 4) and the results of the mismatch analysis (Figure 4) confirmed that B. schroederi had been subjected to a sudden demographic expansion in the past , and this was consistent with the hypothesis of Zhou et al. (2013). Furthermore, across the three habitats of giant pandas examined, the prevalence and intensity of B. schroederi infections is relatively high (50–100%) ,. No commercial vaccine is available against B. schroederi, and its control relies largely on anthelmintic drugs. However, some studies have reported the occurrence of drug resistance in B. schroederi populations in China . Given high levels of gene flow among B. schroederi populations and the potential of rare resistance alleles to spread, an investigation of drug resistance in B. schroederi is urgently needed.
This investigation initially combined two different mitochondrial markers, atp6 and cox1, to explore population genetic variability and structure of 57 B. schroederi isolates sampled from the three main giant panda habitats in China. The low levels of genetic diversity and high levels of gene flow indicated that B. schroederi was not genetically differentiated and had experienced lineage re-arrangement. Moreover, the negative Fu’s Fs and Tajima’s D values coupled with mismatch analysis consistently showed that the B. schroederi population could have been subjected to a sudden demographic expansion in the past. These results, together with previous studies, further suggest that B. schroederi populations did not follow a pattern of isolation by distance, revealing the possible existence of physical connections before the populations became geographically separated. This finding should contribute to a deeper understanding of population genetics and evolutionary biology of B. schroederi, a nematode parasite of the giant panda.
- atp6 :
ATPase subunit 6
- cox1 :
Cytochrome c oxidase subunit I
- cox3 :
Cytochrome c oxidase subunit III
- nad4 :
Nicotinamide dehydrogenase subunit IV
Analysis of molecular variance
Ling DI, Zwerling AA, Pai M: GenoType MTBDR assays for the diagnosis of multidrug-resistant tuberculosis: a meta-analysis. Eur Respir J. 2008, 32: 1165-1174. 10.1183/09031936.00061808.
Small PM, Hopewell PC, Singh SP, Paz A, Parsonnet J, Ruston DC, Schecter GF, Daley CL, Schoolnik GK: The epidemiology of tuberculosis in San Francisco–a population-based study using conventional and molecular methods. N Engl J Med. 1994, 330: 1703-1709. 10.1056/NEJM199406163302402.
Brudey K, Driscoll JR, Rigouts L, Prodinger WM, Gori A, Al-Hajoj SA, Allix C, Aristimuño L, Arora J, Baumanis V, Binder L, Cafrune P, Cataldi A, Cheong S, Diel R, Ellermeier C, Evans JT, Fauville-Dufaux M, Ferdinand S, Garcia de Viedma D, Garzelli C, Gazzola L, Gomes HM, Guttierez MC, Hawkey PM, van Helden PD, Kadival GV, Kreiswirth BN, Kremer K, Kubin M, et al: Mycobacterium tuberculosis complex genetic diversity: mining the fourth international spoligotyping database (SpolDB4) for classification, population genetics and epidemiology. BMC Microbiol 2006, 6:23.,
Gasser RB, Jabbar A, Mohandas N, Schnyder M, Deplazes P, Littlewood DTJ, Jex AR: Mitochondrial genome of Angiostrongylus vasorum: Comparison with congeners and implications for studying the population genetics and epidemiology of this parasite. Infect Genet Evol. 2012, 12: 1884-1891. 10.1016/j.meegid.2012.07.022.
Mas-Coma S, Valero MA, Bargues MD: Fasciola, lymnaeids and human fascioliasis, with a global overview on disease transmission, epidemiology, evolutionary genetics, molecular epidemiology and control. Adv Parasitol. 2009, 69: 41-146. 10.1016/S0065-308X(09)69002-3.
Betson M, Nejsum P, Bendall RP, Deb RM, Stothard JR: Molecular epidemiology of ascariasis: a global perspective on the transmission dynamics of Ascaris in people and pigs. J Infect Dis. 2014, 210: 932-941. 10.1093/infdis/jiu193.
Barrett LG, Thrall PH, Burdon JJ, Linde CC: Life history determines genetic structure and evolutionary potential of host-parasite interactions. Trends Ecol Evol. 2008, 23: 678-685. 10.1016/j.tree.2008.06.017.
Blouin MS, Yowell CA, Courtney CH, Dame JB: Host movement and the genetic structure of populations of parasitic nematodes. Genetics. 1995, 141: 1007-1014.
Frankham R, Briscoe D, Ballou J: Introduction to Conservation Genetics. 2002, Cambridge University Press, Cambridge
Simon C, Frati F, Beckenbach A, Crespi B, Liu H, Flook P: Evolution, weighting, and phylogenetic utility of mitochondrial gene sequences and a compilation of conserved polymerase chain reaction primers. Ann Entomol Soc Am. 1994, 87: 651-701.
Sunnucks P: Efficient genetic markers for population biology. Trends Ecol Evol. 2000, 15: 199-203. 10.1016/S0169-5347(00)01825-5.
Farias IP, Ortí G, Sampaio I, Schneider H, Meyer A: The cytochrome b gene as a phylogenetic marker: the limits of resolution for analyzing relationships among cichlid fishes. J Mol Evol. 2001, 53: 89-103. 10.1007/s002390010197.
Franssen F, Xie K, Sprong H, van der Giessen J: Molecular analysis of Baylisascaris columnaris revealed mitochondrial and nuclear polymorphisms. Parasit Vectors 2013, 6:124.,
Höglund J, Morrison D, Mattsson J, Engström A: Population genetics of the bovine/cattle lungworm (Dictyocaulus viviparus) based on mtDNA and AFLP marker techniques. Parasitology. 2006, 133: 89-99. 10.1017/S0031182006009991.
Troell K, Engström A, Morrison DA, Mattsson JG, Höglund J: Global patterns reveal strong population structure in Haemonchus contortus, a nematode parasite of domesticated ruminants. Int J Parasitol. 2006, 36: 1305-1316. 10.1016/j.ijpara.2006.06.015.
Yin F, Gasser RB, Li F, Bao M, Huang W, Zou F, Zhao G, Wang C, Yang X, Zhou Y, Zhao J, Fang R, Hu M: Genetic variability within and among Haemonchus contortus isolates from goats and sheep in China. Parasit Vectors 2013, 6:279.,
Iñiguez AM, Leles D, Jaeger LH, Carvalho-Costa FA, Araújo A, Amazonas Research Group: Genetic characterisation and molecular epidemiology of Ascaris spp. from humans and pigs in Brazil. Trans R Soc Trop Med Hyg. 2012, 106: 604-612. 10.1016/j.trstmh.2012.06.009.
Zhou X, Xie Y, Zhang ZH, Wang CD, Sun Y, Gu XB, Wang SX, Peng XR, Yang GY: Analysis of the genetic diversity of the nematode parasite Baylisascaris schroederi from wild giant pandas in different mountain ranges in China. Parasit Vectors 2013, 6:233.,
Zhang ZH, Wei FW: Giant Panda: Ex-Situ Conservation Theory and Practice. 2006, Science Press, Beijing
State Forestry Administration: The Third National Survey Report on Giant Panda in China. 2006, Chinese Science Press House, Beijing
Xie Y, Zhang Z, Wang C, Lan J, Li Y, Chen Z, Fu Y, Nie H, Yan N, Gu X, Wang S, Peng X, Yang G: Complete mitochondrial genomes of Baylisascaris schroederi, Baylisascaris ailuri and Baylisascaris transfuga from giant panda, red panda and polar bear. Gene. 2011, 482: 59-67. 10.1016/j.gene.2011.05.004.
Yang GY, Zhang ZH: Parasitic Diseases of Wildlife. 2013, Science Press, Beijing
Loeffler K, Montali R, Rideout B: Diseases and pathology of giant pandas. Giant Pandas: Biology, Veterinary Medicine and Management. Edited by: Wildt DE, Zhang AJ, Zhang HM, Janssen DL, Ellis S. 2006, Cambridge University Press, Cambridge, 377-409. 10.1017/CBO9780511542244.017.
Wang T, He G, Yang G, Fei Y, Zhang Z, Wang C, Yang Z, Lan J, Luo L, Liu L: Cloning, expression and evaluation of the efficacy of a recombinant Baylisascaris schroederi Bs-Ag3 antigen in mice. Vaccine. 2008, 26: 6919-6924. 10.1016/j.vaccine.2008.09.079.
Lin Q, Li H, Gao M, Wang X, Ren W, Cong M, Tan X, Chen C, Yu S, Zhao G: Characterization of Baylisascaris schroederi from Qinling subspecies of giant panda in China by the first internal transcribed spacer (ITS-1) of nuclear ribosomal DNA. Parasitol Res. 2012, 110: 1297-1303. 10.1007/s00436-011-2618-7.
McIntosh A: A new nematode, Ascaris schroederi, from a giant panda, Ailuropoda melanoleuca. Zoologica. 1939, 24: 355-357.
Sprent JF: Notes on Ascaris and Toxascaris, with a definition of Baylisascaris gen. nov. Parasitology. 1968, 58: 185-198. 10.1017/S0031182000073534.
Xie Y, Chen S, Yan Y, Zhang Z, Li D, Yu H, Wang C, Nong X, Zhou X, Gu X, Wang S, Peng X, Yang G: Potential of recombinant inorganic pyrophosphatase antigen as a new vaccine candidate against Baylisascaris schroederi in mice. Vet Res 2013, 44:90.,
Zhao GH, Li HM, Ryan UM, Cong MM, Hu B, Gao M, Ren WX, Wang XY, Zhang SP, Lin Q, Zhu XQ, Yu SK: Phylogenetic study of Baylisascaris schroederi isolated from Qinling subspecies of giant panda in China based on combined nuclear 5.8S and the second internal transcribed spacer (ITS-2) ribosomal DNA sequences. Parasitol Int. 2012, 61: 497-500. 10.1016/j.parint.2012.02.009.
Zhao GH, Xu MJ, Zhu XQ: Identification and characterization of microRNAs in Baylisascaris schroederi of the giant panda. Parasit Vectors 2013, 6:216.,
Hu M, Gasser RB: Mitochondrial genomes of parasitic nematodes–progress and perspectives. Trends Parasitol. 2006, 22: 78-84. 10.1016/j.pt.2005.12.003.
Jex AR, Littlewood DTJ, Gasser RB: Toward next-generation sequencing of mitochondrial genomes–focus on parasitic worms of animals and biotechnological implications. Biotechnol Adv. 2010, 28: 151-159. 10.1016/j.biotechadv.2009.11.002.
Wang Y, Liu GH, Li JY, Xu MJ, Ye YG, Zhou DH, Song HQ, Lin RQ, Zhu XQ: Genetic variability among Trichuris ovis isolates from different hosts in Guangdong Province, China revealed by sequences of three mitochondrial genes. Mitochondrial DNA. 2013, 24: 50-54. 10.3109/19401736.2012.710210.
Thompson JD, Higgins DG, Gibson TJ: CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res. 1994, 22: 4673-4680. 10.1093/nar/22.22.4673.
Castresana J: Selection of conserved blocks from multiple alignments for their use in phylogenetic analysis. Mol Biol Evol. 2000, 17: 540-552. 10.1093/oxfordjournals.molbev.a026334.
Falush D, Stephens M, Pritchard JK: Inference of population structure using multilocus genotype data: linked loci and correlated allele frequencies. Genetics. 2003, 164: 1567-1587.
Posada D: jModelTest: Phylogenetic model averaging. Mol Biol Evol. 2008, 25: 1253-1256. 10.1093/molbev/msn083.
Xia X, Lemey P: Assessing substitution saturation with DAMBE. Phylogen Handbook. 2009, 2009: 615-630. 10.1017/CBO9780511819049.022.
Librado P, Rozas J: DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics. 2009, 25: 1451-1452. 10.1093/bioinformatics/btp187.
Bandelt HJ, Forster P, Röhl A: Median-joining networks for inferring intraspecific phylogenies. Mol Biol Evol. 1999, 16: 37-48. 10.1093/oxfordjournals.molbev.a026036.
Excoffier L, Laval G, Schneider S: Arlequin (version 3.0): an integrated software package for population genetics data analysis. Evol Bioinform Online. 2007, 1: 47-50.
Nei M: Estimation of average heterozygosity and genetic distance from a small number of individuals. Genetics. 1978, 89: 583-590.
Wright S: The genetical structure of populations. Ann Eugen. 1951, 15: 323-354. 10.1111/j.1469-1809.1949.tb02451.x.
Hedrick PW: A standardized genetic differentiation measure. Evolution. 2005, 59: 1633-1638. 10.1111/j.0014-3820.2005.tb01814.x.
Fu YX: Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics. 1997, 147: 915-925.
Schneider S, Excoffier L: Estimation of past demographic parameters from the distribution of pairwise differences when the mutation rates vary among sites: application to human mitochondrial DNA. Genetics. 1999, 152: 1079-1089.
Tajima F: Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989, 123: 585-595.
Poulin R: Network analysis shining light on parasite ecology and diversity. Trends Parasitol. 2010, 26: 492-498. 10.1016/j.pt.2010.05.008.
Archie EA, Ezenwa VO: Population genetic structure and history of a generalist parasite infecting multiple sympatric host species. Int J Parasitol. 2011, 41: 89-98. 10.1016/j.ijpara.2010.07.014.
Höglund J, Morrison DA, Engström A, Nejsum P, Jansson DS: Population genetic structure of Ascaridia galli re-emerging in non-caged laying hens. Parasit Vectors 2012, 5:97.,
Tydén E, Morrison D, Engström A, Nielsen M, Eydal M, Höglund J: Population genetics of Parascaris equorum based on DNA fingerprinting. Infect Genet Evol. 2013, 13: 236-241. 10.1016/j.meegid.2012.09.022.
Hu M, Chilton NB, Gasser RB: The mitochondrial genomics of parasitic nematodes of socio-economic importance: recent progress, and implications for population genetics and systematics. Adv Parasitol. 2004, 56: 133-212. 10.1016/S0065-308X(03)56003-1.
Ngui R, Mahdy MA, Chua KH, Traub R, Lim YA: Genetic characterization of the partial mitochondrial cytochrome oxidase c subunit I (cox 1) gene of the zoonotic parasitic nematode, Ancylostoma ceylanicum from humans, dogs and cats. Acta Trop. 2013, 128: 154-157. 10.1016/j.actatropica.2013.06.003.
Blouin M: Mitochondrial DNA diversity in nematodes. J Helminthol. 1998, 72: 285-289. 10.1017/S0022149X00016618.
Zhao S, Zheng P, Dong S, Zhan X, Wu Q, Guo X, Hu Y, He W, Zhang S, Fan W: Whole-genome sequencing of giant pandas provides insights into demographic history and local adaptation. Nat Genet. 2013, 45: 67-71. 10.1038/ng.2494.
Uddin W: Host Parasite Co-evolution: Case of nematode Parasites (Syphacia obvelata, Trichuris muris) and Host (Mus musculus) in Mouse Hybrid Zone. In: Current methods on speciation research: population genetic approaches, Průhonice, 2012.
Cutter AD, Payseur BA: Genomic signatures of selection at linked sites: unifying the disparity among species. Nat Rev Genet. 2013, 14: 262-274. 10.1038/nrg3425.
Anderson JL, Morran LT, Phillips PC: Outcrossing and the maintenance of males within C. elegans populations. J Hered. 2010, 101: S62-S74. 10.1093/jhered/esq003.
Wildt DE, Zhang AJ, Zhang HM, Janssen DL, Ellis S: Giant Pandas: Biology, Veterinary Medicine and Management. 2006, Cambridge University Press, Cambridge
Nejsum P, Frydenberg J, Roepstorff A, Parker ED: Population structure in Ascaris suum (Nematoda) among domestic swine in Denmark as measured by whole genome DNA fingerprinting. Hereditas. 2005, 142: 7-14. 10.1111/j.1601-5223.2005.01864.x.
Zhou C, Yuan K, Tang X, Hu N, Peng W: Molecular genetic evidence for polyandry in Ascaris suum. Parasitol Res. 2011, 108: 703-708. 10.1007/s00436-010-2116-3.
Posada D, Crandall KA: Intraspecific gene genealogies: trees grafting into networks. Trends Ecol Evol. 2001, 16: 37-45. 10.1016/S0169-5347(00)02026-7.
Anderson RM: Complex dynamic behaviours in the interaction between parasite populations and the host’s immune system. Int J Parasitol. 1998, 28: 551-566. 10.1016/S0020-7519(97)00207-5.
Louhi KR, Karvonen A, Rellstab C, Jokela J: Is the population genetic structure of complex life cycle parasites determined by the geographic range of the most motile host?. Infect Genet Evol. 2010, 10: 1271-1277. 10.1016/j.meegid.2010.08.013.
Bowen BW, Grant WS: Phylogeography of the sardines (Sardinops spp): assessing biogeographic models and population histories in temperate upwelling zones. Evolution. 1997, 51: 1601-1610. 10.2307/2411212.
This work was supported by grants from the Science & Technology Ministry, China (No. 200910188) and the Research Fund for the Chengdu Giant Panda Breeding (No. CPF-2012-13). We would like to thank Zuqing Chen and Yu Wang (College of Veterinary Medicine, Sichuan Agricultural University, China) for their technical assistance; Chengdong Wang and his staff (Ya’an Bifengxia Research Base of Giant Panda Breeding, China) for materials.
The authors declare that they have no competing interests.
XY and YGY conceived the concept of the project. ZZH, WCD, GXB, WT, PXR and YGY collected samples. XY, ZX and SY performed the lab work. XY and ZX carried out the statistical analyses. XY and ZX wrote the initial draft of the manuscript. YGY and XY revised the manuscript. All authors read and approved the final manuscript.
Yue Xie, Xuan Zhou contributed equally to this work.