Regulating Strategies of Transcription and Alternative Splicing for Cold Tolerance Harpadon nehereus Fish

In recent years, Harpadon nehereus gradually become a dominant species with great potential for exploitation in the East China Sea, and it is worth investigating whether H. nehereus would tolerate cold stress to continue to expand into the colder northern waters. The molecular regulation level is favorable evidence to explore the cold tolerance of H. nehereus, a total of 6,650, 1,936, and 2,772 differentially expressed genes (DEGs) in transcription regulation, and 4,409, 1,250, and 2,303 differential alternative splicing genes (DASGs) in alternative splicing regulation were identified in H. nehereus at 13, 15, and 17°C, respectively, importantly, 47 genes were identified as the key candidate genes for cold tolerance in H. nehereus. In transcription regulation, up-regulated DEGs were enriched in metabolic process terms and ribosome, spliceosome pathway, etc., while down-regulated DEGs were enriched in signal transduction terms, focal adhesion, proteoglycans in cancer pathway, etc., at 13, 15, and 17°C, respectively. In alternative splicing regulation, spliceosome, mRNA surveillance pathway, etc., were significantly enriched in DASGs. In a word, H. nehereus adapts to cold environments mainly through transcription and translation, transmembrane transport, protein modification, etc., while cold stress may also induce some diseases in H. nehereus.


INTRODUCTION
Global warming is causing a series of climate changes that threaten the survival and distribution of fish, such as ocean warming and frequent extreme winter cold events (Fossheim et al., 2015;Johnson et al., 2018). Water temperature is the master factor to determine nearly all life activities of most fishes, including growth, development, reproduction, behavior, and geographical distribution (Lopez-Olmeda and Sanchez-Vazquez, 2011;Poloczanska et al., 2016). Like other ectotherms, changes in fish distribution areas may vary with daily and seasonal water temperature variations (Morgan et al., 2013). Currently, the impact of ocean warming on fish is a hot topic of interest for researchers, which often overlooks the survival and distribution pressure on fish caused by extremely winter cold temperatures triggered by climate warming. The research on the cold tolerance of fish was more concerned with cultured economics (Xu et al., 2018;Sun et al., 2019;Jiao et al., 2020) and model fish species (Valdez et al., 2005;Ren et al., 2021), while studies on natural marine fish were relatively scarce. In addition, the study found that fish can adapt to variations in water temperature through physiological plasticity or microevolution (Whitehead et al., 2013), and the adaptive evolution of temperature tolerance determines their ecological status and distribution range, so the study of cold tolerance of fish in natural marine areas has significance to estimate the variations in their distribution range.
Lizardfish Harpadon nehereus, also known as Bombay duck, is mostly discontinuously distributed in the benthic along the coasts and estuaries of East Africa, Asia, and Indonesia, e.g., the Ganges Estuary, the Bay of Bengal, and the Cambay Bay (Bapat, 1970). In China, H. nehereus is widely distributed in the south of the Yellow Sea, the East China Sea (Liu and Cheng, 2015), and the South China Sea (Zhu et al., 2014). The morphological characteristics of H. nehereus including small eyes, wide mouth, incomplete skull, membranous skin, and gelatinous body with high moisture (>90%) (Roy et al., 2014) give it the resemblance to deep-sea fishes , whereas H. nehereus inhabits offshore sandy mud bottom in most of a year. Historically, H. nehereus have been overlooked for its low value and less yield in China. In recent decades, with the marine environmental change and the traditional fisheries decline, the population size of H. nehereus in the East China Sea has rapidly expanded (Lin, 2009). Resource investigations of fish communities indicate that H. nehereus has become dominant species in newly expanded areas of its current range along the southern coast of China (Lin, 2009;Sun et al., 2015). The growth of H. nehereus populations may be accompanied by a lack of survival resources and interspecific competition, and it may have forced them to expand to colder northern waters for survival, deserving detailed studies for more evidence besides biological habits (Firdaus et al., 2018), and nutritional processing (Taqwa et al., 2020). This study not only provides molecular insights into the cold tolerance of H. nehereus, but also enriches the genetic information.
The molecular response mechanism of fish to cold stress is subject to complex regulation at multiple levels, e.g., epigenetic modifications, transcription, translation, alternative splicing, and post-translational modifications. Exploring the cold tolerance of fish at the molecular level is the basis for the study of fish adaptation to low-temperature environments. With the development of high-throughput sequencing technology and bioinformatics, RNA-Seq has developed rapidly into a powerful tool for better understanding the candidate genes and potential biological pathways in responding to different environmental stresses (Sun et al., 2020;Sun et al., 2021), as well as predicting new transcripts, simple sequence repeats (SSR), and alternative splicing (AS) (Lan et al., 2014). Although the flexibility of AS contributes to environmental adaption and phenotypic plasticity in organisms (Mastrangelo et al., 2012), studies on cold tolerance of fish were mostly at the transcriptional level, such as Nibea albiflora (Xu et al., 2018), Epinephelus coioides (Sun et al., 2019), Larimichthys crocea (Mininni et al., 2014), Larimichthys polyactis (Liu F. et al., 2020), and Nibea albiflora (Jiao et al., 2020), and lack of studies on alternative splicing. Therefore, the study of cold tolerance in H. nehereus from both transcription and alternative splicing levels is certainly an innovative way that could search for more evidence from multiple regulation levels.
In this study, we performed the transcription and alternative splicing regulation to explore the regulating strategies of cold tolerance in H. nehereus. This study is the first to reveal the molecular mechanism of cold tolerance in H. nehereus, and fill a partial gap in transcription and alternative splicing research on cold tolerance in fish.

