Transcriptomic analyses on muscle tissues of Litopenaeus vannamei provide the first profile insight into the response to low temperature stress

PLOS ONE, Jun 2017

The Pacific white shrimp (Litopenaeus vannamei) is an important cultured crustacean species worldwide. However, little is known about the molecular mechanism of this species involved in the response to cold stress. In this study, four separate RNA-Seq libraries of L. vannamei were generated from 13°C stress and control temperature. Total 29,662 of Unigenes and overall of 19,619 annotated genes were obtained. Three comparisons were carried out among the four libraries, in which 72 of the top 20% of differentially-expressed genes were obtained, 15 GO and 5 KEGG temperature-sensitive pathways were fished out. Catalytic activity (GO: 0003824) and Metabolic pathways (ko01100) were the most annotated GO and KEGG pathways in response to cold stress, respectively. In addition, Calcium, MAPK cascade, Transcription factor and Serine/threonine-protein kinase signal pathway were picked out and clustered. Serine/threonine-protein kinase signal pathway might play more important roles in cold adaptation, while other three signal pathway were not widely transcribed. Our results had summarized the differentially-expressed genes and suggested the major important signaling pathways and related genes. These findings provide the first profile insight into the molecular basis of L. vannamei response to cold stress.

A PDF file should load here. If you do not see its contents the file may be temporarily unavailable at the journal website or you do not have a PDF plug-in installed and enabled in your browser.

Alternatively, you can download the file locally and open with any standalone PDF reader:

http://journals.plos.org/plosone/article/file?id=10.1371/journal.pone.0178604&type=printable

Transcriptomic analyses on muscle tissues of Litopenaeus vannamei provide the first profile insight into the response to low temperature stress

