Introduction
Mitochondrial DNA (mtDNA) sequences have been used extensively for studies of population or ecological genetics, and also for phylogenetic and evolutionary analyses at various taxonomic levels of different organisms, including nematodes, trematodes, cestodes and protozoans (M.W. Li et al., Reference Li, Lin, Song, Sani, Wu and Zhu2008; Zhao et al., Reference Zhao, Mo, Zou, Weng, Lin, Xia and Zhu2009; Shekhovtsov et al., Reference Shekhovtsov, Katokhin, Kolchanov and Mordvinov2010; Lin et al., Reference Lin, Qiu, Liu, Wu, Weng, Xie, Hou, Pan, Yuan, Zou, Hu and Zhu2011; Liu et al., Reference Liu, Lin, Li, Liu, Liu, Yuan, Song, Zhao, Zhang and Zhu2011; Dai et al., Reference Dai, Liu, Song, Lin, Yuan, Li, Huang, Liu and Zhu2012), due to the maternal inheritance, rapid evolutionary rate and lack of recombination of mtDNA. However, there is limited information on the genetic variation in populations of some important parasite groups from China, such as the liver fluke Clonorchis sinensis infecting animals and humans.
Clonorchis sinensis is an important liver fluke parasite of the family Opisthorchiidae and it parasitizes the bile duct of many mammalian hosts, including humans. Clonorchiasis caused by C. sinensis is endemic in South-East Asia, including China, Vietnam, Japan and Korea (Lun et al., Reference Lun, Gasser, Lai, Li, Zhu, Yu and Fang2005). It is estimated that approximately 35 million people globally are infected with C. sinensis, of whom about 15 million are in China (Lun et al., Reference Lun, Gasser, Lai, Li, Zhu, Yu and Fang2005). Clonorchis sinensis has been considered to be carcinogenic to humans (Shin et al., Reference Shin, Oh, Masuyer, Curado, Bouvard, Fang, Wiangnon, Sripa and Hong2010).
Genetic variation is common in parasite populations, and it is a valuable resource for studying population biology, epidemiology and genetic structure of parasites (Gasser & Newton, Reference Gasser and Newton2000; Li et al., Reference Li, Lin, Song, Sani, Wu and Zhu2008; Cerutti et al., Reference Cerutti, Citterio, Bazzocchi, Epis, D'Amelio, Ferrari and Lanfranchi2010). Recently, sequences of nuclear ribosomal DNA (rDNA) and mtDNA have been used to study genetic variations in helminths of human and animal health significance (M.W. Li et al., Reference Li, Lin, Song, Sani, Wu and Zhu2008; Zhao et al., Reference Zhao, Mo, Zou, Weng, Lin, Xia and Zhu2009; J. Li et al., Reference Li, Zhao, Zou, Mo, Yuan, Ai, Li, Weng, Lin and Zhu2010; Ai et al., Reference Ai, Chen, Alasaad, Elsheikha, Li, Li, Lin, Zou, Zhu and Chen2011). Although intra- and inter-specific variations have recently been studied using nuclear rDNA and mtDNA sequences among C. sinensis from different geographic regions (Park, Reference Park2007), there is a paucity of information on genetic variation among populations of C. sinensis from China, where C. sinensis infection remains an important human and animal health problem.
The objectives of the present study were to examine sequence variability in four mtDNA regions, namely cytochrome c oxidase subunits 1 (cox1) and 2 (cox2), and NADH dehydrogenase subunits 1 and 2 (nad1 and nad2), among C. sinensis isolates from different endemic regions in China. Based on the combined sequences of these four mtDNA regions, phylogenetic relationships of C. sinensis in China with other zoonotic trematodes were also reconstructed.
Materials and methods
Parasites and isolation of genomic DNA
All adult worms of C. sinensis were obtained from naturally infected cats and dogs from different geographic locations in China. These C. sinensis were washed in physiological saline, identified preliminarily to species based on host preference, morphological characters and predilection sites (Kaewkes, Reference Kaewkes2003), fixed in 70% (v/v) ethanol and stored at − 20°C until use. Their sample codes, gender and GenBank accession numbers are listed in table 1.
Table 1 Geographic origins of Clonorchis sinensis samples used in the present study, as well as their GenBank accession numbers for partial sequences of mitochondrial cox1 and cox2 genes (pcox1 and pcox2), and nad1 and nad2 genes (pnad1 and pnad2).
Enzymatic amplification and sequencing
A portion of the cox1, cox2, nad1 and nad2 were amplified by polymerase chain reaction (PCR) with primers listed in table 2. PCR reactions (25 ml) were performed in 10 mm Tris–HCl (pH 8.4), 50 mm KCl, 4 mm MgCl2, 200 mm each of deoxynucleoside triphosphate (dNTP), 50 pmol of each primer and 2 U Taq polymerase (Takara, Dalian, China) and 1 μl of DNA sample in a thermocycler (Biometra, Göttingen, Germany) under the following conditions: after an initial denaturation at 94°C for 5 min, then 94°C for 30 s (denaturation); 50°C (for pcox1) or 52°C (for pcox2) or 58°C (for pnad1) or 61°C (for pnad2) for 30 s (annealing); 72°C for 30 s (extension) for 35 cycles, followed by a final extension at 72°C for 10 min. These optimized cycling conditions for the specific and efficient amplification of individual mtDNA fragments were obtained after varying annealing temperatures. Samples without genomic DNA (no-DNA controls) and host genomic DNA (host-DNA controls) were included in each amplification run, and in no case were amplicons detected in the no-DNA and host-DNA controls (not shown). Each amplicon (5 μl) was examined by agarose gel (1%) electrophoresis to validate amplification efficiency. PCR products were sent to Sangon Company (Shanghai, China) for sequencing from both directions.
Table 2 Sequences of primers used to amplify a portion of the mitochondrial cox1 and cox2 genes (pcox1 and pcox2), and nad1 and nad2 genes (pnad1 and pnad2) from Clonorchis sinensis isolates in China. These primers were designed on mt genomes of C. sinensis.
Sequence analysis and reconstruction of phylogenetic relationships
Sequences of the four mt genes were separately aligned using the computer program Clustal X 1.83 (Thompson et al., Reference Thompson, Gibson, Plewniak, Jeanmougin and Higgins1997). Pairwise comparisons were made of the level of sequence differences (D) among C. sinensis isolates using the formula D = 1 − (M/L), where M is the number of alignment positions at which the two sequences have a base in common, and L is the total number of alignment positions over which the two sequences are compared (Chilton et al., Reference Chilton, Gasser and Beveridge1995).
Representative samples whose four gene sequences were available in this study were used for phylogenetic analyses. Three inference methods, namely, Bayesian analysis (Bayes), neighbour joining (NJ) and maximum likelihood (ML), were used for phylogenetic reconstructions. Bayesian analyses were conducted with four independent Markov chains run for 1,000,000 metropolis-coupled MCMC generations, sampling a tree every 1000 generations in MrBayes 3.1.1 (Ronquist & Huelsenbeck, Reference Ronquist and Huelsenbeck2003). The first 250 trees were omitted as burn-in and the remaining trees were used to calculate Bayesian posterior probabilities (PP). NJ analyses were carried out using the PAUP 4.0 Beta 10 program (Swofford, Reference Swofford2002), and ML analyses were performed using PhyML 3.0 (Guindon & Gascuel, Reference Guindon and Gascuel2003), and the general time reversible (GTR) model with its parameter for the concatenated dataset was determined for the ML analysis using JModeltest (Posada, Reference Posada2008) based on the Akaike information criterion (AIC). Bootstrap support for ML trees was calculated using 100 bootstrap replicates. To study the phylogenetic relationships with other zoonotic trematode species, C. sinensis Russia isolate (FJ381664), C. sinensis Korea isolate (JF729304), C. sinensis China isolate (JF729303), Schistosoma japonicum (NC_002544), S. mansoni (NC_002545), S. haematobium (DQ157222), Opisthorchis felineus (EU921260), O. viverrini (JF739555), Fasciola hepatica (NC_002546) and Paragonimus westermani (NC_002354) were included in the present study, with Ascaris suum (NC_001327) as the outgroup. Phylograms were drawn using the Tree View program version 1.65 (Page, Reference Page1996).
Results and discussion
Genomic DNA was prepared from 31 individual adult C. sinensis isolates from three endemic provinces in mainland China. Amplicons of pcox1, pcox2, pnad1 and pnad2 (~420, 460, 480 and 580 bp, respectively) were amplified individually and subjected to agarose gel electrophoresis. For each mtDNA region, no size variation was detected on agarose gel among any of the amplicons examined (not shown).
To examine sequence variation in the four mtDNA regions among C. sinensis isolates from different endemic provinces, amplicons of pcox1, pcox2, pnad1 and pnad2 representing different endemic types were subjected to direct sequencing. The sequences of pcox1, pcox2, pnad1 and pnad2 were 381 bp, 420 bp, 432 bp and 531 bp in length, respectively. The A+T contents of the sequences were 57.74–58.53% (pcox1), 59.05–59.76% (pcox2), 60.65–61.34% (pnad1) and 63.16–63.74% (pnad2), respectively. The intra-specific sequence variations within populations of all C. sinensis isolates were 0–1.6% for pcox1, 0–1.4% for pcox2, 0–0.9% for pnad1 and 0–1.2% for pnad2. These results were consistent with that of previous studies for this parasite (Park, Reference Park2007; Lai et al., Reference Lai, Wang, Chen, Cai, Wu, Zhu and Lun2008). Similarly, genetic variability has also been detected among populations of the other members of the Opisthorchiidae, such as Opisthorchis viverrini and O. felineus, by analysis of mtDNA sequences (Katokhin et al., Reference Katokhin, Shekhovtsov, Konkow, Yurlova, Serbina, Vodianitskai, Fedorov, Loktev, Muratov, Ohyama, Makhnev, Pel'tek and Mordvinov2008; Saijuntha et al., Reference Saijuntha, Sithithaworn, Wongkham, Laha, Chilton, Petney, Barton and Andrews2008).
For pcox1, pcox2, pnad1 and pnad2, intra-specific nucleotide variation was related mainly to changes at the first and third codon positions, consistent with results of other helminths, including trematodes and cestodes (M.W. Li et al., Reference Li, Lin, Song, Sani, Wu and Zhu2008; Zhao et al., Reference Zhao, Mo, Zou, Weng, Lin, Xia and Zhu2009; Dai et al., Reference Dai, Liu, Song, Lin, Yuan, Li, Huang, Liu and Zhu2012). For pcox1 and pcox2, there was only one change at the second codon position. For pnad1, there was no change in the second codon (table 3). Intra-specific nucleotide variations represented transitions (A ↔ G or C ↔ T: n = 13 for pcox1, n = 15 for pcox2, n = 13 for pnad1 and n = 20 for pnad2) and transversions (A ↔ C, A ↔ T, and/or T ↔ G: n = 3 for pcox1, n = 2 for pcox2, n = 2 for pnad1, n = 7 for pnad2) (table 3).
Table 3 Number and position of codon variations in a portion of the mitochondrial cox1 and cox2 genes (pcox1 and pcox2), and nad1 and nad2 genes (pnad1 and pnad2) among Clonorchis sinensis isolates in China.
a The first codon position of each sequence was determined in relation to the complete mitochondrial DNA sequence of C. sinensis (Shekhovtsov et al., Reference Shekhovtsov, Katokhin, Kolchanov and Mordvinov2010).
The combined sequences of pcox1, pcox2, pnad1 and pnad2 were aligned over a consensus length of 1764 bp. Topologies of all trees inferred by different methods (NJ, MP and ML) with different building strategies and/or different distance models were identical or similar, with only the small difference of bootstrap values (fig. 1). In this tree, the Schistosomatidae, Fasciolidae, Paragonimidae and Opisthorchiidae form monophyletic groups, respectively. Within the clade consisting of the Fasciolidae, Paragonimidae and Opisthorchiidae, all the isolates of C. sinensis were sister to O. felineus. These results were consistent with that of a previous study using internal transcribed spacer (ITS) and cox1 sequences (Kang et al., Reference Kang, Sultana, Loktev, Wongratanacheewin, Sohn, Eom and Park2008), but were not consistent with that of a previous study (Cai et al., Reference Cai, Liu, Song, Wu, Zou, Yan, Yuan, Lin and Zhu2012) using concatenated amino acid sequences of the 12 protein-coding genes, which showed that C. sinensis was more closely related to O. viverrini than to O. felineus. Here, we will not discuss the phylogenetic relationships between O. viverrini and C. sinensis further, since phylogenetic relationships are not robustly supported. However, our results provided some useful information for future studies.
Fig. 1 Phylogenetic relationship of Clonorchis sinensis isolates from mainland China with other zoonotic trematodes inferred by Bayesian (Bayes), maximum likelihood (ML) and neighbour-joining (NJ) analyses using the combined dataset (cox1+ cox2+ nad1+ nad2), with Ascaris suum as outgroup. Bootstrap values (in percentage) above 50% from 1000 pseudo–replicates are shown for the neighbour-joining (the first value), maximum likelihood (the second value) and Bayesian analyses (the third value). Weak = node resolved by method but very weak ( < 50%).
In conclusion, genetic variations among C. sinensis isolates from different geographic regions in China were revealed by sequence analyses of mt cox1, cox2, nad1 and nad2, and demonstrated that mitochondrial DNA sequences provide reliable genetic markers for phylogenetic studies of zoonotic trematodes. For the four mt DNA genes, sequence variation was higher in pcox1 than in pcox2, pnad2 and pnad1.
Acknowledgements
This work was supported in part by the State Key Laboratory of Veterinary Etiological Biology, Lanzhou Veterinary Research Institute, Chinese Academy of Agricultural Sciences (Grant nos SKLVEB2009KFKT014, SKLVEB2010KFKT010 and SKLVEB2011KFKT004), the Yunnan Provincial Program for Introducing High-level Scientists (Grant no. 2009CI125) and the Beef and Yak Production System Programme, Ministry of Agriculture (CARS-38).
