Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Global Transcriptome Profiles of 'Meyer' Zoysiagrass in Response to Cold Stress

  • Shanjun Wei,

    Affiliation College of Life & Environmental Sciences, Minzu University of China, Beijing, PR China

  • Zhenlin Du,

    Affiliation Beijing Institute of Genomics, Chinese Academy of Sciences, Beijing, PR China

  • Fei Gao,

    Affiliation College of Life & Environmental Sciences, Minzu University of China, Beijing, PR China

  • Xiang Ke,

    Affiliation College of Life & Environmental Sciences, Minzu University of China, Beijing, PR China

  • Jing Li,

    Affiliation College of Life & Environmental Sciences, Minzu University of China, Beijing, PR China

  • Jianxiu Liu ,

    queenzhou@263.net (YZ); turfunit@aliyun.com (JL)

    Affiliation Jiangsu Province and the Chinese Academy of Sciences, No.1 Qianhu Houcun, Zhongshanmenwai, Nanjing, Jiangsu Province, PR China

  • Yijun Zhou

    queenzhou@263.net (YZ); turfunit@aliyun.com (JL)

    Affiliation College of Life & Environmental Sciences, Minzu University of China, Beijing, PR China

Abstract

A long green period is essential for a turfgrass species with high ornamental value and a wide area of use. Zoysiagrasses (Zoysia spp. Willd.) are perennial turfgrass species popular in tropical, subtropical and temperate zones, possessing many properties necessary to be economically useful turfgrass. They do not have a long green period because of cold sensitivity. A main focus in zoysiagrass research is to develop cold tolerant cultivars. Understanding the cold response in zoysiagrass is a fundamental area of research. In the present study, ‘Meyer’ zoysiagrass (Zoysia japonica), a widely cultivated variety in the genus, is used. We employed RNA-Seq to investigate genome-wide gene expression profiles in leaves under cold stress (4°C). Using the Illumina sequencing platform, we obtained approximately 206 million high-quality paired-end reads from three libraries (0 h, 2 h, and 72 h cold treatment at 4°C). After de novo assembly and quantitative assessment, 46,412 unigenes were generated with an average length of 998 bp and an N50 of 1,522 bp. A total of 25,644 (55.2%) unigenes were annotated by alignment with public protein databases including NR, SwissProt, KEGG and KOG. Differentially expressed genes (DEGs) were investigated using the RPKM method. A total of 756 DEGs were identified between 0h and 2h-cold treatment, with 522 up-regulated and 234 down-regulated; and 5327 DEGs were identified between 0h and 72h-cold treatment, with 2453 up-regulated and 2874 down-regulated. The expression profile of 15 DEGs selected randomly was confirmed with qRT-PCR. The results suggest that cold stress can induce desiccation and oxidative stress, inhibit photosynthesis and substance transport. In response to the stress, genes involved in proline synthesis, in starch hydrolysis, in methionine and ascorbic acid metabolism, in SOD activity, and in DREBs response pathway were up-regulated. GA metabolism, ABA and JA stimulus response were affected under cold exposure. This is the first transcriptome sequencing of Z. japonica, providing a large set of sequence data as well as gene expression profiles under cold stress. It will improve our current understanding of the cold response of zoysiagrass and be beneficial in breeding research.

Introduction

Zoysiagrasses (Zoysia spp. Willd.) are widely used, environmentally friendly warm-season perennial grass (Poaceae) species naturally distributed on seacoasts and grasslands around the Pacific. They spread by both rhizomes and stolons, and can tolerate wide variations in temperature, sunlight, and water. There are 16 species within the genus [1], most of which are sexually compatible [2]. Of these species, Z.japonica Steud is excellent at repelling weeds, holds up well under traffic, http://en.wikipedia.org/wiki/Zoysia-citenote-4resists disease, and requires little fertilization, so it is grown extensively for forage and for sustainable turf and golf courses in transitional and warm climatic regions throughout the world [3]. Lawns built up with zoysiagrasses can provide a higher quality playing surface, require fewer cultural inputs, have lower water and nutritional requirements, and possess good heat and drought tolerance, with generally good disease and insect resistance compared to other turfgrass species [4]. In transitional and temperate regions, however, the green period of zoysiagrass is largely limited by low temperature, which is the primary reason why it is not widely planted [5]. Breeding cultivars with cold tolerance and a long green period is one of the main focuses in zoysiagrass research. Understanding the cold response of this species at the molecular level should lead to the characterization of candidate genes important for genetic improvement.

Significant progress has been made in the past decade with model plants in elucidating transcriptional networks in the response to low temperature [67]. Generally, cold stress-inducible genes in plants can be classified into two groups. The first group consists of genes encoding proteins to protect cells from abiotic stresses, including genes for antioxidant enzymes, osmolyte biosynthesis, and molecular chaperones. The second group is made up mainly of genes that regulate stress signal transduction and gene expression, such as transcription factors and protein kinases. These findings are of use in understanding the cold response of zoysiagrass. Several studies have investigated cold-induced physiological changes, showing changes in the concentrations of sugars, organic acids, proline, soluble proteins and polar lipids under cold stress [810]. In a mapping group consisting of 96 zoysiagrass accessions, three SSR (simple sequence repeat) loci and one SRAP (sequence-related amplified polymorphism) locus were found to be associated with cold tolerance, and three SSR loci and two SRAP loci were found to be associated with the green period [11]. At the protein level, 70 protein spots were found differentially accumulated in stolons after a 14d-cold exposure [12]. A genome-wide investigation of the gene expression profiles would be helpful for our understanding of the cold response in plants as the gene regulation network is very complex [13].

Next-generation sequencing (NGS)-based RNA sequencing for transcriptome methods (RNA-seq) allows for transcript identification as well as potentially discovering genes involved in specific biological processes. This is especially suitable for non-model organisms whose genomic sequences are unknown [14]. In recent years, RNA-seq has emerged as a powerful method for discovering and identifying genes involved in abiotic stress resistance, dramatically improving our understanding of the interaction between plant and environment [1516]. No previous genomic information has been reported in Z.japonica, and fewer than 500 ESTs have been deposited in Genbank for this species. In this study with cultivar ‘Meyer’, we used Illumina sequencing technology to determine the first comprehensive transcriptome for Z.japonica. By comparing the gene expression profiles of cold-treated and control zoysiagrass, we identified numerous differentially and specifically expressed transcripts of cold-regulated genes. These results will be a valuable resource for genetic and genomic studies for breeding cold tolerance of zoysiagrass.

Material and Methods

Plant material and cold treatment

‘Meyer’ zoysiagrass (Zoysia japonica) was collected from the nursery garden of turfgrass resources, Jiangsu Province and the Chinese Academy of Sciences. Dr Jianxiu Liu, one of the corresponding authors of this study, is in charge of the nursery garden and provided the plant material. This collection does not involve endangered or protected species. Grass patches of 10cm×10cm×10cm in size were grown in pots (20cm) in a greenhouse for 6 months, with a 14h light period and 20°C/25°C (Light/Dark) temperature conditions. For cold treatment, grass pots were transferred to a growth chamber set at 4°C for 2h (CT_2h) and 72h (CT_72h). Grass pots transferred to a chamber set at 25°C were used as the control (CK). The light period was 14h, and light intensity was 100 μmol m-2s-1. The second expanded leaf from the top was collected and stored at -80°C after freezing in liquid nitrogen.

Physiological response assay

Relative water content [17], levels of total soluble carbohydrate, proline [8], and malondialdehyde (MDA) [18] were measured using previously described protocols. Activities of detoxification enzymes including superoxide dismutase (SOD), catalase (CAT) and peroxidase (POD) were determined according to Li &Yi (2012) [19]. One unit of SOD activity was defined as the amount of enzyme inhibiting NBT reduction by 50%. One unit of POD and CAT activity was defined as the amount of enzyme causing a change of 0.01 in absorbance per minute. Experiments were conducted with five biological replicates, and results were expressed as the mean ±standard error (SE). All data were subjected to one-way analysis of variance (ANOVA) and the LSD multiple comparison test (p < 0.05) using the SPSS statistical package.

Total RNA extraction, RNA-seq library construction and sequencing

Total RNA was extracted using the Universal Plant Total RNA Extraction Kit RP3301 (Bio Take Corporation, Beijing, China) following the manufacturer’s instructions, with three biological replicates. RNA integrity was confirmed using an Agilent 2100 Bio Analyzer. Samples with RIN ≥ 7 and 28S:18S ratio ≥ 1.5:1, total amount 0.1∽4 ug, were considered acceptable. RNA samples from the same treatments were mixed and used for the determinations. Pair-end libraries were constructed as previously described [20] and sequenced on a HiSeq2000 system according to the manufacturer’s instructions (Capitalbio Company, Beijing, China). The raw sequence reads were deposited in the NCBI (National Center for Biotechnology Information) SRA database.

RNA-seq data processing, de novo assembly and annotation

