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

Transcriptome analysis of nitric oxide-responsive genes in upland cotton (Gossypium hirsutum)

  • Juan Huang,

    Roles Formal analysis, Methodology, Software, Writing – original draft

    Affiliations Research Center of Buckwheat Industry Technology, Guizhou Normal University, Guiyang, Guizhou, P. R. China, State Key Laboratory of Cotton Biology, Institute of Cotton Research, Chinese Academy of Agricultural Sciences, Anyang, Henan, P. R. China

  • Hengling Wei,

    Roles Formal analysis, Writing – review & editing

    Affiliation State Key Laboratory of Cotton Biology, Institute of Cotton Research, Chinese Academy of Agricultural Sciences, Anyang, Henan, P. R. China

  • Libei Li,

    Roles Formal analysis, Writing – review & editing

    Affiliation State Key Laboratory of Cotton Biology, Institute of Cotton Research, Chinese Academy of Agricultural Sciences, Anyang, Henan, P. R. China

  • Shuxun Yu

    Roles Conceptualization, Funding acquisition, Supervision, Validation

    yu@cricaas.com.cn

    Affiliation State Key Laboratory of Cotton Biology, Institute of Cotton Research, Chinese Academy of Agricultural Sciences, Anyang, Henan, P. R. China

Abstract

Nitric oxide (NO) is an important signaling molecule with diverse physiological functions in plants. It is therefore important to characterize the downstream genes and signal transduction networks modulated by NO. Here, we identified 1,932 differentially expressed genes (DEGs) responding to NO in upland cotton using high throughput tag sequencing. The results of quantitative real-time polymerase chain reaction (qRT-PCR) analysis of 25 DEGs showed good consistency. Gene Ontology (GO) and KEGG pathway were analyzed to gain a better understanding of these DEGs. We identified 157 DEGs belonging to 36 transcription factor (TF) families and 72 DEGs related to eight plant hormones, among which several TF families and hormones were involved in stress responses. Hydrogen peroxide and malondialdehyde (MDA) contents were increased, as well related genes after treatment with sodium nitroprusside (SNP) (an NO donor), suggesting a role for NO in the plant stress response. Finally, we compared of the current and previous data indicating a massive number of NO-responsive genes at the large-scale transcriptome level. This study evaluated the landscape of NO-responsive genes in cotton and identified the involvement of NO in the stress response. Some of the identified DEGs represent good candidates for further functional analysis in cotton.

Introduction

Nitric oxide (NO) is a small gaseous molecule with a free radical nature. It participates in a wide range of important biological processes [1]. Since the significant role of NO in plants was realized in the mid-1900s [2], the biosynthesis of NO in plants has been extensively studied, and two NO biosynthetic pathways of plants have been described. The main pathway is dependent on nitrate reductase, which catalyzes the reduction of nitrite to NO [3, 4]. The second pathway is an arginine-dependent reaction catalyzed by NOA1 (previously referred to as NOS1) [5]. However, NOA1 was subsequently found to play a primary role in chloroplast function, and its effects on NO accumulation are likely indirect [6]. Therefore, nitrate reductase appears to be the primary enzyme for the biosynthesis of NO in plants.

NO acts as a key signaling molecule in many biological processes, including seed dormancy and germination, root organogenesis, the transition to flowering, the regulation of pollen tubes, fruit ripening, plant senescence, and responses to biotic and abiotic stresses [4, 710]. In addition, NO mediates the biological effects of small molecules, such as reactive oxygen species and iron [11, 12]. The NO-donor sodium nitroprusside (SNP) induces ferritin at both the mRNA and protein levels [11]. After iron treatment, NO quickly accumulates in plastids, and the expression of AtFer1 increases [12]. NO also interacts with signal pathways involving hormones and Ca2+ [1315].

NO play a role in plant responses to various abiotic stresses, such as drought, heat, and salt [8, 16]. The application of two NO donors, SNP and diethylenetriamine NO adduct (DETA/NO), leads to a significant reversal in the dry weight and water loss in wheat seedling roots [8]. Treatment of Arabidopsis seedlings with heat results in a reduction in endogenous NO levels, and NO functions in signaling and acts upstream of AtCaM3 in thermo tolerance [16]. In addition, NO interacts with endogenous signals, such as abscisic acid (ABA), ethylene (ET), and auxin (AUX), to regulate diverse biological processes [13, 15].

Efforts have been made to identify NO-responsive genes in plants using large-scale transcriptomic methods. In Arabidopsis thaliana, four transcriptomic analyses identified large amounts of differentially expressed genes (DEGs) and common transcription factor-binding sites (TFBSs) in the promoter regions of DEGs. These DEGs are related to signal transduction, hormone biosynthesis and signaling, transcription factors (TFs), disease resistance, stress responses, photosynthesis, cellular transport, and basic metabolism [1720]. NO-responsive genes have also been identified in other plants, such as Nicotiana tabacum [21], Medicago truncatula [22], Litchi chinensis [23], Betula platyphylla [24], and Larix olgensis [25]. These analyses provide useful information about how NO influences downstream genes and networks under cell death, pathogenic and symbiotic conditions, in addition to triterpene biosynthesis. Upland cotton (Gossypium hirsutum) is an economically import species, as it produces fiber that keeps humans warm. A proteomics analysis was previously performed to identify NO-responsive proteins in cotton [26], but the NO-responsive genes in this species remain unknown. In the present study, we used a high-throughput tag-sequencing method to identify DEGs responding to NO after treatment of cotton with SNP. We obtained sets of up-regulated and down-regulated genes. Several NO-responsive genes were identified and are discussed in association with the diverse functions of NO in plants.

Materials and methods

Plant growth and treatment with NO donor

