Early response to heat stress in Chinese tongue sole (Cynoglossus semilaevis): performance of different sexes, candidate genes and networks
BMC Genomics volume 21, Article number: 745 (2020)
Temperature is known to affect living organisms and alter the expression of responsive genes, which affects a series of life processes, such as development, reproduction and metabolism. Several genes and gene families have been involved in high temperature responses, such as heat shock protein (hsp) family, Jumonji family and genes related to cortisol synthesis. Gonad is a vital organ related to the existence of a species. However, the comprehensive understanding of gonadal responses to environmental temperature is limited.
To explore the effects of environmental temperature on genes and gene networks in gonads, we performed acute heat treatment (48 h) on Chinese tongue sole (Cynoglossus semilaevis). Gonadal transcriptome analysis was conducted on females, pseudomales and males exposed to high (28 °C) and normal (22 °C) temperatures. A total of 1226.24 million clean reads were obtained from 18 libraries. Principal component analysis (PCA) and differentially expressed gene (DEG) analysis revealed different performance of sex responses to heat stress. There were 4565, 790 and 1117 specific genes altered their expression level in females, pseudomales and males, respectively. Of these, genes related to hsp gene family, cortisol synthesis and metabolism and epigenetic regulation were involved in early heat response. Furthermore, a total of 1048 DEGs were shared among females, pesudomales and males, which may represent the inherent difference between high and normal temperatures. Genes, such as eef1akmt3, eef1akmt4, pnmt and hsp family members, were found.
Our results depicted for the first time the gonadal gene expression under acute high temperature treatment in Chinese tongue sole. The findings may provide a clue for understanding the responses of genes and networks to environmental temperature.
Water temperature is a major environmental factor that affects the development, metabolism and reproduction of aquatic ectotherms [1, 2]. Due to global warming, the water temperature of the ocean could rise between 1.1–6.4 °C by the end of the twenty-first century . The increase of temperature may change the responsory genes expression, which can affect the regulatory networks. The response to high temperature has the difference between females and males . The gonad is a vital organ related to the existence of a species. Therefore, the particular response mechanism in gonad under high temperature is needed to be clarified.
Many genes can respond rapidly to high temperature. For example, heat shock proteins (hsps) are molecular chaperones that function in protein folding, localization, secretion, and degradation . Heat stress induces the expression of several hsp gene families. Of these, hsp70 and hsp90 are the most prominent ones. Up-regulated hsps have been observed in many fish species, such as red garra (Garra rufa), snakehead murrel (Channa striatus), and killifish (Fundulus heteroclitus) [6,7,8]. Cortisol is another factor in which production is increased by stress. Many studies on different fish species have found significantly higher cortisol concentration under high temperature and indicated that cortisol was involved in high temperature caused masculinization, such as pejerrey (Odontesthes bonariensis), olive flounder (Paralichthys olivaceus), and medaka (Oryzias latipes) [9,10,11]. Epigenetic regulators can link the external factors to internal gene regulation. Among them, the Jumonji family is one of the important factors. The members of this family, such as jarid2 and jmjd3 (kdm6b), can respond to temperature rapidly and affect sex-related gene expression [12, 13]. Although these genes and gene families have been involved in high temperature perception, a comprehensive understanding of gonad responses to environmental temperature is still lacking.
Recently, transcriptome analysis has become an excellent approach for investigating the responses of organisms to environmental changes, including heat [14,15,16]. Most of the reports are mainly focused on short-term high temperature stimulation in tissues, such as muscle and liver. In addition, there are several reports of the effects of long-term high temperature exposure on fish [17,18,19]. However, to the best of our knowledge, research is still lacking regarding the transcriptome data of fish gonads under short-term high temperature stimulation.
Chinese tongue sole (Cynoglossus semilaevis) is an economically important marine flatfish widely distributed in Chinese coastal waters, which possesses a ZW/ZZ sex-determination system. In our previous study, the genetic female: male ratio was close to 1:1 in this species. Approximately 14% of genetic females can spontaneously sex revered into phenotypic males (also known as pseudomales) [20, 21]. The genome of C. semilaevis has been well-sequenced, and the sex-specific simple sequence repeat (SSR) markers have been developed [21, 22]. C. semilaevis is an excellent model to evaluate the response to high temperature among different sexes. In this study, we conducted acute heat treatment on C. semilaevis and performed transcriptome sequencing on gonads of females, males and pseudomales. Our results may provide a clue for understanding gene and network responses to environment temperature.
Gonadal transcriptome of C. semilaevis
To identify genes involved in response to heat stress in C. semilaevis, transcriptome sequencing was performed in the gonads of 3 females, 3 pseudomales, and 3 males of the control (CT_F1–3, CT_P1–3 and CT_M1–3; 22 °C treatment, CT group) and heat stress (HS_F1–3, HS_P1–3 and HS_M1–3; 28 °C treatment, HS group) groups. A total of 1354.19 million raw reads were generated from 18 libraries, which yielded 1226.24 million clean reads after quality control. An average of 88.87% of the clean reads was mapped to the high-quality C. semilaevis reference genome [National Center for Biotechnology Information (NCBI), Cse_v1.0, BioProject no. PRJNA73987] (Additional file 1). A total of 24,230 genes were detected in the transcriptome.
Performance of sexes under heat stress
As shown in Fig. 1a, CT and HS groups were separated into two distinct clusters. Under normal temperature, females, pseudomales and males were apart from each other. After heat stress, females and males exhibited a separation, while pseudomales were close to males. The grouping was supported by differentially expressed gene (DEG) number analysis of three sexes between heat stress and normal temperature. The comparison of HS_F and CT_F showed 7267 DEGs (405 up-regulated and 6862 down-regulated genes). Comparison of HS_P and CT_P observed 2687 DEGs (1111 up-regulated and 1576 were down-regulated genes). In the comparison between HS_M and CT_M, a total of 3702 DEGs were identified (424 up-regulated and 3278 down-regulated genes) (Fig. 1b).
Specific DEGs in each sex under heat stress
To detect different responses of sexes encountering heat stress, female, pseudomale and male-specific DEGs were analyzed. Females, pseudomales and males exhibited 4565, 790, and 1117 specific DEGs under heat stress (Fig. 1c).
In this study, some members of DEGs related to the hsp gene family, epigenetic regulators, cortisol biosynthesis and sex steroid receptors showed sex-biased characteristics under high temperature. Heat shock protein family A (Hsp70) member 4 like (hspa4l), heat shock factor binding protein 1 (hsbp1) and heat shock 70 kDa protein 14-like (hspa14l) were differentially expressed in heat stress treated females. Heat shock protein family A (Hsp70) member 5 (hspa5), heat shock protein family A (Hsp70) member 12A (hspa12a) and heat shock cognate 70 kDa protein (hsc70) genes were found in pseudomales. DnaJ heat shock protein family (Hsp40) member A2 (dnaja2) and heat shock protein beta-7-like (hspβ7l) genes were found in males. Epigenetic regulation-related lysine demethylase 6B (kdm6b) was specifically down-regulated in females under heat stress. S-adenosylmethionine synthase-like (saml) was specifically up-regulated in pseudomales. DNA (cytosine-5)-methyltransferase 1-like (dnmt1l) was down-regulated in males. In addition, cortisol biosynthesis-related genes, including hydroxysteroid 11-beta dehydrogenase 2 (hsd11b2) and hydroxysteroid 17-beta dehydrogenase 1 (hsd17b1) were down-regulated in females under heat stress. On the other hand, estrogen related receptor alpha (esrrα) and androgen receptor (ar) were down-regulated in females, and estrogen-related receptor gamma (errγl) was down-regulated in males in HT vs. CT comparison (Fig. 2, Additional file 2).
Gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) analysis of sex-specific DEGs
There were 100, 55 and 85 GO terms significantly enriched in females, pseudomales and males, respectively (p < 0.05). In HS_F vs. CT_F comparison, the most enriched GO terms were immune system process, integral component of membrane, transferase activity-transferring glycosyl groups and positive regulation of actin filament polymerization. Besides, the GO term of regulation of histone H3-K4 methylation was significantly enriched (p < 0.05). In HS_P vs. CT_P comparison, the GO terms, including cell morphogenesis, RNA polymerase II transcription factor activity, actin binding and signal transducer activity, were mostly enriched. In HS_M vs. CT_M comparison, the most enriched GO terms were chromatin, nucleosome, tissue homeostasis and retina homeostasis (Additional file 3).
There were 103, 6 and 12 KEGG pathways significantly enriched in females, pseudomales and males, respectively (p < 0.05). In females, enriched KEGG pathways mainly referred to the immune system, such as NF-kappa B signaling pathway, cytokine-cytokine receptor interaction, natural killer cell mediated cytotoxicity and chemokine signaling pathway. In pseudomales, the most enriched KEGG pathways were cell adhesion molecules (CAMs), glycosaminoglycan biosynthesis-chondroitin sulfate/dermatan sulfate, amino sugar and nucleotide sugar metabolism and axon guidance. In males, the most enriched KEGG pathways included systemic lupus erythematosus, neuroactive ligand-receptor interaction, cellular senescence and pancreatic secretion (Additional file 3).
Shared DEGs among three sexes under heat stress
The Venn diagram showed that there were 1048 shared DEGs among female, pesudomale and male individuals from the HS vs. CT comparison (Fig. 1c). EEF1A lysine methyltransferase 3 (eef1akmt3), EEF1A lysine methyltransferase 4 (eef1akmt4), phenylethanolamine N-methyltransferase (pnmt), and hsp family members, such as heat shock protein family B (small) member 1 (hspb1), heat shock protein family B (small) member 9 (hspb9), heat shock 70 kDa protein 1 (hspa1), heat shock transcription factor 4 (hsf4), and heat shock transcription factor 2 binding protein (hsf2bp) were found (Additional file 2). GO annotations of the shared DEGs were performed. In the biological process, cellular process, metabolic process and biological regulation were abundant. Cell, membrane and membrane part were abundant in the cellular component. In molecular function, binding, catalytic activity and transporter activity were enriched (Fig. 3a). A total of 26 KEGG pathways were significantly enriched (p < 0.05), including cell cycle, DNA replication and oocyte meiosis (Fig. 3b).
Validation of RNA-seq by real-time qualitative polymerase chain reaction (qRT-PCR)
To validate the transcriptome results, a total of nine significantly differentially expressed DEGs in at least one comparison group were selected for qRT-PCR analysis. Of these, three genes were the members of the hsp family, including hspa1 (LOC103393458), hsf4 and hsc70 (LOC103395053). The eef1akmt4 gene was a candidate regulator of HSP70. The saml gene (LOC103387328) was involved in S-adenosylmethionine synthesis and epigenetic regulation. The hsd11b2 and hsd17b1 genes were involved in cortisol synthesis. In addition, esrrα and ar genes were encoded for sex hormone receptors. For hsf4, RNA-Seq analysis showed significantly down-regulation in females, pseudomales and males under heat-stress, while qRT-PCR indicated no significant differences in females between high and control temperatures. For esrrα, RNA-Seq showed significantly down-regulation in females under heat-stress, while qRT-PCR indicated that there were significant differences in both females and males. Nevertheless, the down-regulation patterns of hsf4 in females and esrrα in males were consistent between qRT-PCR and RNA-Seq. For the other seven genes, RNA-Seq and qRT-PCR showed similar significant up- and down-regulation patterns (Fig. 4). These results indicated that gene expression patterns determined by qRT-PCR were consistent with those determined by RNA-seq analysis, which supports the reliability of our transcriptome data.
Water temperature is an important external factor, which can significantly influence various physiological processes in fish. Females and males may have a sex-specific adaptation to the changing environments, which is a form of sexual dimorphism . To access the gonadal performance in different sexes under heat stress, high-temperature stimulation in female, pseudomale and male C. semilaevis were performed in the present study, and the gonadal transcriptome data were analyzed.
In C. semilaevis, the response to acute high temperature was different among sexes. Females exhibited more DEGs than pseudomales and males. Similarly, Ribas et al. reported that females exhibited more DEGs than males in heat-treated zebrafish . Interestingly, differences among sexes are so obvious in only 48-h acute heat treatment. In C. semilaevis, ovarian differentiation can be observed earlier than the testicular differentiation. Histological observations show that gonadal differentiation begins at ~ 62 days post hatching (dph) . At this stage, primordial germ cells (PGCs) begin mitosis rapidly to form clusters of oocytes. The appearance of the ovarian cavity is observed at ~ 100 dph. On the other hand, in deduced testis, PGCs start mitosis at ~ 80 dph to develop into spermatogonia, and testicular differentiation begins at ~ 100 dah . In this study, 90 dph females showed more variable responses to heat stress, indicating that the early gonadal differentiation is vulnerable to high temperature, and different cell types may have different performance under heat stress.
The 1048 DEGs shared in females, pesudomales and males may represent the inherent difference between high and normal temperatures. Of these, several members of hsp gene family were found, which is consistent with previous studies [26, 27]. The Hsps, mainly Hsp70s and their chaperons, can re-fold thermally damaged proteins and prevent their cytotoxic aggregation [28, 29]. Potential regulators of hsp gene family were differentially expressed under high temperature. For example, RNA-seq showed that eef1akmt3 was up-regulated in females and pseudomales, but down-regulated in males. The eef1akmt4 was down-regulated in all three sexes, both in RNA-seq and qRT-PCR validation. The translation elongation factor, eEF1A participates in heat shock response (HSR) in mammalian cells. It can rapidly activate HSP70 transcription by recruiting HSF1 to its promoter, stabilize HSP70 mRNA and facilitate its transport from the nucleus to ribosomes . The eef1akmt3 and eef1akmt4 are methyltransferase-related genes, which specifically target Lys-165 and Lys-36 in eEF1A [31, 32]. The differentially expressed eef1akmt3 and eef1akmt4 genes may indicate the change in the eEF1A methylation state and the eEF1A mRNA expression level.
In addition, some DEGs are sex-specific under heat stress. In transcriptome data of females, kdm6b was down-regulated after heat stress. kdm6b is an epigenetic regulator demethylating H3K27me3 at dmrt1 promoter, which plays a causal role in male sex determination in red-eared slider turtle (Trachemys scripta elegans). The kdm6b can rapidly respond to temperature change, and it can be down-regulated in the high temperature . In turtles, high temperature increased the production of females, while it increased the production of males in C. semilaevis. Our results may indicate kdm6b as a conserved responder to high temperatures among fish and reptiles. Its role in promoting masculinization of C. semilaevis under high temperatures deserves further study. Another candidate epigenetic regulator, saml, was specifically up-regulated in heat-treated pseudomales both in RNA-seq and qRT-PCR analyses. The S-adenosylmethionine synthase enzyme is sensitive to abiotic stress and catalyzes the S-adenosylmethionine synthesis, which provides methyl groups for the methylation of DNA, RNA and proteins . These results may indicate the role of epigenetic regulation during the physiological response to high temperature.
Two cortisol related genes, hsd11b2 and hsd17b1 were specifically down-regulated in heat-treated females, both in RNA-seq and qRT-PCR analyses. Cortisol is a stress-related hormone, which plays an important role in high-temperature induced masculinization [9,10,11]. There is crosstalk between cortisol and androgen synthesis . The esrrα and ar genes were down-regulated in heat-treated females, while errγl was specifically down-regulated in males. These changes in sex hormones indicate the role of cortisol and androgen pathways response to high temperature. However, expression of sex-related genes, such as cyp19a1a, foxl2, dmrt1 and amh was not significantly different under heat stress. This may show that 48-h heat stress can only affect the upstream hierarchy but cannot change the expression level of sex-related genes.
Our study revealed different performances among sexes under heat stress in C. semilaevis. Besides, several genes were identified to participate in high temperature perception by comparative transcriptome analysis, and some of them showed a sex-biased characteristic. These results may provide a clue for understanding the responses of genes and networks to environmental temperature and expand our current understanding of the effects of gonadal gene expression in teleosts. However, future studies are required to determine the role of these temperature-sensitive genes on sex determination and differentiation.
Heat stress treatment and sample collection
Three-month-old sex determined C. semilaevis were used in this study. The number of fish was determined in order to sample a sufficient number of pseudomales. The genetic female: male ratio was close to 1:1, and approximately 14% of genetic females could spontaneously sex revered into pseudomales [20, 21]. At least five pseudomales of both control and heat-stress groups with three biological replicates for each group were randomly selected for transcriptome sequencing. Thus, a total number of 140 fish were used in the experiment. Fish were obtained from Laizhou Mingbo Aquatic Co., Ltd., Yantai, China (body length: 5.457 ± 0.663 cm), transported to the laboratory and reared in 120 L tanks at 22 °C. After the acclimation period, fish were randomly divided into two groups. The heat stress group (HS, n = 70) was exposed to high temperature (28 °C) with a gradual increase of temperature (2 °C per hour) up to 48 h. The control group (CT, n = 70) was maintained at 22 °C. The temperature treatment used in this study would not induce pseudomale formation. After treatment, fish were anesthetized with 0.05% MS-222 (Sigma, Shanghai, China) via immersion bath and decapitation . The gonad of each fish was dissected and immediately stored in liquid nitrogen for RNA extraction, and the caudal fin was stored in ethanol for genetic sex determination.
Sex identification of C. semilaevis
The sex of each fish was identified after the heat treatment. Genomic DNA (gDNA) was extracted from the caudal fin using the phenol-chloroform method. The genetic sex was identified by PCR using sex-F and sex-R primers, which amplify two bands of 169 and 134 bp in females, and one band of 169 bp in males . The pseudomales were distinguished from females by qRT-PCR against the sex-determining gene, dmrt1  (Additional file 4).
Total RNA extraction, cDNA library construction and sequencing
Total RNA from each gonad was extracted using Trizol (Invitrogen, Carlsbad, USA) according to the manufacturer’s instructions. The quantity and quality of the RNA were determined using an Agilent 2100 bioanalyzer (Thermo Fisher Scientific, Santa Clara, USA). High-quality RNA (RIN > 7) was used for mRNA library construction via the conventional protocol. Three biological replicates for each group were used. The libraries of heat-stressed females (HS_F1–3), males (HS_M1–3) and pseudomales (HS_P1–3) were sequenced on BGISEQ-500 platform by PTM Biolabs (Hangzhou, China). The raw reads were deposited in NCBI’s Sequence Read Archive (SRA) database (BioProject no. PRJNA605682). The transcriptome data of CT group were reported in our previous study  and downloaded from SRA database (BioProject no. PRJNA576366).
Quality control, mapping and annotation of sequencing reads
The quality of the raw reads was assessed using SOAPnuke v1.4.0 . Adapters and low-quality bases were trimmed by Trimmomatic v0.36 . All filtered clean reads from CT and HS groups were mapped to the C. semilaevis reference genome (NCBI, Cse_v1.0, BioProject no. PRJNA73987) using HISAT2 v2.1.0  and aligned with the reference transcript sequence by Bowtie2 v2.2.5 .
Identification of DEGs
Gene expression levels were estimated using RNA-Seq by Expectation Maximization, RSEM v1.2.8 , and Fragments Per Kilobase Million (FPKM) was calculated to represent the expression level of each gene. The DEGseq and DEseq2 were used to determine the DEGs [43, 44]. Genes with an adjusted p-value of less than 0.001 and fold change greater than 2 were defined as DEGs. PCA was performed using OmicShare tools (www.omicshare.com/tools). The heatmap of sex-biased genes was drawn by R Package pheatmap v1.0.12 (https://CRAN.R-project.org/package=pheatmap). The average expression level of three biological replicates was used to represent the corresponding group and conduct the PCA, DEG, GO enrichment and KEGG enrichment analyses.
GO and KEGG pathway enrichment analyses
GO (http://www.geneontology.org/) and KEGG (http://www.kegg.jp/) enrichment analyses of annotated DEGs were performed by Phyper (http://en.wikipedia.org/wiki/Hypergeometric_distribution) based on Hypergeometric test. The threshold of p < 0.05 is considered to be significant.
RNA-seq data validation by qRT-PCR
Nine candidate genes related to high-temperature response (hspa1, hsf4, hsc70, eef1akmt4, saml, hsd11b2, hsd17b1, esrrα and ar) were selected to validate the RNA-seq data by qRT-PCR. The primers were designed based on their sequences from the NCBI database (Additional file 5). β-actin gene was used as the internal control. One microgram of total RNA for high-throughput transcriptome sequencing was reverse transcribed into cDNA with the PrimeScript™ RT reagent Kit with gDNA Eraser (Takara, Japan). Then, qRT-PCR was performed using QuantiNova™ SYBR Green PCR Kit (Qiagen, Germany) in 20-μl reactions, containing 10 μl 2 × SYBR Green PCR Master Mix, 2 μl QN ROX Reference Dye, 0.7 μM forward primer, 0.7 μM reverse primer and 1 μl cDNA. The cycling program was carried out at 95 °C for 2 min, followed by 40 cycles of 95 °C for 5 s and 60 °C for 10 s; this was followed by a melting curve analysis in an ABI StepOnePlus Real-Time PCR system (Applied Biosystems, USA). Reactions were performed in triplicate. The relative expression fold changes of these genes were analyzed using the 2−ΔΔCt method.
The results of the qRT-PCR analysis were expressed as means ± S.E.M.. The values were compared by multiple t-test using GraphPad Prism 7.0 (GraphPad, USA). Statistically significant differences were defined as p < 0.05.
Availability of data and materials
The transcriptome data used in this study were uploaded to the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA) with accession numbers PRJNA576366 and PRJNA605682. The high-quality C. semilaevis reference genome was downloaded from NCBI (Cse_v1.0, BioProject no. PRJNA73987). The sequences of genes selected for RNA-seq validating were acquired from the NCBI database with their corresponding GeneIDs [hspa1 (103393458), hsf4 (103379546), hsc70 (103395053), eef1akmt4 (103395991), saml (103387328), hsd11b2 (103378825), hsd17b1 (103393165), esrrα (103394741), ar (103395329), β-actin (103393304), dmrt1 (103397807)].
Days after hatching
Differentially expressed gene
Fragments per kilobase million
National Center for Biotechnology Information
Principal component analysis
Heat stress group
Heat shock transcription factor 1
Heat shock protein
Heat shock response
Kyoto Encyclopedia of Genes and Genomes
Primordial germ cells
Real-time quantitative-polymerase chain reaction
Sequence Read Archive
Sylvester JR. Possible effects of thermal effluents on fish: a review. Environ Pollut. 1972;3:205–15.
Dadras H, Dzyuba B, Cosson J, Golpour A, Siddique MAM, Linhart O. Effect of water temperature on the physiology of fish spermatozoon function: a brief review. Aquac Res. 2017;48:729–40.
Mimura N. Sea-level rise caused by climate change and its implications for society. Proc Jpn Acad Ser B Phys Biol Sci. 2013;89:281–301.
Hsu SK, Jakšić AM, Nolte V, Lirakis M, Kofler R, Barghi N, et al. Rapid sex-specific adaptation to high temperature in Drosophila. Elife. 2020;9:e53237.
Chandel G, Dubey M, Meena R. Differential expression of heat shock proteins and heat stress transcription factor genes in rice exposed to different levels of heat stress. J Plant Biochem Biotechnol. 2013;22:277–85.
Oksala NK, Ekmekçi FG, Ozsoy E, Kirankaya S, Kokkola T, Emecen G, et al. Natural thermal adaptation increases heat shock protein levels and decreases oxidative stress. Redox Biol. 2014;3:25–8.
Purohit GK, Mahanty A, Suar M, Sharma AP, Mohanty BP, Mohanty S. Investigating hsp gene expression in liver of Channa striatus under heat stress for understanding the upper thermal acclimation. Biomed Res Int. 2014;2014:381719.
Healy TM, Tymchuk WE, Osborne EJ, Schulte PM. Heat shock response of killifish (Fundulus heteroclitus): candidate gene and heterologous microarray approaches. Physiol Genomics. 2010;41:171–84.
Hattori RS, Fernandino JI, Kishii A, Kimura H, Kinno T, Oura M, et al. Cortisol-induced masculinization: does thermal stress affect gonadal fate in pejerrey, a teleost fish with temperature-dependent sex determination? PLoS One. 2009;4:e6548.
Yamaguchi T, Yoshinaga N, Yazawa T, Gen K, Kitano T. Cortisol is involved in temperature-dependent sex determination in the Japanese flounder. Endocrinology. 2010;151:3900–8.
Hayashi Y, Kobira H, Yamaguchi T, Shiraishi E, Yazawa T, Hirai T, et al. High temperature causes masculinization of genetically female medaka by elevation of cortisol. Mol Reprod Dev. 2010;77:679–86.
Deveson IW, Holleley CE, Blackburn J, JAM G, Mattick JS, Waters PD, et al. Differential intron retention in Jumonji chromatin modifier genes is implicated in reptile temperature-dependent sex determination. Sci Adv. 2017;3:e1700731.
Ge CT, Ye J, Weber C, Sun W, Zhang HY, Zhou YJ, et al. The histone demethylase KDM6B regulates temperature-dependent sex determination in a turtle species. Science. 2018;360:645–8.
Evans TG. Considerations for the use of transcriptomics in identifying the 'genes that matter' for environmental adaptation. J Exp Biol. 2015;218:1925–35.
López-Maury L, Marguerat S, Bähler J. Tuning gene expression to changing environments: from rapid responses to evolutionary adaptation. Nat Rev Genet. 2008;9:583–93.
Mahat DB, Salamanca HH, Duarte FM, Danko CG, Lis JT. Mammalian heat shock response and mechanisms underlying its genome-wide transcriptional regulation. Mol Cell. 2016;62:63–78.
Diaz N, Piferrer F. Lasting effects of early exposure to temperature on the gonadal transcriptome at the time of sex differentiation in the European sea bass, a fish with mixed genetic and environmental sex determination. BMC Genomics. 2015;16:679.
Sun LX, Teng J, Zhao Y, Li N, Wang H, Ji XS. Gonad Transcriptome analysis of high-temperature-treated females and high-temperature-induced sex-reversed neomales in Nile tilapia. Int J Mol Sci. 2018;19:E689.
Fernandino JI, Popesku JT, Paul-Prasanth B, Xiong H, Hattori RS, Oura M, et al. Analysis of sexually dimorphic expression of genes at early gonadogenesis of pejerrey Odontesthes bonariensis using a heterologous microarray. Sex Dev. 2011;5:89–101.
Shao CW, Li QY, Chen SL, Zhang P, Lian JM, Hu QM, et al. Epigenetic modification and inheritance in sexual reversal of fish. Genome Res. 2014;24:604–15.
Chen SL, Zhang GJ, Shao CW, Huang QF, 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:253–60.
Liu Y, Chen SL, Gao FT, Meng L, Hu QM, Song WT, et al. SCAR-transformation of sex-specific SSR marker and its application in half-smooth tongue sole (Cynoglossus semiliaevis). J Agr Biotechnol. 2014;65:S22.
Ribas L, Liew WC, Diaz N, Sreenivasan R, Orban L, Piferrer F. Heat-induced masculinization in domesticated zebrafish is family-specific and yields a set of different gonadal transcriptomes. Proc Natl Acad Sci U S A. 2017;114:E941–50.
Zhu Y, Hu Q, Xu W, Li H, Guo H, Meng L, et al. Identification and analysis of the β-catenin1 gene in half-smooth tongue sole (Cynoglossus semilaevis). PLoS One. 2017;12:e0176122.
Ma XK, Liu XZ, Wen HS, Xu YJ, Zhang LJ. Histological observation on gonadal sex differentiation in Cynoglossus semilaevis Günther. Mar Fish Res. 2006;27:55–61 (in Chinese).
Poonlaphdecha S, Pepey E, Canonne M, de Verdal H, Baroiller JF, D'Cotta H. Temperature induced-masculinisation in the Nile tilapia causes rapid up-regulation of both dmrt1 and amh expressions. Gen Comp Endocrinol. 2013;193:234–42.
Lee KH, Yamaguchi A, Rashid H, Kadomura K, Yasumoto S, Matsuyama M. Germ cell degeneration in high-temperature treated pufferfish, Takifugu rubripes. Sex Dev. 2009;3:225–32.
Logan CA, Buckley BA. Transcriptomic responses to environmental temperature in eurythermal and stenothermal fishes. J Exp Biol. 2015;218:1915–24.
Rosenzweig R, Nillegoda NB, Mayer MP, Bukau B. The Hsp70 chaperone network. Nat Rev Mol Cell Biol. 2019;20:665–80.
Vera M, Pani B, Griffiths LA, Muchardt C, Abbott CM, Singer RH, et al. The translation elongation factor eEF1A1 couples transcription to translation during heat shock response. Elife. 2014;3:e03164.
Malecki J, Aileni VK, Ho AYY, Schwarz J, Moen A, Sørensen V, et al. The novel lysine specific methyltransferase METTL21B affects mRNA translation through inducible and dynamic methylation of Lys-165 in human eukaryotic elongation factor 1 alpha (eEF1A). Nucleic Acids Res. 2017;45:4370–89.
Jakobsson ME, Małecki J, Falnes P. Regulation of eukaryotic elongation factor 1 alpha (eEF1A) by dynamic lysine methylation. RNA Biol. 2018;15:314–9.
Bauerle MR, Schwalm EL, Booker SJ. Mechanistic diversity of radical S-adenosylmethionine (SAM)-dependent methylation. J Biol Chem. 2015;290:3995–4002.
Fernandino JJ, Hattori RS, Kishii A, Strussmann CA, Somoza GM. The cortisol and androgen pathways cross talk in high temperature-induced masculinization: the 11beta-hydroxysteroid dehydrogenase as a key enzyme. Endocrinology. 2012;153:6003–11.
Xu YJ, Zhang YX, Wang B, Liu XZ, Liu Q, Song XS, et al. Leptin and leptin receptor genes in tongue sole (Cynoglossus semilaevis): molecular cloning, tissue distribution and differential regulation of these genes by sex steroids. Comp Biochem Phys A. 2018;224:11–22.
Cui ZK, Liu Y, Wang WW, Wang Q, Zhang N, Lin F, et al. Genome editing reveals dmrt1 as an essential male sex-determining gene in Chinese tongue sole (Cynoglossus semilaevis). Sci Rep. 2017;7:42213.
Wang Q, Liu KQ, Feng B, Zhang ZH, Wang RK, Tang LL, et al. Transcriptome of gonads from high temperature induced sex reversal during sex determination and differentiation in Chinese tongue sole, Cynoglossus semilaevis. Front Genet. 2019;10:1128.
Chen Y, Chen Y, Shi C, Huang Z, Zhang Y, Li S, et al. SOAPnuke: a MapReduce acceleration-supported software for integrated quality control and preprocessing of high-throughput sequencing data. Gigascience. 2018;7:1–6.
Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20.
Kim D, Langmead B, Salzberg SL. HISAT: a fast spliced aligner with low memory requirements. Nat Methods. 2015;12:357–60.
Langmead B, Salzberg SL. Fast gapped-read alignment with bowtie 2. Nat Methods. 2012;9:357–9.
Li B, Dewey CN. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics. 2011;12:323.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.
Wang L, Feng Z, Wang X, Wang X, Zhang X. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010;26:136–8.
This work was supported by the National Nature Science Foundation of China [grant numbers 31722058, 31802275 and 31472269]; the National Key R&D Program of China [grant number 2018YFD0900301]; the AoShan Talents Cultivation Program Supported by Qingdao National Laboratory for Marine Science and Technology [grant number 2017ASTCP-ES06]; the Taishan Scholar Project Fund of Shandong of China to C.S.; the National Ten-Thousands Talents Special Support Program to C.S.; the Central Public-interest Scientific Institution Basal Research Fund, CAFS [grant number No.2020TD19]; and the International Scientific Partnership Program ISPP at King Saud University [grant number No. 0050]. The funding bodies did not participate in the design of the study, the collection, analysis, and interpretation of data and in the writing the manuscript.
Ethics approval and consent to participate
The collection and handling of fish and experimental procedures were performed in accordance with the Guidelines for Experimental Animals of the Ministry of Science and Technology (Beijing, China), and approved by the Institutional Animal Care and Use Committee, IACUC of Yellow Sea Fisheries Research Institute, CAFS (Qingdao, China) (No. YSFRI-2020005).
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.
Statistics of transcriptome sequencing. Data in females, pseudomales and males from each group (control group: CT_F, CT_P and CT_M; heat-stress group: HS_F, HS_P and HS_M) were listed.
Shared and sex-specific genes involved in high temperature response. Interested differentially expressed genes (DEGs) shared in females, pseudomales and males, and DEGs specifically changed in females, pseudomales or males between control (CT) and heat-stress (HS) groups were listed.
Enriched Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways of sex-specific DEGs under heat stress. The differentially expressed genes (DEGs) specifically changed in females, pseudomales, or males in comparison of heat-stress (HS) vs. control (CT) groups were subjected to GO enrichment analysis and resulted significant GO terms (p < 0.05) were listed in Sheet 1. KEGG pathway enrichment analysis of these sex-specific DEGs was listed in Sheet 2 (p < 0.05).
Sex identification of C. semilaevis. (A, B) Genetic sex identification of C. semilaevis using sex-linked SSR marker. (A) Genetic sex identification of C. semilaevis in the control (CT) group; (B) Genetic sex identification of C. semilaevis in the heat-stress (HS) group. Genetic females have two sex-linked SSR bands, while genetic males have only one band. (C, D) Relative expression level of male-determining gene dmrt1 in females, pseudomales and males. (C) Relative expression level of dmrt1 in CT group; (D) Relative expression level of dmrt1 in HS group.
The primers used in this study. The column of Gene ID lists the accession number of each gene in the NCBI database.
About this article
Cite this article
Wang, Q., Hao, X., Liu, K. et al. Early response to heat stress in Chinese tongue sole (Cynoglossus semilaevis): performance of different sexes, candidate genes and networks. BMC Genomics 21, 745 (2020). https://doi.org/10.1186/s12864-020-07157-x
- Transcriptome analysis
- Heat stress
- Heat shock protein
- Epigenetic regulator