June Transcriptomic analyses on muscle tissues of Litopenaeus vannamei provide the first profile insight into the response to low temperature stress Wen Huang 0 1 Chunhua Ren 0 1 Hongmei Li 0 1 2 Da Huo 0 1 2 Yanhong Wang 0 1 Xiao Jiang 0 1 Yushun Tian 0 1 2 Peng Luo 0 1 Ting Chen 0 1 Chaoqun Hu 0 1 0 CAS Key Laboratory of Tropical Marine Bio-resources and Ecology (LMB), South China Sea Institute of Oceanology, Chinese Academy of Sciences , Guangzhou, Guangdong , China , 2 Key Laboratory of Applied Marine Biology of Guangdong Province and Chinese Academy of Sciences (LAMB), South China Sea Institute of Oceanology, Chinese Academy of Sciences , Guangzhou, Guangdong, China, 3 South China Sea Bio-Resource Exploitation and Utilization Collaborative Innovation Center, Guangzhou, Guangdong , China 1 Editor: Sebastian D. Fugmann, Chang Gung University , TAIWAN 2 University of the Chinese Academy of Sciences , Beijing , China The Pacific white shrimp (Litopenaeus vannamei) is an important cultured crustacean species worldwide. However, little is known about the molecular mechanism of this species involved in the response to cold stress. In this study, four separate RNA-Seq libraries of L. vannamei were generated from 13ÊC stress and control temperature. Total 29,662 of Unigenes and overall of 19,619 annotated genes were obtained. Three comparisons were carried out among the four libraries, in which 72 of the top 20% of differentially-expressed genes were obtained, 15 GO and 5 KEGG temperature-sensitive pathways were fished out. Catalytic activity (GO: 0003824) and Metabolic pathways (ko01100) were the most annotated GO and KEGG pathways in response to cold stress, respectively. In addition, Calcium, MAPK cascade, Transcription factor and Serine/threonine-protein kinase signal pathway were picked out and clustered. Serine/threonine-protein kinase signal pathway might play more important roles in cold adaptation, while other three signal pathway were not widely transcribed. Our results had summarized the differentially-expressed genes and suggested the major important signaling pathways and related genes. These findings provide the first profile insight into the molecular basis of L. vannamei response to cold stress. - Data Availability Statement: All relevant data are within the paper and its Supporting Information files. Introduction The Pacific white shrimp, Litopenaeus vannamei, is one of the most important cultural shrimp species, and its production (by weight) has reached nearly 71% of the total economic penaeid shrimp production worldwide [ 1, 2, 3, 4 ]. Naturally, L. vannamei can inhabit estuaries, lagoons or marine areas, but they are sensitive to low temperatures [1]. With high economic value, the farming areas of L. vannamei have been exploded from the southern to the northern coast, and (A201601A03), the National High Technology Research & Development Program of China (863 Program) (2012AA10A404-4), and the Science and Technology Service Network Initiative of the Chinese Academy of Sciences (KFJ-SW-STS-146). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. to the inland water in China. Therefore, low temperatures have become one of the major constrained factors to t L. vannamei culture, since that the culturing water temperature of the Northern China in winter is usually lower than that needed for L. vannamei survival. Water temperature is one of the most common abiotic stress factors for aquatic ectotherms [ 5 ]. It influences nearly all biochemical and physiological progresses of aquatic organisms [ 6 ], and low temperatures may cause severe growth inhibition and mortality [ 6, 7 ]. In fishes, the understanding of biochemical and molecular mechanisms under temperature stress have been well developed. In species such as rainbow trout (Oncorhynchus mykiss) [8], channel catfish (Ictalurus punctatus) [ 9 ] and common carp (Cyprinus carpio L.) [ 10 ], studies showed that about 10 percent of the investigated genes were found to be up-regulated under temperature stress, and they were involved in an extensive gene expression network [ 11 ]. Digital gene expression in Nile tilapia (Oreochromis niloticus) revealed the differentially expressed-genes when comparing the fish cultured in 10ÊC and 30ÊC [ 12 ]. These studies provided useful insights and promoted the molecular breeding work in fish species. However, in L. vannamei, the gene expression profile in response to temperature stress has not yet been illustrated. RNA-Seq technology is a crucial method for quantifying the transcriptional expression in non-model species or in those organisms that lack whole-genome sequencing data [ 13, 14 ]. Studies on channel catfish (I. punctatus) have suggested that both the de novo synthesis of cold-induced proteins and the modification of existing proteins were required for adaptation to low temperatures [9]. Transcriptional regulation research in the coral reef fish Pomacentrus moluccensis revealed the existence of differentially expressed genes in this species [ 15 ]. Moreover, RNA-Seq has already been applied as a useful method on other stress response studies in L. vannamei, such as high salinity [ 16 ] and low salinity [ 17, 18, 19, 20 ]. However, the specific expression profile report on L. vannamei under the stress of cold temperature is still lacked. Molecular signaling pathways can help to improve the understanding of the activity and coordination of the cell-expression network in response to extracellular stress. Calcium is a ubiquitous second messenger in eukaryotic signal transduction cascades [ 21, 22 ]. It impacts nearly every aspect of cellular life [23] and fulfills its key role in transferring extracellular stress signals into chromatin and in regulating protein phosphorylation [ 24 ]. Structural comparisons of mammalian Ca2+/calmodulin-dependent protein kinase II (CaMKII) and plant calciumdependent protein kinases (CDPKs) of calcium signaling showed that their kinase domain containing the highly conserved subdomains of typical eukaryotic Ser/Thr protein kinases [ 22 ]. Mitogen-activated protein kinases (MAPKs), regulated by MAPK kinase kinase-MAPK kinase-MAPK (MKKK-MKK-MAPK) phosphorelay system, are Ser/Thr protein kinases and are key components for vital signal transduction pathways that converting extracellular stimuli into a wide range of cellular responses in the organisms from yeast to humans [ 25,26,27 ]. Consequently, in order to understand the molecular basis of L. vannamei under cold temperatures, investigations conducted on differentially expressed genes in signaling were necessary. Shrimps are poikilothermic species, for which changing of water temperatures can conduct to each part of their bodies and influence their cellular activities. Muscle compose nearly 60% body weight of L. vannamei, and it is the main tissue for converting biochemical energy in shrimp [ 6 ]. Furthermore, we found that the muscle tissues of L. vannamei turned whitish immediately when they were chilled in low temperature water, but regained pellucidity after being returned to the control temperature, indicating that the muscle tissues of L. vannamei may process overt physiology activities in response to acute cold stress. In the present study, a transcriptomic analysis has been carried out on the muscle of L. vannamei. Gene expression profiles under the stress of low temperature has been characterized by comparing various libraries from corresponding temperature treatments, including acute cold stress for 2 hours, endurance cold stress for 48 hours and recovering to normal temperature for 2 hours. The 2 / 21 identified candidate genes will be useful and helpful in revealing a deeper insight into the molecular adaptation basis under low temperature stress in L. vannamei. Materials and methods Ethics statement Shrimp care and experiments were carried out according to the Care and Use of Agricultural Animals in Agricultural Research and Teaching and approved by the Science and Technology Bureau of China. Approval from the Department of Wildlife Administration was not required for the experiments conducted in this paper. All experiments in this paper were performed with permits obtained from the Government of the People's Republic of China and endorsed by the Animal Experimentation Ethics Committee of Chinese Academy of Sciences. Treatments and cold challenges of the shrimp materials Healthy shrimp with an average weight of 7.09±3.22 g were treated in gradient temperatures (9, 10, 11, 12, 13, 14, 15, 18, 21, 24, 27, 30 and 33ÊC) for 2 hours, 48 hours and then recovered at a control temperature (28ÊC) for 2 hours. Survival was evaluated by observing the shrimp's ability to move spontaneously or after gentle prodding [ 28 ]. The temperature of semi-death group was set as the stress cold temperature. Thirteen degree was set as the stress cold temperature. Four groups were placed either in control temperature (CT), 13ÊC for 2 hours (T_2), 13ÊC for 48 hours (T_48) and then, the T_48 group was placed to recover temperature for 2 hours (RC). Muscle tissues from the surviving individuals were put in centrifuge tubes and dipped into liquid nitrogen immediately. Three biological replicates in equal amounts of each group (2 μg of RNA for each individual, and total 6 μg of RNA for each group) were pooled for subsequent RNA-Seq experiments. RNA isolation and cDNA library construction Four separate libraries (CT, T_2, T_48 and RC) were sampled. Total RNA was extracted by TRIzol Reagent (Invitrogen, USA) according to the manufacturer's instructions. The quality of the RNA products was verified as described by Ren et al. (2014) [ 14 ]; RNase-free DNase I (Invitrogen, USA) was used to degrade any possible DNA, and oligo (dT)-coated magnetic beads were mixed with the total RNA to concentrate the ployA tailed mRNA. Fragmentation was performed by incubation in an NEB Next First Strand Synthesis Reaction Buffer (NEB, USA), and the second strand was generated with the buffer, dNTPs, RNase H and DNA polymerase-I. Adapters were ligated to the synthesized cDNA fragments after an end repair step. RNA-Seq and annotation The library was constructed using an Illumina HiSeq 4000 sequencing platform located at the BGI Company (Shenzhen, China; http://www.genomics.cn/index). The clean read data were deposited to the US National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA, http://www.ncbi.nlm.nih.gov/Traces/sra) with the accession No. SRP095377. After assembling, functional annotation by gene ontology terms (GO, http://www. geneontology.org) was analyzed using the Blast2GO program [ 29 ]. The Kyoto Encyclopedia of Genes and Genomes Pathway (KEGG; http://www.genome.jp/kegg), the Non-redundant NCBI collection of nucleotide (NT) and protein (NR) sequence database (ftp://ftp.ncbi.nih. gov/blast/db/FASTA/), the Swiss-Prot database (http://web.expasy.org/docs/swiss-prot_ guideline.html) and the Cluster of Orthologous Groups (COG) protein database (http://www. ncbi.nlm.nih.gov/COG/) were used to predict and classify the possible functions [ 30 ]. 3 / 21 Quantification of differentially expressed genes Three comparisons were taken to the four libraries, including CT-VS-T_2 (T_2/CT), T_2-VS-T_48 (T_48/T_2), and T_48-VS-RC (RC/T_48). The differentially expressed genes (DEGs) in comparison of T_2/CT represented the transcripts in response to acute cold stress (placed from control temperature to acute 13ÊC). The DEGs in comparison of T_48/T_2 represented the transcripts in response to endurance cold stress (endured the cold stress for 48 hours). The DEGs in comparison of RC/T_48 represented the transcripts for recover adaptation (placed from stress cold to control temperature). The clean reads per kilobase per million (RPKM) value was used to estimate transcript abundance on the basis of eliminating the influence of different gene lengths and sequencing discrepancies [ 31 ]. Therefore, the RPKM values could be directly used to compare the differences of gene expression among samples. The algorithm developed by Audic and Claverie (1997) [ 32 ] was used to compare the differences in various gene expression levels. A value of the FDR (false discovery rate) of less than 0.001 and a value of the Log2 ratio of less than 1 were set as the criteria values. Identification of major response genes and pathways The coexisted DEGs in different comparisons of total genes, GO annotated genes and KEGG annotated genes were filtered. Four important kinds of signaling genes regarding Calcium and MAPK cascades, Transcription factor and Serine/threonine-protein kinase signaling genes were searched from the total annotated genes [ 14 ]. To understand the expression model of target genes, hierarchical clustering analysis was performed with the software Heatmap Illustrator Heml 1.0.3.3. Real-time PCR validation To verify the expression profiles of genes in our RNA-Seq results, the represented 15 DEGs in comparison of T_2/CT were selected for validation of the Illumina sequences by real-time PCR analysis. All primers used were listed in Table 1. The treatments of L. vannamei, sampling of muscle, isolation of total RNA and treatment of DNase were performed as described above. The first cDNA strand was synthesized by PrimeScript1 Reverse Transcriptase Kit (Takara, China) with an oligo (dT) primer, a qPCR analysis was conducted in a RotorGene RG3000 real-time PCR system (Corbett research, Australia), PCR reactions were conducted using a SYBR Premix Ex TaqTM II (Takara, China), and the reactions were exposed to an initial denaturation (94ÊC for 3 min) followed by 40 cycles of 94ÊC for 25 s, 60ÊC for 15 s, and 72ÊC for 15 s. The relative transcript abundance was obtained by normalizing with expression of the L. vannamei β-actin gene based on the 2-ΔΔCT method [ 33 ]. Results Cold stress temperature A batch of L. vannamei were treated with different temperatures, and the survival rates were summarized (Fig 1B). The results showed that L. vannamei could not survive when treated with 11ÊC or lower temperatures, the survival rates were slightly affected in the water at 14ÊC, 15ÊC or higher temperatures (Fig 1B). The survival rates at 12ÊC and 13ÊC were 9.30% and 47.22%, respectively. We used 13ÊC as the cold stress temperature since a nearly half survival rate (47.22%) of shrimp was observed at this temperature. The same batch of shrimp from the same cultural environment was used for the following experiments. 4 / 21 GeneID Unigene7411_All Unigene9779_All Unigene3723_All Unigene9994_All Unigene11894_All Unigene25375_All Unigene5734_All CL1921.Contig2_All CL1624.Contig5_All CL199.Contig2_All Unigene15069_All Unigene5290_All Unigene14415_All Unigene26017_All CL261.Contig3_All GeneBank: AF300705 Product length 105 bp 150 bp Overview of the libraries After filtering, the quality metrics of the reads are showed in Table 2: the total raw reads of the samples were from 87.69 to 90.94 Mb, the number of clean reads per library ranged from 65.23 to 66.02 Mb, and the clean reads ratio ranged from 71.95% to 75.28%. The mean value of clean reads Q20 and Q30 percentages were 97.41% and 94.34%, respectively. The clean reads were assembled and clustered into Unigenes, and the quality metrics of Unigenes are presented in Table 3. The total number of All-Unigenes was 29,662, the total length was 29,632,338 bp and the mean length was 999 bp. The annotated genes in NR, NT, SwissProt, KEGG, COG, GO database were 17,103 (57.66%), 11,976 (40.37%), 15,552 (52.43%), 13,936 (46.98%), 8,558 (28.85%) and 2,504 (8.44%), respectively (Fig 2A). Overall 19,619 Unigenes were annotated (Fig 2A). The main Venn diagram of NR, KEGG, COG and SwissProt annotation are showed in Fig 2B, the largest interior gene numbers were 8,264. The percentage of annotated species with NR database are displayed in Fig 2C, the top three annotated species were Zootermopsis nevadensis (the largest blue part, 12.34%), Daphnia pulex (the second large tenne part, 5.31%) and Tribolium castaneum (the third large gray part, 2.95%), respectively. 5 / 21 Fig 1. The cold challenge and survival rate of L. vannamei. (A): Treatment processes of the shrimp. (B): Gradient temperature challenge and survival rates. Control: treated in control temperature of 28℃; T_2: treated in related temperature for 2 hours; T_48: treated in related temperature for 48 hours; Recover: treated in control temperature of 28℃ for 2 hours recovering. Red frame indicated the stress temperature. GO and KEGG classification of the DEGs GO annotation was conducted for the three comparisons. In the comparison group of T_2/ CT, 870 of 1471 DEGs could be assigned by GO classification, and the equivalent numbers for other comparison groups were 1147 of 1943 DEGs (T_48/T_2) and 690 of 1477 DEGs (RC/ T_48); the detailed GO classification data (Pathway, PathwayID, and GeneID) are listed in S1 Table. The most annotated GO pathway in each comparison was Catalytic activity (GO: 0003824). KEGG annotation was conducted for the three comparisons. In the comparison of T_2/CT, 1202 of 1471 DEGs were assigned to 83 pathways. In the comparison of T_48/T_2, 1556 of aQ20: the rate of bases whitch quality is greater than 20 b Q30: The rate of bases which quality is greater than 30 6 / 21 a N50: a weighted median statistic that 50% of the TotalLength is contained in Unigenes great than or equal to this value b N70: a weighted median statistic that 70% of the TotalLength is contained in Unigenes great than or equal to this value c N90: a weighted median statistic that 90% of the TotalLength is contained in Unigenes great than or equal to this value. d GC (%): the percentage of G and C bases in all Unigenes 1943 DEGs were assigned to 91 pathways. In the comparison of RC/T_48, 1208 of 1477 DEGs were assigned to 77 pathways. The details of the KEGG classification (Pathway, PathwayID, and GeneID) of the three comparisons are listed in S2 Table. The top two annotated KEGG pathways in all the three comparisons were Metabolic pathways (ko01100) and RNA transport pathways (ko03013). The RNA transport pathway maps (there were no Metabolic pathway maps in KEGG database) were traced and downloaded from the KEGG database (S1 Fig), the annotated DEGs in the three comparisons were summarized (S3 Table), the major discrepant DEGs in RNA transport pathway were the components of exon-junction complex (EJC). The annotated up- and down-regulated genes of EJC complex were marked and displayed in Fig 3. The EJC complex mainly participated in mRNA translation (S1 Fig, Fig 3A), the EJC outer Fig 2. Annotation of the transcriptome. (A): The summary results with main annotation database; the X axis were labeled with various items; Y axis represented the proportion genes annotated by the functional database; 29,662 was the total numbers of Unigenes, only the percentage numbers were marked in other blue dots. (B): The Venn diagram between NR, KEGG, COG and Swissprot; the annotated gene numbers were exhibit in the related areas. (C): The pie chart of percentage of annotated species with NR database; the main species names were displayed under the pie chart, species (with the corresponding color) were displayed with descending order. 7 / 21 Fig 3. The main components and DEGs in RNA transport pathways. (A): mRNA translation processes. (B): the annotated different expressed genes in EJC complex of T_2/CT. (C): the annotated different expressed genes in EJC complex of T_48/T_2. (D): the annotated different expressed genes in EJC complex of RC/T_48. The up-regulated genes were marked with red frame, down-regulated genes were marked with green frame. shell and Transiently interacting factor genes were the most two kinds of responding components in EJC complex (S3 Table, Fig 3B, 3C and 3D). Represented GO and KEGG classification columnar figures in different comparisons are showed in Fig 4. GO classes contained three categories, including Biological processes, Cellular components and Molecular functions. The top two annotated classes of Biological processes in the comparison of T_2/CT were Cellular processes (97 annotated genes) and Metabolic processes (96 annotated genes), the major classes of Cellular components were Cell (57 annotated genes) and Cell parts (57 annotated genes), and the major classes of Molecular functions were Binding (94 annotated genes) and Catalytic activity (117 annotated genes). The top two annotated KEGG classes in comparison of T_2/CT were Global and overview maps (107 annotated genes of Metabolism) and Translation (76 annotated genes of Genetic Information Processing). The most annotated GO classes in comparison of T_48/T_2 were Catalytic activity (151 annotated genes) in Molecular function, Cell (65 annotated genes) and Cell junction (65 annotated genes) in Cellular component, Metabolic process (128 annotated genes) in Biological process; the most top two annotated KEGG classes were Global and overview maps (123 annotated genes of Metabolism) and Translation (106 annotated genes of Genetic Information Processing). The most annotated GO classes in comparison of RC/T_48 were Catalytic activity (96 8 / 21 Fig 4. GO and KEGG classification of the DEGs in each comparison. (A): Comparison of T_2/CT. (B): Comparison of T_48/ T_2. (C): Comparison of RC/T_48. The Y-axis represents the number of DEGs in each item. 9 / 21 annotated genes) in Molecular function, Cell (36 annotated genes) and Cell part (36 annotated genes) in Cellular component, Cellular process (84 annotated genes) in Biological process; the most top two KEGG classes were Global and overview maps (90 annotated genes of Metabolism) and Translation (82 annotated genes of Genetic Information Processing). Quantification of transcripts The transcript abundance of each Unigene from the four related libraries is listed in S4 Table. DEGs were identified through a pairwise comparison between various transcripts by setting a threshold FDR (false discovery rate) value of 0.001 and a |log2FoldChange| value of 1 based on an algorithm described by Audic and Claverie (1997) and Ren et al. (2014). From the three comparisons (T_2/CT, T_48/T_2 and RC/T_48), a large number of DEGs were identified, and the filtered results (P value < 0.05, FDR 0.001, and estimated absolute |log2FoldChange| 1) are exhibited in S5±S7 Tables. The columnar and Volcano-plot results of the DEGs are showed in Fig 5. The X-axis of Volcano-plot represented the comparison value of log2FoldChange, the YFig 5. Column and Volcano-plot of the DEGs in comparisons. The columnar results of the DEGs are listed and the corresponding Volcano-plot figures are shown. Plots in the Volcano-plot figures represent the DEGs. The Y-axis in the columnar figures represents the DEG numbers in comparisons. In the Volcano-plot figures, the X-axis represents the comparison value of log2FoldChange, and the Y-axis represents the value of |-log10FDR|. The red color refers to up-regulated DEGs, blue color refers to down-regulated DEGs, and the black plot refers to no significantly different DEGs between the pairs of libraries (FDR value 0.001, |log2FoldChange| value 1). 10 / 21 axis represented value of |-log10FDR|. From the pictures, CT-vs-T_2 (T_2/CT) had 904 up- and 567 down-regulated genes, T_2-vs-T_48 (T_48/T_2) had 634 up- and 1309 down-regulated genes and T_48-vs-RC (RC/T_48) had 997 up- and 480 down-regulated genes. The down-regulated transcripts were far more than the up-regulated genes in comparison of T_48/T_2, indicating that more genes were suppressed under the endurance cold stress. Identification of major response genes The Venn diagrams in Fig 6A were used to exhibit the coexisted DEGs in different comparisons. It showed that 234 genes were coexisting in the three comparisons of Total DEGs, and 35 of coexisted GO annotated DEGs and 48 of coexisted KEGG annotated DEGs were also identified. Top 20% of differently-expressed coexisted genes in the three comparisons were selected and taken to Hierarchical clustering analysis, 72 genes were obtained, the detailed GeneID and other information are showed in S8 Table, the heat map figures are showed in Fig 6B. Take the Unigene of CL532.Contig1_All as an example, significantly up-regulated (log2FoldChange value = 9.33) in comparison of T_2/CT, down-regulated (log2FoldChange value = -9.33) in comparison of T_48/T_2 and up-regulated (log2FoldChange value = 8.60) in comparison of RC/T_48. The results indicated that this gene was sensitive to acute cold temperature and recover temperature, and it was down-regulated under the endurance cold stress environment. This gene might play important role in adaptation of changing temperature but was greatly inhibited in long time cold stress environment. The Unigene of Unigene3735_All was taken as another example, for which was significantly down-regulated (Log2FoldChange Value = -7.63) in comparison of T_2/CT, up-regulated (Log2FoldChange Value = 8.12) in comparison of T_48/T_2 and down-regulated (log2FoldChange Value = -8.12) in comparison of RC/T_48. The results indicated this gene was transcriptional inhibited under the environment of changing temperature but was highly expressed under the environment of endurance cold stress, suggesting that it might play important role in response to long time cold stress. Identification of major pathways Total 35 coexisted GO-annotated DEGs and 48 coexisted KEGG-annotated DEGs were selected, the detailed GeneID information are showed in S9 and S10 Tables. The major (Absolute Value of log2FoldChange 5.0) Unigene in GO and KEGG Pathways were filtered, and 15 GO pathways and 5 KEGG pathways were recognized as the temperature-sensitive pathways in response to cold stress (Table 4). The Unigene of CL596.Contig2_All, described as Eukaryotic translation initiation factor 5 gene (gi|646705355|), was annotated by Translation initiation factor 5 pathway (KEGG: K03262) and Cellular macromolecule metabolic (GO: 0044260) and primary metabolic process (GO: 0044238), implying that this gene might play important role in several molecular functions in response to cold stress. Major signaling genes in response to cold stress Four important kinds of signaling genes regarding Calcium signaling, MAPK cascades, Transcription factors and Serine/threonine-protein kinase were fished out (S11±S14 Tables) and clustered into heat map. Generally, the Calcium signaling genes, MAPK cascades genes and Transcription factor genes were not widely transcribed in response to cold stress (Fig 6C, 6D and 6E). Four Unigenes of Calcium signaling genes, 4 Unigenes of MAPK cascades genes and 5 Unigenes of Transcription factor genes were obtained when filtered by the top two levels of log2FoldChange (absolute value of log2FoldChange 6.0 in Calcium signaling genes, absolute value of log2FoldChange 5.0 in MAPK cascades genes and absolute value of log2FoldChange 3.5 in 11 / 21 Fig 6. Venn diagrams and heat maps of the major coexisted DEGs. (A): The Venn diagrams of Total DEGs, GO annotated DEGs and KEGG annotated DEGs in the three comparisons; DEG numbers were labeled in the corresponding areas. (B): 12 / 21 Hierarchical clustering analysis of the top 20% of Total coexisting DEGs in the three comparisons. (C): Hierarchical clustering analysis of the annotated Calcium signaling genes in cold stress. (D): Hierarchical clustering analysis of the annotated MAPK cascades genes in cold stress. (E): Hierarchical clustering analysis of the annotated Transcription factor genes in cold stress. (F): Hierarchical clustering analysis of the annotated Serine/threonine-protein kinase genes in cold stress. Transcription factor genes, Table 5). Whereas, conspicuous patches of Serine/threonine-protein kinase genes were observed (Fig 6F), 16 Unigenes were obtained when filtered by the top two levels of log2FoldChange (absolute value of log2FoldChange 5.0). The results suggested that Serine/threonine-protein kinase might participate in the adaptation of acute or endurance cold stress. Real-time PCR validation The Real-time PCR outcomes were further summarized. The qPCR results of the 15 genes in the treatments of CT and T_2 are showed in Fig 7A. Grouped comparison results between aT_2/CT: The value of Log2FoldChange in comparison of T_2/CT bT_48/T_2: The value of Log2FoldChange in comparison of T_48/T_2 cRC/T_48: The value of Log2FoldChange in comparison of RC/T_48 PathwayID GO:0044260; GO:0044238 GO:0016706; GO:0043167 GO:0016310; GO:0004674; GO:0000166 GO:0006184; GO:0007188; GO:0031683; GO:0005525; GO:0004871; GO:0003924 GO:0097159; GO:1901363 K03262 K12042 K00472 K14005 K12881 13 / 21 Mapk Ser/Thr TFs aT_2/CT: The value of Log2FoldChange in comparison of T_2/CT b T_48/T_2: The value of Log2FoldChange in comparison of T_48/T_2. Ca: Calcium signaling pathway genes Mapk: Mitogen-activated protein kinases cascades genes; Ser/Thr: Serine/threonine-protein kinase genes; TFs: Transcription factor genes qPCR and RNA-Seq were also displayed, and the results showed all the 15 candidate genes in qPCR verification consisted with the results of the RNA-Seq technology (Fig 7B). 14 / 21 Fig 7. Validation for the 15 transcriptomic DEGs in the comparison of T_2/CT. (A): The Real-time PCR results of the 15 genes in treatments of CT and T_2; Y axis represented the Log2 value of the related expression. (B): Grouped results of the 15 genes; the blue columnar mean the Real-time PCR results, red columnar mean the RNA-Seq results; X axis indicated the gene names, Y axis represented the value of log2FoldChange value in the comparison of T_2/CT, the FoldChange value were obtained by2±ΔΔCT method. Bars in this figure represented the standard deviation (SD). Discussion Temperature is one of the most important environmental variances for aquatic ectotherms, and it has been investigated for its key role in affecting overall biological processes in various species [ 7, 34 ]. In larval zebrafish (Danio rerio), cold stress induced-transcripts were generated from 2 and 48 hours treatments under 16ÊC [11]. In tilapia (O. niloticus), an analysis of gene expressions was carried out under cold stress by comparing different temperatures, from 30ÊC followed by a decrement of 1ÊC/Day to 10ÊC [ 12 ]. In Manila clam (Ruditapes philippinarum), a transcriptome study was conducted in response to cold stress with an extreme low temperature of -1ÊC [ 35 ]. Stresses could be classified into mild, chronic or acute, which represented a 15 / 21 dramatic shift in the environmental conditions [ 36 ]. In our study, shrimp were treated in acute cold temperatures for 2 hours, 48 hours and then in a recovery temperature for 2 hours, and a nearly 50% survival rate (47.22%) was observed in the 13ÊC group. The stress temperature in our result was similar to that reported by Fan et al. (2013) [ 37 ], in which 13ÊC for 24 hours was regarded as a suitable cold stress temperature to investigate the altered proteins and genes in L. vannamei. RNA-Seq technology has been demonstrated as a straightforward method to reveal the general molecular regulation processes in species underlying the response to stressful environments [ 14 ]. In the previous studies for L. vannamei, RNA-Seq technology has been employed to investigate the different gene expression levels involved in various environment factors and developmental processes, such as a transcriptome analysis under acute ammonia stress [ 38 ], some digital gene expression analyses in response to low salinity stress [ 17, 18, 19, 20, 39 ], and the transcriptional profiling during ontogenesis [40]. Parts of the molecular mechanisms have already been characterized in response to temperature in the shrimp. The protein content of CAT, GST, Ferritin and HSP60 have been found to be influenced by temperature stress in L. vannamei [ 41 ]. Lipid peroxidation, catalase activity and glutathione-S-transferase in Palaemon elegans and Palaemon serratus have been proved as the biomarkers for monitoring aquatic environments [ 42 ]. In addition, the decreasing water temperatures has been shown to induce oxidative stress in L. vannamei [ 43, 44 ]. Furthermore, 20 proteomic spots and several genes have demonstrated to play important roles in temperature regulation of L. vannamei [37]. However, the profile information of the regulation of gene expression under low temperatures in L. vannamei was still lacked. In the present study, transcriptomic analyses were conducted under cold temperature stress on the muscle of L. vannamei. To our knowledge, this is the first report of transcriptional profiling under acute low temperature stress in L. vannamei. Transcriptome analyses may provide insights into the molecular basis of cold tolerance in many species. Transcriptomic characterization under temperature stress in zebrafish (D. rerio) revealed the transcriptional and post-transcriptional regulation of the cold-acclimated processes [ 5 ], the multiple temperature-regulated biological processes and pathways [ 11 ], and the over-expressed biological processes (such as circadian rhythm, energy metabolism, lipid transport and metabolism) [19]. An RNA-Seq analysis in Melanotaenia duboulayi provided candidate genes in response to changing temperatures in freshwater fish [ 45 ]. The study of transcriptomic analysis on the marine fish ectoparasitic ciliate Cryptocaryon irritans under a low temperature revealed the required genes for cell survival and a deeper dormancy state [ 46 ]. In our present study, a transcriptomic analysis for L. vannamei under acute low temperature was carried out. In total of the four libraries, 904 up- and 567 down-regulated genes in comparison of T_2/CT, 634 up- and 1309 down-regulated genes in comparisons of T_48/T_2, and 997 up- and 480 down-regulated genes in comparisons of RC/T_48 were detected. The top 20% of coexisted DEGs were fished out. In our result, 15 GO pathways and 5 KEGG pathways were annotated as temperature-sensitive pathways. Our study promoted the understanding of profiling gene expression of L. vannamei in response to cold stress. Acute stress places cells in danger, and rapid adaptation is crucial for maximizing cell survival [ 36 ]. To adapt to either biotic or abiotic stress, cellular signaling pathways may coordinate multiple integrated and regulated layers of different steps of mRNA biogenesis in inner cells [ 24, 36 ]. Thus, the important signaling pathways regarding cellular stress response adaptation have been investigated in many studies. Calcium ions (Ca2+) impact nearly every aspect of cellular life [23]; Teets et al. (2013) demonstrated that goldenrod gall fly (Eurosta solidaginis) tissues used calcium signaling to instantly detect decreases in temperature and to trigger downstream cold-hardening mechanisms [ 47 ]. In the transcriptomic analysis reported by Ren et al. (2014), the gene members of the calcium signaling, including cation/calcium exchangers, 16 / 21 calcium-binding proteins, calmodulin-like proteins, etc., were found to be involved in the response to cold stress in Chrysanthemum morifolium [ 14 ]. The MAPKs influenced the cellular signal transduction in organisms from yeast to human. When yeast under osmostress, the Ste11 MAP3Ks were activated, followed by activation of Pbs2 MAP2K that combined to Hog1 and then phosphorylation of the specific osmostress transcription factors for cell adaption [ 36 ]. The p38 family in mammalian species contains 4 isoforms (MAPK11-14) with many overlapping functions in cells, and downstream targets of p38 MAPKs include several kinases that are involved in the control of gene expression and nuclear proteins, such as transcription factors (CREB, ATF1, NFκB, STAT1 and STAT3) and regulators of chromatin remodeling (H3 and HMG14) [ 36 ]. In this paper, the expression model for the genes in four signal pathways were investigated under cold stress. The major differentially-expressed signaling genes were suggested involve in the adaptation of cold stress. Serine/threonine-protein kinase genes might play more important roles in response to cold stress, while Calcium, MAPK cascades and Transcription factor signaling genes were not widely transcribed in response to cold stress. In conclusion, the overall transcriptional expression and quantification analyses under the cold stress temperature of 13ÊC were conducted for Litopenaeus vannamei. Top 20% of differentially expressed genes were filtered. The most annotated and temperature-sensitive of GO and KEGG pathways were found out. The major differentially-expressed genes in signaling were suggested. The results showed that Ser/Thr kinase signal pathway might play more important roles in participating the cold adaptation, while Calcium, MAPK cascades and Transcription factors signaling genes were not widely transcribed. Our study provided the first insight into the molecular basis and supplied an important reference for the mechanism of L. vannamei responded to cold stress. Supporting information S1 Fig. The RNA transport pathway map from KEGG database (map03013). (TIF) S1 Table. GO classification of DEGs in each comparison. (XLSX) S2 Table. KEGG pathway classification of the DEGs in each comparison. (XLSX) S3 Table. Summary of DEGs of the main components of RNA transport pathway. (XLSX) S4 Table. The transcription level of each unigene derived from the number of relevant reads of the four libraries. The ªGeneLengthº column gives the length of exon sequence. "-" represented absence of this unigene in the related group of transcript. RPKM: clean reads per kilobase per million for each unigene in each libraries, respectively. (XLSX) S5 Table. Different transcribed genes between the comparison of T_2/CT. a The ªGeneLengthº column gives the length of exon sequence. b CT: reads per kb per million reads (RPKM) for each unigene in transcription of Control. c T_2: reads per kb per million reads (RPKM) for each unigene in transcription of T_2. d Log2FoldChange(T_2/CT): the ratio between the RPKM in CT to the RPKM in T_2. The criteria applied for assigning significance were: Pvalue < 0.05, FDR 0.001, and estimated absolute |log2FoldChange| 1. Genes listed in descending order of absolute |log2FoldChange|. NA: No annotation were found in the 17 / 21 related database. (XLSX) S6 Table. Different transcribed genes between the comparison of T_48/T_2. a The ªGeneLengthº column gives the length of exon sequence. b T_2: reads per kb per million reads (RPKM) for each unigene in transcription of T_2. c T_48: reads per kb per million reads (RPKM) for each unigene in transcription of T_48. d Log2FoldChange(T_48/T_2): the ratio between the RPKM in T_2 to the RPKM in T_48. The criteria applied for assigning significance were: Pvalue < 0.05, FDR 0.001, and estimated absolute |log2FoldChange| 1. Genes listed in descending order of absolute |log2FoldChange|. NA: No annotation were found in the related database. (XLSX) S7 Table. Different transcribed genes between the comparison of RC/T_48. a The ªGeneLengthº column gives the length of exon sequence. b T_48: reads per kb per million reads (RPKM) for each unigene in transcription of T_48. c RC: reads per kb per million reads (RPKM) for each unigene in transcription of RC. d Log2FoldChange(RC/T_48): the ratio between the RPKM in T_48 to the RPKM in RC. The criteria applied for assigning significance were: Pvalue < 0.05, FDR 0.001, and estimated absolute |log2FoldChange| 1. Genes listed in descending order of absolute |log2FoldChange|. NA: No annotation were found in the related database. (XLSX) S8 Table. Top 20% of coexisted DEGs in the three comparisons. a T_2/CT: The value of Log2FoldChange in comparison of T_2/CT; b T_48/T_2: The value of Log2FoldChange in comparison of T_48/T_2; c RC/T_48: The value of Log2FoldChange in comparison of RC/ T_48. (XLSX) S9 Table. The coexisted DEGs in GO annotation. a T_2/CT: The value of Log2FoldChange in comparison of T_2/CT; b T_48/T_2: The value of Log2FoldChange in comparison of T_48/ T_2; c RC/T_48: The value of Log2FoldChange in comparison of RC/T_48. (XLSX) S10 Table. The coexisted DEGs of KEGG annotation. a T_2/CT: The value of Log2FoldChange in comparison of T_2/CT; b T_48/T_2: The value of Log2FoldChange in comparison of T_48/T_2; c RC/T_48: The value of Log2FoldChange in comparison of RC/T_48. (XLSX) S11 Table. The differential gene expression of Calcium signaling genes in each comparison. (XLSX) S12 Table. The differential gene expression of MAPK cascades genes in each comparison. (XLSX) S13 Table. The differential gene expression of Transcription factors (TFs) in each comparison. (XLSX) S14 Table. The differential gene expression of Serine/threonine-protein kinase genes in each comparison. (XLSX) 18 / 21 Acknowledgments The authors are grateful to Lin Huang (Boya College, Sun Yat-Sen University, Guangzhou 510275, China) for helping us with English Copy-editing. Author Contributions Conceptualization: WH TC CH. Data curation: WH YW XJ PL TC. Formal analysis: WH CR TC. Funding acquisition: WH RC TC CH. Investigation: WH HL DH YT. Methodology: WH CR TC. Project administration: CR CH. Resources: CR CH. Software: WH YW XJ PL TC. Supervision: CR CH. Validation: WH RC YW XJ PL TC. Visualization: WH TC. Writing ± original draft: WH. Writing ± review & editing: RC TC CH. 19 / 21 20 / 21 1. FAO. Fisheries and Aquaculture Department: The State of World Fisheries and Aquaculture . Rome, Italy: Food and Agriculture Organization of the United Nations; 2012 : 37 . 2. Kim S , Pang Z , Seo H , Cho Y , Samocha T , Jang I . Effect of bioflocs on growth and immune activity of Pacific white shrimp, Litopenaeus vannamei postlarvae . Aquac Res . 2014 ; 45 : 362 ± 371 . 3. Chen T , Ren C , Wang Y , Luo P , Jiang X , Huang W , et al. Molecular cloning, inducible expression and antibacterial analysis of a novel i-type lysozyme (lyz-i2) in Pacific white shrimp, Litopenaeus vannamei . Fish Shellfish Immun . 2016 ; 54 : 197 ± 203 . 4. Ghaffari N , Sanchez-Flores A , Doan R , Garcia-Orozco K , Chen P , Ochoa-Leyva A , et al. Novel transcriptome assembly and improved annotation of the whiteleg shrimp (Litopenaeus vannamei), a dominant crustacean in global seafood mariculture . Sci Rep . 2014 ; 4 : 7081 . https://doi.org/10.1038/ srep07081 PMID: 25420880 5. Long Y , Song G , Yan J , He X , Li Q , Cui Z. Transcriptomic characterization of cold acclimation in larval zebrafish . BMC Genomics . 2013 ; 14 : 612 . https://doi.org/10.1186/ 1471 -2164-14-612 PMID: 24024969 6. Wang Q , Tan X , Jiao S , You F , Zhang P. Analyzing Cold Tolerance Mechanism in Transgenic Zebrafish (Danio rerio) . PLoS ONE . 2014 ; 9 ( 7 ):e102492. https://doi.org/10.1371/journal.pone. 0102492 PMID: 25058652 7. Long Y , Yan J , Song G , Li X , Li X , Li Q , et al. Transcriptional events co-regulated by hypoxia and cold stresses in Zebrafish larvae . BMC Genomics . 2015 ; 16 : 385 . https://doi.org/10.1186/s12864-015-1560- y PMID: 25975375 8. Vornanen M , Hassinen M , Koskinen H , Krasnov A . Steady-state effects of temperature acclimation on the transcriptome of the rainbow trout heart . Am J Physiol Regul Integr Comp Physiol . 2005 ; 289 : R1177 ± 1184 . https://doi.org/10.1152/ajpregu.00157. 2005 PMID: 15932967 9. Ju Z , Dunham RA , Liu Z. Differential gene expression in the brain of channel catfish (Ictalurus punctatus) in response to cold acclimation . Mol Genet Genomics . 2002 ; 268 : 87 ± 95 . https://doi.org/10.1007/ s00438-002 -0727-9 PMID: 12242503 10. Gracey A , Fraser E , Li W , Fang Y , Taylor R , Rogers J , et al. Coping with cold: An integrative, multitissue analysis of the transcriptome of a poikilothermic vertebrate . PNAS . 2004 ; 101 : 16970 ± 16975 . https:// doi.org/10.1073/pnas.0403627101 PMID: 15550548 11. Long Y , Li L , Li Q , He X , Cui Z. Transcriptomic Characterization of Temperature Stress Responses in Larval Zebrafish . PLoS ONE . 2012 ; 7 ( 5 ):e37209. https://doi.org/10.1371/journal.pone. 0037209 PMID: 22666345 12. Yang C , Jiang M , Wen H , Tian J , Liu W , Wu F , et al. Analysis of differential gene expression under lowtemperature stress in Nile tilapia (Oreochromis niloticus) using digital gene expression . Gene . 2015 ; 564 : 134 ± 140 . https://doi.org/10.1016/j.gene. 2015 . 01 .038 PMID: 25617524 13. Hornett E , Wheat C . Quantitative RNA-Seq analysis in non-model species: assessing transcriptome assemblies as a scaffold and the utility of evolutionary divergent genomic reference species . BMC Genomics . 2012 ; 13 ( 1 ): 361 . 14. Ren L , Sun J , Chen S , Gao J , Dong B , Liu Y , et al. A transcriptomic analysis of Chrysanthemum nankingense provides insights into the basis of low temperature tolerance . BMC Genomics . 2014 ; 15 : 844 . https://doi.org/10.1186/ 1471 -2164-15-844 PMID: 25277256 15. Kassahn K , Crozier R , Ward A , Stone G , Caley M. From transcriptome to biological function: environmental stress in an ectothermic vertebrate, the coral reef fish Pomacentrus moluccensis . BMC Genomics . 2007 ; 8 : 358 . https://doi.org/10.1186/ 1471 -2164-8-358 PMID: 17916261 16. Kumlu M , Kumlu M , Turkmen S . Combined effects of temperature and salinity on critical thermal minima of pacific white shrimp Litopenaeus vannamei (Crustacea: Penaeidae) . J Therm Biol . 2010 ; 35 : 302 ± 304 . 17. Gao W , Tan B , Mai K , Chi S , Liu H , Dong X , et al. Profiling of differentially expressed genes in hepatopancreas of white shrimp (Litopenaeus vannamei) exposed to long-term low salinity stress . Aquaculture . 2012 ; 364 ( 365 ): 186 ± 191 . 18. Chen K , Li E , Li T , Xu C , Wang X , Lin H , et al. Transcriptome and Molecular Pathway Analysis of the Hepatopancreas in the Pacific White Shrimp Litopenaeus vannamei under Chronic Low-Salinity Stress . PLoS ONE . 2015 ; 10 ( 7 ):e0131503. https://doi.org/10.1371/journal.pone.0131503 PMID: 26147449 19 . Wang X , Wang S , Li C , Chen K , Qin J , 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 ± 1048 . 20. Zhao Q , Pan L , Ren Q , Hu D. Digital gene expression analysis in hemocytes of the white shrimp Litopenaeus vannamei in response to low salinity stress . Fish Shellfish Immun 2015 ; 42 : 400 ± 407 . 21. Berridge M , Lipp P , Bootman M. The versatility and universality of calcium signaling . Nat Rev Mol Cell Biol . 2000 ; 1 ( 1 ): 11 ± 21 . https://doi.org/10.1038/35036035 PMID: 11413485 22. Cheng S , Willmann M , Chen H , Sheen J. Calcium Signaling through Protein Kinases. The Arabidopsis Calcium-Dependent Protein Kinase Gene Family . Plant Physiol . 2002 ; 129 : 469 ± 485 . https://doi.org/ 10.1104/pp.005645 PMID: 12068094 23. Clapham D. Calcium Signaling . Cell . 2007 ; 131 ( 14 ): 1047 ± 1058 . 24. Dodd A , Kudla J , Sanders D. The Language of Calcium Signaling . Annu Rev Plant Biol . 2010 ; 61 : 593 ± 620 . https://doi.org/10.1146/annurev-arplant- 070109 -104628 PMID: 20192754 25. Cuevas B , Abell A , Johnson G . Role of mitogen-activated protein kinase kinase kinases in signal integration . Oncogene . 2007 ; 26 : 3159 ± 3171 . https://doi.org/10.1038/sj.onc. 1210409 PMID: 17496913 26. Cargnello M , Roux P. Activation and Function of the MAPKs and Their Substrates, the MAPK-Activated Protein Kinases . Microbiol Mol Biol R . 2011 ; 75 ( 1 ): 50 ± 83 . 27. Morrison D. MAP Kinase Pathways . CSH Perspect Biol . 2012 ; 4 : a011254 . 28. Teets N , Elnitsky M , Benoit J , Lopez-Martinez G , Denlinger D , Lee R . Rapid cold-hardening in larvae of the Antarctic midge Belgica antarctica: cellular cold-sensing and a role for calcium . Am J Physiol Regul Integr Comp Physiol . 2008 ; 294 : 1938 ± 1946 . 29. Conesa A , GoÈtz S , GarcÂõa-GoÂmez JM , Terol J , TaloÂn M , Robles M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research . Bioinformatics. 2005 ; 21 ( 18 ): 3674 ± 3676 . https://doi.org/10.1093/bioinformatics/bti610 PMID: 16081474 30. 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 . 2008 ; 36 ( suppl 1 ):D480± D484 . 31. Mortazavi A , Williams B , McCue K , Schaeffer L , Wold B . Mapping and quantifying mammalian transcriptomes by RNA-Seq . Nat Methods . 2008 ; 5 ( 7 ): 621 ± 628 . https://doi.org/10.1038/nmeth.1226 PMID: 18516045 32. Audic S , Claverie J. The significance of digital gene expression profiles . Genome Res . 1997 ; 7 ( 10 ): 986 ± 995 . PMID: 9331369 33. Livak K , Schmittgen T. Analysis of relative gene expression data using real-time quantitative PCR and the 2− ΔΔCT Method . Methods . 2001 ; 25 ( 4 ): 402 ± 408 . https://doi.org/10.1006/meth. 2001 .1262 PMID: 11846609 34. Abele D , Puntarulo S. Formation of reactive species and induction of antioxidant defence systems in polar and temperate marine invertebrates and fish . Comp Biochem Phys A . 2004 ; 138 ( 4 ): 405 ± 15 . 35. Nie H , Jiang L , Huo Z , Liu L , Yang F , Yan X . Transcriptomic responses to low temperature stress in the Manila clam, Ruditapes philippinarum . Fish Shellfish Immunol . 2016 ; 55 : 358 ± 66 . https://doi.org/10. 1016/j.fsi. 2016 . 06 .008 PMID: 27288255 36. Nadal E , Ammerer G , Posas F . Controlling gene expression in response to stress . Nat Rev Genet . 2011 ; 12 : 833 ± 845 . https://doi.org/10.1038/nrg3055 PMID: 22048664 37. Fan L , Wang A , Wu Y . Comparative proteomic identification of the hemocyte response to cold stress in white shrimp, Litopenaeus vannamei . J Proteomics . 2013 ; 80 : 196 ± 206 . https://doi.org/10.1016/j.jprot. 2012 . 12 .017 PMID: 23396037 38. Lu X , Kong J , Luan S , Dai P , Meng X , Cao B , et al. Transcriptome Analysis of the Hepatopancreas in the Pacific White Shrimp (Litopenaeus vannamei) under Acute Ammonia Stress . PloS ONE . 2016 ; 11 ( 10 ):e0164396. https://doi.org/10.1371/journal.pone. 0164396 PMID: 27760162 39. Hu D , Pan L , Zhao Q , Ren Q . Transcriptomic response to low salinity stress in gills of the Pacific white shrimp, Litopenaeus vannamei . Mar Genom . 2015 ; 24 : 297 ± 304 . 40. Quispe R , Justino E , Vieira F , Jaramillo M , Rosa R , Perazzolo L. Transcriptional profiling of immunerelated genes in Pacific white shrimp (Litopenaeus vannamei) during ontogenesis . Fish Shellfish Immunol . 2016 ; 58 : 103 ± 107 . https://doi.org/10.1016/j.fsi. 2016 . 09 .024 PMID: 27637731 41. Zhou J , Wang L , Xin Y , Wang W , He W , Wang A , et al. Effect of temperature on antioxidant enzyme gene expression and stress protein response in white shrimp, Litopenaeus vannamei . J Therm Biol . 2010 ; 35 : 284 ± 289 . 42. Vinagre C , Madeira D , MendoncËa V , Dias M , Roma J , Diniz M. Effect of temperature in multiple biomarkers of oxidative stress in coastal shrimp . J Therm Biol . 2014 ; 41 : 38 ± 42 . https://doi.org/10.1016/j. jtherbio. 2014 . 02 .005 PMID: 24679970 43. Qiu J , Wang W , Wang, Liu Y , Wang A . Oxidative stress, DNA damage and osmolality in the Pacific white shrimp, Litopenaeus vannamei exposed to acute low temperature stress . Comp Biochem Phys C . 2011 ; 154 : 36 ± 41 . 44. Li B , Xian J , Guo H , Wang A , Miao Y , Ye J , et al. Effect of temperature decrease on hemocyte apoptosis of the white shrimp Litopenaeus vannamei . Aquacult Int . 2014 ; 22 : 761 ± 774 . 45. Smith S , Bernatchez L , Beheregaray L . RNA-Seq analysis reveals extensive transcriptional plasticity to temperature stress in a freshwater fish species . BMC Genomics . 2013 ; 14 : 375 . https://doi.org/10.1186/ 1471 -2164-14-375 PMID: 23738713 46. Yin F , Sun P , Wang J , Gao Q . Transcriptome analysis of dormant tomonts of the marine fish ectoparasitic ciliate Cryptocaryon irritans under low temperature . Parasite Vector . 2016 ; 9 : 280 . 47. Teets N , Yi S , Lee R , Denlinger D. Calcium signaling mediates cold sensing in insect tissues . PNAS . 2013 ; 110 ( 22 ): 9154 ± 9159 . https://doi.org/10.1073/pnas.1306705110 PMID: 23671084


This is a preview of a remote PDF: http://journals.plos.org/plosone/article/file?id=10.1371/journal.pone.0178604&type=printable

Wen Huang, Chunhua Ren, Hongmei Li, Da Huo, Yanhong Wang, Xiao Jiang, Yushun Tian, Peng Luo, Ting Chen, Chaoqun Hu. Transcriptomic analyses on muscle tissues of Litopenaeus vannamei provide the first profile insight into the response to low temperature stress, PLOS ONE, 2017, DOI: 10.1371/journal.pone.0178604