Comparative Transcriptome Analysis Reveals Molecular Defensive Mechanism of Arachis hypogaea in Response to Salt Stress

Abiotic stresses comprise all nonliving factors, such as soil salinity, drought, extreme temperatures, and metal toxicity, posing a serious threat to agriculture and affecting the plant production around the world. Peanut (Arachis hypogaea L.) is one of the most important crops for vegetable oil, proteins, minerals, and vitamins in the world. Therefore, it is of importance to understand the molecular mechanism of peanut against salt stress. Six transcriptome sequencing libraries including 24-hour salt treatments and control samples were constructed from the young leaves of peanut. A comprehensive analysis between two groups detected 3,425 differentially expressed genes (DEGs) including 2,013 upregulated genes and 1,412 downregulated genes. Of these DEGs, 141 transcription factors (TFs) mainly consisting of MYB, AP2/ERF, WRKY, bHLH, and HSF were identified in response to salinity stress. Further, GO categories of the DEGs highly related to regulation of cell growth, cell periphery, sustained external encapsulating structure, cell wall organization or biogenesis, antioxidant activity, and peroxidase activity were significantly enriched for upregulated DEGs. The function of downregulated DEGs was mainly enriched in regulation of metabolic processes, oxidoreductase activity, and catalytic activity. Fourteen DEGs with response to salt tolerance were validated by real-time PCR. Taken together, the identification of DEGs' response to salt tolerance of cultivated peanut will provide a solid foundation for improving salt-tolerant peanut genetic manipulation in the future.