Experimental Animals
From November 2020 to January 2021, H. nehereus (average body length: 159.71 ± 24.16 mm, weight: 50.53 ± 12.30 g) were collected from the East China Sea, they were harvested from 13, 15, 17, and 24 • C natural marine areas. During the sampling process, we tried to ensure that temperature was the only variable in different sampling areas and other environmental indicators were similar, the location and environment information of sampling sites were presented in Supplementary Table 1. We immediately anesthetized the freshly caught H. nehereus with ice in the sampling process, followed by rapid dissection of muscle (M), and put them into liquid nitrogen for preservation. After disembarkation, they were transferred to −80 • C refrigerator for storage and subsequent RNA extraction. Each group had three biological duplications and a total of 12 samples were taken for RNA-Seq.

RNA Extraction, Library Preparation, and Sequencing
Total RNA was isolated from tissues sampled by using TRIZOL Reagent (Invitrogen, United States) and treated with RNase-free DNase I (Takara, China) at 37 • C for 1 h to remove residual genomic DNA. RNA purity was checked using NanoPhotometer R spectrophotometer (IMPLEN, CA, United States). RNA concentration was measured using Qubit R RNA Assay Kit in Qubit R 2.0 Flurometer (Life Technologies, CA, United States). Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, United States) was used to assess the quality of extracted RNA. Only high-quality RNA samples were used to construct cDNA libraries, and the libraries were sequenced by paired-end sequencing on the Illumina Novaseq6000 sequencing platform (Illumina, United States).

Identification and Enrichment Analysis of Differentially Expressed Genes
Firstly, clean reads from raw data were obtained using Fastp software by removing reads containing adapter, reads containing ploy-N, and low-quality reads (Chen et al., 2018). Subsequently, the clean reads were mapped to the H. nehereus reference genome using STAR software (Dobin et al., 2013). Based on the numbers of clean reads mapped to the genome sequence, the differential expressed level of genes was calculated using featureCounts (Liao et al., 2014). The expression of each transcript was normalized using the transcripts per million (TPM) algorithm, and then read-count data was subjected to differential expressed analysis using the DESeq2 R package (Love et al., 2014). The differentially expressed genes (DEGs) were identified by | log 2 FC| ≥ 1 (which corresponds to Fold Chang ≥ 2 or ≤ −2) and p-adjusted (padj) ≤ 0.01. Moreover, DEGs were divided into up-regulated DEGs and down-regulated DEGs, and Gene Ontology (GO) (Ashburner et al., 2000) and Kyoto Encyclopedia of Genes and Genomes (KEGG) (Kanehisa and Goto, 2000) pathway enrichment analyses (p-adj < 0.05) were performed for up-and down-regulated DEGs, respectively.

