Liver transcriptome analysis reveals extensive transcriptional plasticity during acclimation to low salinity in Cynoglossus semilaevis
BMC Genomics volume 19, Article number: 464 (2018)
Salinity is an important abiotic stress that influences the physiological and metabolic activity, reproduction, growth and development of marine fish. It has been suggested that half-smooth tongue sole (Cynoglossus semilaevis), a euryhaline fish species, uses a large amount of energy to maintain osmotic pressure balance when exposed to fluctuations in salinity. To delineate the molecular response of C. semilaevis to different levels of salinity, we performed RNA-seq analysis of the liver to identify the genes and molecular and biological processes involved in responding to salinity changes.
The present study yielded 330.4 million clean reads, of which 83.9% were successfully mapped to the reference genome of C. semilaevis. One hundred twenty-eight differentially expressed genes (DEGs), including 43 up-regulated genes and 85 down-regulated genes, were identified. These DEGs were highly represented in metabolic pathways, steroid biosynthesis, terpenoid backbone biosynthesis, butanoate metabolism, glycerolipid metabolism and the 2-oxocarboxylic acid metabolism pathway. In addition, genes involved in metabolism, osmoregulation and ion transport, signal transduction, immune response and stress response, and cytoskeleton remodeling were affected during acclimation to low salinity. Genes acat2, fdps, hmgcr, hmgcs1, mvk, pmvk, ebp, lss, dhcr7, and dhcr24 were up-regulated and abat, ddc, acy1 were down-regulated in metabolic pathways. Genes aqp10 and slc6a6 were down-regulated in osmoregulation and ion transport. Genes abat, fdps, hmgcs1, mvk, pmvk and dhcr7 were first reported to be associated with salinity adaptation in teleosts.
Our results revealed that metabolic pathways, especially lipid metabolism were important for salinity adaptation. The candidate genes identified from this study provide a basis for further studies to investigate the molecular mechanism of salinity adaptation and transcriptional plasticity in marine fish.
Salinity is an important environmental factor for marine fish aquaculture . Fluctuations in salinity have a significant impact on fish reproduction, growth, development, and physiological and metabolic activities . However, the molecular mechanism underlying salinity adaptation is not clear. Euryhaline fish species can adapt to a wide range of salinities and cope with both chronic and rapid osmotic stress , thus provide an excellent model to study osmoregulation during acclimation to various aquatic environments. The half-smooth tongue sole (Cynoglossus semilaevis) is a euryhaline fish species that can survive in a wide range of salinities range from 14 to 37 ppt (parts per thousand) . The fish is mainly distributed in the East Asia and has emerged as an important commercial fish in aquaculture in China. The optimal salinity is 30 ppt in aquaculture [4, 5].
Osmoregulation is one of the most energetically costly metabolic activities in teleosts . A large amount of energy is consumed by fish to maintain their osmotic homeostasis during acclimation to either freshwater or hyper-saline water [2, 7]. The liver, as the essential metabolic organ in fish, has been recognized as the major source of carbohydrate metabolites for osmoregulatory organs [2, 8]. In addition, the liver participates in certain important physiological processes such as antioxidant activity, glycogen synthesis and bile secretion in teleost fish. Long-term salinity stress inevitably leads to metabolic changes in the liver. However, compared with the widely studied osmoregulatory organs such as the gills, kidney and intestine, the liver response to osmotic stress is hardly known in teleosts [2, 9].
Previous studies on salinity stress response were mainly focused on the physiological and biochemical aspects , and molecular studies have been limited to the cloning and expression level detection of few osmoregulation genes [10,11,12,13]. Identification of the candidate genes involved in salinity adaptation is the first step in elucidating the molecular basis underlying osmoregulation . Previous studies have identified several genes involved in the acute osmotic-stress response including Na+/K+-ATPase (NKA) [12, 13], vacuolar-type H+-ATPase (VHA) , cytochrome c oxidase (COX)  and heat shock proteins (HSPs) [10, 11, 14]. In addition, many genes associated with metabolic processes have also been associated with responses to long-term salinity acclimation [2, 7]. With the advancement of next-generation sequencing (NGS) technology, RNA-seq has become a powerful approach to identify genes and pathways involved in the responses to stress [2, 10, 15, 16]. In the last five years, RNA-seq has been used to study salinity regulation in model and non-model fish species, including spotted sea bass (Lateolabrax maculatus) , striped catfish (Pangasianodon hypophthalmus) [10, 17], Nile tilapia (Oreochromis niloticus) , and medaka (Oryzias melastigm; Oryzias latipes) .
Whole-genome sequencing of the half-smooth tongue sole has been recently completed, thus providing a valuable reference to unravel the molecular mechanisms underlying salinity adaptation using RNA-seq . In the present study, we conducted an RNA-seq analysis to characterize the salinity-induced changes in gene expression in the liver of half-smooth tongue sole. By assessing the transcriptional variations caused by different salinities, it was possible to identify the genes and the molecular and biological processes involved in salinity adaptation. With these results, investigators can develop functional markers to monitor for contemporary responses to climate change  or screen broadly across a species range to predict the potential for adaptation . This information would help to make the fish aquaculture industry successful . Investigating the long-term stress adaptation mechanism also provides insights into understanding plastic and evolutionary responses to various aquatic environments [21,22,23].
All animal experiments were conducted in accordance with the guidelines and approval of Institutional Animal Care and Use Committee of Ocean University of China. The field studies did not involve endangered or protected species.
Salinity challenge experiment and fish sampling
Half-smooth tongue sole were obtained from a local fish farm named shuangying aquatic breeding LLC in Lijin, China. The fish were reared in commercial fish ponds (5 m*5 m*1 m) for 10 months, under controlled conditions (20 ± 0.5 °C; 14:10 h light/dark cycle; O2 ≥ 4 ng/ml; 30 ppt) and fed a compound feed. Individuals were randomly divided into six experimental tanks (N = 40 in each tank) and acclimatized at an optimal salinity of 30 ppt for a period of 7 days before the start of the salinity trial. For the trials, we only selected female half-smooth tongue of similar lengths, since different patterns of expression have been detected for certain genes between male and female C. semilaevis [24,25,26,27]. Six experimental tanks were equally divided into two groups including the low salinity group (LS_L group) and high-salinity group (HS_L group, normal living salinity, acted as a control salinity group). Each group had three replicates. The HS_L group was maintained at 30 ppt during the experiment, whereas the salinity in the LS_L group was decreased by 5 ppt per day over a period of three days towards a target of 15 ppt. The experiment lasted for 60 days. During sampling, the experimental fish were anesthetized using 0.1% tricaine methanesulfonate (MS-222), and livers from three female individuals in each tank were dissected, immediately frozen in liquid nitrogen, and stored at − 80 °C until RNA extraction. Each biological replicate represented three pooled female fish livers to make representative samples for deep sequencing analysis. All animal handling was carried out in accordance with the ethical guidelines and protocols of Ocean University of China’s Animal Care Committee.
Total RNA was extracted from livers using RNAiso reagent (TakaRa, Japan) according to the manufacturer’s instruction. Three fish livers from each tank were pooled for RNA extraction. RNA quality was determined using an Agilent Bioanalyzer system (Agilent Technologies, US), and only samples with an RNA Integrity Number (RIN) greater than 7.0 were used for RNA library construction. Six cDNA libraries were prepared from the RNA samples of two groups according to the manufacturer’s recommendation (New England Biolabs, US). The libraries were then sequenced by the Novogene company (China) using the Illumina HiSeq 2500 platform (125 bp, paired-end reads). Raw data (raw reads) in fastq format were first processed through in-house Perl scripts. In this step, the Q20, Q30 and GC content of the clean data were calculated, and the sequencing quality was assessed. At the same time, clean data (clean reads) were obtained by removing the reads containing adaptors, reads with more than 10% poly-N and reads of low quality (the number of bases with sQ < = 5 accounts for more than 50% of the total read length) from the raw data. All downstream analyses were based on high quality clean data.
Reference genome sequence and read mapping
The clean RNA-seq data were aligned to the reference genome (Accession number: PRJNA251742) using TopHat (version 2.0.9) with mismatch 2 and other parameters set as default. Sequence reads that were mapped to multiple genes or positions were removed. HTSeq v0.6.1 with -m union was used to count the number of reads mapped to each gene. For normalization, the count for each gene was divided by the number of fragments per kilobase of transcript sequence per million base pairs sequenced in each sample (FPKM; fragments per kilobase per million reads) .
Differential expression analysis, gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) enrichment analysis
After performing a quality check of the raw count data, the differences in gene expression between the LS_L group and HS_L group were analyzed using DESeq2 . Read count data were used as input for the program. Differential expression was tested at an adjusted P-value < 0.05 using the Benjamini-Hochberg procedure. Annotation of the DEGs was achieved through BLASTN similarity searches against the whole-genome sequence of half-smooth tongue sole. To characterize the DEGs with Gene Ontology (GO) terms, which have three main categories (biological process, cellular component and molecular function), GO enrichment analysis was conducted using the GOseq R package  with correction of gene length bias. Significantly enriched GO terms were ranked using corrected P-values less than 0.05. The results provided a broad overview of groups of genes cataloged in the three ontology vocabularies. KOBAS software was used to retrieve the statistical enrichment of differential expression genes in KEGG pathways (http://www.genome.jp/kegg/) . FDR ≤ 0.05 indicated relevant pathways within which regulated genes were significantly enriched.
Data validation by quantitative real-time PCR (qPCR)
To verify the accuracy of transcriptomic sequencing data, eleven DEGs were selected for quantitative real-time PCR analysis. Total RNA from the liver was extracted individually from the LS_L group and HS_L group. RNA samples were analyzed in biological triplicate and technical triplicate for qPCR. Quantitative PCR was performed using StepOnePlus™ (American) and the SYBR Premix Ex TaqTM (TliRNaseH Plus) Kit (TaKaRa, Japan, Code No. RR420A) according to the manufacturer’s protocol. Gene-specific primers of eleven target genes and one housekeeping gene (18 s)  were designed using Oligo 6.0 software and synthesized by BGI (The Beijing Genomics Institute). The primer sequences are listed in Additional file 1: Table S1. The 20-μl PCR mixture consisted of 10 μl SYBR®Premix Ex Taq (TliRNaseH Plus), 0.4 μl PCR Forward Primer, 0.4 μl PCR Reverse Primer, 0.4 μl ROX Reference Dye (50×), 2 μl DNA template, and RNase-free water to a total volume of 20 μl. The qPCRs were performed with the following conditions: denaturation at 95 °C for 30 s, 40 cycles of denaturation at 95 °C for 5 s, annealing at 60 °C for 30 s, and extension at 72 °C for 30 s. Relative gene expression was calculated using the 2-ΔΔCt method . Excel software was used to calculate the correlation coefficient between quantitative expression by qPCR and by transcriptome analysis.
Availability of supporting data
The sequencing data from this study have been submitted to the NCBI Sequence Read Archive (SRA) (https://www.ncbi.nlm.nih.gov/sra) under the accession number SRP071827. The raw and processed data have been submitted to the NCBI Gene Expression Omnibus (GEO) (https://www.ncbi.nlm.nih.gov/gds) under the accession number GSE111312.
Illumina sequencing and reads mapping
A total of 342.36 million raw reads were obtained, including 172.19 million reads in the LS_L group and 170.17 million in the HS_L group. The correlation coefficients between replicated samples were greater than 0.88 (Additional file 2: Figure S1). After preprocessing and removal of low-quality sequences, a total of 330.39 million clean reads were obtained, including 166.07 million reads in the LS_L group and 164.32 million in the HS_L group. Ten GB of clean bases were generated for each group. More than 94% of bases had a base accuracy of 99% and more than 89% of bases had a base accuracy of 99.9%. These clean reads were mapped against the annotated genome of the half-smooth tongue sole. A total of 282.10 million reads were successfully mapped. A total of 271.09 million reads were uniquely mapped (Table 1).
Differential expression analysis between livers from the LS_L and HS_L groups
To identify the differentially expressed genes, transcriptome data of half-smooth tongue sole livers from the LS_L group and HS_L group were analyzed. One hundred twenty-eight candidate genes were significantly differentially expressed between the LS_L group and HS_L group using the criteria of FDR-adjusted P-value< 0.05 and |Log2(fold change)| > 1. Compared with the HS_L group, 43 DEGs were up-regulated and 85 DEGs were down-regulated in the LS_L group (Fig. 1). Detailed DEG information is shown in Additional file 3: Table S2. To illustrate the differential expression of genes detected in the liver among different salinities, a heat map of FPKM-normalized transcript counts for the DEGs in each pairwise comparison was generated (Additional file 4: Figure S2).
Gene ontology analysis
All annotated genes were divided into three categories: molecular function, cellular component and biological process. Of the 34,751 unigenes, 27,828 were assigned to 881 GO categories (Additional file 5: Table S3). GO enrichment analysis of the DEGs showed the top processes according to the number of genes and the enrichment level (P-value < 0.05). According to biological process, sterol metabolic process (GO: 0016125) and organic hydroxy compound metabolic process (GO: 1901615) were highly ranked in GO enrichment. According to molecular function, iron ion binding (GO: 0005506) was markedly enriched. There was no significant enrichment of GO terms in the cellular component category (Fig. 2). More GO terms were identified based on the up-regulated genes. These terms included sterol metabolic process (GO: 0016125), alcohol metabolic process (GO: 0006066), organic hydroxy compound metabolic process (GO: 1901615), steroid metabolic process (GO: 0008202) and fatty acid biosynthetic process (GO: 0006633), which were all significantly enriched in the biological process category. In the molecular function category, two more significant enriched GO terms were found, namely cholesterol delta-isomerase activity (GO: 0047750) and intramolecular oxidoreductase activity (GO: 0016863) (Additional file 6: Figure S3).
KEGG pathway analysis identifies molecular interaction networks within cells and helps to elucidate the potential biological functions of analyzed genes. As shown in Figs. 3, 6 pathways were found to be highly enriched by KEGG analysis. In which, metabolic pathways were the most frequently represented pathways in response to long-term hypotonic stress. The DEGs enriched in these pathways are shown in Table 2. The metabolic pathways contained 32 DEGs, whereby 15 of these genes were up-regulated in low salinity conditions, and 17 were down-regulated during acclimation to low salinity. Some pathways shared the same DEGs. Genes highlighted in the steroid biosynthesis pathway, terpenoid backbone biosynthesis pathway, glycerolipid metabolism pathway, and 2-oxocarboxylic acid metabolism pathway also functioned in the metabolic pathways. In the butanoate metabolism pathway, all of the highlighted genes except for the aacs gene were also involved in metabolic pathways.
The classification of gene function analysis
Based on the combination of enrichment analysis, annotation and manual literature searches, gene function annotation resources showed various candidate DEGs potentially associated with salinity adaptation and osmoregulation. These DEGs are involved in metabolism, osmoregulation and ion transport, signal transduction, immune response and stress response, and cytoskeleton remodeling (Table 3). To illustrate the differential expression of genes detected in the liver among different salinities, a heat map of FPKM-normalized transcript counts for the DEGs in two groups was generated (Fig. 4).
Validation of RNA-seq results by qPCR
To verify the accuracy of the RNA-seq data, eleven genes were chosen for validation by qPCR. These genes included hpx, LOC103377668, slc40a1, LOC103381466, pmvk, pdia6, g0 s2, ctsk, itpka, anxa3 and sar1b. As shown in Fig. 5, the correlation coefficient of the relative log2(fold changes) was 0.91 between RNA-seq and qPCR after salinity acclimation, suggesting the correctness of the bioinformatics analysis for the transcriptomic sequencing data.
With the advancement of transcriptome sequencing, RNA-Seq has emerged as a powerful tool to study the molecular response in fish impacted by environmental changes . In the present study, we conducted a transcriptome analysis to reveal the changes of gene expression profiles in the liver of half-smooth tongue sole when responding to salinity changes. A total of 128 DEGs were identified during acclimation to low salinity conditions. These DEGs are enriched in several biological processes, including lipid metabolism, steroid biosynthesis, osmoregulation and ion transport, and signal transduction. Data from our studies are consistent with previous findings that alteration of liver metabolism plays a critical role in salinity adaptation in L. maculatus . This is probably due to the fact that additional energy is needed for ion transport during acclimation to different osmotic environments in euryhaline fishes, and liver is the main source of energy for osmoregulatory organs to help with the synthesis and function of enzymes and transporters .
Metabolism was one of the most enriched pathways during acclimation to a low salinity environment
It has been reported that euryhaline fish can utilize lipids as an energy source when encountering osmotic stress [9, 21, 34]. In our study, we found that 32 genes in the metabolic pathways were identified, in which genes ebp, cyp51, fdps, pmvk, hmgcs1, lss, dhcr24, acat2, dhcr7, and abat were involved in lipid metabolism under the low salinity environment. The cyp51 gene encodes lanosterol 14-alpha demethylase, a key enzyme involved in lipid metabolism and steroid biosynthesis [35, 36]. Interestingly, the expression of cyp51-like gene was down-regulated in Oreochromis mossambicus, a euryhaline freshwater fish, but up-regulated in O. niloticus, a stenohaline freshwater fish when responding to increased salinity environment. This indicated the different role of cyp51-like in different species during acclimation to different salinity conditions . It has been shown that reduced expression of cyp51 gene might contribute to the decreased lipid levels in both XBP1- and IRE1α-deficient mice, whose lipid metabolism were significantly down-regulated . In our study, we found that cyp51 was up-regulated during acclimation to low salinity conditions, indicating the increased lipid metabolism and lipid levels in C. semilaevis. The novel liver protein acetyl-CoA acetyltransferase-2 (ACAT2) is involved in beta-oxidation and lipid metabolism . In L. maculatus, acetyl-CoA acetyltransferase (acat) was down-regulated during adaptation to high salinity conditions . The down-regulation of expression of acat1 suggested that the amount of energy produced from the TCA cycle might be low in FW milkfish under hypothermia . In our study, the expression of acat2 was up-regulated when salinity was decreased, suggesting that the amount of energy might be high.
Terpenoids are the most abundant lipids and precursors in the synthesis of complex secondary metabolites such as sterols and steroids . Owing to the increasing amount of energy, the demand for terpenoids is highly urgent. The terpenoid backbone biosynthesis pathway was first found to be involved in salinity adaptation. The production of terpenoids requires six key enzymes and the expression of these genes (acat2, hmgcs1, hmgcr, mvk, pmvk and fdps) were shown to be up-regulated in our study. ACAT2, HMGCS1, HMGCR, MVK, PMVK are all crucial enzymes involved in mevalonate pathway, which is an important process of terpenoid backbone biosynthesis and is considered to be the only pathway of producing the precursors isopentenyl-PP (IPP) and dimethylallyl-PP (DMAPP) . Next, the transformation of DMAPP and IPP was catalyzed into farnesyl-PP (FPP), which finally produces terpenoids with the help of farnesyl diphosphate synthase (FDPS) [42, 43]. These highly activated enzymes involved in lipid metabolism and terpenoids biosynthesis indicated the increased amount of energy and lipids in the present study. It can be deduced that lipid metabolism and biosynthesis were more likely needed during hypotonic acclimation.
Steroid biosynthesis played a pivotal role in response to salinity stress in aquatic animals . Our study identified some biosynthesis-related genes such as ebp, lss, dhcr7, dhcr24, and cyp51, and these genes also played roles in metabolic pathways. Sterols are important compounds in many biological membranes; they not only act as cell membrane components, but also possess transport capability [44, 45]. Emopamil binding protein (ebp), also known as sterol isomerase, is an essential enzyme in the sterol biosynthesis pathway in eukaryotes . In the present study, the ebp gene expression was up-regulated to active sterol biosynthesis in response to salinity stress, indicating an increase in sterol levels in the liver. Lanosterol is an upstream precursor of sterol biosynthesis for all animal and fungal steroids (especially cholesterol) . The lss gene encodes lanosterol synthase and a null mutation for lss decreased cholesterol levels in the Shumiya cataract rat . Increased lss expression in our case may indicate improved cholesterol levels in the liver. After the synthesis of lanosterol, the production of 14-demethylation of lanosterol catalyzed by CYP51 eventually yields cholesterol after a series of complicated reactions. In our study, cyp51 was up-regulated during acclimation to low salinity conditions, indicating active cholesterol synthesis in C. semilaevis. In a previous study of constant salinity change, lss and dhcr24 in the steroid biosynthesis pathway were found to be up-regulated in Nile tilapia (O. niloticus) , which is consistent with the data from our present study.
Cholesterol biosynthesis is crucial in the relationship between osmoregulation and steroid because the cholesterol induced by salinity fluctuation is related to the physical properties of cell membranes [15, 49]. An increasing trend was found for high-density lipoprotein cholesterol as salinity decreased in Eriocheir sinensis . Genes 24-dehydrocholesterol reductase (dhcr24) and 7-dehydrocholesterol reductase (dhcr7) are the final enzymes required for cholesterol biosynthesis . A stronger evolutionary correlation was found between DHCR24 and Na+,K+-ATPase than between DHCR24 and any other membrane protein investigated, indicating that cholesterol evolved together with Na+,K+-ATPase in multicellular animals to support Na+,K+-ATPase activity . Previous studies suggested that the cholesterol production induced by salinity might influence the physical properties of cell membranes . Transcriptome profiling in O. niloticuspivotal showed the vital roles of sterol metabolism in response to salinity stress . Hyposmotic stress may cause damage to cell membranes. Given that five up-regulated DEGs in our study contributed to cholesterol synthesis, we can infer that biosynthesis of cholesterol played an important role in the relationship between osmoregulation and steroid production in order to maintain osmotic homeostasis in a long-term low salinity environment. Cholesterol function was more likely related to the physical properties of cell membranes. In addition, cholesterol is the precursor of cortisol, a hormone which responds to various stress conditions and can be widely used as a biomarker in aquaculture. Previous study has comfirmed that the levels of cortisol were elevated under stress to maintain the homeostasis in fish .
It has been reported that the propanoate metabolism pathway was activated in the gill of the pacific white shrimp Litopenaeus vannamei when responding to acute low salinity stress . Butanoate metabolism was identified as a new pathway for long-term low salinity responses in our study. The 4-aminobutyrate aminotransferase (abat) gene, also known as GABA transaminase, was down-regulated in the metabolic pathways and butanoate metabolism pathway. In the plant Arabidopsis, GABA transaminase (GABA-T), the first step of GABA catabolism, was identified as the most responsive to NaCl levels . This is the first report of abat in teleosts, and further research is needed to reveal the underlying regulation mechanism.
Osmoregulation and ion transport were involved in acclimation to long-term low salinity conditions
Osmoregulation is the active regulation of osmotic pressure by means of maintaining a balance of intracellular solute concentration when adapting to the surrounding conditions. Aquaporins (AQPs) represent a family of transmembrane channel proteins that allow osmotic-driven transport of water and small solutes across biological membranes and are found in all living organisms . Our studies here showed that aquaporin 10 (aqp10) was down-regulated in low salinity water. Aquaporin 10 is an aquaglyceroporin that is located in the membrane of the mammalian small intestine, where it is believed to participate in transport of both water and solutes . Previous studies have demonstrated that other members of the AQP family play a role in osmoregulation. AQP1 was often used as a reference gene because it was identified to change in teleosts after salinity challenge . A decrease in AQP-3 can help prevent loss of glycerol from the cell . The previous study also showed a decrease in AQP-8a and AQP-10b expression abundance in the intestine of Pangasianodon hypophthalmus when salinity levels increased . Collecitvely, these provided a pivotal hint to clarify water transport activity in teleosts.
The transport of ions across the plasma membrane by transporters is essential to osmoregulation. During accumulation to different salinities, cells can depend on transporters and accumulate intracellular osmolytes, including taurine, betaine, sorbitol and myo-inositol, to balance any osmotic change [60, 61]. The sodium/chloride-dependent taurine transporter (SLC6A6; TauT) has specificity for taurine . Functional analysis of slc6a6 in pacific oyster Crassostrea gigas revealed that slc6a6 required higher NaCl concentration . In our study, half-smooth tongue sole slc6a6 mRNA expression was down-regulated under hyposmotic stress. Taurine transporter mRNA has also been found to be highly expressed in the gill of seawater-acclimated Japanese eel (Anguilla japonica) . It is speculated that the down-regulation of slc6a6 in response to hyposmotic stress was induced by a substantial decrease in tissue taurine content following the decrease in internal osmolality.
Signal transduction was indispensable during acclimation to salinity changes
The initial basis for euryhaline fish to adapt to salinity stress is the efficient mechanisms of osmosensing and osmotic stress signaling . They are the upstream mechanisms that regulate effector protein expression and orchestrate adaptive responses . Signal transduction is activated by ligand-receptor binding and then is propagated through several transducer proteins via phosphorylation or dephosphorylation events . In our results, several genes involved in signal transduction were found differentially expressed between the LS_L group and HS_L group, such as abat, rem2, dpf3 and npr1. Natriuretic peptides play a central role in fish osmoregulation (especially in seawater). Natriuretic peptides inhibited salt loading, rather than stimulated salt loss in marine fish . The gene expression of its receptor NPR1 was found down-regulated in the present result, indicating the potential role of npr1 in the signal transduction. These genes found in our study may have important functions in osmo-regulated signaling pathways, which can be used to regulate effector protein expression and orchestrate adaptive responses.
In the present study, we only analyzed coding RNA (mRNA), but not the non-coding RNAs (e.g., lncRNA or microRNA). Future research will be carried out to obtain a comprehensive understanding of the whole transcriptomic regulatory mechanism in salinity adaptation including lncRNA, microRNA and circRNA sequencing.
The results of this study highlight the response of the liver transcriptome to long-term environmental salinity changes using next-generation sequencing technology. Data from this study indicate that the liver plays a critical role in the long-term salinity adaptation of aquatic animals. A list of candidate genes with differential expression patterns were identified in response to C. semilaevis salinity tolerance. These genes were enriched in diverse biological processes including metabolism, osmoregulation and ion transport, signal transduction, immune response and stress response, and cytoskeleton remodeling. Lipid metabolism and biosynthesis were the most remarkable affected genetic pathways during acclimation to low salinity. Some genes such as abat, fdps, hmgcs1, mvk and pmvk and dhcr7 involved in these processes were first reported to be associated with salinity adaptation. This study not only offered a number of candidate genes for salinity adaptation in half-smooth tongue sole but will also facilitate research into extensive transcriptional plasticity to identify underlying physiological adaptations to different salinities for other euryhaline teleosts.
- C. semilaevis :
Cytochrome c oxidase
Differentially expressed genes
Fragments per kilobase per million reads
Heat shock proteins
Kyoto Encyclopedia of Genes and Genomes
Parts per thousand
RNA Integrity Number
Callaway R, Shinn AP, Grenfell SE, Bron JE, Burnell G, Cook EJ, et al. Review of climate change impacts on marine aquaculture in the UK and Ireland. Aquat Conserv Mar Freshwat Ecosyst. 2012;22(3):389–421.
Zhang X, Wen H, Wang H, Ren Y, Zhao J, Li Y. RNA-Seq analysis of salinity stress–responsive transcriptome in the liver of spotted sea bass (Lateolabrax maculatus). PLoS One. 2017;12(3):e0173238.
Lam SH, Lui EY, Li Z, Cai S, Sung WK, Mathavan S, et al. Differential transcriptomic analyses revealed genes and signaling pathways involved in iono-osmoregulation and cellular remodeling in the gills of euryhaline Mozambique tilapia, Oreochromis mossambicus. BMC Genomics. 2014;15(1):921.
Li S, He F, Wen H, Li J, Si Y, Liu M, et al. Low salinity affects cellularity, DNA methylation, and mRNA expression of igf1 in the liver of half smooth tongue sole (Cynoglossus semilaevis). Fish Physiol Biochem. 2017;43(6):1587–602.
Xu H, Cao L, Zhang Y, Johnson RB, Wei Y, Zheng K, et al. Dietary arachidonic acid differentially regulates the gonadal steroidogenesis in the marine teleost, tongue sole (Cynoglossus semilaevis), depending on fish gender and maturation stage. Aquaculture. 2017;468:378–85.
Boeuf G, Payan P. How should salinity influence fish growth? Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology. 2001;130(4):411–23.
Hwang P, Lee T. New insights into fish ion regulation and mitochondrion-rich cells. Comp Biochem Physiol A Mol Integr Physiol. 2007;148(3):479–97.
Tseng Y, Hwang P. Some insights into energy metabolism for osmoregulation in fish. Comparative Biochemistry & Physiology Toxicology & Pharmacology Cbp. 2008;148(4):419–29.
Lai KP, Li JW, Wang SY, Chiu JMY, Tse A, Lau K, et al. Tissue-specific transcriptome assemblies of the marine medaka Oryzias melastigma and comparative analysis with the freshwater medaka Oryzias latipes. BMC Genomics. 2015;16(1):135.
Thanh NM, Jung H, Lyons RE, Chand V, Tuan NV, Thu VTM, et al. A transcriptomic analysis of striped catfish (Pangasianodon hypophthalmus) in response to salinity adaptation: De novo assembly, gene annotation and marker discovery. Comparative Biochemistry & Physiology Part D Genomics & Proteomics. 2014;10C(1):52–63.
Wang SJ, Wu MJ, Chen XJ, Jiang Y, Yan YB. DsHsp90 is involved in the early response of Dunaliella salina to environmental stress. Int J Mol Sci. 2012;13(7):7963–79.
Huang C, Chao P, Lin H. Na+/K+-ATPase and vacuolar-type H+-ATPase in the gills of the aquatic air-breathing fish Trichogaster microlepis in response to salinity variation. Comp Biochem Physiol A Mol Integr Physiol. 2010;155(3):309–18.
Tine M, McKenzie DJ, Bonhomme F, Durand J-D. Salinity-related variation in gene expression in wild populations of the black-chinned tilapia from various west African coastal marine, estuarine and freshwater habitats. Estuar Coast Shelf Sci. 2011;91(1):102–9.
Wang J, Wei Y, Li X, Cao H, Xu M, Dai J. The identification of heat shock protein genes in goldfish (Carassius auratus) and their expression in a complex environment in Gaobeidian Lake, Beijing, China. Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology. 2007;145(3):350–62.
Xu Z, Gan L, Li T, Xu C, Chen K, Wang X, et al. Transcriptome profiling and molecular pathway analysis of genes in association with salinity adaptation in Nile tilapia Oreochromis niloticus. PLoS One. 2015;10(8):e0136506.
Zhao X, Yu H, Kong L, Li Q. Transcriptomic responses to salinity stress in the Pacific oyster Crassostrea gigas. PLoS One. 2012;7(9):e46244.
Nguyen TV, Jung H, Nguyen TM, Hurwood D, Mather P. Evaluation of potential candidate genes involved in salinity tolerance in striped catfish (Pangasianodon hypophthalmus) using an RNA-Seq approach. Mar Genomics. 2016;25:75–88.
Chen S, Zhang G, Shao C, Huang Q, Liu G, Zhang P, et al. Whole-genome sequence of a flatfish provides insights into ZW sex chromosome evolution and adaptation to a benthic lifestyle. Nat Genet. 2014;46(3):253–60.
Hoffmann AA, Willi Y. Detecting genetic responses to environmental change. Nat Rev Genet. 2008;9(6):421.
Hoffmann AA, Sgrò CM. Climate change and evolutionary adaptation. Nature. 2011;470(7335):479.
Gibbons TC, Metzger DC, Healy TM, Schulte PM. Gene expression plasticity in response to salinity acclimation in threespine stickleback ecotypes from different salinity habitats. Mol Ecol. 2017;26(10):2711–25.
Smith S, Bernatchez L, Beheregaray LB. RNA-seq analysis reveals extensive transcriptional plasticity to temperature stress in a freshwater fish species. BMC Genomics. 2013;14(1):375.
Oomen RA, Hutchings JA. Transcriptomic responses to environmental change in fishes: insights from RNA sequencing. FACETS. 2017;2(2):610–41.
Ma Q, Liu SF, Zhuang ZM, Sun ZZ, Liu CL, Tang QS. The co-existence of two growth hormone receptors and their differential expression profiles between female and male tongue sole (Cynoglossus semilaevis). Gene. 2012;511(2):341–52.
Ma Q, Liu SF, Zhuang ZM, Lin L, Sun ZZ, Liu CL, et al. Genomic structure, polymorphism and expression analysis of growth hormone-releasing hormone and pituitary adenylate cyclase activating polypeptide genes in the half-smooth tongue sole (Cynoglossus semilaevis). Genetics & Molecular Research Gmr. 2011;10(4):3828–46.
Chen SL, Li J, Deng SP, Tian YS, Wang QY, Zhuang ZM, et al. Isolation of female-specific AFLP markers and molecular identification of genetic sex in half-smooth tongue sole (Cynoglossus semilaevis). Mar Biotechnol. 2007;9(2):273–80.
Deng SP, Chen SL, Xu JY, Liu BW. Molecular cloning, characterization and expression analysis of gonadal P450 aromatase in the half-smooth tongue-sole, Cynoglossus semilaevis. Aquaculture. 2009;287(1–2):211–8.
Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011;29(7):644.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.
Matthew D Young MJW, Smyth GK, Oshlack A. Method Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010;11:R14.
Kanehisa M, Araki M, Goto S, Hattori M, Hirakawa M, Itoh M, et al. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 2007;36:D480–D4.
Si Y, He F, Wen H, Li J, Zhao J, Ren Y, et al. Genetic polymorphisms and DNA methylation in exon 1 CpG-rich regions of PACAP gene and its effect on mRNA expression and growth traits in half smooth tongue sole (Cynoglossus semilaevis). Fish Physiol Biochem. 2016;42(2):407–21.
Livak J, Schmittgen D. Analysis of relative gene expression data using real time quantitative PCR and the 2−ΔΔCT method. Methods. 2001;25(4):402–8.
Luvizotto-santos R, Bianchini A. Lipids as energy source during salinity acclimation in the euryhaline crab Chasmagnathus granulata Dana, 1851 (crustacea-grapsidae). Journal of Experimental Zoology Part A: Ecological and Integrative Physiology. 2003;295(2):200–5.
So JS, Hur KY, Tarrio M, Ruda V, Frank-Kamenetsky M, Fitzgerald K, et al. Silencing of lipid metabolism genes through IRE1α-mediated mRNA decay lowers plasma lipids in mice. Cell Metab. 2012;16(4):487–99.
Pilloff D, Dabovic K, Romanowski MJ, Bonanno JB, Doherty M, Burley SK, et al. The kinetic mechanism of phosphomevalonate kinase. J Biol Chem. 2003;278(7):4510–5.
Ronkin D, Seroussi E, Nitzan T, Doron-Faigenboim A, Cnaani A. Intestinal transcriptome analysis revealed differential salinity adaptation between two tilapiine species. Comparative Biochemistry and Physiology Part D: Genomics and Proteomics. 2015;13:35–43.
Sodhi SS, Ghosh M, Song KD, Sharma N, Kim JH, Kim NE, et al. An approach to identify SNPs in the gene encoding acetyl-CoA acetyltransferase-2 (ACAT-2) and their proposed role in metabolic processes in pig. PLoS One. 2014;9(7):e102432.
Hu Y-C, Kang C-K, Tang C-H, Lee T-H. Transcriptomic analysis of metabolic pathways in milkfish that respond to salinity and temperature changes. PLoS One. 2015;10(8):e0134959.
Labaronne E, Pinteur C, Vega N, Pesenti S, Julien B, Meugnier-Fouilloux E, et al. Low-dose pollutant mixture triggers metabolic disturbances in female mice leading to common and specific features as compared to a high-fat diet. J Nutr Biochem. 2017;45:83–93.
Lynen F. Biosynthetic pathways from acetate to natural products. Pure Appl Chem. 1967;14(1):137–68.
Banthorpe DV, Bucknall GA, Doonan HJ, Doonan S, Rowan MG. Biosynthesis of geraniol and nerol in cell-free extracts of Tanacetum vulgare. Phytochemistry. 1976;15(1):91–100.
Sagami H, Ogura K, Seto S, Kurokawa T. A new prenyltransferase from Micrococcus lysodeikticus. Biochem Biophys Res Commun. 1978;85(2):572–8.
Bossche HV, Lauwers W, Willemsens G, Marichal P, Cornelissen F, Cools W. Molecular basis for the antimycotic and antibacterial activity of N-substituted imidazoles and triazoles: the inhibition of isoprenoid biosynthesis. Pest Manag Sci. 1984;15(2):188–98.
Hartmann M-A. Plant sterols and the membrane environment. Trends Plant Sci. 1998;3(5):170–5.
Silve S, Dupuy PH, Labit-Lebouteiller C, Kaghad M, Chalon P, Rahier A, et al. Emopamil-binding protein, a mammalian protein that binds a series of structurally diverse neuroprotective agents, exhibits delta8-delta7 sterol isomerase activity in yeast. J Biol Chem. 1996;271(37):22434–40.
Yoshida Y, Aoyama Y, Noshiro M, Gotoh O. Sterol 14-demethylase P450 (CYP51) provides a breakthrough for the discussion on the evolution of cytochrome P450 gene superfamily. Biochem Biophys Res Commun. 2000;273(3):799–804.
Mori M, Li G, Abe I, Nakayama J, Guo Z, Sawashita J, et al. Lanosterol synthase mutations cause cholesterol deficiency–associated cataracts in the Shumiya cataract rat. J Clin Invest. 2006;116(2):395–404.
Parasassi T, Giusti AM, Raimondi M, Ravagnan G, Sapora O, Gratton E. Cholesterol protects the phospholipid bilayer from oxidative damage. Free Radic Biol Med. 1995;19(4):511–6.
Long X, Wu X, Zhao L, Ye H, Cheng Y, Zeng C. Effects of salinity on gonadal development, osmoregulation and metabolism of adult male Chinese mitten crab, Eriocheir sinensis. PLoS One. 2017;12(6):e0179036.
Prabhu AV, Sharpe LJ, Brown AJ. The sterol-based transcriptional control of human 7-dehydrocholesterol reductase (DHCR7): evidence of a cooperative regulatory program in cholesterol synthesis. Biochimica et Biophysica Acta (BBA)-Molecular and Cell Biology of Lipids. 2014;1841(10):1431–9.
Lambropoulos N, Garcia A, Clarke RJ. Stimulation of Na+, K+-ATPase activity as a possible driving force in cholesterol evolution. J Membr Biol. 2016;249(3):251–9.
Mommsen TP, Vijayan MM, Moon TW. Cortisol in teleosts: dynamics, mechanisms of action. and metabolic regulation Rev Fish Biol Fish. 1999;9(3):211–68.
Wang X, Wang S, Li C, Chen K, Qin JG, Chen L, et al. Molecular pathway and gene responses of the pacific white shrimp Litopenaeus vannamei to acute low salinity stress. J Shellfish Res. 2015;34(3):1037–48.
Renault H, Roussel V, El Amrani A, Arzel M, Renault D, Bouchereau A, et al. The Arabidopsis pop2-1 mutant reveals the involvement of GABA transaminase in salt stress tolerance. BMC Plant Biol. 2010;10(1):20.
Calvanese L, Pellegrini-Calace M, Oliva R. In silico study of human aquaporin AQP11 and AQP12 channels. Protein Sci. 2013;22(4):455–66.
Mobasheri A, Shakibaei M, Marples D. Immunohistochemical localization of aquaporin 10 in the apical membranes of the human ileum: a potential pathway for luminal water and small solute absorption. Histochem Cell Biol. 2004;121(6):463–71.
Wong MKS, Ozaki H, Suzuki Y, Iwasaki W, Takei Y. Discovery of osmotic sensitive transcription factors in fish intestine via a transcriptomic approach. BMC Genomics. 2014;15(1):1134.
Campbell EM, Ball A, Hoppler S, Bowman AS. Invertebrate aquaporins: a review. J Comp Physiol B. 2008;178(8):935–55.
Alfieri RR, Petronini PG. Hyperosmotic stress response: comparison with other cellular stresses. Pflügers Archiv-European Journal of Physiology. 2007;454(2):173–85.
Ho SN. Intracellular water homeostasis and the mammalian cellular osmotic stress response. J Cell Physiol. 2006;206(1):9–15.
Hosoi M, Shinzato C, Takagi M, HOSOI-TANABE S, Sawada H, Terasawa E, et al. Taurine transporter from the giant Pacific oyster Crassostrea gigas: function and expression in response to hyper-and hypo-osmotic stress. Fish Sci. 2007;73(2):385–94.
Chow S, Ching L, Wong A, Wong CK. Cloning and regulation of expression of the Na+–cl−–taurine transporter in gill cells of freshwater Japanese eels. J Exp Biol. 2009;212(20):3205–10.
Fiol DF, Kültz D. Osmotic stress sensing and signaling in fishes. FEBS J. 2007;274(22):5790–8.
Evans T. Co-ordination of osmotic stress responses through osmosensing and signal transduction events in fishes. J Fish Biol. 2010;76(8):1903–25.
Evans DH. Cell signaling and ion transport across the fish gill epithelium. J Exp Zool A Ecol Genet Physiol. 2002;293(3):336–47.
This work was supported by The National Natural Science Funds (41476110), State 863 High-Technology R&D Project of China (2012AA10A403) and Fundamental Research Funds for the Central Universities (201712026).
Availability of data and materials
RNA sequence data were deposited in the DDBJ/EBI/NCBI databases with the accession number SRP071827. Raw/processed counts in addition to sequence data were submitted to the NCBI Gene Expression Omnibus (GEO) database with the accession number GSE111312.
Ethics approval and consent to participate
The fish in this study were obtained from a fish farm named shuangying aquatic breeding LLC in Lijin, China. All animal experiments were conducted in accordance with the guidelines and approval of Institutional Animal Care and Use Committee of Ocean University of China. The studies did not involve endangered or protected species.
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Table S1. Primers used in qPCR validation. (XLSX 11 kb)
Figure S1. Pearson’s correlation coefficients among the LS_L group and HS_L group. (PDF 5 kb)
Table S2. List of DEGs from C. semilaevis in response to long-term hypotonic stress. (XLS 50 kb)
Figure S2. Heat map analysis of all DEGs in the liver from the LS_L group and HS_L group. (PDF 13 kb)
Table S3. List of Gene Ontology (GO) enrichment terms in the liver of half-smooth tongue sole. (XLS 221 kb)
Figure S3. Gene Ontology (GO) terms based on up-regulated DEGs of half-smooth tongue sole during acclimation to low salinity (* indicates significantly enriched GO terms). (PDF 5 kb)
About this article
Cite this article
Si, Y., Wen, H., Li, Y. et al. Liver transcriptome analysis reveals extensive transcriptional plasticity during acclimation to low salinity in Cynoglossus semilaevis. BMC Genomics 19, 464 (2018). https://doi.org/10.1186/s12864-018-4825-4