Seeds of G. hirsutum ecotype CCRI10 were sown in nutritional soil in a greenhouse. The daily photoperiod was set to 16 h/8 h (light/dark), and the corresponding temperatures were 30/22°C (light/dark). As NO is gaseous, we selected a widely accepted NO donor, SNP, for NO treatment [18, 23, 25]. Twenty days after germination, plants with robust and consistent growth were irrigated with 100 μM SNP or 250 μM SNP under continuous light. The control plants were treated with double-distilled water. Cotyledons were harvested at 0, 1, 3, 6, 9, 12, and 24 h after treatment, immediately frozen in liquid nitrogen, and stored at –80°C. Each treatment consisted of three biological replicates, with at least 20 plants pooled per group. To measure NO content, 1 g sample was ground in 7 mL of phosphatic buffer solution. The homogenate was then centrifuged at 11,000 g for 10 min, and the supernatant solution was collected. NO content was measured and calculated using the Griess Reagent System kit (Promega Corporation, Madison, WI, USA) with three technical replicates according to the manufacturer’s instructions. Briefly, 50 μL of the supernatant solution and 50 μL of the sulfanilamide solution were mixed, followed by incubation for 10 min at room temperature, protected from light. Then 50 μl of the NED Solution was added, and incubated 10 minutes at room temperature, protected from light. Finally, the absorbance was measured at 540 nm.

Library construction and sequencing

The materials at 3 h after treatment were used for library construction, with three biological replicates for each treatment, coming to a total of nine libraries. Total RNA was extracted using the EASYspin Plant Total RNA Extraction Reagent (Aidlab, Beijing, China) in accordance with the manufacturer’s instructions and was then stored at −80°C. RNA samples with an A260/A280 ≥ 1.8, A260/A230 ≥ 1.8, and RNA integrity number ≥ 6.5 were considered acceptable for library construction. Libraries were constructed using the following steps: 1) Poly (A)-containing mRNA was purified from 6 μg of total RNA using oligo (dT) magnetic adsorption beads. 2) cDNA was synthesized from the bead-bound RNA using oligo (dT) primers. 3) The restriction enzyme NlaIII was used to digest cDNA, leaving a 4-base overhang (CATG) on the digested cDNA. 4) Unnecessary fragments were washed away, and the Illumina adaptor 1 was ligated to the 5' end of the digested cDNA, which was bead-bounded. 5) The cDNA was digested with MmeI, which recognizes the junction between the Illumina adaptor 1 and the NlaIII site but cuts 17 bp downstream of CATG. 6) Illumina adaptor 2 was ligated to the 3' end of the sequence tag, generating tags with different adaptors at each end. 7) After linear PCR amplification, fragments were purified via polyacrylamide gel electrophoresis. The quantification and qualification of the libraries were detected using an Agilent 2100 Bioanalyzer and an ABI StepOnePlus Real-Time PCR system (Applied Biosystems, USA). The libraries were sequenced using the Illumina HiSeq™ 2000 system. The sequence data of this study have been deposited under SRA accession PRJNA391847.

Data processing