RNA-seq reads of each sample were first processed by removing the adapter sequence, then merged and assembled to transcripts using Trinity with default parameters (v.r2013-08-14) [21]. Alternative splicing transcripts were assigned to a unigene cluster. The longest transcript in a cluster was assigned as a unigene to form a set of non-redundant transcript sequences. Unigene sequences were aligned to a series of protein databases by BLAST (blastx, E-value < = 1e-5) (v2.2.6) following a priority order of Uniprot (Swiss-Prot/TrEMBL), Nr (non-redundant protein sequences in NCBI), KEGG (Kyoto Encyclopedia of Genes and Genomes database) and KOG (eukaryotic orthologous groups). The best hit of alignment was used to infer biological function of a unigene. The representative Gene Ontology annotation was performed by the software Blast2GO [22] which assigned homologous sequences aligned by BLAST with Uniprot and NCBI nr database to GO terms. We predicted transcription factors from protein sequences obtained by CDS predictions using hmmsearch to search the domains of plant transcription factors (http://plntfdb.bio.uni-potsdam.de/v3.0/). The transcription factor-related unigenes (E-value < = 1e-5) were classified according to the gene family information.

Gene expression quantification and differential expression analysis

We used the Trinity platform for gene expression quantification and differential expression analysis, which included RSEM (v1.2.6) [23] for transcript abundance estimation, and normalization of expression values as FPKM (Fragments per kilobase of transcript permillion fragments mapped), and edgeR [24] for identifying differentially expressed genes (absolute value of log2 (ratio) ≥ 2and adjusted p-value < 0.001). GO terms and pathways enriched in the set of differentially expressed genes were calculated by the hypergeometrictest [25].

Quantitative real-time PCR analysis

For quantitative real-time PCR (qRT-PCR) analysis of transcripts, 1 μg DNase I-treated total RNA was used to synthesize cDNA by M-MLV (Promega) using poly(dT) oligonucleotides. Actin-1 was used as an internal control. qRT-PCR was performed using One Step SYBR PrimeScript RT-PCR Kit (TAKARA, Dalian, China) according to the manufacturer’s instructions. Products were verified by melting curve analysis. Quantification was achieved by normalizing the number of target transcript copies to the reference Actin-1 gene using the comparative ΔΔCt method [26]. All reactions were performed with at least three biological replicates. Primers used in all quantitative RT-PCR experiments are listed in S1 Table.

Results

Physiological changes induced by low temperature

To investigate physiological responses in leaves, we measured relative water content (RWC) and levels of total soluble carbohydrate, proline and malondialdehyde (MDA). The activities of antioxidant enzymes, including SOD, CAT and POD, were also determined. After 2h at 4°C, RWC declined significantly, and the content of total soluble carbohydrate increased (P < 0.05). After 72h at 4°C, RWC remained constant, whereas soluble carbohydrate, proline and MDA increased (P < 0.05) considerably (Fig 1A–1D). The activities of SOD and CAT in leaves exposed to 2h did not change significantly (P ≥ 0.05), while the activity of POD increased by 14.67% (P < 0.05). After a 72h-cold exposure, SOD activity increased by 10.78%, while activities of CAT and POD decreased by 35.99% (P <0.05) and 13.85% (P <0.05), respectively (Fig 1E–1G).

thumbnail
Fig 1. Physiological change in leaves under cold stress.

‘Meyer’ was placed in growth chambers at 4°C for the cold treatment and at 25°C for the control. The second expanded leaf from the top was used for determinations. Experiments were conducted with five biological replicates, and results were expressed as mean ± standard error (SE). *P < 0.05.

https://doi.org/10.1371/journal.pone.0131153.g001

Illumina sequencing and reads assembly

To measure the transcriptome in the cold, we selected ‘Meyer’ leaves from normal condition (CK), 2h cold treatment (CT_2h), and 72h treatment (CT_72h) for RNA-Seq analysis. Three gene pair-end libraries, constructed from cDNA synthesized from leaf RNA, were sequenced with the Illumina HiSeq2000. Total clean nucleotides generated from each sample exceeded 5.5 Gb. We obtained 206,232,142 raw reads in total. Low-quality reads which contained adapter and unknown or low quality bases were discarded, and a total of 206,164,694 clean reads were obtained. Of these clean reads, the average Q30 percentage (sequencing error rate < 1%) was 88.47%. All clean reads were deposited in the NCBI and can be accessed in the Short Read Archive (SRA) under the study accession number SRP051412 (PRJNA270938).

A total of 89,868 transcripts with length ≥ 300 nt were generated using Trinity software (Table 1). The mean transcript size was 1262 nt, and N50 is 1854 nt. Clustering analysis was used to assign all transcripts to 46,412 unigene clusters. A set of non-redundant unigenes was obtained by keeping the longest transcript in a cluster as representative. Among these unigenes, the mean length was 998 nt, with N50 being 1522 nt. A total of 28103 (60.55%) unigenes were longer than 500 nt, and 5818 (12.53%) unigenes were longer than 2000 nt (Fig 2), implying that a batch of transcripts with complete CDS were obtained. Reads were then mapped back to unigenes, 82.76%, 85.68%, and 83.60% of reads in libraries of CK, CT_2h, and CT_72h respectively, were mapped successfully, indicating that the transcript data was suitable for further analysis.

thumbnail
Fig 2. Length distributions of the unigenes.

The X-axis lists the length of unigenes, and the Y-axis lists shows the number of unigenes of that size.

https://doi.org/10.1371/journal.pone.0131153.g002

Among the unigene clusters, 33,837 only contained one transcript sequence, while 12575 clusters contained more than one sequence (Fig 3). These multiple transcript unigene clusters can represent transcription variants, allelic variants, closely related paralogues, misassembled transcripts, or transcripts that were fragmented due to low coverage. The largest unigene cluster contains 188 transcript sequences and was blasted to a probable trehalose-phosphate phosphatase (TPP) gene family in Setaria italic. TPP removes the phosphate from trehalose-6-phosphate to produce free trehalose. There is also a large TPP gene family in the genome of rice and Arabidopsis, some of which are involved in the response to stress [2728].

thumbnail
Fig 3. Histogram of Unigene clusters.

The X-axis lists how many transcripts are in the unigene clusters, and the Y-axis lists show the number of unigenes clusters of that size.

https://doi.org/10.1371/journal.pone.0131153.g003

Functional annotation and classification

We first annotated the unigenes by homologous searches against public databases in the priority order of Uniprot, Nr, KEGG and KOG. Using a cutoff value of E<10−5, a total of 25,644 (55.2%) unigenes had significant BLASTx hits in at least one of the four databases. Among them, 25438(54.8%), 25468(54.9%), 19031(43.1%), 13187(28.6%) unigenes were blasted to Uniprot, Nr, KEGG and KOG, respectively. A total of 12162 (26.2%) unigenes had hit in all four databases, while 20777(44.8%) unigene were not identified (Table 2). To study the sequence conservation of ‘Meyer’ in other plant species, we analyzed the species distribution of unigene datasets by aligning sequences against the Uniprot database. The species distributions for the best match from each sequence are shown in Fig 4. Specifically, 98.0% of these unigenes had the greatest number of hits to monocotyledonous species. The majority of unigenes had a significant level of sequence identity to foxtail millet (Setariaitalica) genes (52.07%). The second closest reference species was Sorghum bicolor (21.92%). Oryza sativa (9.99%), Zea mays (6.30%), Oryza brachyantha (3.99%), Brachypodium distachyon (3.62%), and Hordeum vulgare (0.15%) are reference species in decreasing order (Fig 4).

We used the Gene Ontology (GO) classification system to assess the possible functions of the unigenes. A total of 16,206 (16206/46421) unigenes were successfully assigned to at least one GO term (Fig 5). The unigenes were then classified into three categories: biological processes, cellular component, and molecular function. The category of biological processes consisted of 2021 GO terms, which were assigned to 82.94% of the unigenes; “metabolic process”, “cellular process” and “single-organism process” were the top three. The cellular component category consisted of 391 GO terms, which were assigned to 93.99% of the unigenes; “cell”, “cell part” and “organelle” were the major components. The category of molecular function consisted of 1359 GO terms, which were assigned to 73.66% of the unigenes; “catalytic activity”, “binding” and “transporter activity” were the top three functions.

We then aligned all the unigenes to the KOG database and 11,308 unigenes were identified. These unigenes were grouped into 25 functional categories (Fig 6). The top three categories were “signal transduction mechanism”, “general function prediction only”, and “posttranslational modification, protein turnover, chaperone”; the smallest category was “cell motility”.

Transcription factors are important upstream regulatory proteins and play critical roles in plant responses to abiotic stresses. In the present study we identified a total of 1979 transcription factors from the ‘Meyer’ unigenes that were classified into 55 different families (Table 3). The largest group of transcription factors was bHLH (203, 10.26%), followed by MYB-related family (144, 7.27%), ERF (132, 6.67%), NAC (130, 6.57%) and WRKY (127, 6.42%). These five families contributed up to 37.19% of the transcription factors identified in ‘Meyer’.

thumbnail
Table 3. Distribution of transcription factors by gene family.

https://doi.org/10.1371/journal.pone.0131153.t003

Identification of differentially expressed unigenes

To identify a significant change in gene expression induced by low temperature, we used the FPKM method (Fragments per kilobase of transcript permillion fragments mapped) allowing us to calculate the expression levels of genes. Compared with CK, unigenes with absolute values of log2 (ratio) ≥ 2 and p-value ≤ 0.001 were considered as DEGs (differentially expressed genes). There were 756 DEGs between CK and CT_2h, including 522 up-regulated unigenes and 234 down-regulated ones (Fig 7, S2 Table). There were 5327 DEGs between CK and CT_72h, including 2453 up-regulated unigenes and 2874 down-regulated ones (Fig 7, S3 Table). To further validate our RNA-seq expression profile data, we performed qRT-PCR assays with primers designed according to fifteen randomly-selected unigenes (S1 Table). qRT-PCR data were normalized using the ‘housekeeping’ gene Actin-1.The results showed that gene expression profiles were consistent using either RNA-seq or qRT-PCR analysis, although the exact fold change varied between the two techniques (Table 4). The correlation coefficient between RNA-Seq data and qRT-PCR data is 0.837 (p < 0.01). In addition, PCR products of these primers showed the expected fragment size by agarose gel electrophoresis, which was consistent with the reliability of the assembled sequences.

thumbnail
Fig 7. Number of DEGs between normal condition and cold treatments.

https://doi.org/10.1371/journal.pone.0131153.g007

GO and KEGG enrichment analysis of DEGs

To investigate the function of the DEGs, we then mapped all DEGs to the Gene Ontology database (http://www.genontology.org/) and calculated the gene numbers from each GO term. Using a hypergeometric test, we identified significantly enriched GO terms (p-value < 0.05) in DEGs compared to normal condition. After 2h at low temperature, 39 GO terms were enriched and can be categorized into 6 types (S4 Table): (1) Cell wall and external encapsulating structure organization; (2) Water and carbohydrate transport; (3) Electron carrier activity; (4) Tetrapyrrole binding; (5) Sequence-specific DNA binding transcription factor activity; (6) Auxin mediated signaling pathway. After a 72h-cold exposure, 112 GO terms were enriched (S5 Table). GO term types listed in (2) to (4) above were enriched more intensively. Other enriched GO terms in the CT_72h treatment were mainly related to functions such as: (1) Photosystem and photosynthesis activity; (2) Oxidation-reduction process; (3) Carbohydrate metabolism; (4) Nitrogenous compound metabolism and transport; (5) Lipid metabolism; (6) GA catabolism; (7) Response to ABA and JA stimulus.

KEGG enrichment analysis identified some biochemical pathways and signal transduction pathways. Compared with CK, 5 pathways showed significant change (p-value < 0.05) in the 2h-cold treatment, including biosynthesis of unsaturated fatty acids, biosynthesis of flavonoid and isoflavonoid, MAPK signaling pathway, and terpenoid backbone biosynthesis (S6 Table). In the 72h-cold treatment, 32 pathways showed significant changes (p-value < 0.05) (S7 Table). These pathways included: (1) Photosynthesis, including “antenna proteins”, “chlorophyll metabolism”, and “carbon fixation”; (2) Carbohydrate metabolism, including “starch and sucrose metabolism” and “pentose and glucuronate interconversions”; (3) Lipid metabolism, including “Linoleic acid metabolism” and “alpha-Linolenic acid metabolism”; (4) Antioxidants metabolism, including “ascorbate and aldarate metabolism” and “cysteine and methionine metabolism”; (5) Biosynthesis of secondary metabolites, including “phenylpropanoid biosynthesis”, “isoquinoline alkaloid biosynthesis”, and “flavonoid biosynthesis”; (6) Plant hormone signal transduction. The results of KEGG enrichment and GO term enrichment in the 72h-cold stress were in general similar.

Taken together, after exposure to 4°C for 2h, electron carrier activity and the biological processes of cell wall organization, water homeostasis and carbohydrate transport were affected. After 72 h of cold exposure electron carrier activity, redox balance, photosynthetic activity, carbohydrate and lipid metabolism, biosynthesis of some secondary metabolite, and organic substance transport were affected significantly. The metabolism of GA and responses to ABA and JA stimulus were also involved in cold stress response at the 72 h time point.

Discussion

Temperature is an important environmental factor that affects the distribution, growth and development of plants. Cold stress, including chilling (<20°C) and/or freezing (<0°C) temperatures, prevents the expression of full genetic potential of plants owing to its direct inhibition of metabolic reactions and, indirectly, through cold-induced osmotic, oxidative and other stresses [67, 29]. To understand the cold response in ‘Meyer’ zoysiagrass at the molecular level, we compared the transcriptomes at normal growth temperature and at 4°C. The characteristics of the cold response in ‘Meyer’ are discussed below.

Photosynthesis and nitrogen assimilation

Photosynthesis is greatly inhibited by low temperature in a number of plant species, of which the main impact is photosystem and its action [29]. Similar inhibition occurred in ‘Meyer’ under cold stress, indicated by the down-regulated expression of unigenes involved in porphyrin and chlorophyll synthesis, photosystem structure, photosynthetic electronic transport, and carbon dioxide fixation (S8 Table).

Chlorophylls are complex molecules exquisitely suited to the light absorption, energy transfer, and electron transfer functions carried out in photosynthesis. During the 72h-4°C stress, most of the DEGs involved in chlorophyll synthesis were down-regulated (Fig 8). At the same time, one unigene of ferrochelatase (Hem H), and one of pheophorbide a oxygenase (PAO) were up-regulated. Hem H catalyzes the insertion of ferrous iron into protoporphyrin IX to form protoheme, which departs from the chlorophyll biosynthetic pathway. PAO, the key enzyme in chlorophyll catabolism, catalyzes porphyrinmacrocycle cleavage of pheophorbide a (pheide a) to a primary fluorescent catabolite (pFCC) [30]; its expression was up-regulated during senescence [31]. Taken together, these gene expression changes suggest that chlorophyll biosynthesis decreases and chlorophyll catabolism increases during low temperature treatment.

thumbnail
Fig 8. Chlorophyll metabolism in ‘Meyer’ under cold stress.

gltX, glutamyl-tRNA synthetase; GluTR, Glutamyl-tRNA reductase; PPOX, protoporphyrinogen oxidase; Chl H, magnesium chelatase; EC 2.1.1.11, magnesium protoporphyrin IX methyltransferase; EC 1.14.13.81, magnesium-protoporphyrin IX monomethy lester (oxidative) cyclase; EC 1.3.1.33, protochlorophyllide reductase; EC 3.1.1.14, chlorophyllase; PAO, Pheophorbide a oxygenase; Hem H, ferrochelatase; COX15, cytochrome c oxidase 15 (Heme A synthase). In brackets, red numbers represent the number of up-regulated unigenes, blue numbers represent the number of down-regulated unigenes.

https://doi.org/10.1371/journal.pone.0131153.g008

An intact photosystem is necessary for effective light-absorption and energy transformation. In this study, 10 unigenes encoding chlorophyll a-b binding proteins, 8 unigenes encoding proteins in the reaction center of photosystem Ⅰ (PSⅠ), and 3 unigenes encoding oxygen evolving complex (PsbQ, PsbO and PsbR) were down-regulated. It can be inferred that light harvesting, the water-splitting reaction, and electron transport in PSⅠdecrease in leaves under low temperature exposure. The photosynthetic electron transport activity was also affected, as demonstrated by the down-regulated unigenes of plastocyanin and ferredoxins-1. These changes might lead to abnormal photosynthetic electron transport.

In the carbon fixation process, gene expressions of four enzymes, including carbonic anhydrases (CA), phosphoenolpyruvate carboxylase (PEPC), ribulosebisphosphate carboxylase (Rubisco) and fructose-1,6-bisphosphatase (FBP) were altered considerably during cold treatment. CA catalyzes the reaction of carbon dioxide with water to form bicarbonate, which is then transported to photosynthesizing cells. In this study, the expression of three unigenes encoding CA was down-regulated, suggesting that the exchange rate of carbon dioxide decrease during the 72h cold treatment. This result was consistent with previous research demonstrating that the exchange rate of carbon dioxide in zoysiagrass fell by 80–95% when canopy air temperature was 10°C-12°C, compared to the control (optimal temperature) [32]. PEPC, a key enzyme during C4 photosynthesis, has been suggested to be related to the cold sensitivity of C4 plants [33]. Findings in Amaranthus hypochondriacus (a NAD-ME type C4 species) have shown that PEPC activity fell to 70% on exposure to a low temperature of 15°C, but the protein content and phosphorylation of the enzyme did not change significantly [34] compared to the control. The present study showed that the expression of six PEPC unigenes changed significantly at 4°C, three were up-regulated and three were down-regulated. A previous study indicated that the activity of PEPC in zoysiagrass decreased drastically under 10°C/7°C (L/D) chilling stress [35]. The depressed enzyme activity might be related to the expression changes of PEPC isogenes under cold stress. In addition, one unigene of Rubisco small chain and one unigene of chloroplast located FBP were down-regulated. The expression changes of these carbon-fixation related unigenes support previous data showing that photosynthesis in zoysiagrass decreases greatly under chilling stress compared with that at 25°C [32].

In contrast to carbon fixation, nitrogen assimilation appears to be stimulated under cold stress. The expression changes of related unigenes are shown in S9 Table. First, the expressions of three nitrate reductase (NR) unigenes were up-regulated with fold changes of 6 to 11. Second, two unigenes, homologous with ferredoxin-6 and FNR root isozyme genes in Zea mays, respectively, were up-regulated. In Zea mays the two genes are nitrate-inducible and might function to transfer electrons from the NADPH to Fd-dependent enzymes such as nitrite reductase (NiR) and glutamine synthetase in root plastids [3637]. Third, eight DEGs with GO term “glutamine family amino acid biosynthetic process” (GO:0009084), including two glutamine synthetase unigenes, were all up-regulated. NR catalyzes the reaction: NO3+ NAD(P)H + H++ 2 e→NO2+ NAD(P)++ H2O. Nitrite can be converted to ammonia nitrogen by NiR, and the latter can then be assimilated as amide nitrogen through the glutamine synthetase pathway. It has been shown that in wheat seedlings the activity of glutamine synthetase and nitrate reductase increases during cold acclimation [38]. The overexpression of a glutamine synthetase gene in rice was found to enhance salt and chilling tolerance and increase photorespiration [39]. The up-regulated expression of these genes might increase nitrogen assimilation and photorespiration. It could be inferred that, owing to the decrease of carbon fixation under cold stress, nitrogen assimilation and photorespiration activities were enhanced so as to divert surplus assimilatory power—ATP and NADPH.

Oxidative stress

Reactive oxygen species (ROS) are produced in plants as by products during many metabolic reactions, such as photosynthesis and respiration. Oxidative stress occurs when there is a serious imbalance between the production of ROS and antioxidant defense, causing rapid cell damage by triggering a chain reaction. In the present study, after a 72h-4°C exposure, a total of 270 DEGs were assigned to the GO term “oxidation-reduction process” (GO:0055114), including 177 down-regulated and 93 up-regulated with the fold changes of -10.54 to 13.16, compared to the control. The expression changes of these genes could affect some oxidation-reduction activities and trigger oxidative stress. GO enrichment analysis of DEGs indicated that GO terms “cellular response to reactive oxygen species” (GO:0034614) and “cellular response to oxidative stress” (GO:0034599) were enriched in the 72h-cold treatment. Oxidative stress can dramatically promote nonenzymatic lipid peroxidation [40], of which malondialdehyde (MDA) is one of the important end products. Physiological measurement showed that MDA levels in leaves increased significantly during the 72h treatment.

Desiccation stress and inhibition in substance transport

Chilling stress can inhibit water uptake [41]. This study showed that the relative water content in leaves of ‘Meyer’ declined significantly (p < 0.05) after being exposed to 4°C for 2h, and remained at this level during the 72h exposure. Genes characterized by the GO term “response to desiccation” (GO:0009269) were enriched under cold stress for 2h and 72h, demonstrated by the down-regulated expression of two and four putative aquaporin unigenes, respectively (S10 Table). These results indicated that ‘Meyer’ became water stressed after cold treatment. Aquaporins are intrinsic membrane proteins that mediate the transport of water, small neutral solutes and CO2. It is well documented that regulation of aquaporins plays a role in early responses and acclimation to water stress [42]. In Arabidopsis, 12 aquaporin genes are down-regulated during cold acclimation and recover to original expression levels during de-acclimation [43]. The over-expression of aquaporin genes of Rhododendron catawbiense in Arabidopsis plants leads to a significant decrease in cold tolerance and the increase of dehydration rate of leaves [43]. The down-regulation of aquaporin unigenes thus might help to prevent water loss in ‘Meyer’ leaves under cold stress.

Desiccation stress might affect the activities of transporters. GO enrichment analysis of DEGs indicated that genes with GO terms “dicarboxylic acid transport” (GO:0006835) and “carbohydrate transport” (GO:0008643) were enriched during 2h-cold stress. After 72h-cold stress, genes with GO terms including “nucleobase transmembrane transporter activity” (GO:0015205), “nitrogen compound transport” (GO:0071705), “amino acid transmembrane transport” (GO:0003333), “peptide transport” (GO:0015833), and “carboxylic acid transport” (GO:0046942), were enriched. A total of49 DEGs, 37 of which were down-regulated, were found to have these GO terms, suggesting that transporter activities are inhibited during cold stress.

The affect of cold stress on “Carbohydrate transport” genes was the most noteworthy. After a 2h-cold exposure, three DEGs with this GO term, one encoding a bidirectional sugar transporter (SWEET1a) and two encoding monosaccharide transporter 3 (MST3) respectively, were up-regulated. SWEET1a mediates both low-affinity uptake and efflux of sugar across the plasma membrane. MST3 in rice is involved in the accumulation of monosaccharides required for cell wall synthesis during cell thickening [44]. It is possible that early sugar transport adjustment is involved in cell wall thickening, which is supported by the enrichment of DEGs with the GO term “cell wall organization or biogenesis” in the 2h-cold treatment. During the 72h-cold exposure, 14 DEG sassigned to the GO term “carbohydrate transport” were differentially regulated (S11 Table), including the nine down-regulated and the five up-regulated genes. Unigenes encoding sugar transporter proteins (STP), sucrose transport protein (SUT), hexose carrier, and polyol transporter were down-regulated. SUTs are responsible for the transport of sucrose into the cell with the concomitant uptake of protons (symport system) [45]. STPs have been shown to catalyze proton symport of hexoses across the plasma membrane into the cell, as demonstrated by heterologous expression in yeast [46]. The down-regulated expression of SUT and STP will restrict the phloem loading activity and, thus can affect the long distance transport of sugar. Unigenes encoding putative isoform plastidic glucose transporter 2, vacuolar monosaccharide symporter, and CMP-sialic acid transporter 3-like isoform X1 were up-regulated in the cold. In Arabidopsis thaliana, plastidic glucose transporter 4 is located in the chloroplast inner membrane and may be involved in the efflux of glucose into the cytosol [47]; the tonoplast monosaccharide transporter 1 has been shown to be involved in the transport of monosaccharides into the vacuole during cold stress [48]; and a CMP-sialic acid transporter is responsible for transporting CMP-sialic acid to the Golgi lumen [49]. It could be inferred that intracellular sugar transport was stimulated in ‘Meyer’ during the exposure to low temperature. Abnormal sugar transport in the cold could possibly lead to abnormal distribution of photosynthetic products in cells. Rogers et al [32] found that, when zoysiagrass was chilled at 10°C the chloroplasts accumulated starch and were substantially inflated.

Resistance to desiccation stress

In response to drought stress brought about by low temperature treatment, osmotic adjustment occurs in ‘Meyer’, including the adjustment of potassium ion homeostasis and the accumulation of proline and soluble sugars. GO and KEGG enrichment analysis indicated that genes with the GO terms “potassium ion homeostasis” (GO:0055075), “proline biosynthetic process” (GO:0006561), and “starch and sucrose metabolic” KEGG pathway (ko00500) were enriched in ‘Meyer’ during the 72h-4°C stress.

Potassium (K+) is the most important and abundant cation in living plant cells and plays crucial roles in osmoregulation. In the enriched GO term “potassium ion homeostasis” two unigenes encoding calcineurin B-likes (CBLs) and ABC transporter C family member 5-like (abcC5-like) respectively were up-regulated, while a unigene encoding calcium-activated outward-rectifying potassium channel 1(KCa1) was down-regulated significantly (S10 Table). CBLs are commonly related to abiotic stress tolerance in plants [50], participating in the regulation of K+ uptake [51]; abcC5 is involved in regulation of K+ and Na+ cell content, and is required for regulation of stomatal opening by auxin, abscisic acid (ABA) and external Ca2+; KCa1, also known as two-pore potassium channel 1(TPK1), is apotassium-selective tonoplast ion channel functioning in K+ transport across the vacuolar membrane and plays a role in intracellular K+ homeostasis. The up-regulation of CBL and abcC5 accompanied with the down- regulation of KCa1 might change the distribution of K+ in apoplast, cytoplast and vacuolar to adapt to cold stress.

Proline is an important osmolyte in plant cells. It is also considered to be a scavenger of free radicals, an energy sink and a stress-related signal. The synthesis of proline is stimulated under multiple abiotic stresses [52]. In the present study, 6 DEGs with the GO term “proline biosynthetic process” (GO:0006561) were found to be up-regulated during the 72h-cold (S10 Table) but not during the 2h-cold exposure, which corresponded to the accumulation of proline in leaves.

Soluble sugars play an important role in osmotic adjustment. Sucrose and glucose accumulate and comprise the majority of the total soluble sugar in the stolon and rhizome of zoysiagrass during cold acclimation [8]. This study showed that the level of total soluble sugar went up significantly at the 72h-cold time point. At the transcriptional level, 102 DEGs were assigned to the “starch and sucrose metabolism” KEGG pathway (ko00500), including 34 up-regulated and 68 down-regulated. Among these, unigenes encoding starch branching enzyme (SBE), starch phosphorylase and α-amylase were up-regulated with fold changes of 3 to 12 in the 72h-cold treatment. SBE influence the structure of starch by catalyzing the formation of a-1,6-branch points with varied frequency and branch chain length. Overexpression of SBEIIb in an SBEIIb-deficient mutant of rice resulted in the accumulation of excessive branched, water-soluble polysaccharides instead of amylopectin [53]. The changes of these genes suggest that amylohydrolysis and starch-branching may contribute to the accumulation of soluble sugars in leaves. Sucrose synthase (SUS), a key enzyme involved in sucrose metabolism, catalyzes the reversible conversion of sucrose and UDP to UDP-glucose and fructose. Induced expression and increased activity of SUSs occurs in several species, such as Arabidopsis [54], wheat [55] and chickpea [56]. In the present study, three of the four SUS unigenes identified were induced with fold changes of 2.50, 2.35 and 1.29, respectively. During the reproductive stage of chickpea, a significantly higher specific activity of SUS was observed in seeds of most of the cold-tolerant lines but not in the cold susceptible lines [56]. It may be worthwhile to determine whether SUSs can be used as candidate genes to evaluate cold tolerance of different germplasm in Zoysia genus.

Resistance to oxidative stress

Antioxidants play an important role in combating oxidative stress. Ascorbic acid (AsA) functions as a major redox buffer and as a cofactor for enzymes involved in regenerating other antioxidants. In the present study, genes with “ascorbate and aldarate metabolism” KEGG pathway were enriched at the 72h time point, in which two encoding GME (GDP-Mannose 3',5'-epimerase) were induced 2.3-fold and 6.2-fold, respectively. GME catalyzes an important step in the AsA biosynthetic pathway in higher plants. Overexpression of GMEs in tomato leads to AsA accumulation and enhanced tolerance to oxidative, cold, and salt stress [57]. The up-regulated expression of GME may enhance tolerance to oxidative stress in ‘Meyer’ exposed to low temperature.

Methionine (Met), a sulfur-containing amino acid, is highly susceptible to oxidation [58]. In arsenic-exposed rats methionine can result in the reversal of oxidative stress with a significant (P < 0.05) decline in tissue arsenic burden [59]. In the present study, genes with GO terms “Cysteine and methionine metabolism” KEGG pathway (ko00270) were enriched in ‘Meyer’ after a 72h-cold exposure, among which four encoding cobalamin-independent methionine synthase (MetE), homocysteine S-methyltransferase (HMT), S-adenosyl-l-methionine synthase (SAMS), and methionine γ-lyase (MGL) respectively, were up-regulated. Both MetE and HMT catalyze the synthesis of methionine, whose role in the stress response has been largely unreported. Methionine is converted to S-adenosylmethionine (SAM) which serves as a precursor for metabolites including glycinebetaine, methylated polyols, polyamines and ethylene. MGL is involved in cellular Met homeostasis and Ile synthesis. Transgenic Arabidopsis lines overexpressing a SAMS gene from potato exhibited higher salt and drought stress tolerance compared to controls [60]. MGL of Arabidopsis was up-regulated by simultaneous water deficit and nematode stress in root tissue, and conferred resistance to nematodes when overexpressed [61]. The up-regulation of the four genes could increase the level of Met and its derivatives, which may alleviate oxidative damage and enhance tolerance to cold in ‘Meyer’. The induced expression of MetE and HMT genes by low temperature provides a new insight into the cold response in plants.

Antioxidant enzymes such as superoxide dismutase (SOD), catalase (CAT) and peroxidase (POD) are important in reducing levels of ROS. Several of these enzymes have been shown to be induced under cold stress at either the transcript or protein level [62]. In the present study, the activity of POD increased during the 2h-cold stress, but declined significantly at the 72h-cold time point; the activity of SOD increased significantly during the 72h-cold stress; the activity of CAT declined during the cold stress. Compared with the control, two POD unigenes were transcriptionally up-regulated in the 2h-cold treatment, one SOD unigene was up-regulated but 12 POD unigenes were down-regulated in the 72h-cold treatment. These findings were consistent with the changes of POD and SOD activity. Three unigenes for CAT were identified and designated as comp35050_c4, comp30190_c0 and comp30190_c4. Under cold stress, the expression of comp35050_c4 increased, while that of the other two decreased. The decrease of CAT activity under cold stress suggests that the protein product of comp30190_c0 and comp30190_c4 could have great activity than that of comp35050_c4. The lower activities of CAT and POD under cold stress may be one of reasons why zoysiagrasses are cold sensitive. Increasing antioxidant enzymes activities should be considered in cold tolerance breeding research.

Expression changes of DREBs and COR/dehydrin genes

The DREBs (Dehydration Responsive Element Binding proteins/C-repeat Binding Factors) proteins namely, DREB1 and DREB2, are thought to be the major transcription factors (TFs) that control stress-inducible gene expression in the ABA-independent pathway [63]. Overexpressing DREB genes in plants can significantly enhance plant tolerance to multiple abiotic stresses, including drought, freezing and high salinity [63]. Moreover, findings in some plants indicate that the different stress tolerance ability is related to diverse gene structures and expression patterns of DREB1 and DREB2. For example, in winter cultivars of wheat (Tritic eaeaestivum and T. monococcum) the constitutive and inducible expression of DREB1/CBFs within five pooideae-specific groups may play a predominant role in their superior freezing tolerance ability [6465]; Among Arabidopsis natural populations, the relaxed selection on the CBF/DREB1s in some accessions from the southern part of Europe has compromised the ability of these genes to act as efficient transcriptional activators during the cold acclimation process [66]; In Zea mays, it was reported that natural variation in drought tolerance at the seedling stage was related to the variation of ZmDREB2.7 [67]. In this study, eight and five unigenes assigned to DREB1 and DREB2 subgroups, respectively, were identified (Table 5). Compared with the control, six DREB1 and two DREB2 unigenes were up-regulated significantly at the 2h-cold time point, and one DREB2 unigene was up-regulated at the 72h-cold time point. The induced expression of these DREB1 and DREB2 genes will be involved in increasing cold resistance in zoysiagrass. The relationship between the DREBs and cold tolerance in zoysiagrass is worthy of further research.

thumbnail
Table 5. Expression changes of unigenes of DREB1, DREB2 and COR/dehydrin.

https://doi.org/10.1371/journal.pone.0131153.t005

A number of cold-regulated (COR) genes have been characterized in plants and some are identified as regulons of DREBs [7]. The dehydrins, belonging to the late-embryogenesis abundant (LEA) class of proteins [68], are among the COR genes that have been extensively studied [69] and shown to have major roles in conferring low temperature tolerance. Studies in wheat (T.aestivum L.) have shown that the expression of the COR genes Wcs120, Wcor410, and Wcor14 was highest for low temperature-tolerant and lowest for more cold sensitive genotypes [70]. In Arabidopsis, the expression levels of COR genes are higher in more cold-tolerant leaves than in cold-sensitive pollen [71]. In this study, ten unigenes annotated as dehydrin/COR were detected, six of which were up-regulated several fold during 72h-cold stress, indicating a potential relationship to the cold tolerance ability in ‘Meyer’.

The relationship of the COR413 family to cold tolerance is of interest. The COR413 family consists of two distinct groups: plasma membrane targeted (COR413-PM) and thylakoid membrane targeted (COR413-TM), respectively [72]. The expressions of one gene of group COR413-PM and one gene of group COR413-TM are closely associated with the acquisition of freezing tolerance in several plant species such as wheat, rye, and Arabidopsis. On the other hand, group Cor413-pm and-tm transcripts are not induced in the low temperature-sensitive species rice and maize [72]. In the present study, a unigene encoding COR413-pm 1 was induced by low temperature in ‘Meyer’, which is the reverse of the expression profiles in the low temperature-sensitive species mentioned above. Determining whether this induced expression is related to cold tolerance in zoysiagrass requires further research.

Adjustment in phytohormone metabolism and signaling pathway

Phytohormones play important roles in the response of plants to stress. In this study, genes with GO terms “gibberellin metabolic process” (GO:0009685), “response to abscisic acid stimulus”(GO:0009737), and “response to jasmonic acid stimulus” (GO:0009753) were enriched after the 72h-cold exposure. The expression changes of DEGs in these GO terms are shown in S12 Table.

Gibberellins (GAs) play a central role in plant development and the metabolism of GAs in plants is highly regulated [73]. From a regulatory standpoint, two biosynthetic enzymes—GA 20-oxidase (GA20ox) and GA 3-oxidase (GA3ox)—and an enzyme involved in gibberellin catabolism, GA 2-oxidase (GA2ox), are most notable. The gene transcription of these enzymes is highly regulated. In this study, three unigenes of GA20ox were down-regulated; three of GA2ox were up-regulated (Fig 9). These changes suggest that gibberellin catabolic pathways were activated and lead to the decreased endogenous level of bioactive GA just like other plant species [74].

thumbnail
Fig 9. Gibberellin metabolism in ‘Meyer’ under cold stress.

GA20ox, GA 20-oxidase; GA3ox, GA3-oxidase; GA2ox, GA 2-oxidase. Numbers in brackets represent numbers of unique transcripts regulated by cold stress. Red and blue indicate up- and down-regulated by low temperature, respectively.

https://doi.org/10.1371/journal.pone.0131153.g009

ABA serves as an endogenous messenger in abiotic stress responses in plants, and has been called a ‘stress hormone’. Genes with the GO term of “response to abscisic acid stimulus” was enriched in ‘Meyer’ during the 72h-4°C exposure. Three regulators of ABA signal transduction pathways, including type 2C protein phosphatases (PP2C), sucrose nonfermenting 1-related protein kinase2 (SnRK2/ SAPK4), and proline-rich receptor-like protein kinase (PERK4), were up-regulated significantly in the cold. Without ABA, PP2C dephosphorylates SnRK2 kinases and blocks signal transduction. The ABA-bound receptor inactivates PP2C, which causes the SnRK2 kinases to be in a phosphorylated state, allowing them to activate bZIP group transcription factors (TFs) that then promote downstream ABA-responsive gene transcription [7577]. PERK4, localized to the plasma membrane, has been shown to be an ABA- and Ca2+-activated protein kinase, is required for abscisic acid responses in Arabidopsis [78]. In this study, the expression of three PP2C unigenes, one SnRK2/ SAPK4 unigene and one PERK4 unigene were induced with fold changes of 2.53 to 7.10, suggesting that the ABA signaling pathway was stimulated in ‘Meyer’ exposed to cold.

Jasmonates, including jasmonic acid (JA) and its related metabolites, are essential components of the signaling pathway triggering the expression of plant defense genes in response to various sources of environmental stress [79]. Levels of endogenous jasmonates have been reported to increase upon pathogen infection, insect wounding, and various abiotic stresses [8081]. In this study, genes with the GO term “response to jasmonic acid stimulus” (GO:0009753) was enriched during 72h cold stress, among which the regulation of allene oxide synthase (AOS) unigenes and indole-3-acetic acid-amidosynthetase GH3.12 unigenes were noteworthy. AOS is a major control point in Arabidopsis thaliana octadecanoid signaling. This enzyme acts on a number of unsaturated fatty-acid hydroperoxides forming the corresponding allene oxides. Several AOS enzymes have been cloned from various plant species. Most of them prefer 13-hydroperoxides as substrates and are involved in JA production [8283]. JA can also induce the expression of AOS [84]. OsGH3.12 belongs to group I of GH3 genes in plants [85], which can conjugate JA with isoleucine to produce bioactive JA-Ile [86]. In ‘Meyer’ under 4°C stress for72h, two AOS genes were down-regulated with fold changes of -2 to -5, while AOS2-like, AOS3-like and AOS4-like were up-regulated 4 to 8-fold, and three OsGH3.12 homologs were down-regulated-7 to -9-fold. The expression changes of these unigenes suggest that metabolism and signaling pathways of JA are significantly regulated by cold stress in ‘Meyer’, which has also been reported in rice [81].

Cross-talk between the different plant hormones results in synergistic or antagonistic interactions that play crucial roles in the response to abiotic stress in many plants [87]

As shown in this study, the expression changes of multiple unigenes that are involved in metabolism and signaling pathway of different plant hormones might indicate a coordination of the regulatory network of stress signaling and cold response in zoysiagrass.

Conclusion

In conclusion, this study provided the first large-scale transcriptome data set in warm-season turfgrass in response to low-temperature stress. More than 206 million high-quality 100-bp paired-end reads were generated and assembled into 46,412 unigenes in ‘Meyer’ zoysiagrass. Furthermore, 756 and 5327 unigenes were identified as DEGs in 2h- and 72h-4°C treatments, respectively; in addition, their classification, functions, and metabolic pathways are shown here for the first time. The down-regulated expression of genes involved in electron carrier activity, chlorophyll biosynthesis and carbon-fixation activity, water transport, organic substance transport, and antioxidant enzyme activity are consistent with the cold sensitivity of ‘Meyer’. The induced expression of genes involved in osmotic adjustment, antioxidants, nitrogen assimilation, phytohormone metabolism and signaling, and the DREB response pathway may play a vital role in cold resistance. The high-resolution expression patterns presented here further our understanding of the molecular response to low temperature in zoysiagrass, and would facilitate the breeding of cold tolerant varieties.

Supporting Information

S2 Table. The DEGs between 2h-cold stressand CK.

https://doi.org/10.1371/journal.pone.0131153.s002

(XLSX)

S3 Table. The DEGs between 72h-cold stress and CK.

https://doi.org/10.1371/journal.pone.0131153.s003

(XLSX)

S4 Table. The enriched GO terms under 2h-cold stress.

https://doi.org/10.1371/journal.pone.0131153.s004

(XLSX)

S5 Table. The enriched GO terms under 72h-cold stress.

https://doi.org/10.1371/journal.pone.0131153.s005

(XLSX)

S6 Table. The enriched KEGG pathwaysunder 2h-cold stress.

https://doi.org/10.1371/journal.pone.0131153.s006

(XLSX)

S7 Table. The enriched KEGG pathways under 72h-cold stress.

https://doi.org/10.1371/journal.pone.0131153.s007

(XLSX)

S8 Table. The expression changes of DEGs related to photosystem and photosynthesis.

https://doi.org/10.1371/journal.pone.0131153.s008

(XLSX)

S9 Table. The expression changes of DEGs related to nitrate assimilation.

https://doi.org/10.1371/journal.pone.0131153.s009

(XLSX)

S10 Table. The expression changes of DEGs involved in desiccation stress and osmotic adjustment.

https://doi.org/10.1371/journal.pone.0131153.s010

(XLSX)

S11 Table. The expression changes of DEGs in GO term of carbohydrate transport.

https://doi.org/10.1371/journal.pone.0131153.s011

(XLSX)

S12 Table. The expression changes of DEGs related to plant hormone metabolism or stimulus.

https://doi.org/10.1371/journal.pone.0131153.s012

(XLSX)

Acknowledgments

The authors would like to thank Dr. Yongqiang Qian and Professor Zhenyuan Sun for helping grow plant materials.

Author Contributions

Conceived and designed the experiments: SW FG. Performed the experiments: SW XK J. Li. Analyzed the data: SW ZD. Contributed reagents/materials/analysis tools: SW ZD J. Liu. Wrote the paper: SW YZ.

References

  1. 1. Cai H, Inoue M, Yuyama N, Takahashi W, Hirata M, Sasaki T. Isolation, characterization and mapping of simple sequence repeat markers in zoysiagrass (Zoysia spp.). Theor Appl Genet. 2005 Dec;112(1):158–166. pmid:16235048
  2. 2. Yaneshita M, Kaneko S, Sasakuma T. Allotetraploidy of Zoysia species with 2n = 40 based on a RFLP genetic map. Theor Appl Genet. 1999; 98(5): 751–756. 1999 April; 98(5): 751–756.
  3. 3. Engelke MC, Anderson S. Zoysiagrasses (Zoysia spp.) In: Turfgrass Biology, Genetics, and Breeding (Eds Casler MD, Duncan RR). John Wiley & Sons, Inc, New Jersey. 2003.pp: 271–286.
  4. 4. Mao PS, Wang XG, Wang YH, Han JG. Effect of storage temperature and duration on the vigor of zoysiagrass (Zoysia japonica Steud.) seed harvested at different maturity stages. Grassland Sci. 2009Mar; 55(1):1–5.
  5. 5. Beard JB (1973) Turfgrass: science and culture. Englewood Cliffs: Prentice Hall. 261pp.
  6. 6. Thomashow MF. Plant cold acclimation: freezing tolerance genes and regulatory mechanisms. Annu Rev Plant Physiol Plant Mol Biol. 1999 Jun; 50:571–599. pmid:15012220
  7. 7. Carvallo MA, Pino MT, Jeknic Z, Zou C, Doherty CJ, Shiu SH, et al. A comparison of the low temperature transcriptomes and CBF regulons of three plant species that differ in freezing tolerance: Solanum commersonii, Solanum tuberosum, and Arabidopsis thaliana. J Exp Bot. 2011Jul; 62(11):3807–3819. pmid:21511909
  8. 8. Patton AJ, Cunningham SM, Volenec JJ, Peicher ZJ. Differences in freeze tolerance of zoysiagrasses: II. Carbohydrate and proline accumulation. Crop Sci. 2007 Sep; 47: 2170–2181.
  9. 9. Patton AJ, Cunningham SM, Volenec JJ, Reicher ZJ. Differences in freeze tolerance of zoysiagrasses: I. Role of proteins. Crop Sci. 2007 Sep; 47(5): 2162–2169.
  10. 10. Zhang Q, Fry J, Rajashekar C, Bremer D, Engelke M. Membrane polar lipid changes in zoysiagrass rhizomes and their potential role in freezing tolerance. J Am Soc Hort Sci. 2009 May; 134(3): 322–328.
  11. 11. Guo HL, Xuan JP, Liu JX, Zhang YM, Zheng YQ. Association of molecular markers with cold tolerance and green period in zoysiagrass (Zoysia Willd.). Breed Sci.2012 Dec; 62(4):320–327. pmid:23341745
  12. 12. Xuan J, Song Y, Zhang H, Liu J, Guo Z, Hua Y. Comparative Proteomic Analysis of the Stolon Cold Stress Response between the C4 Perennial Grass Species Zoysia japonica and Zoysia metrella. PLoS One. 2013 Sep 26;8(9):e75705. pmid:24086619
  13. 13. Barah P, Jayavelu N D, Rasmussen S, Nielsen HB, Mundy J and Bone AM. Genome-scale cold stress response regulatory networks in ten Arabidopsis thaliana ecotypes. BMC Genomics. 2013 Oct 22;14:722. pmid:24148294
  14. 14. Varshney RK, Nayak SN, May GD, Jackson SA. Next-generation sequencing technologies and their implications for crop genetics and breeding. Trends Biotechnol. 2009 Sep;27(9):522–530. pmid:19679362
  15. 15. Wang J, Yang Y, Liu X, Huang J, Wang Q, Gu J, et al. Transcriptome profiling of the cold response and signaling pathways in Lilium lancifolium. BMC Genomics. 2014 Mar 17;15:203. pmid:24636716
  16. 16. Hu T, Sun X, Zhang X, Nevo E, and Fu J. An RNA sequencing transcriptome analysis of the high-temperature stressed tall fescue reveals novel insights into plant thermotolerance. BMC Genomics. 2014 Dec 19;15:1147. pmid:25527327
  17. 17. Jain D, Chattopadhyay D. Analysis of gene expression in response to water deficit of chickpea (Cicer arietinum L.) varieties differing in drought tolerance. BMC Plant Biology. 2010 Feb 9;10:24. pmid:20144227
  18. 18. Draper H H & Hadley M. Malondialdehyde determination as index of lipid peroxidation. Method Enzymol. 1990; 186: 421–431.
  19. 19. Li L, Yi H. Effect of sulfur dioxide on ROS production, gene expression and antioxidant enzyme activity in Arabidopsisplants. Plant Physiol Biochem. 2012 Sep;58:46–53. pmid:22771435
  20. 20. Hu T, Sun X, Zhang X, Nevo E, and Fu J. An RNA sequencing transcriptome analysis of the high-temperature stressed tall fescue reveals novel insights into plant thermotolerance. BMC Genomics.2014 Dec 19;15:1147. pmid:25527327
  21. 21. Haas BJ, Papanicolaou A, Yassour M, Grabherr M, Blood PD, Bowden J, et al. De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nat Protoc.2013 Aug;8(8):1494–1512. pmid:23845962
  22. 22. Conesa A, Götz S, García-Gómez JM, Terol J, Talón M, Robles M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005 Sep 15;21(18):3674–3676. pmid:16081474
  23. 23. Li B and Dewey CN. RSEM: accurate transcript quantification from RNA-seq data with or without a reference genome. BMC Bioinformatics.2011 Aug 4;12:323. pmid:21816040
  24. 24. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010 Jan 1;26(1):139–140. pmid:19910308
  25. 25. Zheng Q, Wang XJ. GOEAST: a web-based software toolkit for Gene Ontology enrichment analysis. Nucleic Acids Res. 2008 Jul 1;36(Web Server issue):W358–63. pmid:18487275
  26. 26. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2-ΔΔCT method. Methods. 2001 Dec; 25(4): 402–408. pmid:11846609
  27. 27. Pramanik MH, Imai R. Functional Identification of a Trehalose 6-phosphate phosphatase gene that is involved in transient induction of trehalose biosynthesis during chilling Stress in Rice. Plant Mol Biol.2005 Aug;58(6):751–762. pmid:16240171
  28. 28. Ge L, Chao D, Shi M, Zhu M, Gao J, Lin H. Overexpression of the trehalose-6-phosphate phosphatase gene OsTPP1 confers stress tolerance in rice and results in the activation of stress responsive genes. Planta. 2008 Jun;228(1):191–201. pmid:18365248
  29. 29. Allen DJ and Ort DR. Impacts of chilling temperatures on photosynthesis in warm-climate plants. Trends Plant Sci. 2001 Jan;6(1):36–42. pmid:11164376
  30. 30. Pruzinska A, Tanner G, Anders I, Roca M, Hortensteiner S. Chlorophyll breakdown: pheophorbide a oxygenase is a Rieske-type iron-sulfur protein, encoded by the accelerated cell death1 gene.Proc Natl Acad Sci U S A. 2003 Dec 9;100(25):15259–15264. pmid:14657372
  31. 31. Sakuraba Y, Kim YS, Yoo SC, Hortensteiner S, Paek NC. 7-Hydroxymethyl chlorophyll a reductase functions in metabolic channeling of chlorophyll breakdown intermediates during leaf senescence. Biochem Biophys Res Commun. 2013 Jan 4;430(1):32–37. pmid:23200839
  32. 32. Rogers RA, Dunn JH, Nelson CJ. Photosynthesis and cold hardening in Zoysia and bermudagrass. Crop Sci. 1977 Jan; 17(5): 727–732.
  33. 33. Du YC, Nose A, Wasano K. Effects of chilling temperature on photosynthetic rates, photosynthetic enzyme activities and metabolite levels in three sugarcane species. Plant Cell Environ. 1999 Mar; 22(3): 317–324.
  34. 34. Chinthapalli B, Murmu J, Raghavendra AS. Dramatic difference in the responses of phosphoenolpyruvate carboxylase to temperature in leaves of C3 and C4 plants. J Exp Bot. 2003 Feb; 54(383):707–714. pmid:12554714
  35. 35. Matsuba K, imaizumi N, Kaneko S, Ohsugi R. Photosynthetic responses to temperature of phosphoenolpyruvate carboxykinase type C4 species differing in cold sensitivity. Plant Cell Environ. 1997; 20: 268–274.
  36. 36. Ritchie SW, Redinbaugh MG, Shiraishi N, Vrba JM, Campbell WH. Identification of a maize root transcript expressed in the primary response to nitrate: characterization of a cDNA with homology to ferredoxin-NADP+oxidoreductase. Plant Mol Biol.1994 Oct;26(2):679–690 pmid:7948922
  37. 37. Matsumura T, Sakakibara H, Nakano R, Kimata Y, Sugiyama T,Hase T. Anitrate-inducible ferredoxin in maize roots. Genomic organization and differential expression of two nonphotosynthetic ferredoxin isoproteins. Plant Physiol.1997 Jun;114(2):653–660. pmid:9193097
  38. 38. Yaneva I, Mäck G, Vunkova-Radeva R, Tischner R. Changes in nitrate reductase activity and the protective effect of molybdenum during cold stress in winter wheat grown on acid soil. J Plant Physiol. 1996; 149(1–2): 211–216.
  39. 39. Hoshida H, Tanaka Y, Hibino T, Hayashi Y, Tanaka A, Takabe T, et al. Enhanced tolerance to salt stress in transgenic rice that overexpresses chloroplast glutamine synthetase. Plant Mol Biol.2000 May;43(1):103–111. pmid:10949377
  40. 40. Sattler SE, Mène-Saffrané L, Farmer EE, Krischke M, Mueller MJ, and Penna DD. Nonenzymatic lipid peroxidation reprograms gene expression and activates defense markers in arabidopsis tocopherol-deficient Mutants. Plant Cell.2006 Dec;18(12):3706–3720. pmid:17194769
  41. 41. Chinnusamy V, Zhu J, and Zhu JK. Cold stress regulation of gene expression in plants. Trends Plant Sci. 2007 Oct;12(10):444–451. pmid:17855156
  42. 42. Chaumont F and Tyerman SD. Aquaporins: highly regulated channels controlling plant water relations. Plant Physiol.2014 Apr;164(4):1600–1618. pmid:24449709
  43. 43. Peng Y, Arora R, Li G, Wang X, Fessehaie A. Rhododendron catawbiense plasma membrane intrinsic proteins are aquaporins, and their over-expression compromises constitutive freezing tolerance and cold acclimation ability of transgenic Arabidopsis plants. Plant Cell Environ. 2008 Sep;31(9):1275–1289. pmid:18518915
  44. 44. Toyofuku K, Kasahara M, Yamaguchi J. Characterization and Expression of Monosaccharide Transporters (OsMSTs) in Rice. Plant Cell Physiol.2000 Aug;41(8):940–947. pmid:11038054
  45. 45. Srivastava AC, Ganesan S, Ismail IO, Ayre BG. Functional characterization of the Arabidopsis AtSUC2 Sucrose/H+ symporter by tissue-specific complementation reveals an essential role in phloem loading but not in long-distance transport.Plant Physiol.2008 Sep;148(1):200–211. pmid:18650401
  46. 46. Bu¨ttner M. The monosaccharide transporter(-like) gene family in ArabidopsisFEBS Lett.2007 May 25;581(12):2318–2324. pmid:17379213
  47. 47. Weber A, Servaites JC, Geiger DR, Kofler H, Hille D, Gröner F, et al. Identification, purification, and molecular cloning of a putative plastidic glucose translocator. Plant Cell. 2000 May; 12: 787–802. pmid:10810150
  48. 48. Wormit A, Trentmann O, Feifer I, Lohr C, Tjaden J, Meyer S,et al. Molecular identification and physiological characterization of a novel monosaccharide transporter from Arabidopsis involved in vacuolar sugar transport. Plant Cell. 2006 Dec;18(12):3476–3490. pmid:17158605
  49. 49. Bakker H, Routier F, Ashikov A, Neumann D, Bosch D, Gerardy-Schahn R. A CMP-sialic acid transporter cloned from Arabidopsis thaliana. Carbohydr Res.2008 Aug 11;343(12):2148–2152. pmid:18258224
  50. 50. Kim BG, Waadt R, Cheong YH, Pandey GK, Dominguez-Solis JR, Schültke S, et al. The calcium sensor CBL10 mediates salt tolerance by regulating ion homeostasis in Arabidopsis. Plant J.2007 Nov;52(3):473–484. pmid:17825054
  51. 51. Xu J, Li HD, Chen LQ, Wang Y, Liu LL, He L, et al. A protein kinase, interacting with two calcineurin B-like proteins, regulates K+ transporter AKT1 in Arabidopsis. Cell. 2006 Jun 30;125(7):1347–1360. pmid:16814720
  52. 52. Szabados L, Savouré A. Proline: a multifunctional amino acid. Trends Plant Sci.2010 Feb;15(2):89–97. pmid:20036181
  53. 53. Tanaka N, Fujita N, Nishi A, Satoh H, Hosaka Y, Ugaki M, et al. The structure of starch can be manipulated by changing the expression levels of starch branching enzyme IIb in rice endosperm. Plant Biotechnol J. 2004 Nov;2(6):507–516. pmid:17147623
  54. 54. Baud S, Vaultier M and Rochat C. Structure and expression profile of the sucrose synthase multigene family in Arabidopsis.J Exp Bot.2004 Feb;55(396):397–409. pmid:14739263
  55. 55. Crespi MD, Zabaleta EJ, Pontis HG, and Salerno GL. Sucrose synthase expression during cold acclimation in wheat. Plant Physiol. 1991 Jul; 96(3):887–891. pmid:16668270
  56. 56. Kaur S, Gupta AK, Kaur N, Sandhu JS, and Gupta SK. Antioxidative Enzymes and Sucrose Synthase Contribute to Cold Stress Tolerance in Chickpea. J Agron Crop Sci.2009 Oct; 195(5): 393–397.
  57. 57. Zhang C, Liu J, Zhang Y, Cai X, Gong P, Zhang J, et al. Overexpression of SlGMEs leads to ascorbate accumulation with enhanced oxidative stress, cold, and salt tolerance in tomato.Plant Cell Rep. 2011 Mar;30(3):389–398. pmid:20981454
  58. 58. Davies MJ. The oxidative environment and protein damage. Biochim Biophys Acta. 2005 Feb; 1703:93–109. pmid:15680218
  59. 59. Nandi D, Patra RC, Swarup D. Effect of cysteine, methionine, ascorbic acid and thiamine on arsenic-induced oxidative stress. Toxicology.2005 Jul 1;211(1–2):26–35. pmid:15863245
  60. 60. Kim SH, Kim SH, Palaniyandi SA, Yang SH, Suh JW. Expression of potato S-adenosyl-l-methionine synthase (SbSAMS) gene altered developmental characteristics and stress responses in transgenic Arabidopsis plants.Plant Physiol Biochem.2015 Feb;87:84–91. pmid:25559387
  61. 61. Atkinson NJ, Lilley CJ, and Urwin PE. Identification of genes involved in the response of Arabidopsis to simultaneous biotic and abiotic Stresses.Plant Physiol.2013 Aug;162(4):2028–2041. pmid:23800991
  62. 62. Yan SP, Zhang QY, Tang ZC, Su WA, Sun WN. Comparative proteomic analysis provides new insights into chilling stress responses in rice.Mol Cell Proteomics.2006 Mar;5(3):484–496. pmid:16316980
  63. 63. Lata C, Prasad M. Role of DREBs in regulation of abiotic stress responses in plants. J Exp Bot. 2011 Oct;62(14):4731–4748. pmid:21737415
  64. 64. Badawi M, Danyluk J, Boucho B, Houde M, Sarhan F. The CBF gene family in hexaploid wheat and its relationship to the phylogenetic complexity of cereal CBFs. Mol Genet Genomics. 2007 May;277(5):533–554. pmid:17285309
  65. 65. Knox AK, Li C, Vagujfalvi A, Galiba G, Stockinger EJ, Dubcovsky J. Identification of candidate CBF genes for the frost tolerance locus Fr-Am2 in Triticum monococcum.Plant Mol Biol.2008 Jun;67(3):257–270. pmid:18317935
  66. 66. Zhen Y and Ungerer MC. Relaxed Selection on the CBF/DREB1 Regulatory Genes and Reduced Freezing Tolerance in the Southern Range of Arabidopsis thaliana. Mol Biol Evol.2008 Dec;25(12):2547–2555. pmid:18775899
  67. 67. Liu S, Wang X, Wang H, Xin H, Yang X, Yan JB, et al. Genome-Wide Analysis of ZmDREB Genes and Their Association with Natural Variation in Drought Tolerance at Seedling Stage of Zea mays L. PLoS Genet.2013;9(9):e1003790. pmid:24086146
  68. 68. Thomashow MF. Role of cold-responsive genes in plant freezing tolerance. Plant Physiol. 1998 Sep;118(1):1–8. pmid:9733520
  69. 69. Borovskii GB, Stupnikova IV, Antipina AI, Voinikov VK. Accumulation of dehydrins and ABA-inducible proteins in wheat seedlings during low-temperature acclimation. Russ J Plant Physl+. 2002 Feb; 49(2): 229–234.
  70. 70. Ganeshan S, Vitamvas P, Fowler B and Chibbar RN. Quantitative expression analysis of selected COR genes reveals their differential expression in leaf and crown tissues of wheat (Triticum aestivum L.) during an extended low temperature acclimation regimen.J Exp Bot.2008;59(9):2393–2402. pmid:18508811
  71. 71. Lee JY, Lee DH. Use of serial analysis of gene expression technology to reveal changes in gene expression in Arabidopsis pollen undergoing cold stress.Plant Physiol.2003 Jun;132(2):517–529. pmid:12805584
  72. 72. Breton G, Danyluk J, Charron JF, and Sarhan F. Expression profiling and bioinformatic analyses of a novel stress-regulated multispanning transmembrane protein family from cereals and Arabidopsis. Plant Physiol.2003 May;132(1):64–74. pmid:12746512
  73. 73. Yamauchi Y, Ogawa M, Kuwahara A, Hanada A, Kamiya Y,Yamaguchi S. Activation of gibberellin biosynthesis and response pathways by low temperature during imbibition of Arabidopsis thaliana seeds. Plant Cell. 2004 Feb; 16(2):367–378. pmid:14729916
  74. 74. Vettakkorumakankav NN, Falk D, Saxena P, Fletcher RA.A crucial role for gibberellins in stress protection of plants. Plant & Cell Physiol. 1999(5); 40:542–548.
  75. 75. Rodriguez PL, Benning G, Grill E. ABI2, a second protein phosphatase 2C involved in abscisic acid signal transduction in Arabidopsis. FEBS Lett.1998 Jan 16;421(3):185–190. pmid:9468303
  76. 76. Fujii H, Chinnusamy V, Rodrigues A, Rubio S, Antoni R, et al. In vitro reconstitution of an abscisic acid signaling pathway. Nature.2009 Dec 3;462(7273):660–664. pmid:19924127
  77. 77. Umezawa T, Sugiyama N, Mizoguchi M, Hayashi S, Myouga F, Yamaguchi-Shinozaki K. et al. Type 2C protein phosphatases directly regulate abscisic acid-activated protein kinases in Arabidopsis. Proc Natl Acad Sci U S A.2009 Oct 13;106(41):17588–17593. pmid:19805022
  78. 78. Umezawa T, Nakashima K, Miyakawa T, Kuromori T, Tanokura M, Shinozaki K, et al. Molecular basis of the core regulatory network in ABA responses: sensing, signaling and transport. Plant Cell Physiol.2010 Nov;51(11):1821–1839. pmid:20980270
  79. 79. Turner JG, Ellis C, and Devoto A. The jasmonate signal pathway. Plant Cell. Plant Cell.2002;14 Suppl:S153–164.
  80. 80. Pauwels L, Inze D and Goossens A. Jasmonate-inducible gene: what does it mean? Trends Plant Sci.2009 Feb;14(2):87–91. pmid:19162528
  81. 81. Du H, Liu H and Xiong L. Endogenous auxin and jasmonicacid levels are differentially modulated by abiotic stresses in rice. Front Plant Sci.2013 Oct 9;4:397. pmid:24130566
  82. 82. Froehlich JE, Itoh A and Howe GA. Tomato allene oxide synthase and fatty acid hydroperoxide lyase, two cytochrome P450s involved in oxylipin metabolism, are targeted to different membranes of chloroplast envelope.Plant Physiol.2001 Jan;125(1):306–317. pmid:11154338
  83. 83. Stumpe M, Göbel C, Demchenko K, Hoffmann M, Klösgen RB, Pawlowski K, et al. Identification of an allene oxide synthase (CYP74C) that leads to formation of a-ketols from 9-hydroperoxides of linoleic and linolenic acid in below-ground organs of potato.Plant J.2006 Sep;47(6):883–896. pmid:16899083
  84. 84. Ha SB, Lee BC, Lee DE, Kuk YI, Lee AY, Han O, et al. Molecular characterization of the gene encoding rice allene oxide synthase and its expression.Biosci Biotechnol Biochem.2002 Dec;66(12):2719–2722. pmid:12596875
  85. 85. Terol J, Domingo C, Talón M. The GH3 family in plants: genome wide analysis in rice and evolutionary history based on EST analysis. Gene.2006 Apr 26; 371(2):279–290. pmid:16488558
  86. 86. Staswick PE, Tiryaki I. The oxylipin signal jasmonic acid is activated by an enzyme that conjugates it to isoleucine in Arabidopsis.Plant Cell.2004 Aug;16(8):2117–2127. pmid:15258265
  87. 87. Peleg Z, Blumwald E. Hormone balance and abiotic stress tolerance in crop plants. Curr Opin Plant Biol. 2011 Jun;14(3):290–295. pmid:21377404