3.4.1 Validation of RNA-seq results by quantitative real-time RT-PCR (qPCR)
-
A cDNA library was constructed from equal amounts of RNA extracted from P. indica leaves subjected to 4℃. The transcriptome sequencing of 18 samples was finally completed, and a total of 128.88 Gb Clean Data was obtained, and the Clean Data of each sample reached 5.73 Gb, with the percentage of Q30 bases at 93.42% and above (Tab. S2). A total of 90,725 Unigenes were obtained after assembly, among which 31,073 Unigenes were above 1kb in length. After obtaining high-quality sequencing data, the sequence was assembled using Trinity. A total of 90,725 Unigenes were obtained from the assembly, and the N50 of Unigene was 2963, which showed high assembly integrity. The Clean Data of each sample was compared with the assembled Transcript or Unigene libraries, and the statistics of the comparison results were shown in Tab. S2. The statistical results indicate that the transcriptome data are sufficient and valid for the following analysis step.
-
Twelve randomly selected DEGs were subjected to qPCR for transcript abundance to verify that the RNA-seq results were accurate. There were 7 Unigene data up-regulated expression and 5 down-regulated (Fig. 5). The qPCR results of the tested genes Unigene 30052, Unigene 69681 and Unigene79700, Unigene 11281, Unigene40972 were generally consistent with the changes in the RNA-seq data. Unigene08110, Unigene64040, Unigene03061, Unigene82860, Unigene 43750, Unigene 44592, and Unigene 51682. These seven genes were weakly up- or down-regulated in expression within 4h, which may have led to a slight bias in the q-PCR test results, but the overall expression changes were consistent. This result was validated the reliability of the RNA-seq data and supports further bioinformatics analysis.
Note
There are 7 up-regulated and 5 down-regulated expressed genes in the graph; the X-axis represents the time (h) of low-temperature treatment, the left Y-axis displays the corresponding RNA-seq expression data (blue histogram), and the right Y-axis displays the relative gene expression level by qPCR (black dash line). Bars represent ± SE (n = 3) and broken lines represent ± SE (n = 9).
3.4.2 GO and KEGG Enrichment Reveals Cold Stress Response is Closely Linked to Transcriptional Regulation and Hormonal Signaling
The BLAST parameter E-value ≤ 1e-5 and the HMMER parameter E-value ≤ 1e-10 were selected to obtain 41,111 Unigenes with annotated information in the total database (Tab. S3). The GO database categorized the annotated 33,538 Unigenes into 45 functional groups (Fig. 5a), which were further enriched into three GO level II functions, including "cellular components," "molecular functions," and "biological processes." The top three "cellular components" were "cellular anatomical entity (19589)", "intracellular (11698)" and "protein-containing complex (3515)"; The top three "molecular functions" were "structural domain binding (16993)", "catalytic activity (15944)" and "transcription regulator activity (1112)"; the top three "biological process" were "cellular process (18945) ", "metabolic process (16477)" and "biological regulation (5812)". The enrichment degree of Pathway was analyzed in KEGG using the Enrichment Factor and the significance of enrichment was calculated using the hypergeometric test. Analysis shows that in the 0h vs. 24h group, the most significant KEGG enrichment pathways are "Plant hormone signal transmission," "Plant pathway interaction," "MAPK signaling pathway - plant," and "Circadian rhythm - plant" (Fig. 5b). By ranking GO and KEGG, it was not difficult to find that P. indica cold stress response was closely related to transcriptional regulation and hormone signaling, which were also coincided with our previous finding of increased ACC content, ACO and ACS enzyme activities in the phytohormone ethylene synthesis pathway.
3.4.3 Differential expression gene annotation of ethylene response factors PiERF1 and PiDREB1A involved in P. indica cold response
Transcription factor prediction analysis by BMKCloud (BMKCloud) was showed that AP2/ERF transcription factors were induced in the most significant number (108) during the cold response of P. indica within 24h, suggested that transcription factors ERFs play an essential role in the molecular mechanism of the cold response of P. indica (Fig. 7a). Analysis of the differentially expressed ERFs genes in the four treatment groups of 0h vs 4h, 0h vs. 8h, 0h vs 12h and 0h vs 24h revealed that the Unigene_83531 gene (named PiERF1) was significantly up-regulated in response to cold stress, with |log2FC| reaching 2.46 at 24h. Combined with the correlation matrix analysis, PiERF1 was found to have a common relationship with several colds. In combination with the correlation matrix analysis, PiERF1 was positively correlated with several physiological indicators related to cold stress, including ACC content, ACO and ACS activity. It was hypothesized that PiERF1 might be regulated by the ethylene precursor ACC (Fig. 7b) and positively regulated by cold stress.
In the ethylene signaling pathway, included ETR1 (Unigene 12465), EIN3 (Unigene 53149) with the ACC (1-aminocyclopropane-1-carboxylic acid) synthase gene, ACS6 (Unigene 81360), as well as the upstream genes of ACS6, EIN3/EIN1, MPK3 (Unigene 78362, Unigene 64590), both responded positively to cold stress (Fig. 7c). In the DREB-COR pathway, PiDREB1A, a fast cold-responsive gene, was the highest expression at 8h, followed by a gradual decrease. The COR gene (Unigene 28684) was regulated downstream by PiDREB1A, with a gradual up-regulation of its expression at 24h (Fig. 7d), which were regulated the genes of the plant group to resist cold stress. The results were indicated that both key pathways, ethylene signaling and DREB-COR, responded to cold stress, but the crosstalk between them needs further exploration.
-
3.5 Evidence of negative regulation of cold hardiness of P. indica by exogenous spraying of ethylene promoter ACC and inhibitor AVG by ethylene signaling
P. indica leaves were sprayed with 50 µmol/L ACC and AVG respectively (control group sprayed with RO water), and it was found that the ACC group was more sensitive to cold relative to the control group, and the foliage crumpling was evident at 72h of 4℃ stress, and the leaves were dry and severely damaged after 72h of recovery incubation (Fig. 8a). Chlorophyll fluorescence parameters were greatly reduced, and the photochemical efficiency was significantly inhibited (Fig. 8b); at the same time, ACC increased EL and membrane permeability; ROS scavenging mechanism was weakened, resulting in an increase in MDA, which caused P. indica to exhibit cold sensitivity (Fig. 8c).In contrast to this, the AVG group was able to increase P. indica cold tolerance, and the degree of damage was lower relative to that of the CONTROL group after 72h of restoration culture (Fig. 8a) and higher photochemical efficiency (Fig. 8b); it could reduce EL permeability, maintain cellular homeostasis, enhance ROS enzyme activity, and inhibit MDA overaccumulation (Fig. 8c). Excess accumulation of ACC in P. indica induced reduced photochemical efficiency, elevated cell membrane permeability, and increased lipid peroxidation, leading to reduced cold tolerance. On the contrary, AVG inhibited ACC accumulation and suppressed ethylene signaling, resulting in a certain protective effect on P. indica; however, further experiments need to demonstrate how ethylene signaling regulates its downstream genes in response to cold.