Xiaocong Ma1†
Xiuying Yan2†
Ren Ke1Huiquan Shan2
Saif ur Rehman2Tong Feng2Yalin Zheng2Chen Chuang3Weiguan Zhou4
Qingyou Liu2,5*Jinghui Zheng1*- 1Department of Cardiology, Ruikang Hospital Affiliated to Guangxi University of Chinese Medicine, Nanning, China
- 2State Key Laboratory for Conservation and Utilization of Subtropical Agro-bioresources, Guangxi University, Nanning, China
- 3Guangxi Medical University Cancer Hospital, Nanning, China
- 4Thai Natural Hirudin Co, Ltd., Bangkok, Thailand
- 5Guangdong Provincial Key Laboratory of Animal Molecular Design and Precise Breeding, School of Life Science and Engineering, Foshan University, Foshan, China
Hirudo nipponia is the only blood-sucking leech included in Chinese Pharmacopoeia having distinct features of anticoagulation, exorcizing blood stasis, and promoting menstruation. Despite such significant characteristics, very little is known about its molecular genetics and related physiological mechanisms. In this study, the transcriptomes of H. nipponia at three developmental stages (larvae, young, and adults), revealed a total of 1,348 differentially expressed genes (DEGs), 223 differentially expressed lncRNAs, and 88 novel mRNAs. A significant diverse gene expression patterns were observed at different developmental stages which were analyzed by differential gene expression trends, and the overall gene expression trends consist of three overall down-regulated trends, and two overall up-regulated trends. Furthermore, the GO and KEGG enrichment functional annotation analysis revealed that these DEGs were mainly associated with protein hydrolysis, signal transduction, energy metabolism, and lipid metabolism while growth, development, metabolism, and reproduction-related DEGs were also found. Additionally, real-time quantitative PCR results confirmed deep sequencing results based on the relative expression levels of nine randomly selected genes. This is the first transcriptome-based comprehensive study of H. irudo nipponia at different developmental stages which provided considerable deep understanding related to gene expression patterns and their relevant developmental pathways, neurodevelopmental and reproductive characteristics of the leech.
Introduction
Hirudo nipponia (Annelida, Clitellata, and Hirudiniformes) is widely distributed in many regions of the globe (Liu et al., 2018). Physiologically, leeches are narrow, slightly flattened, cylindrical, and small-sized creatures used as raw material in Chinese traditional medicine with diverse functioning such as anticoagulation, exorcizing blood stasis, and promoting menstruation (Saglam, 2019). They are hermaphroditic where the male part matures first and then mates heterozygously, and oviparously (Ceylan, 2020). After hatching, leech seedlings moult once every 10 days, leeches grow for one age after one molt, 30–90 days is the optimum growth period for leeches, and it takes approximately 80 days to grow from seedling to adult size (Manav et al., 2019). However, H. nipponia has a prolonged maturation cycle with low reproduction ability as the H. nipponia begins to lay eggs in the second year. The average number of spawning cocoons of H. nipponia was the lowest compared with Whitmania pigra Whitman and Hirudinaria Poecilobdella manillensis (Ceylan, 2020).
H. nipponia as a traditional Chinese medicine contains hirudin which is the strongest thrombin natural inhibitor found so far (Straub et al., 2011; Liu et al., 2019) and it has unique properties of hypolipidemic, anti-tumor, anti-fibrotic, and anti-inflammatory. Furthermore, it has widely been used clinically to treat a variety of diseases such as hypertension, hyperlipidemia, cerebral infarction and malignancies (Sig et al., 2017; Pickrell et al., 2020; Wang et al., 2021). In addition, leeches have the ability to self-repair damage to the central nervous system and restore function, and are often used as model organisms in the field of neurobiology (Tasiemski and Salzet, 2017; Meiser et al., 2019; Heath-Heckman et al., 2021). In recent years, the number of wild leeches has gradually been decreased due to environmental destruction and increased hunting which results in high market demand but low supply. Thus, H. nipponia breeding at the farm level with selective traits is emerging everywhere all around the world (Chen et al., 2018; Yu et al., 2021).
Many studies have focused on the role of effective components in leech and their cultivation (Iwama et al., 2019; Wang et al., 2021). Recent advancements in high-throughput technologies provide an opportunity to dig out in-depth and comprehensive high quality leech genomes. Previously Liu et al. (2018) conducted a comparative transcriptomic study on three leech species including rod-shaped bull leech, broad-bodied golden leech and the cave mountain leech, and reported differentially expressed genes related to the olfactory transduction pathway among them. Moreover, another transcriptome-based study of the H. nipponia salivary gland identified a large number of anticoagulants, anti-inflammatory, antibacterial and anti-tumour proteins (Lu et al., 2018). But very little is reported related to the omics of H. nipponia, which is the major obstacle in the development and application of H. nipponia gene resources effectively. In this study, RNA-seq was used to sequence the transcriptome of H. nipponia at different developmental stages including larvae (L), young (Y), and adult (A) to provide a scientific basis for later research on the breeding, growth, and development of H. nipponia and related functional genes as well as the lncRNAs which are important for chromosomal modification, transcriptional activation, and interference.
Materials and Methods
Sample Collection
H. nipponia used in the experiment was collected from Qinzhou Zhanhong traditional Chinese medicine company. The samples were divided into a larval stage (within 10 days of hatching and have not sucked the blood), young insect stage (about 3 months after hatching and sucked blood three times) and adult stage (about 2 years after hatching, sucked blood more than six times and were adult and sexually mature) (Tan et al., 2002; Lu et al., 2018). The leeches were directly killed with scissors, then 24 samples were collected from the anterior and posterior segments from each of the larvae and young insects, along with adult oral sucker, gonad, tail sucker, and remaining insects. All the collected samples were placed in liquid nitrogen and stored at −80°C for further analysis. Meanwhile, a total of nine intact individuals, three each of larvae, young and adults, were taken for histological observation, fixed in 4% paraformaldehyde, and stored at 4°C. Three individuals from each group were taken as biological replicates.
Histomorphological Analysis
The histological features of H. nipponia were observed through paraffin embedding, sectioning, hematoxylin, and eosin staining. The muscles’ ultra-morphology was observed by microscope (EVOS FL auto), and the images were collected via micrograph system (EVOS auto).
cDNA Library Construction and Sequence Data Analysis
Total RNA was extracted from different parts of the H. nipponia at different developmental stages (L, Y, and A) using TRIzol reagent (Throme, United Sttaes). The RNA was purified using the RNeasy Mini Kit (Oiagen, Germany). The cDNA library was sent to Beijing Anuo Uda Gene Technology and sequenced using the Illumina HiSeq 4,000 platform (NEB, United Sttaes), which generated 150 bp paired-end reads. Each sample quality was checked using Fastp (v0.20.1) (Chen et al., 2018). After that HISAT2 (v2.2.1) (Kim et al., 2015) was used to map the reads from each sample to the reference genome and Samtools (v1.11) (Li et al., 2009) was employed to convert the sequence SAM files into BAM. StringTie (v2.1.4) (Pertea et al., 2015) was then applied to assemble and merge the transcripts from each sample. Gffcompare v0.11.2 software (Pertea and Pertea, 2020) was used to compare the merged transcripts to the reference annotation file in GTF, and StringTie was used to compare the transcripts to the reference annotationfile in GTF. StringTie was used to estimate transcript abundance using the option "-e-B-P20”. The Abundance results files were ending in ". balltown” and the prepDE.py was used to compare these files.
Analysis of Differentially Expressed Genes
The FPKM (Fragments Per Kilobase of exon model per Million mapped fragments) method was used to calculate the expression of all genes in each sample. The differentially expressed genes (DEGs) between different groups were analyzed using DESeq2 software. Genes with FDR (false discovery rate) < 0.05 and |log2FC|>1 was set as a criterion to find the DEGs between two stages and Ggplot2 software was used to locate the DEGs.
Sample Time Series Analysis
The Tcseq (time-course sequencing data analysis) software was used to analyze the DEGs at three-time points, and after clustering the patterns of DEGs were visualized at different developmental stages. The FDR value <0.05 and |log2FC|>1 was considered as the level of significance to determine the DEGs between the two stages.
Functional Enrichment Analysis of DEGs
The functional annotation analysis of DEGs was done by GO (gene ontology) and KEGG (Kyoto Encyclopedia of Genes and genomes) using a cluster profiler (statistical analysis and visualization of gene and gene cluster functional maps). The value of FDR <0.05 was considered to be significantly enriched.
Identification and Analysis of lncRNAs
Potential lncRNAs were screened by the following highly stringent criteria: 1) transcript length was not less than 200 nt; 2) transcripts were labeled as “i”, “u”, “324” and “x”; 3) the coding potential calculator (CPC2) score was “non-coding”; 4) transcripts containing an open reading frame (ORF) greater than 50 AA were removed; 5) transcripts were identified by alignment with the Swiss lncRNAs act as a cis-regulatory factor and can regulate the expression of adjacent genes. Furthermore, the co-expression networks of candidate lncRNAs and their upstream or downstream 10 kb mRNAs were constructed. Whereas, the ligation and enrichment were performed according to the position-frequency matrix.
Synthesis of cDNA and Real-Time qPCR
We have selected nine genes randomly to validate the expression of RNA-Seq and the genes related to the growth and development of the H. nipponia as suggested by Shan et al. (2022). The polyubiquitinase gene (ubiquitin, UBQ) (Robledo et al., 2014; Peña et al., 2010) was selected as the internal reference gene and sequenced together with the validation gene used the oligo7 designed (Supplementary Table S1), with primers synthesized by Sangon Biotech (Shanghai, China). The RNA was transformed into cDNA using the Hiscript III Reverse Transcription Kit (Vazyme, Nanjing, China), and fluorescent quantitative PCR was performed on an ABI 7500 (Thermo, MA, United Sttaes) using the AceQ qPCR SYBR Green Master Mix system (Vazyme, Nanjing, China). The relative quantitative 2−ΔΔCt method was used to analyze genes expression (Xu et al., 2021).
Results
Developmental Features and Histological Observations of H. nipponia at Different Developmental Stages
The formaldehyde-fixed H. nipponia at different developmental stages such as larvae (L), young (Y), and adults (A) are shown in Figure 1A. The calculated bodyweight of H. nipponia larvae (0.1362 ± 0.0075 g), young insects (1.184 ± 0.031 g), and adults (1.524 ± 0.045 g), were found significantly different in three developmental stages (p < 0.0001) as presented in Figure 1B. Additionally, the results of paraffin sections at different growth stages showed that there was a significantly higher ratio of muscle bundle tissue in the sections of adult worms as compared to larvae and young leech (Figure 1C).
FIGURE 1. Phenotypic data of H. nipponia. (A) Observation on the body shape of H. nipponia in different periods. (B) Weight statistics of H. nipponia in different periods. (C) Histological observation of H. nipponia in different periods.
RNA Sequencing of H. nipponia
We constructed a total of 24 cDNA libraries of H. nipponia at different growth stages (HNLB1, HNLB2, HNLB3, HNLF1, HNLF2, HNLF3, HNYB1, HNYB2, HNYB3, HNYF1, HNYF2, HNYF3, HNO1, HNO2, HNO3, HNG1, HNG2, HNG3 (HNP1, HNP2, HNP3, HNB1, HNB2, and HNB3) as shown in Supplementary Table S1. Three samples were taken from each of the larval, young, and adult stages (HNL1, HNL2, HNL3, HNY1, HNY2, HNY3, HNA1, HNA2, and HNA3). A total of 541,124,257 raw reads were generated from 24 libraries and after filtering the adapters, N ratios greater than 10%, basic reads, and low-quality reads, about 517805628 clean reads were obtained. The GC content of all the samples ranging between 41% and 45% (Supplementary Table S2). These results indicate that the sequencing data are of good quality and meet the requirements for subsequent analysis. The specific data are listed in Supplementary Table S2.
A total of 20,430 mRNAs were detected by RNA-seq, of which 88 were novel. The transcriptional expression of mRNAs was expressed as FPKM values and the FPKM distribution of all mRNAs is shown in Figure 2A, while the gene expression of the different samples is shown in a box plot (Figure 2B). To find the “main” elements and structures in the data, the complex sample composition relationship is reflected in the two eigenvalues of horizontal and vertical coordinates, which were helpful to explore the distance relationship between samples. Moreover, based on a permutation multivariate analysis of variance (PERMANOVA) we obtained p-value = 0.04 (p < 0.05) and then the nine samples were further divided into three sections, showing a tendency to cluster within the same group and a significant separation between the different groups (Figure 2C) and a relational clustering graph (Figure 2D) was constructed to show the correlation between each sample. We observed that the distribution of genes expression data for the samples we acquired was relatively consistent, and the PCA analysis results exhibited that the different subgroups can be more clearly distinguished from one another, and that intra-group correlations are high, providing reliable data for subsequent analysis.
FIGURE 2. The mRNA expression analysis. (A) The density distribution of mRNAs was according to log10 (FPKM). (B) The 9 Samples expression (L1, L2, L3, Y1, Y2, Y3, A1, A2, and A3) violin plot, which was replaced by log10 (FPKM). (C) The PCA distribution of 9 samples. (D) The sample relationship cluster analysis.
Differential Gene Expression Analysis
We used FDR value <0.05 and |log2FC| > 1 as criteria for screening DEGs by pairwise comparisons of three different growth stages of the H. nipponia. We obtained 238 DEGs in LvsY, 976 DEGs in LvsA, and 763 DEGs in YvsA (Figures 3A,B). We plotted the VEEN diagram for DEGs of these different growth periods (Figure 3C). In total, we obtained 1,348 DEGs, and the expression of these genes was shown in the heat map using cluster analysis (Figure 3D). The result showed that the adult stage had more ratio of DEGs and active expression than that of larvae and young, while on the other hand, the differential genes between larvae and young worms were fewer but more stable. Thus, it was observed that the expression of different genes in H. nipponia differed significantly at different developmental periods.
FIGURE 3. Visual analysis of expression differences between contrasts at different growth periods. (A) The barplot of the DEGs in three comparisons Volcano plot displaying DEGs between HNL vs. HNY, HNL vs. HNA, HNY vs. HNA stages. (B) Genetic statistics for differences between growth periods. (C) The ven plots of differential genes in different growth periods. (D) Heat map of genes in different growth periods.
Sample Time Series Analysis of DEGs
In order to fully reveal the growth and development of H. nipponia at different time points, we analyzed the expression trends of differential genes using FDR <0.05 and |log2FC|>1 as a level of significance. The overall gene expression trends included three overall down-regulated trends (cluster1, 5, and 7), two overall up-regulated trends (cluster 4 and 6), and three highest trends of DEG expression in HNY (cluster1, 2, and 8) (Figure 4). The changes in DEGs and their trend enrichment indicated significant differences in gene expression profiles at different developmental stages, effectively revealing the gene expression status of H. nipponia during development.
FIGURE 4. The sample time series analysis of DEGs. [Note: The horizontal axis represents different growth periods, and the vertical axis represents normalized gene expression].
GO and KEGG Enrichment Analysis of DEGs
To further unravel the functions of differentially expressed genes, GO analysis was performed and grouped the functional annotation into three categories: molecular function, cellular component, and biological process (Figure 5A). The DEGs were mainly enriched in 71 biological processes, 23 cellular components, and 130 molecular functions. The molecular functions of DEGs were found to be associated with ion transport, cell adhesion, proteolysis, and G protein-coupled receptor signaling pathway, while cellular components were plasma membrane and hemoglobin complex and the biological processes included protein dephosphorylation, homophilic cell adhesion via plasma membrane adhesion molecules, oxygen transport, response to wounding and cyclic nucleotide biosynthetic process. The KEGG enrichment is mainly related to metabolism, for example, amino acid metabolism such as Phenylalanine, tyrosine and tryptophan biosynthesis, and energy metabolism like glycerophospholipid metabolism, carbon sequestration in photosynthetic organisms, and carbon fixation pathways in prokaryotes as shown in Figure 5B. In addition, GO and KEGG profiles of differentially expressed genes were compared between the two H. nipponia at different periods as shown in Supplementary Figure S1, S2.
FIGURE 5. Differential gene enrichment analysis. (A) GO functional enrichment analysis differential genes. (B) Enrichment analysis of the KEGG pathway for differential genes.
Analysis of Differentially Expressed lncRNAs (DElncRNAs)
Furthermore, we compared two-by-two differences between the three different growth stages of H. nipponia for screening DelncRNAs by using FDR <0.05 and |log2FC|>1 as a criterion and in total, 223 DelncRNAs were obtained (Figure 6A) of these, 44 DelncRNAs were found in LvsY, 102 in LvsA, and 77 in YvsA (Figure 6B). The clustering analysis of these DelncRNAs is presented in a heat map (Figure 6C). Additionally, the GO enrichment analysis showed that DelncRNA target genes were mainly enriched in protein deubiquitination, ribosome, and metallo-carboxypeptidase activities (Figure 6D) whereas, the KEGG pathway analysis exhibited that the DelncRNA target genes were mainly correlated to protein secretion, bacterial secretion system, and secretion system (Figure 6E).
FIGURE 6. Screening and enrichment analysis of DelncRNAs in H. nipponia at different growth stages. (A) DelncRNA statistics at different growth stages. (B) Venn diagram of DElncRNAs. (C) DelncRNA clustering analysis. (D) KEGG analysis for revealing the role of predicted DelncRNA target genes. (E) GO enrichment plot showing GO analysis of predicted DelncRNA target genes.
The Expression of Important Genes Related to Growth and Development
Furthermore, the DEGs were evaluated and several genes associated with leech growth and development including muscle development-related genes such as promyosin, myosin heavy chain, myosin light chain, actin, etc. were determined which showed a progressive increase in gene expression as the leech developed (Figure 7A). However, the genes involved in gonadal development or regulation of egg production in leeches like surface protein (SPA17), IQ motif containing G (IQCG), β-actin (ACTB), 5-hydroxytryptamine (5-HT), and HOX gene family, etc. were highly expressed in the adult stage (Figure 7B). Similarly, the genes associated with neurodevelopment e.g., Tubulin Tyrosine Ligase Like 5 (TTLL5), Filamin-A (FLNA), Protocadherin Gamma Subfamily A10 (PCDHGA10), Solute Carrier Family 5 Member 8 (SLC5A8) and NK2 Homeobox 8 (NKX2-8), etc. were expressed in all growth period but were highly expressed in the larval stage (Figure 7C). On the other hand, in young worms, the anticoagulation-related genes Antistasin (ANTA), a disintegrin and metalloproteinase with thrombospondin motif family (ADAMTS), therostasin (THST), and decorsin (DECO) etc. were highly expressed (Figure 7D).
FIGURE 7. The expression of important genes related to growth and development. (A) Muscle development-related genes. (B) Reproduction-related genes. (C) Neurodevelopment-related genes. (D) Anticoagulation-related genes.
Confirmation of Gene Expression With Quantitative PCR
To verify the expression of genes from transcriptome sequencing data, we randomly selected 9 genes, including 5-Hydroxytryptamine Receptor 7 (5HT1R), ACTB, and UPF3B Regulator of Nonsense Mediated MRNA Decay (REN3B), etc. for RT-qPCR validation. The results show that the expression trends of the genes verified by RT-PCR are in line with the RNA-seq trends (Figure 8), and indicated that the obtained transcriptome sequencing data are highly reliable.
FIGURE 8. The validation of candidate genes. Red means q-PCR and blue means RNA-seq.
Discussion
In this study, we reported transcriptomic analyses from three different developmental stages (L, Y, and A) of H. nipponia, thus significantly expanding the genetic resources available for the leech. Our study also provides a new theoretical basis for studying the molecular genetics-based physiological mechanisms involved in the development of H. nipponia which ultimately benefits to improve the aquaculture privilege of H. nipponia. Based on growth characteristics, the leeches have categorized into three growth or developmental stages such as the larva that have not sucked blood for about 10 days after hatching, young worm included those that have sucked blood three times after hatching for about 2 months, and the adults sucked blood more than six times after 2 years of hatching and are sexually mature. Totally, 1,348 DEGs and 223 DelncRNAs were found in three growth periods, among them, the highest number of differential genes were identified in larval and adult leeches, while the lowest number of differential genes were figured out in young and adult leeches, which may be related to the fact that 3-month-old young leeches are similar in size to adult leeches (Tan et al., 2002; Lu et al., 2018). Furthermore, the GO and KEGG enrichment analyses showed that differential genes expression patterns were mainly associated with protein hydrolysis, enzyme activity, metabolism, and cell growth and development. As leeches grow rapidly during the young worm period, requiring large amounts of food intake and energy through synthesis, digestion, and metabolism of food, our data show that genes associated with digestion and metabolism were highly expressed in the young worm stage. In addition, genes encoding enzymes with protease activity (metalloproteases and thiol proteases) were found upregulated in young worms, consistent with studies showing that the non-reproductive stage (often referred to as the feeding stage) exhibits more aggressive feeding behavior than the adult stage (Zhao et al., 2021).
lncRNAs are long-stranded non-coding RNAs that are greater than 200 nt in length and do not have a coding role in living organisms. Although lncRNAs do not encode proteins but they have a crucial role in life activities. It has been demonstrated that lncRNAs are involved in X chromosome inactivation, chromosomal remodeling, transcriptional activation and interference as well as growth and development, metabolism and differentiation, and tumor production. Compared to other species, the study on leech lncRNAs is not well documented. In this study, during three developmental periods (L, Y, and A) a total of 223 differentially expressed lncRNA molecules were identified. Further analysis of these lncRNA molecules indicated that these lncRNA molecules are likely to play an important role in the growth and development of H. nipponia. Our study data provide new insights into the regulatory mechanisms of lncRNAs in the growth and development of H. nipponia.
Typically, actin, Ef1a, and GAPDH are used as internal reference genes for RT-qPCR testing, however we selected UBQ as internal reference gene. This is because β-actin or GAPDH are differentially expressed genes according to the results of our RNA-seq sequencing datat therefore we cannot use β-actin or GAPDH as internal reference genes for qRT-PCR analysis. Since there is lack of literature on the selection of internal reference genes in H. nipponia, we looked to other species and found that UBQ could be used as a regular internal reference gene (Robledo et al., 2014; Peña et al., 2010) and was not a differentially expressed gene in our RNA-seq data, so it was selected as an internal reference gene.
Previous studies have shown that the degree of muscular growth has changed significantly at different growth stages along with the underlying mechanisms regulating muscle growth might differ (Li et al., 2021). Our results indicated that the genes attributing muscular development were highly expressed and showed a tendency of upregulation as the leech body develops. The muscle-building block protein of annelids, paramyosin (MYSP), plays a key function in the muscular tension contraction process (Sonobe et al., 2016; Silva et al., 2020). In addition to performing basic contractile functions, proto-myosin (TPM) is the actin-binding protein in muscle filaments is also involved in many vital activities in the body, such as cytoplasmic division and signal transduction. It can bind helically to actin filaments and regulate interactions with myosin (Zhang et al., 2018; Jurewicz et al., 2020; Pavadai et al., 2020). Actin is an important member of the actin family that controls almost all forms of cell and organelle motility which also regulates growth, development, bioregulation, and repair (Hammer et al., 2019; Kadzik et al., 2020; Liu et al., 2020). Myosins consist of the myosin heavy chain (MYH) and the myosin light chain. The MYH family is one of the key members in the regulation of muscle development and affects muscle production and repair after injury, mainly through the regulation of cell proliferation and differentiation (Pan et al., 2016; Squire, 2019; Hill et al., 2021).
Leeches are typically hermaphroditic, with the male gonads maturing first than female gonads. Compared to other species of leeches, H. nipponia has low reproductive performance. Though analyzing the DEGs, we screened the genes or proteins associated with gonadal development or oogenesis. For example, SPA17, IQCG, ACTB, 5-HT, and the HOX gene family. The IQCG is also involved in spermatogenesis and is a key regulator of cilia/flagellar motility and is associated with sperm motility (Harris et al., 2014; Gonzalez-Gutierrez et al., 2021). Additionally, the studies have also reported that β-actin is potentially involved in sperm maturation, sperm excretion, and the ability to maintain sperm motility. Early studies on mammalian (guinea pig) sperm showed that actin was involved in sperm motility and random severance of actin filaments inhibited the flagella motility (Yenia, et al., 2007). Moreover, comparative proteomic study on normal and weak spermatozoa revealed that reduced expression of β-actin results in feeble spermatozoa (Martinez-Heredia et al., 2008). Similarly, reports on ovarian development have shown that 5-HT significantly play an important biological role in gonadal maturation and serial egg production (Kondo et al., 2017). The HOX gene is a family of regulatory molecules encoding highly conserved transcription factors whose expression is regulated by sex steroids, which are significantly involved in reproductive tract development, endometrial cycle growth, and embryo implantation (Du and Taylor, 2015; Cohen et al., 2020).
Furthermore, the leech is a fascinating model organism for researching brain regeneration mechanisms because of its capacity to regenerate the central nervous system (CNS) throughout its life cycle (Meriaux et al., 2011; Lu et al., 2018). Second, the leech’s CNS neurons are constantly expanding, implying that the process through which the leech develops and promotes synaptic connectivity may never be turned off. In our study, genes associated with neural development such as TTLL5, FLNA, PCDHGA10, SLC5A8, and NKX2-8 were expressed in all three developmental stages but were highly expressed in the larval stage. These genes are thought to have a major role in the formation of certain synaptic connections as well as the maintenance of neuronal energy status and function. They are regulated in neurons which re-establish axonal and synaptic connections. In addition to hirustasin, other anticoagulant genes including ANTA, ADAMTS, THST, and DECO were also expressed in the leech, and they were highly expressed in young leeches. Hirustasin and THST inhibited trypsin, chymotrypsin, histone G and tissue kinase-releasing enzyme, but not coagulation factor Xa activity (Kwak et al., 2019). In contrast, ANTA is a potent inhibitor of Factor Xa. Additionally, DECO inhibits the interaction of fibrinogen with platelet receptors expressed, on the glycoprotein IIb-IIIa complex. It prevents blood clotting when feeding and/or storing ingested blood (El-Obeid and Al-Badr, 1982). Adult leeches are frequently used as medicine in production applications and their anticoagulant features make it worthwhile to utilize leeches in medicine, like H. nipponia, which is commonly used for the treatment of various ailments (Koeppen et al., 2019; Singh et al., 2020) Our findings indicated that the anticoagulant gene is highly expressed in young leech suggesting the possibility that young leeches are more effective in medicine while further experimental verification of the exact function is needed.
Conclusion
In this study, the transcriptome of H. nipponia was systematically analyzed in three different developmental stages. Trait and histological observations revealed that adult leeches had significantly higher body length and body weight and a significantly higher number of myenteric tissues than larval and young leeches. Transcriptome analysis revealed that the number of DEGs was highest in L vs. A; DEGs enrichment analysis was mainly related to metabolic pathways; genes related to muscle development, genes related to gonadal development or oviposition regulation and genes related to neural development were significantly period-specific, indicating that these changes were associated with leech growth and development. In addition, we found that hirudin and anticoagulation-related genes were highly expressed in the young stage.
Data Availability Statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.
Author Contributions
Conceptualization, QL and JZ; methodology, XY, XM, HS, TF, and YZ; software, CC, TF, HS, and WZ; formal analysis, XM, XY, SR, and YZ; investigation, KR, TF, XY, and CC; resources, QL and JZ; data curation, XM, XY, TF, and HS; writing—original draft preparation, XM, and XY; writing—review and editing, SR, QL, TF, KR, HS, YZ, CC, WZ, and JZ; supervision, QL and JZ; project administration, QL and JZ; funding acquisition, QL. All authors have read and agreed to the published version of the manuscript.
Funding
This work was granted and supported by the National Natural Science Fund (U20A2051), Guangxi Natural Science Foundation (AB18221120), Guangxi Distinguished scholars Program (201835), and Qihuang High-level Talent Team Training Projects of Guangxi University of Chinese Medicine: Application of Systems Biology in Chinese Medicine Research (2021005).
Conflict of Interest
Author WZ is employed by the company Thai Natural Hirudin Co, Ltd.
The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
The reviewer MTUQ declared a shared affiliation with the authors XY, HS, TF, SR, YZ, and QL to the handling editor at the time of review.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphys.2022.873831/full#supplementary-material
References
Ceylan M. (2020). Effects of Maternal Age on Reproductive Performance of the Southern Medicinal Leech, Hirudo Verbana Carena, 1820. Animal Reproduction Sci. 218, 106507. doi:10.1016/j.anireprosci.2020.106507
Chen S., Zhou Y., Chen Y., Gu J. (2018). Fastp: an Ultra-fast All-In-One FASTQ Preprocessor. Balt. Md 34 (17), i884–i890. doi:10.1093/bioinformatics/bty560
Cohen S., Richter-Levin A., Shefi O. (2020). Brief Electrical Stimulation Triggers an Effective Regeneration of Leech CNS. eNeuro 7 (3), 2020. doi:10.1523/ENEURO.0030-19.2020
Du H., Taylor H. S. (2015). The Role ofHoxGenes in Female Reproductive Tract Development, Adult Function, and Fertility. Cold Spring Harb. Perspect. Med. 6 (1), a023002. doi:10.1101/cshperspect.a023002
El-Obeid H. A., Al-Badr A. A. (1982). Mestranol. Anal. Profiles Drug Subst., 375–406. doi:10.1016/s0099-5428(08)60270-1
González-Gutiérrez A. G., Gutiérrez-Mora A., Verdín J., Rodríguez-Garay B. (2021). An F-Actin Mega-Cable Is Associated with the Migration of the Sperm Nucleus during the Fertilization of the Polarity-Inverted Central Cell of Agave Inaequidens. Front. Plant Sci. 12, 774098. doi:10.3389/fpls.2021.774098
Hammer J. A., Wang J. C., Saeed M., Pedrosa A. T. (2019). Origin, Organization, Dynamics, and Function of Actin and Actomyosin Networks at the T Cell Immunological Synapse. Annu. Rev. Immunol. 37, 201–224. doi:10.1146/annurev-immunol-042718-041341
Harris T. P., Schimenti K. J., Munroe R. J., Schimenti J. C. (2014). IQ Motif-Containing G (Iqcg) Is Required for Mouse Spermiogenesis. G3 (Bethesda) 4 (2), 367–372. doi:10.1534/g3.113.009563
Heath-Heckman E., Yoo S., Winchell C., Pellegrino M., Angstadt J., Lammardo V. B., et al. (2021). Transcriptional Profiling of Identified Neurons in Leech. BMC Genomics 22 (1), 215. doi:10.1186/s12864-021-07526-0
Hill C., Brunello E., Fusi L., Ovejero J. G., Irving M. (2021). Myosin-based Regulation of Twitch and Tetanic Contractions in Mammalian Skeletal Muscle. Elife 10, 68211. doi:10.7554/elife.68211
Iwama R., Oceguera-Figueroa A., Giribet G., Kvist S. (2019). The Salivary Transcriptome of Limnobdella Mexicana (Annelida: Clitellata: Praobdellidae) and Orthology Determination of Major Leech Anticoagulants. Parasitology 146 (10), 1338–1346. doi:10.1017/s0031182019000593
Jurewicz E., Robaszkiewicz K., Moraczewska J., Filipek A. (2020). Binding of S100A6 to Actin and the Actin-Tropomyosin Complex. Sci. Rep. 10 (1), 12824. doi:10.1038/s41598-020-69752-y
Kadzik R. S., Homa K. E., Kovar D. R. (2020). F-actin Cytoskeleton Network Self-Organization through Competition and Cooperation. Annu. Rev. Cell. Dev. Biol. 36, 35–60. doi:10.1146/annurev-cellbio-032320-094706
Kim D., Langmead B., Salzberg S. L. (2015). HISAT: a Fast Spliced Aligner with Low Memory Requirements. Nat. Methods. 12 (4), 357–360. doi:10.1038/nmeth.3317
Koeppen D., Aurich M., Pasalar M., Rampp T. (2019). Medicinal Leech Therapy in Venous Congestion and Various Ulcer Forms: Perspectives of Western, Persian and Indian Medicine. J. Tradit. Complement. Med. 10 (2), 104–109. doi:10.1016/j.jtcme.2019.08.003
Kondo M., Yamamoto T., Takahashi S., Taira M. (2017). Comprehensive Analyses Ofhoxgene Expression inXenopus Laevisembryos and Adult Tissues. Dev. Growth Differ. 59 (6), 526–539. doi:10.1111/dgd.12382
Kwak H. J., Park J. S., Medina Jiménez B. I., Park S. C., Cho S. J. (2019). Spatiotemporal Expression of Anticoagulation Factor Antistasin in Freshwater Leeches. Int. J. Mol. Sci. 20 (16), 3994. doi:10.3390/ijms20163994
Li H., Handsaker B., Wysoker A., Fennell T., Ruan J., Homer N., et al. (2009). Genome Project Data Processing Subgroup 1.0.0.0The Sequence Alignment/map Format and SAMtools. Bioinformatics 25 (16), 2078–2079. doi:10.1093/bioinformatics/btp352
Li H., Li Q., Yu H., Du S. (2021). Characterization of Paramyosin Protein Structure and Gene Expression during Myogenesis in Pacific Oyster (Crassostrea gigas). Comp. Biochem. Physiology Part B Biochem. Mol. Biol. 255, 110594. doi:10.1016/j.cbpb.2021.110594
Liu X., Shu S., Korn E. D. (2020). Muscle Myosins Form Folded Monomers, Dimers, and Tetramers during Filament Polymerization In Vitro. Proc. Natl. Acad. Sci. U.S.A. 117 (27), 15666–15672. doi:10.1073/pnas.2001892117
Liu Z., Tong X., Su Y., Wang D., Du X., Zhao F., et al. (2019). In-depth Profiles of Bioactive Large Molecules in Saliva Secretions of Leeches Determined by Combining Salivary Gland Proteome and Transcriptome Data. J. Proteomics 200, 153–160. doi:10.1016/j.jprot.2019.03.009
Liu Z., Zhao F., Tong X., Liu K., Wang B., Yang L., et al. (2018). Comparative Transcriptomic Analysis Reveals the Mechanism of Leech Environmental Adaptation. Gene 664, 70–77. doi:10.1016/j.gene.2018.04.063
Lu Z. H., Shi P., You J. H., Li J. D., Chen S. J. (2018). Analysis of Growth Performance and Fitting Growth Curve of Hirudo Nipponia in Early Developmental Stages. J. Anhui Agric. Sci. 46 (25), 75–77. doi:10.13989/j.cnki.0517-6611.2018.25.021
Lu Z., Shi P., You H., Liu Y., Chen S. (2018). Transcriptomic Analysis of the Salivary Gland of Medicinal Leech Hirudo Nipponia. PLoS One 13 (10), e0205875. doi:10.1371/journal.pone.0205875
Manav M., Ceylan M., Büyükçapar H. M. (2019). Investigation of Reproductive Efficiency, Growth Performance and Survival of the Southern Medicinal Leech, Hirudo Verbana Carena, 1820 Fed with Mammalian and Poultry Blood. Animal Reproduction Sci. 206, 27–37. doi:10.1016/j.anireprosci.2019.05.004
Martínez-Heredia J., de Mateo S., Vidal-Taboada J. M., Ballescà J. L., Oliva R. (2008). Identification of Proteomic Differences in Asthenozoospermic Sperm Samples. Hum. Reprod. 23 (4), 783–791. doi:10.1093/humrep/den024
Meiser S., Ashida G., Kretzberg J. (2019). Non-synaptic Plasticity in Leech Touch Cells. Front. Physiol. 10, 1444. doi:10.3389/fphys.2019.01444
Meriaux C., Arafah K., Tasiemski A., Wisztorski M., Bruand J., Boidin-Wichlacz C., et al. (2011). Multiple Changes in Peptide and Lipid Expression Associated with Regeneration in the Nervous System of the Medicinal Leech. PLoS One 6 (4), e18359. doi:10.1371/journal.pone.0018359
Pan M., Zhang Q., Liu P., Huang J., Wang Y., Chen S. (2016). Inhibition of the Nuclear Export of P65 and IQCG in Leukemogenesis by NUP98-IQCG. Front. Med. 10 (4), 410–419. doi:10.1007/s11684-016-0489-0
Pavadai E., Lehman W., Rynkiewicz M. J. (2020). Protein-Protein Docking Reveals Dynamic Interactions of Tropomyosin on Actin Filaments. Biophysical J. 119 (1), 75–86. doi:10.1016/j.bpj.2020.05.017
Peña A. A., Bols N. C., Marshall S. H. (2010). An Evaluation of Potential Reference Genes for Stability of Expression in Two Salmonid Cell Lines after Infection with Either Piscirickettsia Salmonis or IPNV. BMC Res. Notes 3 (3), 101. doi:10.1186/1756-0500-3-101
Pertea G., Pertea M. (2020). GFF Utilities: GffRead and GffCompare. F1000Res 9, 304. ISCBCommJ-304. doi:10.12688/f1000research.23297.1
Pertea M., Pertea G. M., Antonescu C. M., Chang T.-C., Mendell J. T., Salzberg S. L. (2015). StringTie Enables Improved Reconstruction of a Transcriptome from RNA-Seq Reads. Nat. Biotechnol. 33 (3), 290–295. doi:10.1038/nbt.3122
Pickrell B. B., Daly M. C., Freniere B., Higgins J. P., Safa B., Eberlin K. R. (2020). Leech Therapy Following Digital Replantation and Revascularization. J. Hand Surg. 45 (7), 638–643. doi:10.1016/j.jhsa.2020.03.026
Robledo D., Hernández-Urcera J., Cal R. M., Pardo B. G., Sánchez L., Martínez P., et al. (2014). Analysis of qPCR Reference Gene Stability Determination Methods and a Practical Approach for Efficiency Calculation on a Turbot (Scophthalmus maximus) Gonad Dataset. BMC Genomics 15 (1), 648. doi:10.1186/1471-2164-15-648
Saglam N. (2019). Internal and External Morphological Characteristics of the Medicinal Leech Species Hirudo Sulukii and Hirudo Verbana. Turk. Parazitol. Derg. 43 (4), 204
Shan H., Ren K., Liu J., Rehman S. u., Yan X., Ma X., et al. (2022). Comprehensive Transcriptome Sequencing Analysis of Hirudinaria Manillensis in Different Growth Periods. Front. Physiol. 13, 897458. doi:10.3389/fphys.2022.897458
Sig A. K., Guney M., Uskudar Guclu A., Ozmen E. (2017). Medicinal Leech Therapy-An Overall Perspective. Integr. Med. Res. 6 (4), 337–343. doi:10.1016/j.imr.2017.08.001
Silva A. M. M., Kennedy L. S., Hasan S. C., Cohen A. M., Heeley D. H. (2020). Demonstration of Beta-Tropomyosin (Tpm2) and Duplication of the Alpha-Slow Tropomyosin Gene (TPM3) in Atlantic Salmon Salmo salar. Comp. Biochem. Physiology Part B Biochem. Mol. Biol. 245, 110439. doi:10.1016/j.cbpb.2020.110439
Singh S. K., Rajoria K. (2020). Medical Leech Therapy in Ayurveda and Biomedicine - A Review. J. Ayurveda Integr. Med. 11 (4), 554–564. doi:10.1016/j.jaim.2018.09.003
Sonobe H., Obinata T., Minokawa T., Haruta T., Kawamura Y., Wakatsuki S., et al. (2016). Characterization of Paramyosin and Thin Filaments in the Smooth Muscle of Acorn Worm, a Member of Hemichordates. J. Biochem. 160 (6), 369–379. doi:10.1093/jb/mvw047
Squire J. (2019). Special Issue: The Actin-Myosin Interaction in Muscle: Background and Overview. Int. J. Mol. Sci. 20 (22), 5715. doi:10.3390/ijms20225715
Straub A., Roehrig S., Hillisch A. (2011). Oral, Direct Thrombin and Factor Xa Inhibitors: the Replacement for Warfarin, Leeches, and Pig Intestines? Angew. Chem. Int. Ed. 50 (20), 4574–4590. doi:10.1002/anie.201004575
Tan E. G., Huang L. Y., Guan Y., Wang B. (2002). Studies on Growth and Reproduction of Hirudinaria Manillensis in Guangdong Province. Chin. Traditional Herb. Drugs 33 (9), 837. doi:10.3321/j.issn:0253-2670.2002.09.036
Tasiemski A., Salzet M. (2017). Neuro-immune Lessons from an Annelid: The Medicinal Leech. Dev. Comp. Immunol. 66, 33–42. doi:10.1016/j.dci.2016.06.026
Wang C.-H., Pandey S., Sivalingam K., Shibu M. A., Kuo W.-W., Yu-LanYeh L., et al. (2021). Leech Extract: A Candidate Cardioprotective against Hypertension-Induced Cardiac Hypertrophy and Fibrosis. J. Ethnopharmacol. 264, 113346. doi:10.1016/j.jep.2020.113346
Xu D., Yang H., Zhuo Z., Lu B., Hu J., Yang F. (2021). Characterization and Analysis of the Transcriptome in Opisina Arenosella from Different Developmental Stages Using Single-Molecule Real-Time Transcript Sequencing and RNA-Seq. Int. J. Biol. Macromol. 169, 216–227. doi:10.1016/j.ijbiomac.2020.12.098
Yenia Azamar Y., Uribe S., Mújica A. (2007). F-actin Involvement in guinea Pig Sperm Motility. Mol. Reprod. Dev. 74 (3), 312. doi:10.1002/mrd.20578
Yu M., Zhou M., Cao M., Feng X. L., Zeng D. J., Wu J. Y. (2021). Review on the Development of the Medicinal Leech. J. Aquac. 42 (06), 39. doi:10.3969/j.issn.1004-2091.2021.06.009
Zhang X., Zhang X., Yuan J., Du J., Li F., Xiang J. (2018). Actin Genes and Their Expression in Pacific White Shrimp, Litopenaeus Vannamei. Mol. Genet. Genomics 293 (2), 479–493. doi:10.1007/s00438-017-1397-y
Keywords: Hirudo nipponia, different growth periods, transcription analysis, developmental, reproductive pathways
Citation: Ma X, Yan X, Ke R, Shan H, Rehman Su, Feng T, Zheng Y, Chuang C, Zhou W, Liu Q and Zheng J (2022) Comparative Transcriptome Sequencing Analysis of Hirudo nipponia in Different Growth Periods. Front. Physiol. 13:873831. doi: 10.3389/fphys.2022.873831
Received: 16 February 2022; Accepted: 02 June 2022;
Published: 23 June 2022.
Edited by:
Ya-Nan Zhang, Huaibei Normal University, ChinaCopyright © 2022 Ma, Yan, Ke, Shan, Rehman, Feng, Zheng, Chuang, Zhou, Liu and Zheng. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Qingyou Liu, qyliu-gene@gxu.edu.cn; Jinghui Zheng, zhengjinghui@hotmail.com
†These authors have contributed equally to this work