Identification of Differential Alternative Splicing
In addition to transcription regulation, mRNA data was used to study alternative splicing regulation by evaluating exon expression. Differential exon usage was assessed using the DEXSeq R package. The operational steps are described in detail in the DEXSeq manual. 1 Threshold as p-adjusted (p-adj) ≤ 0.01 and log 2 FC ≥ 2 was set for screening DAS.

Ethics Statement
All animal experiments were conducted following the guidelines and approval of the respective Animal Research and Ethics Committees of the Ocean University of China, and the field studies did not involve endangered or protected species.

Base Data Statistics About Transcriptome
In this study, 12 cDNA libraries were constructed and sequenced, and clean Q20 and Q30 were over 95 and 90%, respectively. The range of clean reads mapped to the genome sequences was 84.41-97.39%, among which uniquely mapped reads over 75% (Supplementary Table 2). All the above results proved that the data were highly accurate, and ensured the reliability of subsequent analysis results.

Analysis of Differential Alternative Splicing
As a result, 18,062, 2,674, and 6,424 exons in 4,409, 1,250, and 2,303 genes were differential alternative splicing (DAS) at 13, 15, and 17 • C, respectively, among which 382 DAS genes (DASGs) were common in three low-temperature environments (Figure 2). Moreover, 47 genes were common to both DEGs and DASGs, among which 19 genes were upregulated expressed in both transcription and alternative splicing regulation ( Table 1). With the combination of GO annotation and manual literature searches, 19 genes were grouped roughly into four categories, including transcription and translation, transmembrane transport, posttranslational modification, and others ( Table 1).
Kyoto encyclopedia of genes and genomes enrichment results showed that DASGs were enriched in 34, 12, and 24 pathways at 13, 15, and 17 • C, respectively. "Spliceosome, " "mRNA surveillance pathway, " and other five pathways were significantly enriched in three low-temperature environments (Figure 6).

DISCUSSION
Overall, numerous DEGs and their relevant pathways were involved in cold stresses, and the number of down-regulated DEGs was much more than up-regulated DEGs in transcription regulation, this may reduce the energy requirement for body storage in low-temperature environments (Wen et al., 2017).  Taking all the results together, H. nehereus adapts to lowtemperature environments mainly through transcription and translation, transmembrane transport, and protein modification.

Transcription and Translation
The temperature has a decisive influence on the structure and function of cellular components such as nucleic acids and proteins (Somero and Hochachka, 1971). When fish faces cold stress, biomolecules such as nucleic acids and proteins do not fold and assemble correctly, resulting in reduced activity (Zhang et al., 2018). In this study, UDEGs were involved in many biological processes and pathways related to transcription and translation such as translation, ribosome biogenesis, ribosome, spliceosome, and nucleocytoplasmic transport (Figures 3A, 4), similar to the results of Danio rerio (Long et al., 2013) and Oryzias latipes (Ikeda et al., 2017). The spliceosome pathway played an important role in regulating the adaptation of H. nehereus to cold stress, so alternative splicing analysis was of great significance in the study of cold tolerance of H. nehereus. DASGs enrichment results showed that transcription and translation related pathways such as spliceosome, ribosome, and mRNA surveillance pathway were also significantly enriched. Up-regulated DEGs and DASGs in ribosome biogenesis and pathways may be a compensatory response to the cell's ability to enhance protein synthesis under low-temperature environments.
In addition, Rpap3, Nploc4, Endouc, Cnot9, Exosc2, Rplp0, Wdr43, and Exoc5 proteins were classified as transcription and translation. Rpap3 functions as co-chaperone for both Hsp70 and Hsp90 (Benbahouche et al., 2014), the increase of Rpap3 appeared to be accompanied by an increase in Hsp70 and Hsp90 proteins, but this did not turn out to be the case, and Hsp70 and Hsp90 were reduced in H. nehereus from cold environments. Rpap3 may play a more role in regulating cell survival and apoptosis (Benbahouche et al., 2014), renewing cells to maintain normal physiological functions of H. nehereus cells  under low-temperature environments. Endou proteins belong to the Eukaryotic Endou ribonuclease family of enzymes. Recent studies demonstrated that Endou plays important roles in cellular processes, including regulating ER structure, RNA degradation, and cell survival (Poe et al., 2014;Schwarz and Blower, 2014). The reason for the increase of Endouc might also be the maintenance of cell normal physiological processes. Noteworthy, autophagy, as a cell cleaner (Galluzzi and Green, 2019), was heavily enriched in the results (Figures 4A, 6A-C and Table 1), this result reinforced our speculation that H. nehereus maintained normal cellular physiological functions by accelerating the removal of damaged cell structures by autophagy and the renewal of damaged cell structures by associated proteins, such as Rpap3 and Endou. In addition, studies showed that Exosc2 knockout zebrafish showed larval lethality 13 days post-fertilization, with microcephaly neurons, myelin deficiency, and loss of spinal motor neurons (Yatsuka et al., 2021). Furthermore, Exosc8, which is in the same family as Exosc2, its mutation causes hypomyelination with spinal muscular atrophy (Boczonadi et al., 2014), combined with the result that thecardiac muscle contraction pathway and Cald1 protein (Zheng et al., 2009) were significantly expressed ( Figure 6A and Table 1), we suspected that the increase of Excos2 was to avoid excessive muscle contraction in H. nehereus under low-temperature environments.

