Gene 564 (2015) 134–140
Contents lists available at ScienceDirect
Gene journal homepage: www.elsevier.com/locate/gene
Research paper
Analysis of differential gene expression under low-temperature stress in Nile tilapia (Oreochromis niloticus) using digital gene expression Changgeng Yang, Ming Jiang, Hua Wen ⁎, Juan Tian, Wei Liu, Fan Wu, Gengwu Gou Key Laboratory of Freshwater Biodiversity Conservation and Utilization of Ministry of Agriculture, Yangtze River Fisheries Research Institute, Chinese Academy of Fishery Sciences, Wuhan 430223, China
a r t i c l e
i n f o
Article history: Received 18 August 2014 Received in revised form 8 December 2014 Accepted 17 January 2015 Available online 21 January 2015 Keywords: Digital gene expression Tilapia Low-temperature stress
a b s t r a c t Tilapia (Oreochromis niloticus) do not survive well at low temperatures. Therefore, improvement of the lowtemperature resistance has become an important issue for aquaculture development of tilapia. The objective of this study was to construct a digital gene expression tag profile to identify genes potentially related to low temperature in tilapia. In this study, tilapia was treated at 30 °C to lethal temperature 10 °C in decrement of 1 °C D−1. Digital gene expression analysis was performed using the Illumina technique to investigate differentially expressed genes in tilapia cultured at different temperatures (30 °C, 26 °C, 20 °C, 16 °C, and 10 °C). A total of 206,861, 188,082, 185,827, 188,067, and 214,171 distinct tags were obtained by sequencing these five libraries, respectively. Compared with the 30 °C library, there were 304, 407, 709, and 772 upregulated genes and 342, 793, 771, and 1466 downregulated genes in 26 °C, 20 °C, 16 °C, and 10 °C libraries, respectively. Trend analysis of these differentially expressed genes identified six statistically significant trends. Functional annotation analysis of the differentially expressed genes identified various functions associated with the response to lowtemperature stress. When tilapia are subjected to low-temperature stress, expression changes were observed in genes associated with nucleic acid synthesis and metabolism, amino acid metabolism and protein synthesis, lipid and carbohydrate content and types, material transport, apoptosis, and immunity. The differentially expressed genes obtained in this study may provide useful insights to help further understand the effects of low temperature on tilapia. © 2015 Elsevier B.V. All rights reserved.
1. Introduction Tilapia (Oreochromis niloticus) are native to Africa and belong to the Cichlidae of the Perciformes, which is the third-ranking aquaculture breed after Cyprinidae and Salmonidae. Tilapia is widely cultivated in the torrid and subtropical zones. China is the largest tilapia-breeding country, with an output of 1,330,000 t in 2010 and 1,550,000 t in 2011 (Fisheries Bureau, 2013). Tilapia's growth temperature is between 16 and 38 °C, with an optimum temperature between 25 and 28 °C (Wohlfarth and Hulata, 1983), indicating that they are unable to thrive at low temperatures (Cnaani et al., 2000). In 2008, more than 300,000 t of tilapia died because of the continuous low temperature (10–13 °C) in China (Danmei et al., 2010). Therefore, improving the low-temperature resistance of tilapia has become an important issue in aquaculture Abbreviations: DGE, digital gene expression; PCR, polymerase chain reaction; GEO, NCBI Gene Expression Omnibus database; TPM, the number of transcripts per million clean tags; FDR, the false discovery rate; DEGs, the number of differentially expressed genes; GO, gene ontology; KO, KEGG Orthology. ⁎ Corresponding author at: Yangtze River Fisheries Research Institute, Chinese Academy of Fishery Sciences, Wuhan 430223, China. E-mail address: [email protected] (H. Wen).
http://dx.doi.org/10.1016/j.gene.2015.01.038 0378-1119/© 2015 Elsevier B.V. All rights reserved.
development. Research on the impact of low temperatures on the performance of tilapia has mainly focused on the effects on appearance, physiology, and biochemistry (Ndong et al., 2007; Jun et al., 2013; Jun et al., 2014). However, there have been few reports on how tilapia react at the molecular level, for example, the changes in gene expression under low-temperature stress (Desale et al., 2010; Jun et al., 2014). The rapid development of high-throughput sequencing technology, which can quickly and economically detect the transcription of thousands of genes in a test sample, has led to its wide use in the expression analysis of various organisms under differing conditions, as well as gene functional studies and other screening research (Lixin et al., 2010; Zhongyu et al., 2012; Yubing et al., 2014). Therefore, in this study, we selected O. niloticus as the research object and treated them with different temperatures, using digital gene expression (DGE) to test the differences in the expressions of genes in the liver of the tilapia during low-temperature stress. Using the published tilapia genome database, the DGE analysis identified a large number of differentially expressed genes (DEGs). The functions of these DEGs were analyzed and relate to various aspects of the life activities of tilapia, and will form the basis for future research on the mechanism of resistance to low temperature in tilapia. These studies will deepen our understanding of the molecular response of tilapia to
C. Yang et al. / Gene 564 (2015) 134–140
low-temperature stress, and will identify low-temperature-sensitive genes or genes associated with low-temperature tolerance. 2. Material and methods 2.1. Ethics statement Nile tilapia is widely cultivated in South China and is not listed as endangered or protected species. All the experimental animal programs involved in this study were approved by the Yangtze River Fisheries Research Institute's animal care and use committee, and followed the experimental basic principles; and the field studies did not involve endangered or protected species. 2.2. Experimental fish and sample preparation Fifty Nile tilapia with body weights of 100 ± 15 g were purchased from a fish hatchery in Guangxi Province of China, and then transported to the experimental base of the Yangtze River Fisheries Research Institute (Wuhan, Hubei Province, China). Fish used in this study were cultured in 405 m3 polyethylene cultivation tanks. The fish were acclimatized for 2 weeks at 28 °C, and then grown for 1 week at 30 °C before the beginning of the experiments. The tilapia was subject to a lethal low temperature of 10 °C (Sifa et al., 2002) by cooling at a rate of 1 °C D−1. Three individuals from each temperature point were sedated and then dissected at 30 °C, 26 °C, 20 °C, 16 °C, and 10 °C. The liver samples were collected at each temperature point and stored at − 80 °C until needed. 2.3. DGE-tag profiling The TRIZOL reagent (Invitrogen, USA) was used to extract total RNA (approximately 6 μg). The extracted RNA was treated with DNase I (Takara, China), detected by gel electrophoresis and then quantified using ultraviolet spectrophotometry. Oligo(dT) beads were used to adsorb and purify mRNA, and then an Oligo(dT) primer was used to synthesize the first- and second-strand cDNAs. The Illumina Gene Expression Sample Prep Kit was used to prepare the sequence tags, according to the manufacturer's instructions. Illumina adaptors 1 and 2 were added to the 5′ and 3′ ends of the tags, respectively, to obtain a tag library with different 21-bp adaptor sequences at either end for polymerase chain reaction (PCR) amplification. After 15 cycles of linear PCR amplification, 6% TBE PAGE gel electrophoresis was used to purify the 105-bp fragments. The DNA was purified and subjected to Illumina sequencing. Raw sequences are transformed into clean tags after data-processing: removal of 3′ adaptor
135
sequences; removal of tags representing only adaptor sequences; removal of low-quality tags; removal of tags that were too long or too short; and removal of tags with a copy number of 1. Finally, the cleaned raw data were then deposited in the NCBI Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi? acc=GSE55412). 2.4. Analysis and mapping of DGE tags All tags were mapped to the reference tilapia genome (ftp://ftp. ensembl.org/pub/release-72/fasta/oreochromis_niloticus/). To monitor mapping events on both strands, both sense and complementary antisense sequences were included in the mapping process. Clean tags that only mapped to one gene were designated as unambiguous clean tags. The number of unambiguous clean tags for each gene was calculated and then normalized to the number of transcripts per million clean tags (TPM) (Hoen et al., 2008; Morrissy et al., 2009). Referring to “the significance of digital gene expression profiles”(Audic and Claverie, 1997), we screened the DEGs between pairs of samples. The P-value corresponds to the differential gene expression test. The false discovery rate (FDR) was used to determine the threshold of the P-value in multiple tests. In our analysis, the FDR was ≤0.001 and the absolute value of the log2 ratio ≥1 was used as the threshold to judge the significance of differences in gene expression. Cluster software (Eisen et al., 1998) was used to perform cluster analysis of gene expression patterns at different temperature treatments. Additionally, Short Time-series Expression Miner (STEM, version 1.2.2b) software was used to perform trend analysis of the gene expression patterns. The trend analysis of gene expression patterns is a gene expression clustering method (a curve shape of the gene expression in multiple stages) aimed at analyzing multiple continuous samples. It can select genes that have certain biological characteristics (such as continuously increasing expression) to form a gene set. Finally, gene ontology (GO) and pathway classification were used to determine the possible functions of all DEGs by mapping them to GO (http://www. geneontology.org/)(Ashburner et al., 2000; Xin and Zhen, 2007) and KEGG Orthology (KO) (http://www.genome.jp/kegg/) databases (Kanehisa et al., 2008). 3. Results 3.1. Analysis of DGE libraries Libraries produced from the total RNA isolated from the livers of tilapia treated at different temperatures (30 °C, 26 °C, 20 °C, 16 °C, and 10 °C) were named as Tpm-T30, Tpm-T26, Tpm-T20, Tpm-T16, and
Table 1 Categorization and abundance of sequence tags. Summary Raw data Clean tags All tags mapping to genes
Unambiguous tags mapping to genes
Unknown tags
Total Distinct tags Total number Distinct tag number Total number Total % of clean tags Distinct tag number % distinct clean tags Total number Total % of clean tags Distinct tag number % distinct clean tags Total number Total % of clean tags Distinct tag number % distinct clean tags
Tpm-T30
Tpm-T26
Tpm-T20
Tpm-T16
Tpm-T10
4,966,613 206,861 4,847,788 90,650 2,418,439 49.89% 32,479 35.83% 2,348,807 48.45% 31,317 34.55% 2,429,349 50.11% 58,171 64.17%
4,913,802 188,082 4,810,662 87,600 2,534,143 52.68% 32,816 37.46% 2,467,463 51.29% 31,697 36.18% 2,276,519 47.32% 54,784 62.54%
4,756,878 185,827 4,656,085 87,566 2,477,298 53.21% 32,888 37.56% 2,404,072 51.63% 31,730 36.24% 2,178,787 46.79% 54,678 62.44%
4,778,837 188,067 4,676,025 87,979 2,603,667 55.68% 32,730 37.20% 2,544,608 54.42% 31,629 35.95% 2,072,358 44.32% 55,249 62.80%
4,946,774 214,171 4,828,918 98,854 2,568,054 53.18% 36,657 37.08% 2,501,470 51.80% 35,476 35.89% 2,260,864 46.82% 62,197 62.92%
136
C. Yang et al. / Gene 564 (2015) 134–140
709 upregulated and 771 downregulated. At 10 °C, there were 2238 DEGs: 772 upregulated and 1466 downregulated (Fig. 1). 3.3. Trend analysis and clustering of gene expression
Fig. 1. Cluster analysis of differential gene expression patterns. Each column represents a temperature treatment (30 °C, 26 °C, 20 °C, 16 °C, and 10 °C), and each row represents a gene. Expression differences are shown in different colors. Red indicates upregulated and green indicates downregulated.
Tpm-T10, respectively. The major characteristics of these five libraries are described in Table 1, including numbers of raw tags, clean tags (which represented 97.61%, 97.90%, 97.88%, 97.85%, and 97.61% of the original tags in these five libraries, respectively), and distinct tags. The clean tags were compared with the reference unigenes and with the reference genome sequence (Table 1). The tags that unambiguously mapped to gene might represent key genes associated with lowtemperature stress.
3.2. Detection of DEGs at different temperatures Differential expression label screened at different temperatures based on the detection method of different genes with digital gene expression profiling was used to identify DEGs at different temperatures. With 30 °C as the control, at 26 °C, there were 646 DEGs, of which 304 were upregulated and 342 were downregulated. When the temperature dropped to 20 °C, 1200 DEGs were identified, of which 407 were upregulated and 793 were downregulated. At 16 °C, there were 1480 DEGs:
STEM version 1.2.2b software performed trend analysis and clustering of gene expression. STEM uses a specialized calculation method to cluster genes according to their changing trend of expression level with, for example, time or dosage. Moreover, each class has a similar gene expression curve. The genes expressed at different temperatures could be classified into 20 clusters, of which six were judged to be statistically significant. These six clusters included 2320 probe sets (Figs. 2 and 3). Clusters 1, 19, 2, 0, 3, and 17 had significant clustering trends (P b 0.05). Cluster 0 represented genes whose expressions consistently decreased with decreasing temperature and cluster 19 represented genes whose expressions consistently increased with decreasing temperature. Cluster 1 represented the genes whose expressions initially decreased and then were restored to their original level as the temperature decreased. Cluster 2 represented the genes whose expressions initially decreased and then increased to higher than their original level with decreasing temperature. Cluster 3 represented genes whose expressions varied but were lower than the control level. Cluster 17 represented genes whose expressions varied but were higher than the control level. 3.4. GO functional enrichment analysis GO is used to link sets of differentially expressed genes to functional categories that they might be associated with. In this study, we used P b 0.05 for the significant difference GO Term as the threshold. GO analysis was conducted for clusters 1, 19, 2, 0, 3, and 17, which all showed significant trends (Fig. 4). Genes in clusters 0, 1, 2, and 19 were significantly associated with ‘cellular component’, which belongs to biological process. This category includes many biological processes, such as cellular process, metabolic process, biological regulation, pigmentation, localization, and establishment of localization. The genes in cluster 3 were associated with biological process including metabolic process, cellular process, biological regulation, pigmentation, and localization. The genes in cluster 17 mainly take part in the biological process such as cellular process, metabolic process, developmental process, and
Fig. 2. Sketch map of the cluster analysis of differentially expressed genes. Filled color clusters show clusters with significant trends (P value b 0.05). The number at the lower left corner of each cluster represents the number of genes in the cluster.
C. Yang et al. / Gene 564 (2015) 134–140
multicellular organismal process. These clusters showed similar results in terms of their ‘Molecular Function’, which included binding, catalytic activity, and transporter activity. 3.5. KEGG pathway annotations In vivo, the protein product of each gene carries out its biological functions through interaction and coordination with other proteins, and pathway-based analysis is helpful to further understand the biological function of these gene products. For cluster 0, the significantly enriched metabolic pathways included DNA replication, glutathione
137
metabolism, histidine metabolism, tryptophan metabolism, and biosynthesis of unsaturated fatty acids. Compared with clusters 0 and 19, clusters 1, 2 and 3, each of which has 16, 15 and 16, have much more metabolic pathways whose P-value was less than 0.05, and they are similar to involved pathway such as metabolic pathways and various amino acid metabolism. But when it comes to different pathways, it mainly includes ribosome, citrate cycle (TCA cycle) and RNA transport. The genes in cluster 17 are associated with the following pathways: protein export, ribosome biogenesis in eukaryotes, porphyrin and chlorophyll metabolism, and protein processing in the endoplasmic reticulum (Table 2).
Fig. 3. Cluster analysis of differentially expressed genes. The STEM software grouped the probe sets into six statistically significant clusters. The notations “T30”, “T26”, “T20”, “T16”, and “T10” on the x-axis represent different temperature treatment times. The numbers at the upper right indicate the number of genes belonging to each cluster and the P-value significance. A: Cluster 0; B: Cluster 19; C: Cluster 1; D: Cluster 2; E: Cluster 3; F: Cluster 17.
138
C. Yang et al. / Gene 564 (2015) 134–140
4. Discussion 4.1. Effects of low-temperature stress on tilapia 4.1.1. Effects on nucleic acid synthesis and metabolism Fish response to low temperature stress is a process controlled by multiple genes. Low temperature affects DNA replication, transcription, and translation, which represents the basic response to lowtemperature stress in fish (Ju et al., 2002; Buckley et al., 2006; Logan and Somero, 2010, 2011). In the study of low temperature adaptation mechanisms of zebrafish, 43.33% of the spliceosome associated genes were up-regulated by cold stress. Ribosome biogenesis in eukaryotes is another highly represented pathway with 30.88% of the associated genes up-regulated by cold stress (Long et al., 2013). In our study, firstly, we found that the expression levels of some genes are affected by low temperature, such as those related to DNA replication, the mRNA surveillance pathway, nucleotide excision repair, ribosome, ribosome biogenesis in eukaryotes, RNA degradation, RNA polymerase, and RNA transport. Among them, the expression of ectonucleoside triphosphate diphosphohydrolase 2, DNA polymerase, and amidophosphoribosyltransferase showed a downward trend
under low temperature stress. Secondly, under low-temperature stress, cellular DNA is vulnerable to damage. For example, in DNA replication, base pairing errors, base tautomerism, base deamination, base modification, and base loss will occur that results may affect the function and genetic properties of cells. But DNA repair is the foundation of resistance to all types of stress damage. For example, high mobility group box 1 and several other key genes involved in the spliceosome, mismatch repair, and base excision repair were upregulated by lowtemperature stress, as were certain genes involved in the regulation of gene transcription and translation. We identified 19 transcription factors whose expression levels consistently increased under low temperature stress. A further six unknown transcription factors were downregulated. These genes will be candidate genes for follow-up studies. 4.1.2. Effects on the synthesis and metabolism of amino acids and proteins During tilapia's response to low temperature stress, certain key genes regulate amino acid metabolism and protein metabolism. However, in the process of low temperature stress, protein catabolism was frequently mentioned and examined. Low temperature not only is able to affect the stability of the protein, but also affect the rate
Fig. 4. Gene classification of the clusters showing significant trends based on gene ontology. A: Cluster 0; B: Cluster 19; C: Cluster 1; D: Cluster 2; E: Cluster 3; F: Cluster 17.
C. Yang et al. / Gene 564 (2015) 134–140
and correctness of protein folding. That means degradation or modification of unfolded and misfolded protein is very important in low temperature stress (Long et al., 2012; Gracey et al., 2004). In this experiment, we also found some DEGs that are related to protein degradation. These genes were mainly associated with the ubiquitin–proteasome complex pathways (UPPs). UPP is a highly efficient, protein degradation pathway in eukaryotic cells, which is involved in many important biological functions in the body, including the cell cycle, muscle atrophy, antigen presentation, and apoptosis. In particular, UPP regulates different cell pathways via the degradation of activating or inhibiting factors (Dorts et al., 2012). In this study, we found that the expressions of proteasome subunit beta 4, RPN5, RPN7, and RPT4 were always upregulated at low temperature; however, the expression of proteasome subunit beta 1 was lower than the control. Ubiquitin synthesis-related genes, UBE1C UBE2D–E, Cul3, Apc3, and Apc7 were upregulated at low temperatures, but UBE2C and Cdh1 were downregulated. 4.1.3. Effects on lipid metabolism It may have some relevance between proper membrane fluidity and the ability of resistance to low temperature. Low temperature can lead to changes in the fatty acid composition of fish cells (Jun et al., 2014). One of the most obvious is the rapidly rising proportion of unsaturated fatty acids, which help maintain the fluidity of the cell membrane (Peggy and Jeffrey, 1982). Hsieh (Hsieh and Kuo, 2005) found that low temperature led to the regulation of cell membrane fluidity, where transducers perceive the changes in the fluidity of membrane, and then transmit information, ultimately leading to altered gene expressions. The expression levels of DEDS, PAFAH, SPT, and GHMT, which are fatty acid composition genes in tilapia, were upregulated, while the expressions of PLA2, NADPH, carbonyl reductase, palmitoyl-CoA hydrolase, NAD, and YDC1 were downregulated. These genes mainly affect the content and proportion of fatty acids in the body by influencing fatty acid metabolism, ether lipid metabolism, biosynthesis of
139
unsaturated fatty acids, fatty acid elongation in mitochondria, cyanoamino acid metabolism, linoleic acid metabolism, sphingolipid metabolism, alpha-linolenic acid metabolism, and arachidonic acid metabolism. 4.1.4. Effects on glucose metabolism We found that the expressions of genes for malate dehydrogenase, lactase, ALDH, carboxy-lyases, α-trehalase, and α-amylase in the TCA cycle, as well as genes associated with galactose metabolism, pentose and glucuronate interconversions, starch and sucrose metabolism, and seleno compound metabolism, were downregulated. This suggested that when tilapia are subjected to low-temperature stress, the balance of carbohydrates is damaged, which would affect the physiological functions of tilapia (Vergauwen et al., 2010; Jun et al., 2014). 4.1.5. Effects on material transport Under low temperature stress, soluble substances such as electrolyte leak out to the extra-membrane, and destroy the balance of extra- and intra-cellular ions. Low temperature can make low the concentration of sodium and chlorine ions in plasma and drop the osmotic pressure in tilapia mossambica (Allanson et al., 1971). In our study, GO analysis identified a class of genes associated with material transport in cells, including the Na+, K+, Ca2+, and Cl− transport genes, such as ATP1A1, ATP1A3, ATP1B3, and ATP1B1. The expression levels of these genes were significantly below the control levels under low-temperature stress. Similarly a study in Oreochromis mossambicus confirmed that low temperature can decrease the osmotic pressure and the concentration of sodium and chlorine ions in plasma (Mohammad et al., 2007). Meanwhile activity of Na+–K+ ATPase was reduced, resulting in a large number of tilapia deaths in winter (Sardella and Brauner, 2007). Furthermore, the rate of cellular fluidphase endocytosis is an important indicator of cell material transport (Padron et al., 2000). In the DGE analysis, we found that certain significant genes associated with endocytosis, such as ArfGAP, VPS45, CHMP1,
Table 2 KEGG pathways significantly associated with clusters of DEGs. Cluster 0
Cluster 19
Cluster 17
Pathway
P-value
Pathway
P-value
Pathway
P-value
DNA replication Glutathione metabolism Histidine metabolism Tryptophan metabolism
0.000276202 0.001827799 0.03075674 0.04954585
Ribosome Ribosome biogenesis in eukaryotes RNA transport RNA polymerase
2.18E−07 2.33E−06 3.60E−05 0.004576703
Protein export Ribosome biogenesis in eukaryotes Porphyrin and chlorophyll metabolism Protein processing in endoplasmic reticulum
0.02293469 1.15E−02 4.12E−02 4.36E−02
Biosynthesis of unsaturated fatty acids
0.04954585
Glycine, serine and threonine metabolism
0.02156604
Pyrimidine metabolism
0.02886629
Cluster 1
Cluster 2
Cluster 3
Pathway
P-value
Pathway
P-value
Pathway
P-value
Spliceosome Metabolic pathways Ribosome biogenesis in eukaryotes RNA transport Isoquinoline alkaloid biosynthesis Oxidative phosphorylation Ribosome Sphingolipid metabolism
3.92E−06 2.12E−05 1.96E−04 8.13E−04 5.61E−03 0.006436675 0.00704044 0.008083055
8.45E−05 0.00103832 0.002070217 0.004052493 0.005312817 0.009517344 0.01977622 0.0203079
Biosynthesis of secondary metabolites Propanoate metabolism Ascorbate and aldarate metabolism Citrate cycle (TCA cycle) Histidine metabolism Glutathione metabolism Beta-alanine metabolism Tryptophan metabolism
0.002175909 0.003021286 0.003034778 0.005013339 0.01091505 0.0142616 0.0164469 0.01798152
Folate biosynthesis RNA polymerase Glycerolipid metabolism Phenylalanine metabolism Glycerophospholipid metabolism Circadian rhythm — plant mRNA surveillance pathway Pyrimidine metabolism
0.01033338 0.01483869 0.01873486 0.02069493 0.02379591 0.02511005 0.02556651 0.043702
Endocytosis Butanoate metabolism Metabolic pathways Ribosome Citrate cycle (TCA cycle) Synthesis and degradation of ketone bodies Glutathione metabolism Protein processing in endoplasmic reticulum SNARE interactions in vesicular transport Phagosome Valine, leucine and isoleucine degradation Biosynthesis of secondary metabolites Oxidative phosphorylation Terpenoid backbone biosynthesis Fructose and mannose metabolism Proteasome
0.02112345 0.02366301 0.02909268 0.03142113 0.03190187 0.04128445 0.04360911 0.04911374
Starch and sucrose metabolism Metabolic pathways Pentose phosphate pathway Fatty acid metabolism RNA transport Valine, leucine and isoleucine degradation Limonene and pinene degradation
0.02122405 0.02164336 0.03426903 0.03426903 0.03952539 0.04280115 0.04625282
140
C. Yang et al. / Gene 564 (2015) 134–140
and CHMP1, were significantly upregulated at low temperatures. Further research on the effect of these changes is needed. 4.1.6. Effects on apoptosis or programmed cell death Yan (WenJing et al., 2005) showed that a low temperature induced apoptosis in Colossoma brachypomum cells, which provides a reference for the study of tilapia cell death under low-temperature stress. In our DGE analysis, we identified six downregulated genes in the NOD-like receptor signaling pathway, cytosolic DNA-sensing pathway and mTOR signaling pathway. By contrast, the expressions of natural killer cellmediated cytotoxicity, cellular tumor antigen p53, BCL2-related ovarian killer, myeloid cell leukemia sequence 1 (BCL2-related), ubiquitinmediated proteolysis, and aspartic-type endopeptidase were all upregulated. All of these genes are associated with cell death, and will be the focus of future research. 4.1.7. Effects on the immune system The previous studies in tilapia showed that low temperature can regulate the levels of catecholamines, cortisol, epinephrine, norepinephrine, and epinephrine agonists in the blood, thereby inhibiting phagocytosis of white blood cells and lowering antibody levels (WenHsiung et al., 2002). Thus, low-temperature stress causes a decline in fish immunity (Ndong et al., 2007). In our experiment, certain genes related to the immune system showed reduced expression in the DGE analysis, including SOD, F-actin, SGT1, and ATG1. These genes can affect the peroxisome, phagosome, and the regulation of autophagy. This would increase the risk of infection by pathogens when tilapia is cultured at low temperatures. This could be one explanation for the death of tilapia at low temperatures (Leslie et al., 1990). 5. Conclusions Low temperature caused a series of physiological and biochemical changes in fish, which were related to the temperature sensing and regulation mechanisms of fish. When tilapia are subjected to lowtemperature stress, expression changes were observed for genes associated with nucleic acid synthesis and metabolism, amino acid metabolism and protein synthesis, lipid and carbohydrate content and types, material transport, apoptosis, and immunity. Acknowledgments This work was supported by the China Agriculture Research System (Grant No. CARS-49) and funded by the Yangtze River Fisheries Research Institute (Grant No. SZ2012-04). References Allanson, B.R., Bok, A., Vanwyk, N.I., 1971. The influence of exposure to low temperature on Tilapia mossambica Peters (Cichlidae). J. Fish Biol. 3, 181–185. Ashburner, M., Ball, C.A., Blake, J.A., Botstein, D., Butler, H., Cherry, J.M., Davis, A.P., Dolinski, K., Dwight, S.S., Eppig, J.T., Harris, M.A., Hill, D.P., Issel-Tarver, L., Kasarskis, A., Lewis, S., Matese, J.C., Richardson, J.E., Ringwald, M., Rubin, G.M., Sherlock, G., 2000. Gene ontology: tool for the unification of biology. The Gene Ontology Consortium. Nat. Genet. 25, 25–29. Audic, S., Claverie, J.M., 1997. The significance of digital gene expression profiles. Genome Res. 7, 986–995. Buckley, B.A., Gracey, A.Y., Somero, G.N., 2006. The cellular response to heat stress in the goby Gillichthys mirabilis: a cDNA microarray and protein-level analysis. J. Exp. Biol. 209, 2660–2677. Cnaani, A., Gall, A.E., Hulata, G., 2000. Cold tolerance of tilapia species and hybrids. Aquac. Int. 8, 289–298. Danmei, M., GuangPing, C., Haiyan, Y., 2010. Lethal reaction of Nile tilapia (Oreochromis niloticus) under low temperature stress. Guangxi Agric. Sci. 41, 726–728. Desale, B.Z., Kevin, M.F., Robert, J.C., 2010. Transcriptional response of delta-9-desaturase gene to acute and chronic cold stress in Nile tilapia, Oreochromis niloticus. J. World Aquacult. Soc. 41, 800–806.
Dorts, J., Bauwin, A., Kestemont, P., Jolly, S., Sanchez, W., Silvestre, F., 2012. Proteasome and antioxidant responses in Cottus gobio during a combined exposure to heat stress and cadmium. Comp. Biochem. Physiol. C Toxicol. Pharmacol. 155, 318–324. Eisen, M.B., Spellman, P.T., Brown, P.O., Botstein, D., 1998. Cluster analysis and display of genome-wide expression patterns. Proc. Natl. Acad. Sci. U. S. A. 95, 14863–14868. Fisheries Bureau, M.o.A, PRC., 2013. China Fisheries Yearbook. Gracey, A.Y., Fraser, E.J., Li, W., Fang, Y., Taylor, R.R., Rogers, J., Brass, A., Cossins, A.R., 2004. Coping with cold: an integrative, multitissue analysis of the transcriptome of a poikilothermic vertebrate. Proc. Natl. Acad. Sci. U. S. A. 101, 16970–16975. Hoen, P.A.C., Ariyurek, Y., Thygesen, H.H., Vreugdenhil, E., Vossen, R.H.A.M., de Menezes, R.X., Boer, J.M., van Ommen, G.J.B., den Dunnen, J.T., 2008. Deep sequencing-based expression analysis shows major advances in robustness, resolution and inter-lab portability over five microarray platforms. Nucleic Acids Res. 36. Hsieh, S.L., Kuo, C.M., 2005. Stearoyl-CoA desaturase expression and fatty acid composition in milkfish (Chanos chanos) and grass carp (Ctenopharyngodon idella) during cold acclimation. Comp. Biochem. Physiol. B Biochem. Mol. Biol. 141, 95–101. Ju, Z., Dunham, R.A., Liu, Z., 2002. Differential gene expression in the brain of channel catfish (Ictalurus punctatus) in response to cold acclimation. Mol. Genet. Genomics 268, 87–95. Jun, Q., Hui, W., Kpundehb, M.D., Jie, H., Pao, Xu., 2013. Effect of water temperature, salinity, and their interaction on growth, plasma osmolality, and gill Na+, K+-ATPase activity in juvenile GIFT tilapia Oreochromis niloticus (L.). J. Therm. Biol. 38, 331–338. Jun, Q., Jie, H., Hong, Y., Hui, W., Kpundehb, M.D., Pao, Xu., Zhixiang, Z., 2014. Temperature modulates hepatic carbohydrate metabolic enzyme activity and gene expression in juvenile GIFT tilapia (Oreochromis niloticus) fed a carbohydrate-enriched diet. J. Therm. Biol. 40, 25–31. Kanehisa, M., Araki, M., Goto, S., Hattori, M., Hirakawa, M., Itoh, M., Katayama, T., Kawashima, S., Okuda, S., Tokimatsu, T., Yamanishi, Y., 2008. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 36, D480–D484. Leslie, L.B., John, B.K., Michael, J.B., 1990. Cold tolerance in maternal mouthbrooding tilapias: phenotypic variation among species and hybrids. Aquaculture 85, 271–280. Lixin, X., Ding, H., Weiren, D., Yiwen, Z., Jianzhong, S., 2010. Deep sequencing-based transcriptome profiling analysis of bacteria-challenged Lateolabrax japonicus reveals insight into the immune-relevant genes in marine fish. BMC Genomics 11, 472. Logan, C.A., Somero, G.N., 2010. Transcriptional responses to thermal acclimation in the eurythermal fish Gillichthys mirabilis (Cooper 1864). Am. J. Physiol. Regul. Integr. Comp. Physiol. 299, R843–R852. Logan, C.A., Somero, G.N., 2011. Effects of thermal acclimation on transcriptional responses to acute heat stress in the eurythermal fish Gillichthys mirabilis (Cooper). Am. J. Physiol. Regul. Integr. Comp. Physiol. 300, R1373–R1383. Long, Y., Li, L., Li, Q., He, X., Cui, Z., 2012. Transcriptomic characterization of temperature stress responses in larval zebrafish. PLoS One 7, e37209. Long, Y., Song, G., Yan, J., He, X., Li, Q., Cui, Z., 2013. Transcriptomic characterization of cold acclimation in larval zebrafish. BMC Genomics 14, 612. Mohammad, A.A., AbdelFattah, M.E., Ahmed, E.G., 2007. Effects of thermal and thermohaline shocks on survival and osmotic concentration of the tilapias Oreochromis mossambicus and Oreochromis aureus × Oreochromis niloticus hybrids. J. World Aquacult. Soc. 27, 456–461. Morrissy, A.S., Morin, R.D., Delaney, A., Zeng, T., McDonald, H., Jones, S., Zhao, Y., Hirst, M., Marra, M.A., 2009. Next-generation tag sequencing for cancer gene expression profiling. Genome Res. 19, 1825–1835. Ndong, D., Chen, Y.Y., Lin, Y.H., Vaseeharan, B., Chen, J.C., 2007. The immune response of tilapia Oreochromis mossambicus and its susceptibility to Streptococcus iniae under stress in low and high temperatures. Fish Shellfish Immunol. 22, 686–694. Padron, D., Bizeau, M.E., Hazel, J.R., 2000. Is fluid-phase endocytosis conserved in hepatocytes of species acclimated and adapted to different temperatures? Am. J. Physiol. Regul. Integr. Comp. Physiol. 278, R529–R536. Peggy, A.S., Jeffrey, R.H., 1982. Time course of changes in fatty acid composition of gills and liver from rainbow trout (Salmo gairdneri) during thermal acclimation. J. Exp. Zool. 221, 159–168. Sardella, B.A., Brauner, C.J., 2007. Cold temperature-induced osmoregulatory failure: the physiological basis for tilapia winter mortality in the Salton Sea? Calif. Fish Game 93, 200–213. Sifa, L., Chenhong, L., Madan, D., Florabelle, G., Rex, D., 2002. Cold tolerance of three strains of Nile tilapia, Oreochromis niloticus, in China. Aquaculture 213, 123–129. Vergauwen, L., Benoot, D., Blust, R., Knapen, D., 2010. Long-term warm or cold acclimation elicits a specific transcriptional response and affects energy metabolism in zebrafish. Comp. Biochem. Physiol. A Mol. Integr. Physiol. 157, 149–157. WenHsiung, C., LianTien, S., ChingLin, T., YenLin, S., ChingFong, C., 2002. Cold-stress induced the modulation of catecholamines, cortisol, immunoglobulin M, and leukocyte phagocytosis in tilapia. Gen. Comp. Endocrinol. 126, 90–100. WenJing, Y., Hanmin, C., XiaoLi, Z., Hong, Y., Ming, Z., 2005. Cold stress causes early apoptosis in caudal fin cells of Colossoma brachypomum. Acta Biol. Exp. Sin. 38, 105–110. Wohlfarth, G.W., Hulata, G.I., 1983. Applied genetics of tilapias. ICLARM Stud. Rev. 1–26. Xin, Z., Zhen, S., 2007. EasyGO: gene ontology-based annotation and functional enrichment analysis tool for agronomical species. BMC Genomics 8, 246. Yubing, L., Meiling, L., Xinrong, L., Bo, C., Xiaofei, M., 2014. Identification of differentially expressed genes in leaf of Reaumuria soongorica under PEG-induced drought stress by digital gene expression profiling. PLoS One 9, e94277. Zhongyu, Z., Yan, Z., Huizhen, X., Xiaolei, L., Dawei, S., Shaomin, P., 2012. A gene expression profile of the developing human retinal pigment epithelium. Mol. Vis. 18, 2961–2975.
Contents lists available at ScienceDirect
Gene journal homepage: www.elsevier.com/locate/gene
Research paper
Analysis of differential gene expression under low-temperature stress in Nile tilapia (Oreochromis niloticus) using digital gene expression Changgeng Yang, Ming Jiang, Hua Wen ⁎, Juan Tian, Wei Liu, Fan Wu, Gengwu Gou Key Laboratory of Freshwater Biodiversity Conservation and Utilization of Ministry of Agriculture, Yangtze River Fisheries Research Institute, Chinese Academy of Fishery Sciences, Wuhan 430223, China
a r t i c l e
i n f o
Article history: Received 18 August 2014 Received in revised form 8 December 2014 Accepted 17 January 2015 Available online 21 January 2015 Keywords: Digital gene expression Tilapia Low-temperature stress
a b s t r a c t Tilapia (Oreochromis niloticus) do not survive well at low temperatures. Therefore, improvement of the lowtemperature resistance has become an important issue for aquaculture development of tilapia. The objective of this study was to construct a digital gene expression tag profile to identify genes potentially related to low temperature in tilapia. In this study, tilapia was treated at 30 °C to lethal temperature 10 °C in decrement of 1 °C D−1. Digital gene expression analysis was performed using the Illumina technique to investigate differentially expressed genes in tilapia cultured at different temperatures (30 °C, 26 °C, 20 °C, 16 °C, and 10 °C). A total of 206,861, 188,082, 185,827, 188,067, and 214,171 distinct tags were obtained by sequencing these five libraries, respectively. Compared with the 30 °C library, there were 304, 407, 709, and 772 upregulated genes and 342, 793, 771, and 1466 downregulated genes in 26 °C, 20 °C, 16 °C, and 10 °C libraries, respectively. Trend analysis of these differentially expressed genes identified six statistically significant trends. Functional annotation analysis of the differentially expressed genes identified various functions associated with the response to lowtemperature stress. When tilapia are subjected to low-temperature stress, expression changes were observed in genes associated with nucleic acid synthesis and metabolism, amino acid metabolism and protein synthesis, lipid and carbohydrate content and types, material transport, apoptosis, and immunity. The differentially expressed genes obtained in this study may provide useful insights to help further understand the effects of low temperature on tilapia. © 2015 Elsevier B.V. All rights reserved.
1. Introduction Tilapia (Oreochromis niloticus) are native to Africa and belong to the Cichlidae of the Perciformes, which is the third-ranking aquaculture breed after Cyprinidae and Salmonidae. Tilapia is widely cultivated in the torrid and subtropical zones. China is the largest tilapia-breeding country, with an output of 1,330,000 t in 2010 and 1,550,000 t in 2011 (Fisheries Bureau, 2013). Tilapia's growth temperature is between 16 and 38 °C, with an optimum temperature between 25 and 28 °C (Wohlfarth and Hulata, 1983), indicating that they are unable to thrive at low temperatures (Cnaani et al., 2000). In 2008, more than 300,000 t of tilapia died because of the continuous low temperature (10–13 °C) in China (Danmei et al., 2010). Therefore, improving the low-temperature resistance of tilapia has become an important issue in aquaculture Abbreviations: DGE, digital gene expression; PCR, polymerase chain reaction; GEO, NCBI Gene Expression Omnibus database; TPM, the number of transcripts per million clean tags; FDR, the false discovery rate; DEGs, the number of differentially expressed genes; GO, gene ontology; KO, KEGG Orthology. ⁎ Corresponding author at: Yangtze River Fisheries Research Institute, Chinese Academy of Fishery Sciences, Wuhan 430223, China. E-mail address: [email protected] (H. Wen).
http://dx.doi.org/10.1016/j.gene.2015.01.038 0378-1119/© 2015 Elsevier B.V. All rights reserved.
development. Research on the impact of low temperatures on the performance of tilapia has mainly focused on the effects on appearance, physiology, and biochemistry (Ndong et al., 2007; Jun et al., 2013; Jun et al., 2014). However, there have been few reports on how tilapia react at the molecular level, for example, the changes in gene expression under low-temperature stress (Desale et al., 2010; Jun et al., 2014). The rapid development of high-throughput sequencing technology, which can quickly and economically detect the transcription of thousands of genes in a test sample, has led to its wide use in the expression analysis of various organisms under differing conditions, as well as gene functional studies and other screening research (Lixin et al., 2010; Zhongyu et al., 2012; Yubing et al., 2014). Therefore, in this study, we selected O. niloticus as the research object and treated them with different temperatures, using digital gene expression (DGE) to test the differences in the expressions of genes in the liver of the tilapia during low-temperature stress. Using the published tilapia genome database, the DGE analysis identified a large number of differentially expressed genes (DEGs). The functions of these DEGs were analyzed and relate to various aspects of the life activities of tilapia, and will form the basis for future research on the mechanism of resistance to low temperature in tilapia. These studies will deepen our understanding of the molecular response of tilapia to
C. Yang et al. / Gene 564 (2015) 134–140
low-temperature stress, and will identify low-temperature-sensitive genes or genes associated with low-temperature tolerance. 2. Material and methods 2.1. Ethics statement Nile tilapia is widely cultivated in South China and is not listed as endangered or protected species. All the experimental animal programs involved in this study were approved by the Yangtze River Fisheries Research Institute's animal care and use committee, and followed the experimental basic principles; and the field studies did not involve endangered or protected species. 2.2. Experimental fish and sample preparation Fifty Nile tilapia with body weights of 100 ± 15 g were purchased from a fish hatchery in Guangxi Province of China, and then transported to the experimental base of the Yangtze River Fisheries Research Institute (Wuhan, Hubei Province, China). Fish used in this study were cultured in 405 m3 polyethylene cultivation tanks. The fish were acclimatized for 2 weeks at 28 °C, and then grown for 1 week at 30 °C before the beginning of the experiments. The tilapia was subject to a lethal low temperature of 10 °C (Sifa et al., 2002) by cooling at a rate of 1 °C D−1. Three individuals from each temperature point were sedated and then dissected at 30 °C, 26 °C, 20 °C, 16 °C, and 10 °C. The liver samples were collected at each temperature point and stored at − 80 °C until needed. 2.3. DGE-tag profiling The TRIZOL reagent (Invitrogen, USA) was used to extract total RNA (approximately 6 μg). The extracted RNA was treated with DNase I (Takara, China), detected by gel electrophoresis and then quantified using ultraviolet spectrophotometry. Oligo(dT) beads were used to adsorb and purify mRNA, and then an Oligo(dT) primer was used to synthesize the first- and second-strand cDNAs. The Illumina Gene Expression Sample Prep Kit was used to prepare the sequence tags, according to the manufacturer's instructions. Illumina adaptors 1 and 2 were added to the 5′ and 3′ ends of the tags, respectively, to obtain a tag library with different 21-bp adaptor sequences at either end for polymerase chain reaction (PCR) amplification. After 15 cycles of linear PCR amplification, 6% TBE PAGE gel electrophoresis was used to purify the 105-bp fragments. The DNA was purified and subjected to Illumina sequencing. Raw sequences are transformed into clean tags after data-processing: removal of 3′ adaptor
135
sequences; removal of tags representing only adaptor sequences; removal of low-quality tags; removal of tags that were too long or too short; and removal of tags with a copy number of 1. Finally, the cleaned raw data were then deposited in the NCBI Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi? acc=GSE55412). 2.4. Analysis and mapping of DGE tags All tags were mapped to the reference tilapia genome (ftp://ftp. ensembl.org/pub/release-72/fasta/oreochromis_niloticus/). To monitor mapping events on both strands, both sense and complementary antisense sequences were included in the mapping process. Clean tags that only mapped to one gene were designated as unambiguous clean tags. The number of unambiguous clean tags for each gene was calculated and then normalized to the number of transcripts per million clean tags (TPM) (Hoen et al., 2008; Morrissy et al., 2009). Referring to “the significance of digital gene expression profiles”(Audic and Claverie, 1997), we screened the DEGs between pairs of samples. The P-value corresponds to the differential gene expression test. The false discovery rate (FDR) was used to determine the threshold of the P-value in multiple tests. In our analysis, the FDR was ≤0.001 and the absolute value of the log2 ratio ≥1 was used as the threshold to judge the significance of differences in gene expression. Cluster software (Eisen et al., 1998) was used to perform cluster analysis of gene expression patterns at different temperature treatments. Additionally, Short Time-series Expression Miner (STEM, version 1.2.2b) software was used to perform trend analysis of the gene expression patterns. The trend analysis of gene expression patterns is a gene expression clustering method (a curve shape of the gene expression in multiple stages) aimed at analyzing multiple continuous samples. It can select genes that have certain biological characteristics (such as continuously increasing expression) to form a gene set. Finally, gene ontology (GO) and pathway classification were used to determine the possible functions of all DEGs by mapping them to GO (http://www. geneontology.org/)(Ashburner et al., 2000; Xin and Zhen, 2007) and KEGG Orthology (KO) (http://www.genome.jp/kegg/) databases (Kanehisa et al., 2008). 3. Results 3.1. Analysis of DGE libraries Libraries produced from the total RNA isolated from the livers of tilapia treated at different temperatures (30 °C, 26 °C, 20 °C, 16 °C, and 10 °C) were named as Tpm-T30, Tpm-T26, Tpm-T20, Tpm-T16, and
Table 1 Categorization and abundance of sequence tags. Summary Raw data Clean tags All tags mapping to genes
Unambiguous tags mapping to genes
Unknown tags
Total Distinct tags Total number Distinct tag number Total number Total % of clean tags Distinct tag number % distinct clean tags Total number Total % of clean tags Distinct tag number % distinct clean tags Total number Total % of clean tags Distinct tag number % distinct clean tags
Tpm-T30
Tpm-T26
Tpm-T20
Tpm-T16
Tpm-T10
4,966,613 206,861 4,847,788 90,650 2,418,439 49.89% 32,479 35.83% 2,348,807 48.45% 31,317 34.55% 2,429,349 50.11% 58,171 64.17%
4,913,802 188,082 4,810,662 87,600 2,534,143 52.68% 32,816 37.46% 2,467,463 51.29% 31,697 36.18% 2,276,519 47.32% 54,784 62.54%
4,756,878 185,827 4,656,085 87,566 2,477,298 53.21% 32,888 37.56% 2,404,072 51.63% 31,730 36.24% 2,178,787 46.79% 54,678 62.44%
4,778,837 188,067 4,676,025 87,979 2,603,667 55.68% 32,730 37.20% 2,544,608 54.42% 31,629 35.95% 2,072,358 44.32% 55,249 62.80%
4,946,774 214,171 4,828,918 98,854 2,568,054 53.18% 36,657 37.08% 2,501,470 51.80% 35,476 35.89% 2,260,864 46.82% 62,197 62.92%
136
C. Yang et al. / Gene 564 (2015) 134–140
709 upregulated and 771 downregulated. At 10 °C, there were 2238 DEGs: 772 upregulated and 1466 downregulated (Fig. 1). 3.3. Trend analysis and clustering of gene expression
Fig. 1. Cluster analysis of differential gene expression patterns. Each column represents a temperature treatment (30 °C, 26 °C, 20 °C, 16 °C, and 10 °C), and each row represents a gene. Expression differences are shown in different colors. Red indicates upregulated and green indicates downregulated.
Tpm-T10, respectively. The major characteristics of these five libraries are described in Table 1, including numbers of raw tags, clean tags (which represented 97.61%, 97.90%, 97.88%, 97.85%, and 97.61% of the original tags in these five libraries, respectively), and distinct tags. The clean tags were compared with the reference unigenes and with the reference genome sequence (Table 1). The tags that unambiguously mapped to gene might represent key genes associated with lowtemperature stress.
3.2. Detection of DEGs at different temperatures Differential expression label screened at different temperatures based on the detection method of different genes with digital gene expression profiling was used to identify DEGs at different temperatures. With 30 °C as the control, at 26 °C, there were 646 DEGs, of which 304 were upregulated and 342 were downregulated. When the temperature dropped to 20 °C, 1200 DEGs were identified, of which 407 were upregulated and 793 were downregulated. At 16 °C, there were 1480 DEGs:
STEM version 1.2.2b software performed trend analysis and clustering of gene expression. STEM uses a specialized calculation method to cluster genes according to their changing trend of expression level with, for example, time or dosage. Moreover, each class has a similar gene expression curve. The genes expressed at different temperatures could be classified into 20 clusters, of which six were judged to be statistically significant. These six clusters included 2320 probe sets (Figs. 2 and 3). Clusters 1, 19, 2, 0, 3, and 17 had significant clustering trends (P b 0.05). Cluster 0 represented genes whose expressions consistently decreased with decreasing temperature and cluster 19 represented genes whose expressions consistently increased with decreasing temperature. Cluster 1 represented the genes whose expressions initially decreased and then were restored to their original level as the temperature decreased. Cluster 2 represented the genes whose expressions initially decreased and then increased to higher than their original level with decreasing temperature. Cluster 3 represented genes whose expressions varied but were lower than the control level. Cluster 17 represented genes whose expressions varied but were higher than the control level. 3.4. GO functional enrichment analysis GO is used to link sets of differentially expressed genes to functional categories that they might be associated with. In this study, we used P b 0.05 for the significant difference GO Term as the threshold. GO analysis was conducted for clusters 1, 19, 2, 0, 3, and 17, which all showed significant trends (Fig. 4). Genes in clusters 0, 1, 2, and 19 were significantly associated with ‘cellular component’, which belongs to biological process. This category includes many biological processes, such as cellular process, metabolic process, biological regulation, pigmentation, localization, and establishment of localization. The genes in cluster 3 were associated with biological process including metabolic process, cellular process, biological regulation, pigmentation, and localization. The genes in cluster 17 mainly take part in the biological process such as cellular process, metabolic process, developmental process, and
Fig. 2. Sketch map of the cluster analysis of differentially expressed genes. Filled color clusters show clusters with significant trends (P value b 0.05). The number at the lower left corner of each cluster represents the number of genes in the cluster.
C. Yang et al. / Gene 564 (2015) 134–140
multicellular organismal process. These clusters showed similar results in terms of their ‘Molecular Function’, which included binding, catalytic activity, and transporter activity. 3.5. KEGG pathway annotations In vivo, the protein product of each gene carries out its biological functions through interaction and coordination with other proteins, and pathway-based analysis is helpful to further understand the biological function of these gene products. For cluster 0, the significantly enriched metabolic pathways included DNA replication, glutathione
137
metabolism, histidine metabolism, tryptophan metabolism, and biosynthesis of unsaturated fatty acids. Compared with clusters 0 and 19, clusters 1, 2 and 3, each of which has 16, 15 and 16, have much more metabolic pathways whose P-value was less than 0.05, and they are similar to involved pathway such as metabolic pathways and various amino acid metabolism. But when it comes to different pathways, it mainly includes ribosome, citrate cycle (TCA cycle) and RNA transport. The genes in cluster 17 are associated with the following pathways: protein export, ribosome biogenesis in eukaryotes, porphyrin and chlorophyll metabolism, and protein processing in the endoplasmic reticulum (Table 2).
Fig. 3. Cluster analysis of differentially expressed genes. The STEM software grouped the probe sets into six statistically significant clusters. The notations “T30”, “T26”, “T20”, “T16”, and “T10” on the x-axis represent different temperature treatment times. The numbers at the upper right indicate the number of genes belonging to each cluster and the P-value significance. A: Cluster 0; B: Cluster 19; C: Cluster 1; D: Cluster 2; E: Cluster 3; F: Cluster 17.
138
C. Yang et al. / Gene 564 (2015) 134–140
4. Discussion 4.1. Effects of low-temperature stress on tilapia 4.1.1. Effects on nucleic acid synthesis and metabolism Fish response to low temperature stress is a process controlled by multiple genes. Low temperature affects DNA replication, transcription, and translation, which represents the basic response to lowtemperature stress in fish (Ju et al., 2002; Buckley et al., 2006; Logan and Somero, 2010, 2011). In the study of low temperature adaptation mechanisms of zebrafish, 43.33% of the spliceosome associated genes were up-regulated by cold stress. Ribosome biogenesis in eukaryotes is another highly represented pathway with 30.88% of the associated genes up-regulated by cold stress (Long et al., 2013). In our study, firstly, we found that the expression levels of some genes are affected by low temperature, such as those related to DNA replication, the mRNA surveillance pathway, nucleotide excision repair, ribosome, ribosome biogenesis in eukaryotes, RNA degradation, RNA polymerase, and RNA transport. Among them, the expression of ectonucleoside triphosphate diphosphohydrolase 2, DNA polymerase, and amidophosphoribosyltransferase showed a downward trend
under low temperature stress. Secondly, under low-temperature stress, cellular DNA is vulnerable to damage. For example, in DNA replication, base pairing errors, base tautomerism, base deamination, base modification, and base loss will occur that results may affect the function and genetic properties of cells. But DNA repair is the foundation of resistance to all types of stress damage. For example, high mobility group box 1 and several other key genes involved in the spliceosome, mismatch repair, and base excision repair were upregulated by lowtemperature stress, as were certain genes involved in the regulation of gene transcription and translation. We identified 19 transcription factors whose expression levels consistently increased under low temperature stress. A further six unknown transcription factors were downregulated. These genes will be candidate genes for follow-up studies. 4.1.2. Effects on the synthesis and metabolism of amino acids and proteins During tilapia's response to low temperature stress, certain key genes regulate amino acid metabolism and protein metabolism. However, in the process of low temperature stress, protein catabolism was frequently mentioned and examined. Low temperature not only is able to affect the stability of the protein, but also affect the rate
Fig. 4. Gene classification of the clusters showing significant trends based on gene ontology. A: Cluster 0; B: Cluster 19; C: Cluster 1; D: Cluster 2; E: Cluster 3; F: Cluster 17.
C. Yang et al. / Gene 564 (2015) 134–140
and correctness of protein folding. That means degradation or modification of unfolded and misfolded protein is very important in low temperature stress (Long et al., 2012; Gracey et al., 2004). In this experiment, we also found some DEGs that are related to protein degradation. These genes were mainly associated with the ubiquitin–proteasome complex pathways (UPPs). UPP is a highly efficient, protein degradation pathway in eukaryotic cells, which is involved in many important biological functions in the body, including the cell cycle, muscle atrophy, antigen presentation, and apoptosis. In particular, UPP regulates different cell pathways via the degradation of activating or inhibiting factors (Dorts et al., 2012). In this study, we found that the expressions of proteasome subunit beta 4, RPN5, RPN7, and RPT4 were always upregulated at low temperature; however, the expression of proteasome subunit beta 1 was lower than the control. Ubiquitin synthesis-related genes, UBE1C UBE2D–E, Cul3, Apc3, and Apc7 were upregulated at low temperatures, but UBE2C and Cdh1 were downregulated. 4.1.3. Effects on lipid metabolism It may have some relevance between proper membrane fluidity and the ability of resistance to low temperature. Low temperature can lead to changes in the fatty acid composition of fish cells (Jun et al., 2014). One of the most obvious is the rapidly rising proportion of unsaturated fatty acids, which help maintain the fluidity of the cell membrane (Peggy and Jeffrey, 1982). Hsieh (Hsieh and Kuo, 2005) found that low temperature led to the regulation of cell membrane fluidity, where transducers perceive the changes in the fluidity of membrane, and then transmit information, ultimately leading to altered gene expressions. The expression levels of DEDS, PAFAH, SPT, and GHMT, which are fatty acid composition genes in tilapia, were upregulated, while the expressions of PLA2, NADPH, carbonyl reductase, palmitoyl-CoA hydrolase, NAD, and YDC1 were downregulated. These genes mainly affect the content and proportion of fatty acids in the body by influencing fatty acid metabolism, ether lipid metabolism, biosynthesis of
139
unsaturated fatty acids, fatty acid elongation in mitochondria, cyanoamino acid metabolism, linoleic acid metabolism, sphingolipid metabolism, alpha-linolenic acid metabolism, and arachidonic acid metabolism. 4.1.4. Effects on glucose metabolism We found that the expressions of genes for malate dehydrogenase, lactase, ALDH, carboxy-lyases, α-trehalase, and α-amylase in the TCA cycle, as well as genes associated with galactose metabolism, pentose and glucuronate interconversions, starch and sucrose metabolism, and seleno compound metabolism, were downregulated. This suggested that when tilapia are subjected to low-temperature stress, the balance of carbohydrates is damaged, which would affect the physiological functions of tilapia (Vergauwen et al., 2010; Jun et al., 2014). 4.1.5. Effects on material transport Under low temperature stress, soluble substances such as electrolyte leak out to the extra-membrane, and destroy the balance of extra- and intra-cellular ions. Low temperature can make low the concentration of sodium and chlorine ions in plasma and drop the osmotic pressure in tilapia mossambica (Allanson et al., 1971). In our study, GO analysis identified a class of genes associated with material transport in cells, including the Na+, K+, Ca2+, and Cl− transport genes, such as ATP1A1, ATP1A3, ATP1B3, and ATP1B1. The expression levels of these genes were significantly below the control levels under low-temperature stress. Similarly a study in Oreochromis mossambicus confirmed that low temperature can decrease the osmotic pressure and the concentration of sodium and chlorine ions in plasma (Mohammad et al., 2007). Meanwhile activity of Na+–K+ ATPase was reduced, resulting in a large number of tilapia deaths in winter (Sardella and Brauner, 2007). Furthermore, the rate of cellular fluidphase endocytosis is an important indicator of cell material transport (Padron et al., 2000). In the DGE analysis, we found that certain significant genes associated with endocytosis, such as ArfGAP, VPS45, CHMP1,
Table 2 KEGG pathways significantly associated with clusters of DEGs. Cluster 0
Cluster 19
Cluster 17
Pathway
P-value
Pathway
P-value
Pathway
P-value
DNA replication Glutathione metabolism Histidine metabolism Tryptophan metabolism
0.000276202 0.001827799 0.03075674 0.04954585
Ribosome Ribosome biogenesis in eukaryotes RNA transport RNA polymerase
2.18E−07 2.33E−06 3.60E−05 0.004576703
Protein export Ribosome biogenesis in eukaryotes Porphyrin and chlorophyll metabolism Protein processing in endoplasmic reticulum
0.02293469 1.15E−02 4.12E−02 4.36E−02
Biosynthesis of unsaturated fatty acids
0.04954585
Glycine, serine and threonine metabolism
0.02156604
Pyrimidine metabolism
0.02886629
Cluster 1
Cluster 2
Cluster 3
Pathway
P-value
Pathway
P-value
Pathway
P-value
Spliceosome Metabolic pathways Ribosome biogenesis in eukaryotes RNA transport Isoquinoline alkaloid biosynthesis Oxidative phosphorylation Ribosome Sphingolipid metabolism
3.92E−06 2.12E−05 1.96E−04 8.13E−04 5.61E−03 0.006436675 0.00704044 0.008083055
8.45E−05 0.00103832 0.002070217 0.004052493 0.005312817 0.009517344 0.01977622 0.0203079
Biosynthesis of secondary metabolites Propanoate metabolism Ascorbate and aldarate metabolism Citrate cycle (TCA cycle) Histidine metabolism Glutathione metabolism Beta-alanine metabolism Tryptophan metabolism
0.002175909 0.003021286 0.003034778 0.005013339 0.01091505 0.0142616 0.0164469 0.01798152
Folate biosynthesis RNA polymerase Glycerolipid metabolism Phenylalanine metabolism Glycerophospholipid metabolism Circadian rhythm — plant mRNA surveillance pathway Pyrimidine metabolism
0.01033338 0.01483869 0.01873486 0.02069493 0.02379591 0.02511005 0.02556651 0.043702
Endocytosis Butanoate metabolism Metabolic pathways Ribosome Citrate cycle (TCA cycle) Synthesis and degradation of ketone bodies Glutathione metabolism Protein processing in endoplasmic reticulum SNARE interactions in vesicular transport Phagosome Valine, leucine and isoleucine degradation Biosynthesis of secondary metabolites Oxidative phosphorylation Terpenoid backbone biosynthesis Fructose and mannose metabolism Proteasome
0.02112345 0.02366301 0.02909268 0.03142113 0.03190187 0.04128445 0.04360911 0.04911374
Starch and sucrose metabolism Metabolic pathways Pentose phosphate pathway Fatty acid metabolism RNA transport Valine, leucine and isoleucine degradation Limonene and pinene degradation
0.02122405 0.02164336 0.03426903 0.03426903 0.03952539 0.04280115 0.04625282
140
C. Yang et al. / Gene 564 (2015) 134–140
and CHMP1, were significantly upregulated at low temperatures. Further research on the effect of these changes is needed. 4.1.6. Effects on apoptosis or programmed cell death Yan (WenJing et al., 2005) showed that a low temperature induced apoptosis in Colossoma brachypomum cells, which provides a reference for the study of tilapia cell death under low-temperature stress. In our DGE analysis, we identified six downregulated genes in the NOD-like receptor signaling pathway, cytosolic DNA-sensing pathway and mTOR signaling pathway. By contrast, the expressions of natural killer cellmediated cytotoxicity, cellular tumor antigen p53, BCL2-related ovarian killer, myeloid cell leukemia sequence 1 (BCL2-related), ubiquitinmediated proteolysis, and aspartic-type endopeptidase were all upregulated. All of these genes are associated with cell death, and will be the focus of future research. 4.1.7. Effects on the immune system The previous studies in tilapia showed that low temperature can regulate the levels of catecholamines, cortisol, epinephrine, norepinephrine, and epinephrine agonists in the blood, thereby inhibiting phagocytosis of white blood cells and lowering antibody levels (WenHsiung et al., 2002). Thus, low-temperature stress causes a decline in fish immunity (Ndong et al., 2007). In our experiment, certain genes related to the immune system showed reduced expression in the DGE analysis, including SOD, F-actin, SGT1, and ATG1. These genes can affect the peroxisome, phagosome, and the regulation of autophagy. This would increase the risk of infection by pathogens when tilapia is cultured at low temperatures. This could be one explanation for the death of tilapia at low temperatures (Leslie et al., 1990). 5. Conclusions Low temperature caused a series of physiological and biochemical changes in fish, which were related to the temperature sensing and regulation mechanisms of fish. When tilapia are subjected to lowtemperature stress, expression changes were observed for genes associated with nucleic acid synthesis and metabolism, amino acid metabolism and protein synthesis, lipid and carbohydrate content and types, material transport, apoptosis, and immunity. Acknowledgments This work was supported by the China Agriculture Research System (Grant No. CARS-49) and funded by the Yangtze River Fisheries Research Institute (Grant No. SZ2012-04). References Allanson, B.R., Bok, A., Vanwyk, N.I., 1971. The influence of exposure to low temperature on Tilapia mossambica Peters (Cichlidae). J. Fish Biol. 3, 181–185. Ashburner, M., Ball, C.A., Blake, J.A., Botstein, D., Butler, H., Cherry, J.M., Davis, A.P., Dolinski, K., Dwight, S.S., Eppig, J.T., Harris, M.A., Hill, D.P., Issel-Tarver, L., Kasarskis, A., Lewis, S., Matese, J.C., Richardson, J.E., Ringwald, M., Rubin, G.M., Sherlock, G., 2000. Gene ontology: tool for the unification of biology. The Gene Ontology Consortium. Nat. Genet. 25, 25–29. Audic, S., Claverie, J.M., 1997. The significance of digital gene expression profiles. Genome Res. 7, 986–995. Buckley, B.A., Gracey, A.Y., Somero, G.N., 2006. The cellular response to heat stress in the goby Gillichthys mirabilis: a cDNA microarray and protein-level analysis. J. Exp. Biol. 209, 2660–2677. Cnaani, A., Gall, A.E., Hulata, G., 2000. Cold tolerance of tilapia species and hybrids. Aquac. Int. 8, 289–298. Danmei, M., GuangPing, C., Haiyan, Y., 2010. Lethal reaction of Nile tilapia (Oreochromis niloticus) under low temperature stress. Guangxi Agric. Sci. 41, 726–728. Desale, B.Z., Kevin, M.F., Robert, J.C., 2010. Transcriptional response of delta-9-desaturase gene to acute and chronic cold stress in Nile tilapia, Oreochromis niloticus. J. World Aquacult. Soc. 41, 800–806.
Dorts, J., Bauwin, A., Kestemont, P., Jolly, S., Sanchez, W., Silvestre, F., 2012. Proteasome and antioxidant responses in Cottus gobio during a combined exposure to heat stress and cadmium. Comp. Biochem. Physiol. C Toxicol. Pharmacol. 155, 318–324. Eisen, M.B., Spellman, P.T., Brown, P.O., Botstein, D., 1998. Cluster analysis and display of genome-wide expression patterns. Proc. Natl. Acad. Sci. U. S. A. 95, 14863–14868. Fisheries Bureau, M.o.A, PRC., 2013. China Fisheries Yearbook. Gracey, A.Y., Fraser, E.J., Li, W., Fang, Y., Taylor, R.R., Rogers, J., Brass, A., Cossins, A.R., 2004. Coping with cold: an integrative, multitissue analysis of the transcriptome of a poikilothermic vertebrate. Proc. Natl. Acad. Sci. U. S. A. 101, 16970–16975. Hoen, P.A.C., Ariyurek, Y., Thygesen, H.H., Vreugdenhil, E., Vossen, R.H.A.M., de Menezes, R.X., Boer, J.M., van Ommen, G.J.B., den Dunnen, J.T., 2008. Deep sequencing-based expression analysis shows major advances in robustness, resolution and inter-lab portability over five microarray platforms. Nucleic Acids Res. 36. Hsieh, S.L., Kuo, C.M., 2005. Stearoyl-CoA desaturase expression and fatty acid composition in milkfish (Chanos chanos) and grass carp (Ctenopharyngodon idella) during cold acclimation. Comp. Biochem. Physiol. B Biochem. Mol. Biol. 141, 95–101. Ju, Z., Dunham, R.A., Liu, Z., 2002. Differential gene expression in the brain of channel catfish (Ictalurus punctatus) in response to cold acclimation. Mol. Genet. Genomics 268, 87–95. Jun, Q., Hui, W., Kpundehb, M.D., Jie, H., Pao, Xu., 2013. Effect of water temperature, salinity, and their interaction on growth, plasma osmolality, and gill Na+, K+-ATPase activity in juvenile GIFT tilapia Oreochromis niloticus (L.). J. Therm. Biol. 38, 331–338. Jun, Q., Jie, H., Hong, Y., Hui, W., Kpundehb, M.D., Pao, Xu., Zhixiang, Z., 2014. Temperature modulates hepatic carbohydrate metabolic enzyme activity and gene expression in juvenile GIFT tilapia (Oreochromis niloticus) fed a carbohydrate-enriched diet. J. Therm. Biol. 40, 25–31. Kanehisa, M., Araki, M., Goto, S., Hattori, M., Hirakawa, M., Itoh, M., Katayama, T., Kawashima, S., Okuda, S., Tokimatsu, T., Yamanishi, Y., 2008. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 36, D480–D484. Leslie, L.B., John, B.K., Michael, J.B., 1990. Cold tolerance in maternal mouthbrooding tilapias: phenotypic variation among species and hybrids. Aquaculture 85, 271–280. Lixin, X., Ding, H., Weiren, D., Yiwen, Z., Jianzhong, S., 2010. Deep sequencing-based transcriptome profiling analysis of bacteria-challenged Lateolabrax japonicus reveals insight into the immune-relevant genes in marine fish. BMC Genomics 11, 472. Logan, C.A., Somero, G.N., 2010. Transcriptional responses to thermal acclimation in the eurythermal fish Gillichthys mirabilis (Cooper 1864). Am. J. Physiol. Regul. Integr. Comp. Physiol. 299, R843–R852. Logan, C.A., Somero, G.N., 2011. Effects of thermal acclimation on transcriptional responses to acute heat stress in the eurythermal fish Gillichthys mirabilis (Cooper). Am. J. Physiol. Regul. Integr. Comp. Physiol. 300, R1373–R1383. Long, Y., Li, L., Li, Q., He, X., Cui, Z., 2012. Transcriptomic characterization of temperature stress responses in larval zebrafish. PLoS One 7, e37209. Long, Y., Song, G., Yan, J., He, X., Li, Q., Cui, Z., 2013. Transcriptomic characterization of cold acclimation in larval zebrafish. BMC Genomics 14, 612. Mohammad, A.A., AbdelFattah, M.E., Ahmed, E.G., 2007. Effects of thermal and thermohaline shocks on survival and osmotic concentration of the tilapias Oreochromis mossambicus and Oreochromis aureus × Oreochromis niloticus hybrids. J. World Aquacult. Soc. 27, 456–461. Morrissy, A.S., Morin, R.D., Delaney, A., Zeng, T., McDonald, H., Jones, S., Zhao, Y., Hirst, M., Marra, M.A., 2009. Next-generation tag sequencing for cancer gene expression profiling. Genome Res. 19, 1825–1835. Ndong, D., Chen, Y.Y., Lin, Y.H., Vaseeharan, B., Chen, J.C., 2007. The immune response of tilapia Oreochromis mossambicus and its susceptibility to Streptococcus iniae under stress in low and high temperatures. Fish Shellfish Immunol. 22, 686–694. Padron, D., Bizeau, M.E., Hazel, J.R., 2000. Is fluid-phase endocytosis conserved in hepatocytes of species acclimated and adapted to different temperatures? Am. J. Physiol. Regul. Integr. Comp. Physiol. 278, R529–R536. Peggy, A.S., Jeffrey, R.H., 1982. Time course of changes in fatty acid composition of gills and liver from rainbow trout (Salmo gairdneri) during thermal acclimation. J. Exp. Zool. 221, 159–168. Sardella, B.A., Brauner, C.J., 2007. Cold temperature-induced osmoregulatory failure: the physiological basis for tilapia winter mortality in the Salton Sea? Calif. Fish Game 93, 200–213. Sifa, L., Chenhong, L., Madan, D., Florabelle, G., Rex, D., 2002. Cold tolerance of three strains of Nile tilapia, Oreochromis niloticus, in China. Aquaculture 213, 123–129. Vergauwen, L., Benoot, D., Blust, R., Knapen, D., 2010. Long-term warm or cold acclimation elicits a specific transcriptional response and affects energy metabolism in zebrafish. Comp. Biochem. Physiol. A Mol. Integr. Physiol. 157, 149–157. WenHsiung, C., LianTien, S., ChingLin, T., YenLin, S., ChingFong, C., 2002. Cold-stress induced the modulation of catecholamines, cortisol, immunoglobulin M, and leukocyte phagocytosis in tilapia. Gen. Comp. Endocrinol. 126, 90–100. WenJing, Y., Hanmin, C., XiaoLi, Z., Hong, Y., Ming, Z., 2005. Cold stress causes early apoptosis in caudal fin cells of Colossoma brachypomum. Acta Biol. Exp. Sin. 38, 105–110. Wohlfarth, G.W., Hulata, G.I., 1983. Applied genetics of tilapias. ICLARM Stud. Rev. 1–26. Xin, Z., Zhen, S., 2007. EasyGO: gene ontology-based annotation and functional enrichment analysis tool for agronomical species. BMC Genomics 8, 246. Yubing, L., Meiling, L., Xinrong, L., Bo, C., Xiaofei, M., 2014. Identification of differentially expressed genes in leaf of Reaumuria soongorica under PEG-induced drought stress by digital gene expression profiling. PLoS One 9, e94277. Zhongyu, Z., Yan, Z., Huizhen, X., Xiaolei, L., Dawei, S., Shaomin, P., 2012. A gene expression profile of the developing human retinal pigment epithelium. Mol. Vis. 18, 2961–2975.