To obtain clean tags, quality control was performed by filtering out 39-bp adaptor sequences, empty reads, low-quality tags containing unknown sequences (N), tags not equal to 21 bp, and tags with only one copy number. The clean tags were mapped to the predicted coding sequences from the genomic data of G. hirsutum [27]. Unambiguous clean tags were obtained by filtering out clean tags that mapped to multiple sequences. The number of unambiguous clean tags for each gene was counted and normalized to TPM [28]. All detected genes were annotated based on the non-redundant database of NCBI (ftp://ftp.ncbi.nih.gov/blast/db/FASTA/nr.gz) and the Gene Ontology database (ftp.geneontology.org/pub/go).

To identify DEGs, the following thresholds were used: log2 (ratio of two TPM values) ≥ 1, the false discovery rate (FDR) ≤ 0.001, and genes with a tag number > 10 (at least in one library) [29]. Venn diagrams for comparing lists were constructed using an online tool (http://bioinfogp.cnb.csic.es/tools/venny/index.html). Gene Ontology (GO) annotation was performed using the Blast2GO program [30]. GO enrichment was analyzed at the agriGO website (http://bioinfo.cau.edu.cn/agriGO/index.php) via Customized annotation. KEGG Orthology-Based Annotation System (KOBAS 2.0) was then employed to identify the enriched metabolic pathways [31]. To identify plant hormone-related genes, a local BLASTP search was performed with the Arabidopsis Hormone Database (AHD2.0, http://ahd.cbi.pku.edu.cn/). To identify TFs, a local BLASTP search was performed with the Database of Arabidopsis Transcription Factors (DATF, http://datf.cbi.pku.edu.cn/) [32], and the matched sequences were verified in the Conserved Domain Database of NCBI (CDD, https://www.ncbi.nlm.nih.gov/Structure/cdd/cdd.shtml). The BLAST threshold was set as 1e-5.

Quantitative real-time polymerase chain reaction (qRT-PCR) analysis

qRT-PCR analysis was performed to verify DGE results. RNA samples used for qRT-PCR were identical to those employed for the DGE experiments. Twenty-five DEGs were selected, and gene-specific primers were designed using Primer Premier 5.0 (S1 Table). Actin was used as the internal reference gene. qRT-PCR was carried out using SYBR Green PCR master mix (Roche Applied Science, Germany) on an ABI 7500 real-time PCR system (Applied Biosystems, USA) with three technical replicates. Amplification reactions were initiated with a denaturing step (95°C for 10 min), followed by 40 cycles of denaturing (95°C for 10 s), annealing (60°C for 30 s) and extension (72°C for 40 s).

Determination of hydrogen peroxide (H2O2) and malondialdehyde (MDA) contents

Each one g sample was ground in 6 mL of phosphatic buffer solution. The homogenate was then centrifuged at 6,000 g for 25 min, and the supernatant solution was collected. The H2O2 content was determined via the method of Chen [33]. Then, 600 μL of the supernatant solution was mixed with 2 mL of 0.1% (m/v) titanium chloride in 20% (v/v) H2SO4. The mixture was next centrifuged at 6,000 g for 15 min, and the supernatant solution was collected. The absorbance was subsequently measured at 410 nm. The MDA content was determined using thiobarbituric acid (TBA) reaction [34]. For this assay, 2.5 mL of the supernatant solution was mixed with 2.5 mL of 0.5% TBA solution. The mixture was then boiled for 15 min and immediately chilled on ice. Then the mixture was centrifuged at 6,000 g for 15 min and the supernatant solution was collected. The absorbance was measured at 532 nm and 600 nm. The MDA content was calculated from the difference in absorbance at 532 and 600 nm using an extinction coefficient of 155 mM-1 cm-1.

Results and discussion

Physiological evaluation of NO donor treatment

To ensure the best effect of the treatment, endogenous NO contents were measured according to the time-course after plants were treated with double-distilled water (control), 100 μM SNP, or 250 μM SNP (Fig 1A). Generally, NO content under the 250 μM SNP treatment was higher than that under 100 μM SNP treatment at different times. Compared with control, NO content of both the 100 μM SNP and 250 μM SNP groups first increased, reached a peak at 3 h, then decreased, and remained steady after 9 h. Previous studies have showed that NR is the primary gene involved in NO biosynthesis, and can also be induced by SNP treatment. In birch cells, NR is up-regulated by SNP treatment [24]. In Larix olgensis Henry, two genes annotated to NR are also reported to be up-regulated by SNP treatment [25]. These results suggest that exogenous NO can influence endogenous NO synthesis and signals. Therefore, NR could be used as a marker gene to evaluate the effect of SNP treatment. Interestingly, a candidate homologous to NR, CotAD_40664, was identified in cotton. So we used it to evaluate the effect of SNP treatment at the transcriptional level (Fig 1B). The result showed that NR was induced by both 100 μM SNP and 250 μM SNP, especially at 3–6 h after treatment. This result was in accordance with previous studies [24, 25], suggesting that SNP treatment could influence endogenous NO synthesis. Considering the NO content and NR expression level observed after SNP treatment, we performed the next experiment on the plant materials at 3 h after treatment.

thumbnail
Fig 1. Evaluation of NO donor treatment.

A, NO contents after treatment with SNP. B, Expression pattern of the NO biosynthesis gene NR after treatment with SNP.

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

High-throughput tag-sequencing and mapping of DGE libraries

To gain knowledge about the NO-responsive genes, high-throughput tag-sequencing was performed on RNAs from cotton plants treated with ddH2O (control), 100 μM SNP (100-SNP), or 250 μM SNP (250-SNP). As a result, 13.4–16.6 million total tags were generated with 0.51–0.58 million distinct tags for each library (S2 Table). After quality control, 13.1–16.3 million total clean tags and 266,789–294,679 distinct clean tags were produced for each library (S2 Table). The distribution of total and distinct clean tag copy numbers showed similar tendencies among the nine libraries (S3 Table). Among the distinct clean tags, ~7% exhibited a copy number more than 100, ~5% exhibited a copy number of 51–100, ~9.5% exhibited a copy number of 21–50, ~10.5% exhibited a copy number of 11–20, ~15% exhibited a copy number of 6–10, and ~53% exhibited a copy number of 2–5. These results are consistent with the rule that most genes are expressed at low levels, whereas minorities of genes are highly expressed under certain conditions. To reveal the molecular events underlying the DGE profiles, the tag sequences of the nine libraries were mapped to the predicted coding sequence based on the genomic data of G. hirsutum [27]. Nearly 25% of the distinct clean tags (66,342–72,910) could be unambiguously mapped to a gene in the reference gene database (S2 Table). As a result, 41,071–43,218 tag-mapped genes were detected (69.51%–73.14% of the total reference genes), among which 30,439–32,564 (51.52%–55.11% of the total reference genes) genes belonged to unambiguous tag-mapped genes in each library (S2 Table).

We further investigated whether the differences in gene coverage between libraries were due to differences in sample size in terms of unambiguous mapped tags (ranging from 77,541 to 86,903 tags), or to real differences in transcriptome complexity as a function of the number of unique mapped tags (S1 Fig). The analysis clearly showed that, at the depth of sequencing achieved, a plateau was reached for each sample, indicating that the sequencing depth was sufficient for the transcriptome coverage. In addition, Pearson’s rank correlation analysis was performed to evaluate the reproducibility among biological replicates. As shown in S4 Table, the correlations between samples among the biological replicates were good, with values ranging from 0.960 to 0.994. However, the correlations between samples out of the biological replicates were varied, ranging from 0.566 to 0.907. These results illustrate that our data were with high reliability and reproducibility for the biological replicates and the sequencing method.

Global analysis of the cotton transcriptome in response to NO

Using DESeq package in R programming language, a total of 1,932 DEGs were detected (S5 Table). Among these DEGs, 1,491 were up-regulated, while 433 were down-regulated, and eight showed a complex expression pattern (all eight DEGs were down-regulated by 100 μM SNP and up-regulated by 250 μM SNP). Compared with control, 1,067 and 1,432 DEGs were detected under the 100 μM SNP and 250 μM SNP treatment, respectively (Fig 2A). Among these DEGs, 707 and 1,266 genes were up-regulated, whereas 360 and 166 genes were down-regulated under the 100 μM SNP and 250 μM SNP treatment, respectively. In addition, qRT-PCR was carried out for 25 DEGs to confirm the DGE data. These DEGs included 10 hormone-related DEGs and 15 randomly selected DEGs. The expression patterns obtained via qRT-PCR were consistent with those obtained through DGE analysis, with an overall Pearson R of 0.8599 (ranging from 0.6559 to 0.9999 for each DEG) (P < 0.0001; Fig 2B; S1 Table).

thumbnail
Fig 2. Analysis of the DEGs in response to NO.

A, Number of DEGs compared with the control. B, Comparison of gene expression for 25 DEGs. Data generated using qPCR (y-axis) and DGE (x-axis) are shown.

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

We used the GO classification system to classify the possible functions of the DEGs. A total of 1,363 (70.5%) DEGs were successfully assigned to at least one GO term associated with biological process, cellular component, or molecular function (S5 Table, S2 Fig). In addition, GO enrichment analysis was performed using agriGO. A total of 14 significantly enriched GO items were identified with an FDR≤0.05, including translation, structural molecule activity, structural constituent of ribosome, ribonucleoprotein complex, ribosome, macromolecular complex, cytoplasm, cytoplasmic part, intracellular part, intracellular non-membrane-bounded organelle, non-membrane-bounded organelle, intracellular, intracellular organelle, and organelle (S6 Table). Among biological process, “metabolic process” (31.47%) and “cellular process” (26.97%) were the two items that exhibited the highest percentage of DEGs. Among cellular component, “cell” (20.29%) and “cell part” (20.29%) were the two items that exhibited the highest percentage of DEGs. Among molecular function, “binding” (38.20%) and “catalytic activity” (29.35%) were the two items that exhibited the highest percentage of DEGs (Fig 3). These results are consistent with previous studies suggesting that binding and catalytic activities operate predominantly at transcriptome level in response to NO [20, 35]. Meanwhile, “transcription regulator activity” also showed a notable percentage of DEGs, demonstrating that TFs might be involved in NO response.

thumbnail
Fig 3. GO enrichment of DEGs in response to NO.

The analysis was performed using agriGO.

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

Genes with similar expression patterns are often functionally correlated. To obtain an overview of the expression patterns of the DEGs under different SNP concentrations, we employed MultiExperiment Viewer to cluster the DEGs based on the K-means clustering method. As shown in S3 Fig, these DEGs were clustered into six main clusters. Cluster 1, cluster 2, and cluster 3, which contained 435, 602, and 334 genes, respectively, were up-regulated by both the 100 μM SNP and 250 μM SNP treatment. Cluster 4, which contained 134 genes, was significantly up-regulated by 100 μM SNP, and slightly up-regulated by 250 μM SNP treatment. Cluster 5, which contained 128 genes, was down-regulated by 100 μM SNP, but was not obviously altered by 250 μM SNP treatment. Cluster 6, which contained 299 genes, was down-regulated by both the 100 μM SNP and 250 μM SNP treatment.

In addition, we used KOBAS 2.0 to identify the metabolic pathways in which the NO-responsive genes function. The results are presented in S7 Table and S8 Table. Under the 100 μM SNP treatment, 307 genes were annotated and mapped to 199 KEGG pathways. Among these, the stilbenoid, diarylheptanoid and gingerol biosynthesis, ribosome, starch and sucrose metabolism, and phenylpropanoid biosynthesis pathways were the most enriched pathways (S7 Table). Under the 250 μM SNP treatment, 472 genes were annotated and mapped to 226 KEGG pathways. Among these, plant hormone signal transduction, carbon fixation in photosynthetic organisms, plant-pathogen interaction, carbon metabolism, and ribosome pathways were the most enriched pathways (S8 Table). These results suggested that these pathways might play an important role in the response to NO. Some pathways, such as plant hormone signal transduction and plant-pathogen interaction, have been reported to be involved in NO-responsive processes [7, 13, 36].

TF family genes are regulated by NO

TFs are key regulators that temporarily and spatially turn on or off the transcription of their target genes through binding certain upstream elements [32]. Previous studies have revealed sets of TF families genes that are regulated by NO [20]. We therefore performed a local BLASTP, using the TF sequences in DATF as database and the translated protein sequences of DEGs as queries. The matched sequences were further confirmed by a Batch CD-search in CDD. As a result, a total of 157 TFs belonging to 36 TF families were identified, among which 85 TFs were differentially expressed under the 100 μM SNP treatment, and 118 TFs were differentially expressed under the 250 μM SNP treatment (S9 Table). The most enriched TF families were bHLH (15), DBP (12), MYB (12), C3H (11), AP2-EREBP (9), NAC (8), WRKY (7), PHD (6), C2C2-CO-like (5), GRAS (5), HSF (5), Orphans (5), and TRAF (5) (Table 1). Furthermore, under the 100 μM SNP treatment, the most enriched TF families were bHLH (7), MYB (6), DBP (4), and C3H (4) among the translated protein sequences of up-regulated DEGs, whereas those the translated protein sequences of down-regulated DEGs were bHLH (3), MYB (3), C3H (2), and C2C2-CO-like (2). Under the 250 μM SNP treatment, the most enriched TF families were bHLH (9), NAC (6), WRKY (6), and MYB (5) among the translated proteins sequences of up-regulated DEGs, and bHLH (2) and MYB (2) among the translated protein sequences of down-regulated DEGs. Orthologs of some identified TF family genes have been reported to be associated with diverse NO-mediated processes of plant growth and development. For example, three of the identified TF-encoding genes (NAC2, NAP, and WRKY6) were shown to be down-regulated during dark-induced leaf senescence in comparison of nos1/noa1 (an NO-deficient mutant) with wild-type Arabidopsis [37]. Meanwhile, three other TF-encoding genes (AtbHLH38, AtbHLH39, and AtMYB72) are strongly induced by Fe deficiency and GSNO treatment (an NO donor) in A. thaliana [38]. Orthologs of these TF-encoding genes showed marked fold changes after NO treatment, suggesting that NO-responsive genes and pathways might be regulated by these vital TFs.

Differential expression of hormone-related genes in response to NO

Increasing evidences suggest that NO is an essential component of a wide variety of hormone-regulated processes, including stomatal movement, root growth, senescence, and stress responses [13]. Therefore, we analyzed the NO-responsive DEGs that related to plant hormones. In total, 72 of the DEGs were functionally related to the eight plant hormones, 37 of which were differentially expressed under the 100 μM SNP treatment, while 52 were differentially expressed under the 250 μM SNP treatment (S10 Table). Among these genes, 17, 17, 15, 8, 6, 5, 3, and 1 DEGs were related to the biosynthesis and signal transduction networks of ET, abscisic acid (ABA), auxin (AUX), salicylic acid (SA), brassinosteroid (BR), JA, gibberellin (GA), and cytokinin (CK), respectively (Table 2). Under the treatment of 100 μM SNP, the most abundant hormones were ET (7), ABA (5), AUX (4), and BR (4) among the up-regulated DEGs, whereas they were ABA (5), AUX (2), and BR (2) among the down-regulated DEGs. Under the treatment of 250 μM SNP, the most abundant hormones were AUX (12), ABA (11), and ET (4) in the up-regulated DEGs, whereas they were ABA (1) and AUX (1) among the down-regulated DEGs. These results highlight the role of NO in regulating hormone-related signal transduction pathways.

thumbnail
Table 2. Summary of the hormone related genes in response to NO.

https://doi.org/10.1371/journal.pone.0192367.t002

It is reported that NO induces the expression of genes involved in ET biosynthesis, thus regulating the ET signal pathway [36, 39]. In A. thaliana, genes of ET-insensitive protein (EIN) are pivotal in ET signal transduction. We found that an ortholog of EIN4, CotAD_12236, was induced by both the 100 μM SNP and 250 μM SNP treatments (Fig 4). NO is required in ABA signaling to induce stomatal closure and inhibite light-induced stomatal opening [40]. We identified two orthologs of the ABA biosynthesis gene ABA2 among our results: CotAD_35694 and CotAD_45658. CotAD_35694 was down-regulated by 100 μM SNP, but up-regulated by 250 μM SNP, whereas CotAD_45658 was up-regulated by both 100 μM SNP and 250 μM SNP. It has also been demonstrated that NO is involved in mediating iron deficiency-induced root morphology changes together with AUX in A. thaliana [15]. We also identified four NO-responsive DEGs associated with the AUX signal pathway. Three of them were orthologous to auxin response factors (CotAD_44201, CotAD_06531, CotAD_08111), while one was orthologous to indole-3-acetic acid inducible 29 (CotAD_00131). All four of these DEGs were up-regulated by 100 μM SNP or 250 μM SNP. NO functions together with JA and SA in modulating plant defense responses, including the responses to pathogens and abiotic stresses (e.g., low temperature, drought, and high salinity) [41]. One DEG related to JA biosynthesis (CotAD_29356, orthologous to SSI2) and one DEG related to JA signal transduction (CotAD_42959 orthologous to WRKY70), were identified. Both were up-regulated by 100 μM SNP and 250 μM SNP. Although we did not find obvious evidence of an interaction of NO with GA, a key DEG in GA biosynthesis, CotAD_23661, which is orthologous to GA20OX2, was also identified, with an up-regulation by 250 μM SNP. The expression of these hormone-related DEGs was verified through qRT-PCR (Fig 4), where results for nine of the ten DEGs were consistent with the DGE data, with a Pearson correlation value above 0.8. All the above findings clearly indicate that NO regulates hormone-related signal transduction pathways to influence plant biological processes.

thumbnail
Fig 4. Expression levels of candidate DEGs related to hormone biosynthesis and signaling pathways.

https://doi.org/10.1371/journal.pone.0192367.g004

Stress-responsive genes are regulated by NO

Based on the identification of TFs regulated by NO, we found several TF families were involved in biotic and abiotic stress responses, such as NAC [42], WRKY [43], MYB [44], and HSF [45]. At the same time, ABA, SA, JA and ET are major hormones known to play major roles in mediating the plant defense response against pathogens and abiotic stresses [46]. The identification of hormone-related DEGs provides additional evidence that NO might play a role in the regulation of plant stresses. We therefore determined the contents of H2O2 and MDA, two widely employed physiological markers that are induced by various stresses [47], after treatment with SNP. Consistently, the H2O2 content showed significant increases after 3 h treatment with 100 μM SNP and 250 μM SNP compared with the control (P<0.05) (Fig 5A). The MDA content was also increased after 3 h 100 μM SNP treatment, and showed a significant increase after 3 h 250 μM SNP treatment (Fig 5B). At the mRNA level, a gene encoding catalase (producing H2O2) (CotAD_32224) was found to be up-regulated by 250 μM SNP (Fig 5C). In addition, nine genes encoding antioxidant defense proteins showed altered regulation after treatment with SNP (Fig 5D–5L). These proteins included peroxidase (CotAD_20183, CotAD_29182, CotAD_62732), glutathione S-transferase (CotAD_48300, CotAD_21477, CotAD_30505), ascorbate peroxidase (CotAD_32863), and superoxide dismutase (CotAD_05767, CotAD_73921). With the exception of CotAD_20183 and CotAD_32863, the remaining seven genes showed up-regulation under the 100 μM SNP or 250 μM SNP treatment. In a previous study, hundreds of DEGs, including CAT and APX, were found to be involved in plant defense against various biotic and abiotic stresses after treatment of Arabidopsis leaves with S-nitrosocysteine (CysNO, another NO donor) [20]. In tobacco, significant cross-talk between NO and H2O2 was found, and target genes of the combined action of NO and H2O2 during the induction of plant cell death were identified [21]. In birch, the H2O2 content was also increased after SNP treatment, with genes encoding FeSOD, MnSOD, and CZSOD up-regulated after SNP treatment compared with the control [24]. Our results were consistent with these studies, suggesting the involvement of NO in plant stress responses and defense.

thumbnail
Fig 5. NO is involved in the stress response.

A, Changes in hydrogen peroxide content after treatment with SNP. B, Changes in malondialdehyde content after treatment with SNP. C-L, Expression level of DEGs related to stress responses and defense process.

https://doi.org/10.1371/journal.pone.0192367.g005

Comparison of the NO-responsive genes of cotton with other plants

The identification of NO-responsive genes at a large-scale transcriptome level has been performed in several plant species thus far, including A. thaliana [1720, 47], N. tabacum [21], M. truncatula [22], L. chinensis [23], B. platyphylla Suk.[24], and L. olgensis Henry [25]. We compared the experimental design and results of these 10 studies with those of our study, including examination of the treatment method and time, species, sampled tissues, sequencing method, number of total genes, and number of DEGs. As shown in Table 3, the treatment methods were divergent between these studies. In eight studies, the authors chose SNP as the NO donor to perform the experiment. Among other NO donors chosen for the experiments, GSNO, CysNO, and NOR3 were used twice, once, and once, respectively. Besides, gaseous NO was used to perform the experiment in one case. All the NO donors as well as gaseous NO are acceptable, however, their responsive genes exhibit differences. For example, in the Ferrarini et al. study, among genes responsive to either GSNO or SNP, only approximately 11% were common to both groups in the leaves and 1.6% in the roots. Thus, it appears that the modulation of gene expression is dependent on the form of NO produced (NO+ released by SNP or NOreleased by GSNO) [22, 48]. In addition, even when the same NO donor was used, the treatment times, species, and sampled tissues differed among these studies. This may due to the specific purposes of the studies, and the most effective treatment times in different species. The number of total genes and DEGs identified is increasing with the development of new sequencing technologies. For example, using the cDNA–AFLP, hundreds of DEGs have been identified [17, 21]; whereas using the RNA-Seq, thousands of DEGs have been identified [24, 25].

thumbnail
Table 3. Overview of the transcriptome analysis of NO responsive genes in plants.

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

As the treatment methods and sampled tissues used greatly influence gene expression results, we made a comparison of the 100 μM SNP-responsive genes in our study to the 100 μM SNP-responsive genes reported in A. thaliana [18]. Both experiments were irrigated with 100 μM SNP, and the overground tissues were harvested 3 h after treatment. Among the up-regulated genes, we identified 15 genes (out of 707 DEGs) from our data that aligned with 13 genes (out of 123 DEGs) from the A. thaliana data (Table 4). These aligned genes included some well-defined genes induced by NO. For example, the gene encoding ferritin (CotAD_54212, orthologous to AT5G01600) was induced, confirming previous reports that NO regulates iron homeostasis through Arabidopsis ferritin 1 [11, 12, 14]. An ethylene-responsive TF, ERF13 (CotAD_74946, orthologous to AT2G44840), belonging to the AP2/EREBP TF family, was also identified. TFs in this family participate in the regulation of disease resistance [49], suggesting that NO may regulate disease resistance pathways through interaction with ERF13. In addition, two genes related to calcium signaling were identified (CotAD_56107, orthologous to AT1G76650, and CotAD_28267, orthologous to AT2G41010). One report indicated that NO functions in signaling and acts upstream of AtCaM3 in heat tolerance in A. thaliana [16], providing evidence supporting our findings that genes in the calcium signaling pathway can be induced by NO. Among the down-regulated genes, we identified only two genes (out of 360 DEGs) from our data that aligned with two genes (out of 50 DEGs) from the Arabidopsis data (Table 4). However, there is little genetic or physiological evidence that these two genes could be induced by NO.

thumbnail
Table 4. List of genes differential expressed both in cotton and in Arabidopsis in response to 100 μM SNP.

https://doi.org/10.1371/journal.pone.0192367.t004

Conclusion

In conclusion, we described transcriptional changes in response to NO in upland cotton, using high-throughput tag-sequencing. A total of 1,932 DEGs were identified, among which 1,067 and 1,432 were regulated by 100 μM SNP and 250 μM SNP, respectively. qRT-PCR on 25 DEGs was performed to validate the data and consistent results were obtained. GO classification and KEGG pathway were analyzed to obtain additional knowledge about the transcriptome data. Furthermore, we identified 157 DEGs belonging to 36 TF families and 72 DEGs related to eight hormones. Several TF families (NAC, WRKY, MYB, and HSF) and hormones (ABA, SA, JA, and ET) are involved in stress responses. In addition, H2O2 and MDA contents were shown to increase, and related genes showed marked up-regulation after treatment with SNP. These results suggest that NO might play a role in plant stress response. At last, comparison of NO-responsive genes recorded in the present study with the results of previous studies led to the identification of dozens of vital DEGs that shared between A. thaliana and cotton. To our knowledge, this is the first study that performs deep-sequencing analysis of the response to NO in cotton. These data provide a foundation for future studies concerning NO and stress response in cotton. In addition, some DEGs identified in this study represent good candidates for further functional analysis in cotton.

Supporting information

S1 Fig.

Saturation analyses of the DGE libraries of control (A-C), 100 μM SNP (D-F), and 250 μM SNP (G-I) treatment.

https://doi.org/10.1371/journal.pone.0192367.s001

(TIF)

S2 Fig.

Biological processes (A), cellular component (B), and molecular function (C) of NO-responsive genes involved according to GO annotation.

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

(TIF)

S3 Fig. K-means clusters of the NO-responsive genes.

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

(TIF)

S1 Table. Validation of DGE data through qRT-PCR analyses.

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

(XLSX)

S2 Table. Tag mapping and normalization of tag counts.

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

(XLSX)

S3 Table. Distribution of total and distinct clean tags in the DGE libraries.

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

(XLSX)

S4 Table. Pearson’s rank correlation of nine libraries.

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

(XLSX)

S5 Table. TPM value and annotation of DEGs in response to NO.

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

(XLSX)

S6 Table. GO enrichment analysis using agriGO.

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

(XLSX)

S7 Table. KEGG pathways enriched under the 100 μM SNP treatment.

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

(XLSX)

S8 Table. KEGG pathways enriched under the 250 μM SNP treatment.

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

(XLSX)

S9 Table. List of transcription factors regulated by NO.

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

(XLSX)

S10 Table. Annotation of hormone-related DEGs.

https://doi.org/10.1371/journal.pone.0192367.s013

(XLSX)

Acknowledgments

This research was financially supported by the National Natural Science Foundation of China (grant 31301366).

References

  1. 1. Besson-Bard A, Astier J, Rasul S, Wawer I, Dubreuil-Maurizi C, Jeandroz S, et al. Current view of nitric oxide-responsive genes in plants. Plant Science. 2009;177(4):302–309.
  2. 2. Delledonne M, Xia Y, Dixon RA, Lamb C. Nitric oxide functions as a signal in plant disease resistance. Nature. 1998;394(6693):585–588. pmid:9707120
  3. 3. Desikan R, Griffiths R, Hancock J, Neill S. A new role for an old enzyme: nitrate reductase-mediated nitric oxide generation is required for abscisic acid-induced stomatal closure in Arabidopsis thaliana. Proceedings of the National Academy of Sciences of the United States of America. 2002;99(25):16314–16318. pmid:12446847
  4. 4. Zhang H, Zhao X, Yang J, Yin H, Wang W, Lu H, et al. Nitric oxide production and its functional link with OIPK in tobacco defense response elicited by chitooligosaccharide. Plant cell reports. 2011;30(6):1153–1162. pmid:21336582
  5. 5. Guo FQ, Okamoto M, Crawford NM. Identification of a plant nitric oxide synthase gene involved in hormonal signaling. Science. 2003;302(5642):100–103. pmid:14526079
  6. 6. Van Ree K, Gehl B, Chehab EW, Tsai YC, Braam J. Nitric oxide accumulation in Arabidopsis is independent of NOA1 in the presence of sucrose. Plant J. 2011;68(2):225–233. pmid:21689173
  7. 7. Zeier J, Delledonne M, Mishina T, Severi E, Sonoda M, Lamb C. Genetic elucidation of nitric oxide signaling in incompatible plant-pathogen interactions. Plant Physiol. 2004;136(1):2875–2886. pmid:15347797
  8. 8. Xie Y, Ling T, Han Y, Liu K, Zheng Q, Huang L, et al. Carbon monoxide enhances salt tolerance by nitric oxide-mediated maintenance of ion homeostasis and up-regulation of antioxidant defence in wheat seedling roots. Plant, cell & environment. 2008;31(12):1864–1881. pmid:18811735
  9. 9. Zhu S-h, Zhou J. Effect of nitric oxide on ethylene production in strawberry fruit during storage. Food Chemistry. 2007;100(4):1517–1522.
  10. 10. Manjunatha G, Gupta KJ, Lokesh V, Mur LA, Neelwarne B. Nitric oxide counters ethylene effects on ripening fruits. Plant signaling & behavior. 2012;7(4):476–483. pmid:22499176
  11. 11. Murgia I, Delledonne M, Soave C. Nitric oxide mediates iron-induced ferritin accumulation in Arabidopsis. Plant J. 2002;30(5):521–528. pmid:12047627
  12. 12. Arnaud N, Murgia I, Boucherez J, Briat JF, Cellier F, Gaymard F. An iron-induced nitric oxide burst precedes ubiquitin-dependent protein degradation for Arabidopsis AtFer1 ferritin gene expression. The Journal of biological chemistry. 2006;281(33):23579–23588. pmid:16782706
  13. 13. Simontacchi M, Garcia-Mata C, Bartoli CG, Santa-Maria GE, Lamattina L. Nitric oxide as a key component in hormone-regulated processes. Plant cell reports. 2013;32(6):853–866. pmid:23584547
  14. 14. Besson-Bard A, Gravot A, Richaud P, Auroy P, Duc C, Gaymard F, et al. Nitric oxide contributes to cadmium toxicity in Arabidopsis by promoting cadmium accumulation in roots and by up-regulating genes related to iron uptake. Plant Physiol. 2009;149(3):1302–1315. pmid:19168643
  15. 15. Chen WW, Yang JL, Qin C, Jin CW, Mo JH, Ye T, et al. Nitric Oxide Acts Downstream of Auxin to Trigger Root Ferric-Chelate Reductase Activity in Response to Iron Deficiency in Arabidopsis. Plant Physiol. 2010;154(2):810–819. pmid:20699398
  16. 16. Xuan Y, Zhou S, Wang L, Cheng Y, Zhao L. Nitric oxide functions as a signal and acts upstream of AtCaM3 in thermotolerance in Arabidopsis seedlings. Plant Physiol. 2010;153(4):1895–1906. pmid:20576787
  17. 17. Polverari A, Molesini B, Pezzotti M, Buonaurio R, Marte M, Delledonne M. Nitric oxide-mediated transcriptional changes in Arabidopsis thaliana. Mol Plant Microbe In. 2003;16(12):1094–1105. pmid:14651343
  18. 18. Parani M, Rudrabhatla S, Myers R, Weirich H, Smith B, Leaman DW, et al. Microarray analysis of nitric oxide responsive transcripts in Arabidopsis. Plant Biotechnol J. 2004;2(4):359–366. pmid:17134397
  19. 19. Palmieri MC. Nitric oxide-responsive genes and promoters in Arabidopsis thaliana: a bioinformatics approach. J Exp Bot. 2007;59(2):177–186. pmid:18272923
  20. 20. Hussain A, Mun BG, Imran QM, Lee SU, Adamu TA, Shahid M, et al. Nitric Oxide Mediated Transcriptome Profiling Reveals Activation of Multiple Regulatory Pathways in Arabidopsis thaliana. Frontiers in plant science. 2016;7:975. pmid:27446194
  21. 21. Zago E, Morsa S, Dat JF, Alard P, Ferrarini A, Inze D, et al. Nitric oxide- and hydrogen peroxide-responsive gene regulation during cell death induction in tobacco. Plant Physiol. 2006;141(2):404–411. pmid:16603664
  22. 22. Ferrarini A, De SM, Baudouin E, Pucciariello C, Polverari A, Puppo A, et al. Expression of Medicago truncatula genes responsive to nitric oxide in pathogenic and symbiotic conditions. Molecular plant-microbe interactions: MPMI. 2008;21(6):781–790. pmid:18624641
  23. 23. Liu WW, Chen HB, Lu XY, Rahman MJ, Zhong S, Zhou BY. Identification of nitric oxide responsive genes in the floral buds of Litchi chinensis. Biologia Plantarum. 2014;59(1):115–122.
  24. 24. Zeng F, Sun F, Li L, Liu K, Zhan Y. Genome-scale transcriptome analysis in response to nitric oxide in birch cells: implications of the triterpene biosynthetic pathway. Plos One. 2014;9(12):e116157-e. pmid:25551661
  25. 25. Hu X, Yang J, Li C. Transcriptomic Response to Nitric Oxide Treatment in Larix olgensis Henry. International journal of molecular sciences. 2015;16(12):28582–28597. pmid:26633380
  26. 26. Meng Y, Liu F, Pang C, Fan S, Song M, Wang D, et al. Label-free quantitative proteomics analysis of cotton leaf response to nitric oxide. Journal of proteome research. 2011;10(12):5416–5432. pmid:22029526
  27. 27. Li F, Fan G, Lu C, Xiao G, Zou C, Kohel RJ, et al. Genome sequence of cultivated Upland cotton (Gossypium hirsutum TM-1) provides insights into genome evolution. Nature Biotechnology. 2015;33(5):524–530. pmid:25893780
  28. 28. Morrissy AS, Morin RD, Delaney A, Zeng T, McDonald H, Jones S, et al. Next-generation tag sequencing for cancer gene expression profiling. Genome research. 2009;19(10):1825–1835. pmid:19541910
  29. 29. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome biology. 2014;15(12):550. pmid:25516281
  30. 30. Conesa A, Gotz S. Blast2GO: A comprehensive suite for functional analysis in plant genomics. International journal of plant genomics. 2008;2008:619832. pmid:18483572
  31. 31. Xie C, Mao X, Huang J, Ding Y, Wu J, Dong S, et al. KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic Acids Res. 2011;39(Web Server issue):W316–322. pmid:21715386
  32. 32. Pérezrodríguez P, Corrêa LGG, Rensing SA, Kersten B, Muellerroeber B. PlnTFDB: updated content and new features of the plant transcription factor database. Nucleic Acids Research. 2010;38(Database issue): D822–D827. pmid:19858103
  33. 33. Chen J, Zhu C, Li LP, Sun ZY, Pan XB. Effects of exogenous salicylic acid on growth and H2O2-metabolizing enzymes in rice seedlings under lead stress. Journal of environmental sciences (China). 2007;19(1):44–49.
  34. 34. Heath RL, Packer L. Photoperoxidation in isolated chloroplasts. I. Kinetics and stoichiometry of fatty acid peroxidation. Archives of biochemistry and biophysics. 1968;125(1):189–198. pmid:5655425
  35. 35. Begaramorales JC, Sánchezcalvo B, Luque F, Leyvapérez MO, Leterrier M, Corpas FJ, et al. Differential transcriptomic analysis by RNA-Seq of GSNO-responsive genes between Arabidopsis roots and leaves. Plant & Cell Physiology. 2014;55(6):1080–1095. pmid:24599390
  36. 36. Xie Z, Li D, Wang L, Sack FD, Grotewold E. Role of the stomatal development regulators FLP/MYB88 in abiotic stress responses. The Plant Journal. 2010;64(5):731–739. pmid:21105921
  37. 37. Niu YH, Guo FQ. Nitric oxide regulates dark-induced leaf senescence through EIN2 in Arabidopsis. Journal of integrative plant biology. 2012;54(8):516–525. pmid:22765302
  38. 38. Garcia MJ, Suarez V, Romera FJ, Alcantara E, Perez-Vicente R. A new model involving ethylene, nitric oxide and Fe to explain the regulation of Fe-acquisition genes in Strategy I plants. Plant physiology and biochemistry: PPB / Societe francaise de physiologie vegetale. 2011;49(5):537–44. pmid:21316254
  39. 39. Ahlfors R, Brosche M, Kollist H, Kangasjarvi J. Nitric oxide modulates ozone-induced cell death, hormone biosynthesis and gene expression in Arabidopsis thaliana. Plant J. 2009;58(1):1–12. pmid:19054359
  40. 40. Garcia-Mata C, Lamattina L. Abscisic acid (ABA) inhibits light-induced stomatal opening through calcium-and nitric oxide-mediated signaling pathways. Nitric Oxide. 2007;17(3):143–151. pmid:17889574
  41. 41. Cheong YH, Chang H- S, Gupta R, Wang X, Zhu T, Luan S. Transcriptional profiling reveals novel interactions between wounding, pathogen, abiotic stress, and hormonal responses in Arabidopsis. Plant Physiol. 2002;129(2):661–677. pmid:12068110
  42. 42. Shao H, Wang H, Tang X. NAC transcription factors in plant multiple abiotic stress responses: progress and prospects. Frontiers in plant science. 2015;6:902. pmid:26579152
  43. 43. Phukan UJ, Jeena GS, Shukla RK. WRKY Transcription Factors: Molecular Regulation and Stress Responses in Plants. Frontiers in plant science. 2016;7:760. pmid:27375634
  44. 44. Roy S. Function of MYB domain transcription factors in abiotic stress and epigenetic control of stress response in plant genome. Plant Signal Behav. 2016;11(1):e1117723. pmid:26636625
  45. 45. Ohama N, Kusakabe K, Mizoi J, Zhao H, Kidokoro S, Koizumi S, et al. The Transcriptional Cascade in the Heat Stress Response of Arabidopsis Is Strictly Regulated at the Level of Transcription Factor Expression. The Plant cell. 2016;28(1):181–201. pmid:26715648
  46. 46. Verma V, Ravindran P, Kumar PP. Plant hormone-mediated regulation of stress responses. BMC plant biology. 2016;16:86. pmid:27079791
  47. 47. Tewari RK, Watanabe D, Watanabe M. Chloroplastic NADPH oxidase-like activity-mediated perpetual hydrogen peroxide generation in the chloroplast induces apoptotic-like death of Brassica napus leaf protoplasts. Planta. 2012;235(1):99–110. pmid:21853253
  48. 48. Stamler JS, Singel DJ, Loscalzo J. Biochemistry of nitric oxide and its redox-activated forms. Science. 1992;258(5090):1898–1902. pmid:1281928
  49. 49. Gutterson N, Reuber TL. Regulation of disease resistance pathways by AP2/ERF transcription factors. Current opinion in plant biology. 2004;7(4):465–471. pmid:15231271