Transmembrane Transport
Low-temperature stimulation can reduce the fluidity of biofilm (including cell membrane and organelle membrane), change the structure of biofilm, and affect the function of membrane and membrane binding protein (Los and Murata, 2004). The phenomenon that cold stress affects the fluidity and structure of fish biofilms seems to be common, e.g., Takifugu fasciatus (Wen et al., 2019), O. latipes (Ikeda et al., 2017), and Pearl gentian grouper (Miao et al., 2021). Therefore, the cold stress may have affected the fluidity and structure of the biofilm in H. nehereus, which formed a protective mechanism to maintain the normal physiological function of the biofilm through the increase of transmembrane transport-related proteins such as Slc7a8 and Slc7a6 ( Table 1). The increase in Slc7a8 and Slc7a6 proteins seems unsurprising, and the same family as them, Slc7a11 (Zhou et al., 2019), Slc25a5 (Ren et al., 2021), Slc16a1 (Long et al., 2015), and Slc2a10 (Scott and Johnston, 2012) were also found to be differentially expressed in other fish under low-temperature environments. Slc7a as an amino acid transporter (Fotiadis et al., 2013), the increase of Slc7a6 and Slc7a8 may also be due to the need for more frequent material exchange of cells inside and outside for H. nehereus in low-temperature environments.
Protein export and endocytosis pathway in KEGG, and protein metabolisms, peptide metabolic process, peptide  biosynthetic process, cellular amide metabolic process, and amide biosynthetic process terms in GO that require transmembrane transport were also significantly enriched (Figures 3A, 6A).
The above results indicated that a large number of stress proteins were synthesized in H. nehereus cells and secreted to act extracellularly, while endocytosis could break down damaged macromolecules for recycling. The secretion of stress proteins and the endocytosis of macromolecules are inseparable from the transmembrane transport function. Additionally, the fluidity of the membrane was enhanced by altering the content of membrane lipid and the proportion of unsaturated fatty acids (Johnston and Roots, 1964), and the fatty acid metabolism pathway was enriched probably to maintain the fluidity of the biofilms in this study ( Figure 6C).

