EDP Sciences logo
Submit your paper
Search
Menu
All issues Volume 30 (2023) Parasite, 30 (2023) 59 Full HTML
Open Access
Issue
Parasite
Volume 30, 2023
Article Number 59
Number of page(s) 11
DOI https://doi.org/10.1051/parasite/2023061
Published online 12 December 2023

© H.-Y. Sun et al., published by EDP Sciences, 2023

Licence Creative CommonsThis is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Introduction

Schistosomiasis is a parasitic disease caused by trematode flukes of the genus Schistosoma and is estimated to affect at least 230 million people worldwide. Schistosoma japonicum (digenean trematode) is mainly endemic in China and partly in the Philippines [19, 52]. Schistosomiasis in humans in mainland China is caused by S. japonicum only [2]. After nearly 70 years of schistosomiasis control in China, great achievements have been obtained, with infections in both humans and livestock reduced to a much lower level across the country [16]. As a consequence, the central government of China has set the goal to eliminate this disease by 2030 [3]. However, many challenges remain, one of which may be attributable to its zoonotic transmission, as over 40 species of wild and domestic animals can serve as reservoir hosts for S. japonicum [19].

Schistosoma japonicum has a complex life cycle, involving an obligatory molluscan intermediate host of the genus Oncomelania and a mammalian definitive host. In fresh water, miracidia hatch from eggs, which are produced by sexual reproduction in definitive hosts, and actively penetrate a snail within which asexual reproduction occurs. The infected snail releases free-swimming cercariae that infect humans or reservoir hosts via skin penetration [6]. Snail-to-mammal transmission depends on temporal coincidence between the cercariae and the vertebrate because vertebrates are generally in contact with the water for a short length of time, and cercariae have a short lifespan. Cercarial infectivity declines rapidly after a few hours post-emergence from a snail, whereas their mortality increases [47]. Therefore, the chronobiological characteristics of cercarial emission naturally play a paramount role in schistosome transmission.

Cercarial emergence rhythms are regarded as adaptive activities of parasites to maximize the probability of encountering the right host [7, 28, 31, 41], and moreover, are probably genetically determined [35, 42]. Several studies have been carried out to show the intraspecific polymorphisms of cercarial emergence patterns [28, 33, 39, 41]. For example, on the island of Guadeloupe, S. mansoni cercariae with an early shedding pattern were observed in an urbanized area, where humans were highly infected; whereas cercariae with a late shedding pattern were observed in a sylvatic site, where rats were heavily infected [41]. In Anhui (China), S. japonicum cercariae from a hilly region displayed a late afternoon emergence pattern consistent with a nocturnal rodent reservoir, while those from a marshland region showed an early pattern compatible with a diurnal cattle reservoir [28, 39]. The variations of cercarial emergence between schistosome isolates were considered to be related to different final hosts. A number of studies, based on COX1 or microsatellites [28, 32, 36, 43], have attempted to elucidate any genetic differences between parasites with respect to their biological traits, definitive hosts, or possible habitat types. Evidence from morphological to molecular level indicated that S. japonicum in mainland China may comprise a strain complex, but the low- and middle-valley (including Anhui province) of the Yangtze River were endemic with the same S. japonicum strain [18]. Recent advances in whole genome sequencing technology (WGS) make it feasible to study the genetic architecture of schistosomes with different traits [17, 21], as genome-wide genetic variation can provide a more comprehensive and unbiased framework to identify genomic regions associated with the phenotypic traits.

We have previously reported that in Shitai (ST), a hilly/mountain area (HM) of Anhui, S. japonicum cercariae had a unique late afternoon emergence pattern, and the definitive host-parasite system was driven by local rodents in spite of high infections in humans [27, 28, 37]. Moreover, our recent research suggested that there was an increased prevalence of S. japonicum infections in rodents in HM regions across China [53]. We therefore hypothesized that the specific late-afternoon emergence pattern of S. japonicum was the result of ongoing selection posed by definitive host rodents. As there are no studies to explain the underlining genetics of this kind of adaptive traits with WGS, in this study we sequenced two S. japonicum isolates of Anhui, each with a different cercarial shedding pattern, and then analyzed their genetic relationships with the purpose of detecting any changes in nucleotide diversity and allele frequencies in response to the unique chronobiological characteristics.

Materials and methods

Material and sequencing

Schistosoma japonicum were sampled from two sites in Anhui province, China. One site was in Shitai county (ST), a hilly area where S. japonicum was observed with a late afternoon cercarial shedding pattern and rodents were the main reservoir; another was in Hexian (HX), a lake/marshland area (LM) where S. japonicum was observed with a diurnal shedding of cercariae emergence [28] and bovines were the main reservoir. Details on the field samples used in this study are given in Table 1 and Figure 1. Briefly, intermediate host Oncomelania hupensis hupensis snails [9] were captured via routine field surveys in 2013 and 2020, and isolation of infected snails was performed by releasing cercariae in the laboratory. Infected snails were identified in Shitai in 2013 and 2020 and in Hexian in 2013 only. Cercariae released from infected snails of each site/year were separately and subsequently used to establish their corresponding S. japonicum isolates in mice, in which mice were individually exposed to cercariae from a single infected snail. Four weeks later, adult worms were perfused from infected mice and frozen in absolute ethanol, as a future genetic resource. The Ethics Committee of Soochow University reviewed and approved this research protocol (No. 81971957). Ten to 15 male or female adult worms, either batch of which derived from single infected snails, were selected from each site/year (i.e., ST2020m and ST2020f, for males and females from ST in 2020, respectively; ST2013m and ST2013f, for males and females from ST in 2013, respectively; HX2013m and HX2013f, for males and females from HX in 2013, respectively), and were then sent to Novogene (Tianjin, China) for whole genome sequencing. Pooled DNA was extracted using the DNeasy Blood & Tissue Kit (QIAGEN, Hilden, Germany), according to the manufacturer’s instructions. Paired-end reads (2 × 150 bp) for each sample were produced on a NovaSeq 6000 system (Illumina, San Diego, CA, USA). The data that support the findings of this study are openly available in the BioProject database at https://dataview.ncbi.nlm.nih.gov/, reference No. PRJNA793934.