Introduction
Abiotic stresses comprise all nonliving factors, such as soil salinity, drought, extreme temperatures, and metal toxicity, giving rise to a serious threat to agriculture and affecting the plant production around the world [1]. With the increasing problems of global changes of climate, serious desertification of arable land, heavy metal population of soil, severe shortage of fresh water, and rapid growth of the population, it is very important to develop stress-resistant crops for sustaining growth and productivity in the abiotic stress conditions [2]. Soil salinity is one of the most destructive environment problems that can cause remarkable decrease of cultivated land area and crop productivity [3,4] As the FAO (Food and Agriculture Organization) Land and Plant Nutrition Management Service reported, more than 6% of the Earth's lands are affected by salt. About 45-million-hectare irrigated lands, which accounts for 19.5% of 230 million hectares, are affected by salinity. Over 1500 million hectares are under dry land agriculture, while 32 million (2.1%) are salt-affected to varying degrees [5].
Soil salinity has two main effects on the plant growth either by forming an osmotic potential or by the ionic toxicity effects of Na + and Cl − ions on plant cells [6]. When salt stress happens, the high concentration of salt lowers the osmotic pressure; thus, the plant cannot take enough water as in the soil condition; the concentration is higher than that in the plant cells [7]. The result of salt stress is that the plant stomata will be closed to conserve water and will stimulate the production of reactive oxygen species (ROS) like hydrogen peroxide and superoxide anion in cells. The ROS further disturbs a series of the plant cell processes by causing damage to lipids, proteins, and nucleic acids [8]. Ionic toxicity is associated with the balance concentration of Na + /K + and Na + /Ca + ratio as accumulation of Na + and Cl − ions and can inhibit cellular metabolism processes by inhibiting the relative enzymatic activities in the cells [9,10]. Na + is the crucial ion toxicity in most of the saline soils. It can facilitate entry into or exit out of plant cells by several ion transporters. One of the key responses to plant salinity is to regulate and balance the cellular ion homeostasis via restricting the cumulation of the toxic Na + [11]. Intracellular Na + can be transported out of plant cells by SOS1 transporter via activation of the Salt Overly Sensitive (SOS) signaling pathway, or into the root xylem via the high-affinity potassium transporters (HKTs) [12][13][14]. Intracellular Na + can also be isolated into vacuole by Na + /H + antiporters localized in the tonoplast membrane [15].
Peanut (Arachis hypogaea L.), which is a Leguminosae family plant that is cultivated in semiarid tropical and subtropical regions, is one of the most economically important crops for vegetable oil, proteins, minerals, and vitamins in the world [16]. The cultivated peanut variety is originated from the two diploid wild progenitor species of Arachis ipaensis and Arachis duranensis. As a settled plant, peanut is unable to escape from abiotic stresses. Soil salinity stress is a major threat to peanut productivity [17,18]. Plants can tolerate salinity stresses through modulating numerous genes and by coordinating the function of multiple genes from different metabolic pathways or regulatory systems [19]. Many studies have been conducted on physiological and biochemical changes in cultivated peanut under salt stress, but little information regarding genome and complete transcriptome of cultivated peanut is reported, and various studies about molecular mechanisms against salt stress focus on a few independent genes; therefore, it is difficult to obtain systematic biological genetic information about peanut against salt stress [18,[20][21][22][23][24].
Here, the next-generation transcriptome sequencing technique was performed to investigate the molecular foundation of wild peanut salinity tolerance. The high-throughput RNA sequencing as a revolutionary tool has been widely used to accurately monitor and quantify gene expression differences across the transcriptome under different treatments. In this study, a paired-end RNA sequencing was performed to examine the differential gene expression of cultivated peanut variety under long-term salinity stress conditions, based on our group's previous screening work [25] and the release of cultivated variety reference genome sequence by the International Peanut Genome Initiative (IPGI) group (https://www.peanutbase.org/peanut_genome). Several salt resistance genes were selected and validated using quantitative real-time reverse transcription polymerase chain reaction (qRT-PCR). The identification of the genes related to salt tolerance in cultivated peanut will provide a better understanding and cognition of the tolerance mechanism and further provide a wide gene resource to improve salt tolerance in the future peanut genetic manipulation.

Plant
Material. The seeds of the drought-resistant cultivar Fenghua3 were incubated in pots with a mixture of ver-miculite, perlite, and soil (1:1:1) at 26°C with photoperiod of 16 h. Seedlings were then grown in an artificial climate incubator for 14 days under restrained conditions (16 h light 26°C/8 h dark 22°C cycles). At the three-leaf stage, the robust overground parts and more uniform growing seedlings were harvested as control. The parts of these plants were then grown in salinity stress conditions with Hoagland's solution and 200 mM (1%) NaCl. The robust growing seedlings were harvested at 6, 12, 18, 24, and 48 hours later. Young leaves were immediately frozen with liquid nitrogen and stored at −80°C for further experiment. Three replicates were used for each time point.

Physiological
Measurements. The catalase (CAT), peroxidase (POD), and superoxide dismutase (SOD) activities and malondialdehyde (MDA) content were measured in each sample under salinity stress using physiological assay kits (Nanjing Jiancheng Bioengineering Institute, Nanjing, China) and a UV-Vis spectrophotometer (Shimadzu, Kyoto, Japan) following the manufacturers' instructions. Three replicates were used for each time point of each sample. Statistical calculations were performed using SPSS 12.0 (SPSS Inc., Chicago, IL, USA) with the level of significance setting at P value ≤ 0.05.

RNA Sequencing.
Total RNA was extracted from 20 mg young leaves using GeneJET Plant RNA Purification Mini Kit (Thermo Fisher Scientific, Waltham, Massachusetts, USA), and quantified by NanoDrop ND-2000 (Thermo Fisher Scientific, Waltham, Massachusetts, USA). RNA integrity was assessed using Agilent Bioanalyzer 2100 (Agilent Technologies, Santa Clara, California, USA). The samples with integrity number greater than 8 was used for library construction. The library of each sample was constructed using the TruSeq RNA Sample Preparation Kit (Illumina, San Diego, CA, USA) and SuperScript II Reverse Transcriptase (Invitrogen, Carlsbad, CA, USA). Paired-end sequencing was performed using an Illumina sequencing HiSeq 4000 platform with PE-150 module (Illumina, San Diego, CA, USA). The quality of all raw reads was firstly assessed by FastQC version 0.11.7 with default parameters (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). Adaptors and low-quality bases were trimmed using Trimmomatic version 0.38 with default parameters [26]. Reads with a Phred quality score above 30 were used for the following transcriptome analysis.

Transcript
Assembly and Annotation. The peanut genomic sequences and gene annotation were downloaded from PeanutBase (https://www.peanutbase.org/peanut_genome). The clean reads from all six samples were mapped to peanut genome assembly using HISAT2 version 2.1.0 in strandspecific mode and other parameters with default settings [27]. The transcripts were assembled using StringTie version 1.3.4d with default parameters [28]. The abundance of all genes was determined using the "scaledTPM" method in txImport version 1.8.0 package with default parameters [29].
2.5. Differential Gene Expression Profiling. Differential gene expression analysis was performed using DEGseq2 version 2 International Journal of Genomics 1.22.2 [30] in a Bioconductor package. The false discovery rate (FDR) and log2FC (log of fold change) was calculated for all expressed genes, and only the genes with the absolute log 2FC ≥ 2 and FDR ≤ 0:01 were considered differentially expressed between two groups. GO and KEGG functional enrichment analyses were conducted using the topGO package version 2.3.4 (http://bioconductor.org/ packages/topGO/) with default parameters in a Bioconductor package and KOBAS version 3.0 with default parameters [31], respectively.
Fourteen DEGs with a response to salinity stress were chosen for experimental validation by qRT-PCR. The primers of these genes were designed using the Beacon Designer 7.0 software. Total RNA was isolated using the same method as mentioned in RNA sequencing. Single-stranded cDNAs were reverse-transcribed using PrimeScript RT reagent kit with gDNA Eraser (Perfect Real Time) (TaKaRa, Otsu, Shiga, Japan). The SYBR Premix Ex Taq™ TM II reagent with SYBR Green I was used for qRT-PCR analysis by ABI 7500 Fast Real-Time PCR System (Thermo Fisher Scientific, Waltham, Massachusetts, USA). Each reaction comprised 10 μL 2X SYBR Premix Ex Taq II, 2 μL of single-stranded cDNAs, and 0.4 μM of forward and reverse primers in a final volume of 20 μL, upon initial denaturation at 95°C for 30 s, 40 cycles of denaturation at 95°C for 5 s and 30 s of annealing at 60°C, and 10 s extension steps at 72°C. qRT-PCR assays were with three biological and technical replicates. The relative gene expression levels were measured using the 2 -ΔΔCT method [32] and normalized against the amount of the Actin gene [33]. All specific primers are listed in Table S1.

Physiological Characteristics.
The value of POD activity in leaf tissues showed no significance for the first 6 hours and then exhibited a marked increase by 137%, 176%, 192%, and 197% under 12 h, 18 h, 24 h, and 48 h with salinity stress, respectively, compared with the control (P ≤ 0:01) (Figure 1(a)). CAT activity was observed to have the same tendency with POD. Furthermore, CAT activity had a significant increase by 123%, 169%, 205%, and 213% under 12 h, 18 h, 24 h and 48 h with salt treated, respectively, as compared with the control (P ≤ 0.05) (Figure 1(b)). As compared with the control, results showed that SOD activity had a significant decrease by 144% in the first 6 h (P ≤ 0:05) and then was significantly enhanced by 118%, 147%, 152%, and 144% under 12 h, 18 h, 24 h, and 48 h with salinity stress, respectively (P ≤ 0:01) (Figure 1(c)). The oxidative degradation was recognized as malondialdehyde (MDA) content, which is the product of lipid peroxidation. The results showed that in all the time points, lipid peroxidation was significantly affected by salinity stress (P ≤ 0:05). As compared with the control, salinity levels at 200 mM caused 129%, 163%, 177%, and 186% increase in 12 h, 18 h, 24 h, and 48 h, respectively ( Figure 1(d)). In conclusion, the results of physiological characteristics revealed that the genes related to salinity stress were changed after 6 h salt treatment, which is consistent with our group's previous report [25]. The young leaves of peanut seedlings growing within salinity stress conditions for 24 h were used to perform transcriptome analysis as the curves of physiological traits reached their peaks after salt stress.

Transcriptome Sequencing and Assembly.
To obtain the global scenario of peanut gene expression under salinity stress, six cDNA libraries (three replicates for each group) were performed for Illumina RNA sequencing. After adaptor removal and quality filtering, a total of 163.33 million clean reads were obtained with a Q30 ratio more than 90.63% ( Table 1). The quality-checked high-quality sequencing reads of each sample were individually aligned to peanut genome assembly using HISAT2. Of the total reads, a range of 76.7%-77.30% clean reads were mapped to peanut genome assembly (Table 1). After read aligning and assembling, 53,653 genes were identified in the transcriptome and used to perform analysis of the differentially expressed genes.

Identification of DEGs Involved in Salinity
Stress. Based on gene expression level, differential gene expression analysis (see Differential Gene Expression Profiling for details) was implemented between salt-treated (HSR) and the control (CKR) samples. The result showed that in salinity stress, a total of 3,425 genes were found to be differentially expressed with FDR ≤ 0:01 and absolute log FC ≥ 2, of which comprised 2,013 upregulated genes and 1,412 downregulated genes ( Figure 2, Table S1). Among the DEGs, 141 transcription factors (TFs) (89 upregulated and 52 downregulated) were identified in response to salinity stress, which were assigned to ten TF families (Table S2). The TF families of MYB, AP2/ERF, WRKY, bHLH, and HSF were included with relatively largest volumes representing 87.94% of all TFs (Table 2).

Functional Enrichment Analysis of DEGs to Salinity
Stress. To determine the GO functions of DEGs involved in salinity stress, topGO tool was used for enrichment analysis of GO terms in the whole transcriptome analysis. The GO functional enrichment analysis showed that the upregulated and downregulated DEGs were significantly enriched in different GO categories (Table S3). For the upregulated DEGs, 18 GO categories were significantly enriched at term level of 3 with FDR ≤ 0:01 ( Figure 3, Table S4). These GO terms were highly related to regulation of extensive biological activities, such as oxidation-reduction process (FDR = 1.  Table S3). In comparison with the functional enrichment of upregulated DEGs, the downregulated DEGs were mainly involved in regulation of metabolic process, oxidoreductase activity, and catalytic activity. Meanwhile, the upregulated genes in the salinity stress were involved mostly in regulation of cell growth and cell periphery, sustained external encapsulating structure, cell wall organization or biogenesis, antioxidant activity, and peroxidase activity ( Figure 3, Table S3), which was consistent with our physiological characteristics.
To understand the gene functions and pathways, the DEGs were analyzed using KOBAS tool with KEGG database. The results showed that three pathways were significantly enriched in upregulated DEGs and 19 pathways enriched for downregulated DEGs. Of these, the pathways of phenylpropanoid biosynthesis (FDR = 1:72E-10), pentose    Table S4). KEGG enrichment analysis indicated that the DEGs involved in those pathways were possibly linked to salt tolerance.

qRT-PCR Validation.
To validate the expression profile, fourteen DEGs with response to salt tolerance, involved in plant hormone signaling, transcription factors, secondary metabolism, and oxidative damage, were chosen for experimental validation by qRT-PCR (Table S5). As shown in Figure 5, the results indicated that the gene expression changes analyzed by qRT-PCR were mostly consistent with those obtained by Illumina RNA sequencing except for Aradu.8H8GJ, representing 92.86% consistency.

Discussion
Transcriptomic approach is an effective tool which provides a global information towards the gene expression patterns of a plant organism under any conditions, such as drought, light, temperature, aflatoxin, and salinity stress. Numerous studies about transcriptomic characteristics of peanut have been conducted to evaluate the expression pattern of induc-ible genes under abiotic and biotic stress conditions [25,[34][35][36][37][38][39]. Transcriptome data are valuable resources for exploring plants under stress conditions. To our knowledge, there have been an increasing number of studies on transcriptome analysis of legume species using high-throughput RNA sequencing approach, among those regarding the effects of salinity stress were reported on Medicago [40][41][42], glycine [43][44][45], and common bean [46][47][48]. Studies focusing on transcriptome research of the peanut under salt stress are rare, but there exists a substantial amount of reports regarding the physiological responses to salt stress in peanut based on a de novo transcriptomic sequence assembly method [49][50][51]. In this study, a high-throughput RNA sequencing of cultivated peanut variety under long-term salinity stress conditions was performed and deeply analyzed, which will contribute to a well annotation of cultivated peanut reference genome by the IPGI group. Furthermore, GO and KEGG functional enrichment analyses were used to better understand the functions of DEGs in salt tolerance conditions. The DEGs were mainly involved in oxidation-reduction process; oxidoreductase activity; cell wall organization or biogenesis; response to stress; cofactor binding; metabolic process, such as starch and sucrose metabolism; circadian rhythm; flavonoid biosynthesis; and nitrogen metabolism, which showed significance in responds to salinity stress (Figures 3 and 4). Salt tolerance is a complicated character that is controlled by multiple different genes especially TFs in plants [52]. During the response to salt tolerance stress, a number of genes were activated, leading to the accumulation of numerous proteins involved in resistance to abiotic stress, which were mostly regulated by specific TFs [53]. Among these TFs, the MYB family has been well characterized for their regulatory roles in the response of plants to abiotic stress [54,55]. Several MYB genes have been identified as key factors in the signaling pathways for Arabidopsis and rice resistance to abiotic stresses [56,57]. MYB family proteins were also demonstrated in plants as a regulator for mediating salt tolerance gene expression under abiotic stress [58,59]. A previous study reported that the MYB family proteins were the second most highly expressed TF family in drought tolerance Arachis plants [34]. In our study, the TF expression pattern analysis indicated that MYB family genes responding to salt stress were in large quantities including 47 DEGs. Eight DEGs involved in the MYB family (MYB48, MYB51, MYB60, MYB62, MYB64, MYB98, MYB114, and MYB118) were further selected and validated by qRT-PCR with root and leaf tissues among all time points under salt tolerance. The result showed that MYB48, MYB60, and MYB98 were consistently increased, and the expressions of MYB51 and MYB118 were decreased among time points in leaf tissues under additional salt (Figure 6(a)). As compared with the gene expression patterns in leaf tissues, all MYB genes exhibited a consistent increase among all time points with salt tolerance (Figure 6(b)). The observed is not difficult to understand as the root is the one initially exposed to salt stress. Among these TFs, one (MYB118) encodes R2R3 MYB proteins and others encode the MYB-related proteins. In contrast to the R2R3-MYB genes, the MYB-related genes  International Journal of Genomics have attracted little attention, and few have been studied functionally. Previous studies found that CCA1-like and CPC-like genes could involve in the maintenance of circadian rhythms [60][61][62] and in control of cellular morphogenesis [63]. A novel potato single MYB-like domain protein (StMYB1R-1) played positive roles in potato drought resistance, according to Shin et al. [64]. Our research provides evidence that the MYB-related proteins may be involved in peanut salinity stress regulation and needs more investigation in the future work.
Salt tolerance significantly affected the core genes related to a series of secondary metabolism, such as flavonoid biosynthesis, phenylpropanoid biosynthesis, starch and sucrose metabolism, circadian rhythm, and nitrogen metabolism. Previous studies reported that the flavonoid pathway was mostly induced with response to a series environmental stress as a protective mechanism to oxidative stresses induced by metal ions, high light, temperature, drought, salt, or C nutrition in plants [65][66][67][68]. Many studies revealed that with plant under drought conditions, the expression levels of two key enzymes participating in flavonoid biosynthesis were always increased [69,70]. In our study, the expression levels of two key enzymes, chalcone synthase (51 DEGs) and chalcone-flavanone isomerase family protein (one DEG) involved in flavonoid biosynthesis, were also observed to significantly increase with salt stress conditions. As the precursors of a wide range of phenolic compounds, phenylpro-panoid compounds were explored to have the important roles in a number of plant protection mechanisms such as plant defense response, abiotic stresses, and signal transduction [71,72]. PAL (phenylalanine ammonia lyase), regarded as a central enzyme in the phenylpropanoid pathway to catalyze the deamination of phenylalanine to provide cinnamic acid, was investigated in Salvia species, and the result showed that salt stress could increase PAL activity and total phenolic accumulation in the early few hours with stress treatment [73]. Our result had a consistence with these findings and showed that four PAL genes (Aradu.IU1HH, Araip.69J63, Araip.GM19P, and Araip.V9S7Z) were increased with at least 4-folds (Table S1). Nitrogen is an essential nutrient for plant growth, development, and reproduction. As a major source of nitrogen provider of plants, nitrate assimilation was specifically affected by salinity with subsequently inactivation of the nitrate reductase [74]. In our study, the expression levels of 11 nitrate transporter coding genes were decreased under salt stress. This indicated that nitrogen metabolism in peanut with salinity stress was decreased. Glutathione transferases (GSTs) are induced with a series of biotic and abiotic stresses in plant, which are important for protecting host plants against oxidative damage [75]. About 25 glutathione S-transferase family coding genes showed differential expression increasing patterns when exposed to salt stress in our data (Table S1), which indicated that the GSTs in peanut could enhance salt International Journal of Genomics tolerance. Moreover, in our study, the most of PODencoding genes showed significantly decrease under salt stress; also, the all CAT-and SOD-encoding genes showed slight decrease with salt stress (Table S1), which were inconsistent with our physiological characteristics. This phenomenon may be caused by several complex biological processes and posttranscriptional regulation by noncoding RNAs and still need to be further studied. Multiple membrane ion transporters and pH-related transporters were proposed to mediate plant cellular signaling under salinity stress [11,13]. The preset studies indicated that the Arabidopsis transporter HKT functioned as a salt tolerance determinant and mediates Na + influx and K + transport in plants [14,76]. Research in Arabidopsis revealed that reduction of the AtHKT1 gene expression could lead to hypersensitivity to salt ion with more Na + accumulated in the leaves while directly stimulating K + loading into the xylem by AtHKT [77]. In our findings, the expression level of Araip.58313 (high-affinity K + transporter 1, homologous Arabidopsis HKT1) was found significantly decreased under salt stress conditions (Table S1). Our result may indicate that the peanut growing in salinity stress could reduce the Na + /K + ratio in the leaves via reducing the inhibition of K + uptake. The SOS signaling pathway comprising SOS3, SOS2, and SOS1 has been recognized as key mechanism to maintain ion homeostasis by controlling the cellular signaling under salinity stress [78]. In our study, the expression level of Araip.TT7Q9 (homologous Arabidopsis SOS1), Aradu.FWA9B (homologous Arabidopsis SOS2), and Araip.IZ8TF (homologous Arabidopsis SOS3) was observed to be significantly changed by comparing with control (Table S1), which indicated that salinity stress may activate the Na + /H + antiporter and facilitate Na + efflux through cellular plasma membrane to regulate vacuolar H + -    International Journal of Genomics ATPase activity. The expression and activity of vacuolar H + -ATPase (V-ATPase) or vacuolar pyrophosphatase (V-PPase) facilitated the Na + sequestration into vacuoles [79]. In our findings, the expression levels of three V-ATPase coding genes (Aradu.0U1AX, Aradu.ZRX72, and Araip.08L1N) and two V-PPase coding genes (Aradu.PJ77P and Araip.CC90S) were also significantly changed with salinity stress conditions (Table S1). These results indicated that V-ATPase and V-PPase in peanut may enhance the Na + sequestration into vacuoles and modulate ion homeostasis in high salinity stress conditions.
In summary, we performed the transcriptome analysis of a cultivated peanut variety under long-term salinity stress conditions based on a well annotation of cultivated peanut reference genome. Functional analysis of DEGs showed that salt tolerance genes were mainly involved in oxidationreduction process; oxidoreductase activity; cell wall organization or biogenesis; response to stress; cofactor binding; metabolic process, such as starch and sucrose metabolism; circadian rhythm; flavonoid biosynthesis; and nitrogen metabolism, which showed significant in response to salinity stress. The identification of differentially expressed salt 9 International Journal of Genomics tolerance genes in cultivated peanut will provide a better understanding and cognition of the tolerance mechanism and will further provide reference for improving salttolerant peanut genetic manipulation in the future.

Data Availability
The RNA sequencing data used to support the findings of this study have been deposited in the Short Read Archive at NCBI with accession number of SRR8177741. The qRT-PCR validated genes were submitted to NCBI with GeneBank accession number as following: Aradu.

Conflicts of Interest
The authors declare that they have no conflict of interest.

Authors' Contributions
HZ, XBZ, and QXS contributed to the study design, data collection, and data analysis; wrote the first draft; and revised the manuscript. CXY and JW contributed to the study design, sample preparation, and data analysis. CLY and CJL contributed to the study design, physiological analysis and data collection, and qRT-PCR validation. SHS and FZL were the supervisors of the project and contributed to the study design, data analysis, and manuscript revision. All authors reviewed and accepted the content of the final manuscript.