Protein Modification
The low temperature could induce protein modifications such as phosphorylation and acetylation, which change proteins' chemical properties and structure and affect their functions. Prmt1, Csn1, Klhl38, Atg9a, and Mlip proteins were classified as protein modification (Table 1), and protein modification related pathways such as proteasome, ubiquitin mediated proteolysis, and protein processing in endoplasmic reticulum pathways were also significantly enriched in low-temperature environments of H. nehereus (Figures 4A, 6A,B). This result was similar to the pathway enriched by other fish under cold stress, like the proteasome pathway in D. rerio (Long et al., 2013) and Pampus argenteus (Zhang et al., 2022), and protein processing in endoplasmic reticulum pathway in D. rerio (Long et al., 2013), Lota lota lota (Yang et al., 2021), and O. latipes (Ikeda et al., 2017). The ubiquitin proteasome system catalyzes the great majority of the protein degradation, including both the rapid degradation of misfolded and regulation proteins (Zhao et al., 2015), therefore, the ubiquitin protease system may clear the misfolded proteins and accelerate the renewal of damaged cells for H. nehereus under low-temperature environments (Han et al., 2020). Although the endoplasmic reticulum is the primary place of protein synthesis, it is also the site of protein degradation (Christianson and Ye, 2014). According to the significant enrichment of proteasome, ubiquitin mediated proteolysis, autophagy, and other protein degradation pathways, we prefer the protein processing in the endoplasmic reticulum pathway to play a role in protein degradation in H. nehereus under low-temperature environments. Many pathways related to protein degradation were significantly enriched, indicating that cold stress may affect the balance of protein synthesis and degradation. Furthermore, UDEGs were enriched into disease-related pathways, for example, transcriptional misregulation in cancer, systemic lupus erythematosus, alcoholism, viral carcinogenesis pathways, and viral myocarditis (Figures 3A, 4), suggesting that cold stress may induce various diseases in H. nehereus; DDEGs were enriched into signal transduction terms, focal adhesion, regulation of actin cytoskeleton, and other pathways ( Figures 3B, 5), indicating that signal transduction, cell adhesion, and cytoskeleton regulation were reduced in H. nehereus under low-temperature environments. Moreover, as stress proteins, the heat shock proteins (HSPs) family does not only play a role in thermal shock but also plays a regulatory role in other stresses, such as heavy metals, osmotic stress, cold stress, etc. (Basu et al., 2002). Strong evidence suggests that HSPs have critical roles in helping fish cope with environmental change (Weber and Bosworth, 2005), but 12 HSP genes were down-regulated differentially expressed in this study, and the differential expression levels of HSP genes seem to be different in different fish coping with cold stress, e.g., the HSP genes were up-regulated expressed in channel catfish, L. lota lota and Puntius tetrazona (Weber and Bosworth, 2005;Liu L. L. et al., 2020;Yang et al., 2021), while down-regulated expressed in O. latipes, D. rerio, and Oratosquilla oratoria (Ikeda et al., 2017;Lou et al., 2019;Han et al., 2020). In a similar situation, the focal adhesion pathway of Gasterosteus aculeatus was down-regulated expressed (Metzger and Schulte, 2018), which was similar to the result of H. nehereus, while L. polyactis was up-regulated expressed (Liu F. et al., 2020). Therefore, the difference in regulation strategies of different fish for adapting to low-temperature environments is worth further exploration.

CONCLUSION
Cold tolerance of fish is important in determining fish distribution and survival in the changing oceans. In this study, a total of 6,650, 1,936, and 2,772 DEGs in transcription regulation and 4,409, 1,250, and 2,303 DASGs in alternative splicing regulation were identified in H. nehereus at 13, 15, and 17 • C, respectively, importantly, 47 genes were identified as the key candidate genes for cold tolerance in H. nehereus. It's worth noting that ribosome and spliceosome pathways associated with transcription and translation; protein export and endocytosis pathways associated with transmembrane transport; proteasome and protein processing in endoplasmic reticulum pathways associated with protein modification; and transcriptional misregulation in cancer and viral carcinogenesis pathway associated with the disease were significantly enriched in H. nehereus under low-temperature environments. This study provided a certain theoretical basis for future research into the molecular strategies of cold tolerance in fishes. Concurrently, this study is also one of the few papers using transcriptome data to study the regulating strategies of fish under natural low-temperature environments, which can reflect the real molecular mechanism of H. nehereus when coping with cold and other complex natural environmental factors.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: NCBI (accession: PRJNA820367), Sequence Read Archive (accession: SRR18494828, SRR18494829, SRR18494830, SRR18494831, SRR18494832, SRR18494833, SRR18494834, SRR18494835, SRR18494836, SRR18494837, SRR18494838, and SRR18494839). Further queries can be directed to the corresponding author(s).

ETHICS STATEMENT
The animal study was reviewed and approved by Animal Research and Ethics Committees of the Ocean University of China.

AUTHOR CONTRIBUTIONS
ZS: methodology, software, formal analysis, data curation, writing−original draft, and visualization. LH: investigation, resources, validation, and funding acquisition. YK: software, formal analysis, and validation. LW: formal analysis, data curation, and validation. BK: conceptualization, validation, writing−review and editing, supervision, project administration, and funding acquisition. All authors contributed to the article and approved the submitted version.

FUNDING
This study was supported by National Natural Science Foundation of China (Nos. U20A2087 and 41976091).