thumbnail Figure 1

Map of geographical locations of research sites.
Table 1

Geographical information of S. japonicum samples.

Read mapping and variant calling

Raw reads were filtered with FASTP [4] with the parameters -q 5 -u 50 -n 15, and aligned to the reference genome (ASM636876v1) [29] with BWA-MEM (version 0.7.17) [24] with default parameters. The aligned bam files were then sorted using SAMtools [25], and PCR duplicate reads were marked with Picard-MarkDuplicates (http://broadinstitute.github.io/picard/). Mapping statistics were estimated with SAMtools. The Genome Analysis Toolkit (GATK) [10, 30] was used for variant calling. By running HaplotypeCaller, genomic variant call format (GVCF) files were generated for each sample. These GVCF files were next combined into a single GVCF file, which was then used to identify single nucleotide polymorphisms (SNPs) and small indels with GenotypeGVCFs. SNPs and small indels were further filtered using the following criteria: (1) SNPs were filtered with QD < 2.0, QUAL < 30.0, MQ < 40.0, FS > 45.0, SOR > 3.0, MQRankSum < −12.5, and ReadPosRankSum < −8.0, and indels with QD < 2.0, QUAL < 30.0, MQ < 40.0, FS > 200.0 and SOR > 10.0, and ReadPosRankSum < −20.0; (2) only variants with two alleles were retained and missing data prohibited by using VCFtools (version 0.1.13) [8]; (3) linkage disequilibrium pruning was performed with PLINK (version 1.9) (https://zzz.bwh.harvard.edu/plink/ld.shtml) using a 1-kb window with a step size of 10 SNPs and an r 2 threshold of 0.1; (4) variants with a minor allele frequency (MAF) of <0.05 were removed. All identified variants were annotated with Snpeff (version 5.0e) [5]. The VCF files including samples of each site/year were obtained by using VCFtools (version 0.1.13), and annotated with Snpeff (version 5.0e). SNPs and indels were categorized based on their positions on the genome (including intergenic regions, exons, introns, splicing sites, untranslated regions, transcripts, and 1-kb upstream and downstream regions) and on their effects (i.e., missense, nonsense, and silent mutations). Singleton SNPs in each sample were also calculated with VCFtools.

Genetic differentiation

Principal component analysis (PCA) was performed on the filtered SNP sets with PLINK [38]. For each site/year, nucleotide diversity (π) was calculated within a 100-kb window sliding in 10-kb steps, and Tajima’s D was calculated within a nonoverlapping 100-kb window with VCFtools [22]. Pairwise genetic differentiation (FST per gene between sampling sites/years) was calculated within a 100-kb window sliding in 10-kb steps using VCFtools [22].

Identification of selected regions in genomes

We conducted tests to detect selection signals along the genomes of ST and HX. The π ratios (πST/πHX) and FST values were calculated by the sliding-window approach (a 100-kb window sliding in a 10-kb step). To avoid spurious selection signals, windows containing fewer than 10 informative sites from both π and FST analyses were discarded. We used an empirical procedure [26] and selected the windows simultaneously with significantly low or high π ratios (the 5% left or right tails of the distribution) and significantly high FST values (the 5% right tail of the distribution). These windows were regarded as the regions with strong selective sweep signals, which should harbor candidate genes that undergo a selective sweep. Candidate genes of selective regions were identified using annotations of the reference [29] with Snpeff (version 5.0e) [5]. In order to understand what functions these genes were involved in, we tested for enrichment of gene ontology (GO) terms by comparing to the reference genome with the clusterProfiler package [49]. The Benjamini-Hochberg correction was applied and significantly enriched GO terms were identified at a corrected p-value of <0.05.

Impact of sampling time on selection

As ST included two different years (i.e., 2013 and 2020) of samples, we further analyzed allele frequency changes for every SNP on the comparisons of ST2013 vs. ST2020 and HX2013 vs. ST2013 by performing Fisher’s exact tests (FET) and Cochran-Mantel-Haenszel tests (CMH) with popoolation2 [22]. For FET test, the six sorted and marked BAM files were merged into three groups according to sampling site/year using SAMtools [25]. For CMH test, the six files were kept separate and formed pairs as following: 2 × (ST2013 vs. ST2020) and 2 × (HX2013 vs. ST2013). Pile-up files were generated with SAMtools’ mpileup (-d 500 –min-MQ 30 –min-BQ 30 –adjust-MQ 50) using the BAM files as input. The pile-up files were converted into synchronized files (popoolation2 mpileup2sync.pl –min-qual 20), and sequences around indels +5 bp were identified (popoolation2 identify-indel-regions.pl –min-count 2 –indel-window 5) and removed (popoolation2 filter-sync-by-gtf.pl). Significance with a genome-wide Bonferroni’s correction (0.05/(number of genome-wide SNPs)) was used in both tests. Pairwise genetic differentiation (FST per gene) was calculated using popoolation2 (popoolation2 fst-sliding.pl –pool-size 1000 –min-count 4 –min-coverage 20 –max-coverage 2% –window-size 1000000 –step-size 1000000) [22], with the annotated synchronized file (popoolation2 create-genewise-sync.pl) as an input (popoolation2 create-gene-wise-sync.pl).

Results

Sequencing and mapping

Six schistosome samples (i.e., ST2013m, ST2013f, ST2020m, ST2020f, HX2013m, and HX2013f) generated a total of 84.49 Gb of raw sequencing data. After filtering sequences of low quality, 0.55 billion reads were left with an average mapping depth of 35× on the reference genome. The mapping rate was between 95.70% and 96.27%, and the median coverage ranged from 93.58% to 95.18% (Table 2).

Table 2

Summary of mapping statistics.

Variant identification

Single nucleotide variants varied among samples, mostly between two schistosome isolates. Heterozygous SNPs ranged from the lowest 2,093,493 in ST2013m to the highest 4,335,470 in HX2013m, while homozygous SNPs were in range of the lowest 727,770 in HX2013m to the highest 4,335,470 in ST2013m. The Ts/Tv (transitions/transversions) varied between 1.807 in ST2013f and ST2020m and 1.837 in HX2013f. The average Ts/Tv (transitions/transversions) of variation was 1.82. Singleton SNPs ranged from the lowest 196,368 in ST2013m to the highest 445,959 in HX2013f. Nonsynonymous SNPs varied between 39982 in HX2013m and 81449 in ST2020f (see Table 3).

Table 3

Statistics of SNPs of six samples.

Genetic differentiation

PCA analysis classified schistosome samples according to their sampling site and year, in which the percentage of the variation in genetic distances explained by each PC is 46.9 and 34.9% for PC1 and PC2, respectively. When samples were grouped based on sampling site and year (i.e., ST2013m and ST2013f into ST2013, ST2020m and ST2020f into ST2020, and HX2013m and HX2013f into HX2013), the pairwise genetic differentiation (FST) between HX2013 and ST2013 or ST2020 were 0.55 and 0.54, respectively, showing a great genetic difference; whereas the index between ST2013 and ST2020 was −0.01, indicating almost no difference between both. The values of nucleotide diversity (π) of ST2013 and ST2020 was almost the same, and so for the Tajima’s D. Moreover, both π and Tajima’s D were much lower when compared to HX2013 (see Table 4 and Fig. 2). We then further grouped samples based on sampling sites, namely ST and HX. The lower nucleotide diversity (π) was in ST (5.61 × 10−3) and higher in HX (8.48 × 10−3) (Fig. 3A). Selective neutrality, tested with Tajima’s D, showed that ST had more windows observed at extremely high or low values. The overall value of Tajima’s D computed along genomes in ST (0.29) was lower than that in HX (0.67) (Fig. 3B). Significant differences for π and Tajima’s D value distributions between both sites were observed (p-value < 0.05; two-sample Kolmogorov-Smirnov test). FST revealed a high degree of genetic differentiation between two sites (mean genome-wide FST = 0.56) (Fig. 3C).

thumbnail Figure 2

Box plots of diversity, Tajima’s D of S. japonicum and their FST. (A) Nucleotide diversity estimated in 100-kb windows sliding in 10-kb steps throughout the genome. (B) Tajima’s D estimated within a nonoverlapping 100-kb window throughout the genome. (C) pairwise FST computed in 100-kb windows sliding in 10-kb steps throughout the genome.
thumbnail Figure 3

Box plots of diversity and relationship of two sample groups. (A) Nucleotide diversity estimated in 100-kb windows sliding in 10-kb steps throughout the genome. (B) Tajima’s D estimated within a nonoverlapping 100-kb window throughout the genome. (C) FST between computed in 100-kb windows sliding in 10-kb steps throughout the genome.
Table 4

π and FST of three groups.

Genome-wide selective sweep analysis

A sliding-window approach was applied to identify selection signals on genomes. Out of 36,389 windows examined 424 were discarded as they contained fewer than 10 informative sites from both π value and FST analysis. Strong selective sweep signals were detected in the remaining windows. The genome regions in ST affected by natural selection were fewer than in HX. In ST, strong signals were detected within 64 windows containing 78 genes; whereas in HX within 318 windows containing 276 genes (Fig. 4; Supplemental Tables S1 and S2). ST had lower levels of polymorphism than HX (median πST/πHX = 0.658), and selection regions in each site showed higher levels of selection statistics (FST and Tajima’s D) (Supplemental Fig. S1). The results suggested that selection may exist in shaping the genome, resulting in changes in phenotypic and/or behavioral traits of schistosomes. However, no significantly enriched genes were identified in the selection regions in HX, but in ST about 12 significantly enriched genes were identified (Supplemental Table S3). The functions of those genes involve positive regulation of cell migration, cis-regulatory region sequence-specific DNA binding, positive regulation of transcription, DNA-templated, positive regulation of epithelial to mesenchymal transition, sequence-specific DNA binding and RNA polymerase Ⅱ transcription factor binding (Fig. 5). These enriched genes showed a high level of genetic differentiation (average FST = 0.91) and were mainly related to specific DNA binding and intercellular signaling transduction (Table 5).

thumbnail Figure 4

Distribution of π ratios (πST/πHX) and FST values calculated in 100-kb windows sliding in 10-kb steps throughout the genome. Data points colored red and blue were identified as selected regions in ST (red dots) and in HX (blue dots), respectively. These points correspond to the 5% left and right tails of the empirical π ratio distribution, where the π ratios are 0.092 and 1.643, respectively (vertical dashed lines), and the 5% right tail of the empirical FST distribution, where FST is 0.965 (horizontal dashed line).
thumbnail Figure 5

Gene ontology enrichment analysis of genes from genome regions with strong selective signals.
Table 5

Information of enriched genes of S. japonicum from ST.

Impact of sampling time on selection

To investigate whether the two different years of ST samples affect the above selective analysis, we further performed both FET and CMH tests to examine SNP differences on comparisons of ST2013 vs. HX2013 and ST2013 vs. ST2020 (Supplemental Fig. S2). For HX2013 vs. ST2013, FET yielded roughly 2.36 million SNPs in two merged groups, whereas CMH test contained 0.37 million SNPs throughout the comparisons. A high degree of concordance between FET and CMH (Pearson’s correlation; r = 0.98) was identified when comparing the SNPs shared between the two tests (n = 0.36 million). As for ST2013 vs. ST2020, about 1.36 million SNPs were used in FET and 0.92 million in CMH test. Correlation analysis between FET and CMH tests revealed a high level of concordance (n = 0.19 million; Pearson’s correlation; r = 0.96). Both FET and CMH tests between HX2013 and ST2013 included higher number of significant SNPs than between ST2013 and ST2020. Moreover, the SNPs of the enriched genes (red dots in Supplemental Fig. S2) showed differences between HX2013 and ST2013, but none between ST2013 and ST2020. Except for the SNPs of five genes (EWB00_003259, EWB00_008464, EWB00_008467, EWB00_010839, and EWB00_010840) which showed significant changes above the genome-wide correction in the FET but not in the CMH test, the SNPs of all previously found selective genes showed significant changes in both FET and CMH tests. The pairwise FST of these genes showed a high level of difference between ST2013 and HX2013, while none between ST2013 and ST2020 (Supplemental Fig. S3).

Discussion

With the available genome maps of S. haematobium [23], S. mansoni [24], and S. japonicum [25, 26], whole genome sequencing was used to reveal the population history of S. mansoni [27], to investigate introgression events of S. bovis and S. haematobium [28], and to evaluate the village-level relatedness and genetic diversity of S. japonicum [29]. In this work, using whole-genome pooled sequencing (pool-seq), we sequenced and analyzed genomes of six S. japonicum samples from two sites, each endemic for S. japonicum with a different cercarial emergence pattern. More than 11.87 million SNPs were identified, and two S. japonicum isolates showed great genetic differentiation and contrasting diversities. Although the genome region “windows” affected by natural selection were fewer in ST than in HX, more involved genes were identified in ST. The selective regions and relative genes on the schistosome genome might be related to its unique chronobiological characteristics of cercarial emergence.

In this study, the results from whole genome sequences of S. japonicum showed significant genetic differentiation between two isolates. Rudge et al. [36] once reported strong genetic differentiation in S. japonicum between two schistosome isolates and suggested that this differentiation may be associated with the contrasting host reservoirs. Indeed, our previous work in Anhui identified that rodents may be the key host for maintaining S. japonicum transmission within HM regions, while bovines are key hosts only within LM regions [27, 28, 37]. As cercarial emergence is a heritable trait, shaped by the behavior of definitive hosts [7, 41, 42], in the view of evolution, S. japonicum cercariae from the HM region with a late afternoon emergence pattern may be the result of adaptation to nocturnal rodent reservoirs, while those from the LM with an early pattern may be the result of adaptation to diurnal cattle reservoirs [28, 39]. The genetic differentiation observed between the two schistosome isolates in this work may possibly underline their different biological traits, caused by definitive hosts.

Zhao et al. [51] demonstrated that, based on mitochondrial DNA, a much larger degree of genetic variation was observed in low-lying LM regions along the middle and lower reaches of the Yangtze River than in the HM regions in the upper reaches of the River. Similarly, in our work, we observed higher genetic diversity in HX than in ST. The greater genetic diversity and higher positive Tajima’s D in HX could result from balancing selection for the parasites, for example at the stage of asexual reproduction within intermediate host snails, as in the marshland annual floods may have brought about the widespread mixing and dispersal of snails across large geographical areas [12, 48]. However, the positive values of Tajima’s D in both HX and ST indicate that some regions in the genomes of the two S. japonicum isolates may be under selective pressure by definitive hosts or other factors.

More selective signals were detected on schistosome genomes from HX than from ST, suggesting that HX might receive more selective pressure. This may be explained by the fact that bovines and humans are the key hosts in LM regions [37], and as a major source of infection, have been given high priority for control interventions [45]; whereas rodents are the key hosts in HM regions and are not able to be targeted due to logistic difficulties [23]. Significantly enriched genes were identified in selection regions of ST only, suggesting that they might be related to the fact that the schistosomes may have long resided and been transmitted within local rodents. The number of enriched genes were up to 12 and their products are mainly related to transcription and intercellular signal transduction [11, 13, 15, 40, 44], such as zinc finger protein [11], transcriptional enhancer factor, mothers against decapentaplegic member 3 isoform 1, and mothers against decapentaplegic member 3 isoform 2 (SMAD) [40, 46]. The signaling pathways involving SMAD have components sharing high identity with mammalian orthologues. High expression of SMAD may, as identified in the liver of infected mice [50], imply that schistosomes can, in addition to utilizing their own pathways, exploit host growth factors as developmental signals [40]. We thus proposed that schistosome from ST may strengthen its transcription activities to increase the opportunity to invade the right host and prepare to exploit key signaling pathways of their host for growth and metabolism. We further found that one of these enriched genes is responsible for the product of histone-lysine N-methyltransferase. This is of interest as histone modifying enzymes have been selected as targets for druggable epigenetic targets in early-stage schistosome drug discovery projects [34]. In addition, the involved products also included suppressor of hairless protein (transcription factor partner of Notch receptor), which plays a possible role in Notch signaling in sensory architecture [44]. Regarding sensory structures, Hoffmann et al. [20] found that some components of the light detection system, for example S. mansoni rhodopsin (SmRHO), were consistent with the responsiveness of cercariae to light. Their findings were reminiscent of cercarial emergence, a biological trait controlled to some extent by illumination [1, 14]. More studies would be needed to explore the role of Notch signaling in sensory architecture on S. japonicum in the hilly area.

In conclusion, we found that the hilly/mountainous S. japonicum isolate with late afternoon emergence was genetically different from the LM isolate with a morning shedding pattern, although both sites are within the same province and with O. h. hupensis as intermediate hosts. Moreover, we identified certain functional region changes along the genome of the hilly isolate that might be related to the unique late shedding pattern of schistosome cercariae.

Acknowledgments

This study was funded by the National Science Foundation of China (to DL, No. 81971957). The authors report that they have no competing interests.

Supplementary materials

Table S1: Genes of select regions in S. japonicum from ST.

Table S2: Genes of select regions in S. japonicum from HX.

Table S3: GO enrichment results in S. japonicum from ST.

Figure S1: Box plots of (A) π ratio (πST/πHX), (B) FST, and (C) Tajima’s D of selection regions throughout genomes.

Figure S2: Fisher’s exact test (A, B) and CMH test (C, D) were employed to investigate any significant SNP frequency changes (−log10[p-value]) between the two groups. −log10(P) values for common SNPs obtained through performing FET and CMH tests were plotted together (on separate axes) in order to identify significant allele frequency changes of variants in both tests (E, F). Tests between HX2013 and ST2013 are shown in (A, C, and E), while tests between ST2013 and ST2020 are shown in (B, D, and F). The red lines represent genome-wide Bonferroni’s correction of p-value. Dots in red represent SNPs of enriched genes above the genome-wide Bonferroni’s correction level.

Figure S3: FST values for entire genes calculated between HX2013 and ST2013 (A), between ST2013 and ST2020 (B), and between ST and HX (C). Dots in red represent SNPs of enriched genes.

Access here

References

  1. Asch HL. 1972. Rhythmic emergence of Schistosoma mansoni cercariae from Biomphalaria glabrata: control by illumination. Experimental Parasitology, 31(3), 350–355. [CrossRef] [PubMed] [Google Scholar]
  2. Attwood SW, Upatham ES, Meng XH, Qiu DC, Southgate VR. 2002. The phylogeography of Asian Schistosoma (Trematoda: Schistosomatidae). Parasitology, 125(Pt 2), 99–112. [CrossRef] [PubMed] [Google Scholar]
  3. Brattig NW, Bergquist R, Qian MB, Zhou XN, Utzinger J. 2020. Helminthiases in the People’s Republic of China: Status and prospects. Acta Tropica, 212, 105670. [CrossRef] [PubMed] [Google Scholar]
  4. Chen S, Zhou Y, Chen Y, Gu J. 2018. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics, 34(17), i884–i890. [CrossRef] [PubMed] [Google Scholar]
  5. Cingolani P, Platts A, le Wang L, Coon M, Nguyen T, Wang L, Land SJ, Lu X, Ruden DM. 2012. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome of Drosophila melanogaster strain w1118; iso-2; iso-3. Fly (Austin), 6(2), 80–92. [CrossRef] [PubMed] [Google Scholar]
  6. Colley DG, Bustinduy AL, Secor WE, King CH. 2014. Human schistosomiasis. Lancet, 383(9936), 2253–2264. [CrossRef] [PubMed] [Google Scholar]
  7. Combes C, Fournier A, Moné H, Théron A. 1994. Behaviours in trematode cercariae that enhance parasite transmission: patterns and processes. Parasitology, 109(Suppl), S3–S13. [CrossRef] [PubMed] [Google Scholar]
  8. Danecek P, Auton A, Abecasis G, Albers CA, Banks E, DePristo MA, Handsaker RE, Lunter G, Marth GT, Sherry ST, McVean G, Durbin R. 2011. The variant call format and VCFtools. Bioinformatics, 27(15), 2156–2158. [CrossRef] [PubMed] [Google Scholar]
  9. Davis GM, Wu WP, Xu XJ. 2006. Ecogenetics of shell sculpture in Oncomelania (Gastropoda) in canals of Hubei, China, and relevance for schistosome transmission. Malacologia, 48(1–2), 253–264. [Google Scholar]
  10. DePristo MA, Banks E, Poplin R, Garimella KV, Maguire JR, Hartl C, Philippakis AA, del Angel G, Rivas MA, Hanna M, McKenna A, Fennell TJ, Kernytsky AM, Sivachenko AY, Cibulskis K, Gabriel SB, Altshuler D, Daly MJ. 2011. A framework for variation discovery and genotyping using next-generation DNA sequencing data. Nature Genetics, 43(5), 491–498. [CrossRef] [PubMed] [Google Scholar]
  11. Drummond MG, Calzavara-Silva CE, D’Astolfo DS, Cardoso FC, Rajão MA, Mourão MM, Gava E, Oliveira SC, Macedo AM, Machado CR, Pena SD, Kitten GT, Franco GR. 2009. Molecular characterization of the Schistosoma mansoni zinc finger protein SmZF1 as a transcription factor. PLoS Neglected Tropical Diseases, 3(11), e547. [CrossRef] [PubMed] [Google Scholar]
  12. Ebert D, Fields PD. 2020. Host-parasite co-evolution and its genomic signature. Nature Reviews Genetics, 21(12), 754–768. [CrossRef] [PubMed] [Google Scholar]
  13. Eleutério de Souza PR, Valadão AF, Calzavara-Silva CE, Franco GR, de Morais MA, Jr., Abath FG. 2001. Cloning and characterization of SmZF1, a gene encoding a Schistosoma mansoni zinc finger protein. Mémorias do Instituto Oswaldo Cruz, 96(Suppl), 123–130. [Google Scholar]
  14. Favre T, Bogéa T, Rotenberg L, Silva H, Pieri O. 1997. Circadian rhythms in the cercarial emergence of Schistosoma mansoni by Biomphalaria tenagophila at outdoors: a comparative study with Biomphalaria glabrata. Biological Rhythm Research, 28(3), 348–357. [CrossRef] [Google Scholar]
  15. Gu JL, Chen SX, Dou TH, Xu MJ, Xu JX, Zhang L, Hu W, Wang SY, Zhou Y. 2012. Hox genes from the parasitic flatworm Schistosoma japonicum . Genomics, 99(1), 59–65. [CrossRef] [PubMed] [Google Scholar]
  16. Guo JY, Xu J, Zhang LJ, Lv S, Cao CL, Li SZ, Zhou XN. 2020. Surveillance on schistosomiasis in five provincial-level administrative divisions of the People’s Republic of China in the post-elimination era. Infectious Diseases of Poverty, 9(1), 136. [CrossRef] [PubMed] [Google Scholar]
  17. Han ZG, Brindley PJ, Wang SY, Chen Z. 2009. Schistosoma genomics: new perspectives on schistosome biology and host-parasite interaction. Annual Review of Genomics and Human Genetics, 10, 211–240. [CrossRef] [PubMed] [Google Scholar]
  18. He Y-X. 1993. Studies on the strain differences of Schistosoma japonicum in the mainland of China XII. Conclusion. Chinese Journal of Parasitology and Parasitic Diseases, 11(2), 93–97. [Google Scholar]
  19. He YX, Salafsky B, Ramaswamy K. 2001. Host–parasite relationships of Schistosoma japonicum in mammalian hosts. Trends in Parasitology, 17(7), 320–324. [CrossRef] [PubMed] [Google Scholar]
  20. Hoffmann KF, Davis EM, Fischer ER, Wynn TA. 2001. The guanine protein coupled receptor rhodopsin is developmentally regulated in the free-living stages of Schistosoma mansoni. Molecular and Biochemical Parasitology, 112(1), 113–123. [CrossRef] [PubMed] [Google Scholar]
  21. Koboldt DC, Steinberg KM, Larson DE, Wilson RK, Mardis ER. 2013. The next-generation sequencing revolution and its impact on genomics. Cell, 155(1), 27–38. [CrossRef] [PubMed] [Google Scholar]
  22. Kofler R, Pandey RV, Schlötterer C. 2011. PoPoolation2: identifying differentiation between populations using sequencing of pooled DNA samples (Pool-Seq). Bioinformatics, 27(24), 3435–3436. [CrossRef] [PubMed] [Google Scholar]
  23. Jiao-Jiao L. 2019. Endemic status and control of animal schistosomiasis in China. Chinese Journal of Schistosomiasis Control, 31(1), 40–46. [Google Scholar]
  24. Li H, Durbin R. 2009. Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics, 25(14), 1754–1760. [CrossRef] [PubMed] [Google Scholar]
  25. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R. 2009. The sequence alignment/map format and SAMtools. Bioinformatics, 25(16), 2078–2079. [CrossRef] [PubMed] [Google Scholar]
  26. Li M, Tian S, Jin L, Zhou G, Li Y, Zhang Y, Wang T, Yeung CK, Chen L, Ma J, Zhang J, Jiang A, Li J, Zhou C, Zhang J, Liu Y, Sun X, Zhao H, Niu Z, Lou P, Xian L, Shen X, Liu S, Zhang S, Zhang M, Zhu L, Shuai S, Bai L, Tang G, Liu H, Jiang Y, Mai M, Xiao J, Wang X, Zhou Q, Wang Z, Stothard P, Xue M, Gao X, Luo Z, Gu Y, Zhu H, Hu X, Zhao Y, Plastow GS, Wang J, Jiang Z, Li K, Li N, Li X, Li R. 2013. Genomic analyses identify distinct patterns of selection in domesticated pigs and Tibetan wild boars. Nature Genetics, 45(12), 1431–1438. [CrossRef] [PubMed] [Google Scholar]
  27. Lu DB, Rudge JW, Wang TP, Donnelly CA, Fang GR, Webster JP. 2010. Transmission of Schistosoma japonicum in marshland and hilly regions of China: parasite population genetic and sibship structure. PLoS Neglected Tropical Diseases, 4(8), e781. [CrossRef] [PubMed] [Google Scholar]
  28. Lu DB, Wang TP, Rudge JW, Donnelly CA, Fang GR, Webster JP. 2009. Evolution in a multi-host parasite: chronobiological circadian rhythm and population genetics of Schistosoma japonicum cercariae indicates contrasting definitive host reservoirs by habitat. International Journal for Parasitology, 39(14), 1581–1588. [CrossRef] [PubMed] [Google Scholar]
  29. Luo F, Yin M, Mo X, Sun C, Wu Q, Zhu B, Xiang M, Wang J, Wang Y, Li J, Zhang T, Xu B, Zheng H, Feng Z, Hu W. 2019. An improved genome assembly of the fluke Schistosoma japonicum . PLoS Neglected Tropical Diseases, 13(8), e0007612. [CrossRef] [PubMed] [Google Scholar]
  30. McKenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, Garimella K, Altshuler D, Gabriel S, Daly M, DePristo MA. 2010. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Research, 20(9), 1297–1303. [CrossRef] [PubMed] [Google Scholar]
  31. Mouahid G, Idris MA, Verneau O, Théron A, Shaban MM, Moné H. 2012. A new chronotype of Schistosoma mansoni: adaptive significance. Tropical Medicine & International Health, 17(6), 727–732. [CrossRef] [PubMed] [Google Scholar]
  32. Mouahid G, Mintsa Nguema R, Al Mashikhi KM, Al Yafae SA, Idris MA, Moné H. 2019. Host-parasite life-histories of the diurnal vs. nocturnal chronotypes of Schistosoma mansoni: adaptive significance. Tropical Medicine & International Health, 24(6), 692–700. [CrossRef] [PubMed] [Google Scholar]
  33. N’Goran E, Brémond P, Sellin E, Sellin B, Théron A. 1997. Intraspecific diversity of Schistosoma haematobium in west Africa: chronobiology of cercarial emergence. Acta Tropica, 66(1), 35–44. [CrossRef] [PubMed] [Google Scholar]
  34. Padalino G, Ferla S, Brancale A, Chalmers IW, Hoffmann KF. 2018. Combining bioinformatics, cheminformatics, functional genomics and whole organism approaches for identifying epigenetic drug targets in Schistosoma mansoni . International Journal for Parasitology: Drugs and Drug Resistance, 8(3), 559–570. [CrossRef] [Google Scholar]
  35. Pages JR, Théron A. 1990. Analysis and comparison of cercarial emergence rhythms of Schistosoma haematobium, S. intercalatum, S. bovis, and their hybrid progeny. International Journal for Parasitology, 20(2), 193–197. [CrossRef] [PubMed] [Google Scholar]
  36. Rudge JW, Lu DB, Fang GR, Wang TP, Basáñez MG, Webster JP. 2009. Parasite genetic differentiation by habitat type and host species: molecular epidemiology of Schistosoma japonicum in hilly and marshland areas of Anhui Province. China. Molecular Ecology, 18(10), 2134–2147. [CrossRef] [PubMed] [Google Scholar]
  37. Rudge JW, Webster JP, Lu DB, Wang TP, Fang GR, Basáñez MG. 2013. Identifying host species driving transmission of schistosomiasis japonica, a multihost parasite system, in China. Proceedings of the National Academy of Sciences, 110(28), 11457–11462. [CrossRef] [PubMed] [Google Scholar]
  38. Slifer SH. 2018. PLINK: key functions for data analysis. Current Protocols in Human Genetics, 97(1), e59. [CrossRef] [PubMed] [Google Scholar]
  39. Su J, Zhou F, Lu DB. 2013. A circular analysis of chronobiology of Schistosoma japonicum cercarial emergence from hilly areas of Anhui, China. Experimental Parasitology, 135(2), 421–425. [CrossRef] [PubMed] [Google Scholar]
  40. The Schistosoma japonicum Genome Sequencing and Functional Analysis Consortium. 2009. The Schistosoma japonicum genome reveals features of host-parasite interplay. Nature, 460(7253), 345–351. [CrossRef] [PubMed] [Google Scholar]
  41. Théron A. 1984. Early and late shedding patterns of Schistosoma mansoni cercariae: ecological significance in transmission to human and murine hosts. Journal of Parasitology, 70(5), 652–655. [CrossRef] [Google Scholar]
  42. Théron A, Combes C. 1988. Genetic analysis of cercarial emergence rhythms of Schistosoma mansoni. Behavioral Genetics, 18(2), 201–209. [CrossRef] [PubMed] [Google Scholar]
  43. Théron A, Combes C. 1995. Asynchrony of infection timing, habitat preference, and sympatric speciation of schistosome parasites. Evolution, 49(2), 372–375. [CrossRef] [Google Scholar]
  44. Verjovski-Almeida S, DeMarco R, Martins EA, Guimarães PE, Ojopi EP, Paquola AC, Piazza JP, Nishiyama MY Jr, Kitajima JP, Adamson RE, Ashton PD, Bonaldo MF, Coulson PS, Dillon GP, Farias LP, Gregorio SP, Ho PL, Leite RA, Malaquias LC, Marques RC, Miyasato PA, Nascimento AL, Ohlweiler FP, Reis EM, Ribeiro MA, Sá RG, Stukart GC, Soares MB, Gargioni C, Kawano T, Rodrigues V, Madeira AM, Wilson RA, Menck CF, Setubal JC, Leite LC, Dias-Neto E. 2003. Transcriptome analysis of the acoelomate human parasite Schistosoma mansoni. Nature Genetics, 35(2), 148–157. [CrossRef] [PubMed] [Google Scholar]
  45. Wang LD, Chen HG, Guo JG, Zeng XJ, Hong XL, Xiong JJ, Wu XH, Wang XH, Wang LY, Xia G, Hao Y, Chin DP, Zhou XN. 2009. A strategy to control transmission of Schistosoma japonicum in China. New England Journal of Medicine, 360(2), 121–128. [CrossRef] [PubMed] [Google Scholar]
  46. Wang X, Xu X, Lu X, Zhang Y, Pan W. 2015. Transcriptome bioinformatical analysis of vertebrate stages of Schistosoma japonicum reveals alternative splicing events. PLoS One, 10(9), e0138470. [CrossRef] [PubMed] [Google Scholar]
  47. Whitfield PJ, Bartlett A, Khammo N, Clothier RH. 2003. Age-dependent survival and infectivity of Schistosoma mansoni cercariae. Parasitology, 127(Pt 1), 29–35. [CrossRef] [PubMed] [Google Scholar]
  48. Wilke T, Davis GM, Cui EC, Xiao-Nung Z, Xiao Peng Z, Yi Z, Spolsky CM. 2000. Oncomelania hupensis (Gastropoda: rissooidea) in eastern China: molecular phylogeny, population structure, and ecology. Acta Tropica, 77(2), 215–227. [CrossRef] [PubMed] [Google Scholar]
  49. Yu G, Wang LG, Han Y, He QY. 2012. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics, 16(5), 284–287. [CrossRef] [PubMed] [Google Scholar]
  50. Zhang BB, Cai WM, Tao J, Zheng M, Liu RH. 2013. Expression of Smad proteins in the process of liver fibrosis in mice infected with Schistosoma japonicum. Chinese Journal of Parasitology & Parasitic Disease, 31(2), 89–94. [Google Scholar]
  51. Zhao QP, Jiang MS, Dong HF, Nie P. 2012. Diversification of Schistosoma japonicum in Mainland China revealed by mitochondrial DNA. PLoS Neglected Tropical Diseases, 6(2), e1503. [CrossRef] [PubMed] [Google Scholar]
  52. Zhou XN, Guo JG, Wu XH, Jiang QW, Zheng J, Dang H, Wang XH, Xu J, Zhu HQ, Wu GL, Li YS, Xu XJ, Chen HG, Wang TP, Zhu YC, Qiu DC, Dong XQ, Zhao GM, Zhang SJ, Zhao NQ, Xia G, Wang LY, Zhang SQ, Lin DD, Chen MG, Hao Y. 2007. Epidemiology of schistosomiasis in the People’s Republic of China, 2004. Emerging Infectious Diseases, 13(10), 1470–1476. [CrossRef] [PubMed] [Google Scholar]
  53. Zou HY, Yu QF, Qiu C, Webster JP, Lu DB. 2020. Meta-analyses of Schistosoma japonicum infections in wild rodents across China over time indicates a potential challenge to the 2030 elimination targets. PLoS Neglected Tropical Diseases, 14(9), e0008652. [CrossRef] [PubMed] [Google Scholar]

Cite this article as: Sun H-Y, Zhang J-Y, Zhang H-X, Xu Q & Lu D-B. 2023. Genetic difference between two Schistosoma japonicum isolates with contrasting cercarial shedding patterns revealed by whole genome sequencing. Parasite 30, 59.

All Tables

Table 1

Geographical information of S. japonicum samples.

In the text
Table 2

Summary of mapping statistics.

In the text
Table 3

Statistics of SNPs of six samples.

In the text
Table 4

π and FST of three groups.

In the text
Table 5

Information of enriched genes of S. japonicum from ST.

In the text

All Figures

thumbnail Figure 1

Map of geographical locations of research sites.
In the text
thumbnail Figure 2

Box plots of diversity, Tajima’s D of S. japonicum and their FST. (A) Nucleotide diversity estimated in 100-kb windows sliding in 10-kb steps throughout the genome. (B) Tajima’s D estimated within a nonoverlapping 100-kb window throughout the genome. (C) pairwise FST computed in 100-kb windows sliding in 10-kb steps throughout the genome.
In the text
thumbnail Figure 3

Box plots of diversity and relationship of two sample groups. (A) Nucleotide diversity estimated in 100-kb windows sliding in 10-kb steps throughout the genome. (B) Tajima’s D estimated within a nonoverlapping 100-kb window throughout the genome. (C) FST between computed in 100-kb windows sliding in 10-kb steps throughout the genome.
In the text
thumbnail Figure 4

Distribution of π ratios (πST/πHX) and FST values calculated in 100-kb windows sliding in 10-kb steps throughout the genome. Data points colored red and blue were identified as selected regions in ST (red dots) and in HX (blue dots), respectively. These points correspond to the 5% left and right tails of the empirical π ratio distribution, where the π ratios are 0.092 and 1.643, respectively (vertical dashed lines), and the 5% right tail of the empirical FST distribution, where FST is 0.965 (horizontal dashed line).
In the text
thumbnail Figure 5

Gene ontology enrichment analysis of genes from genome regions with strong selective signals.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.