Expression profiling of Chrysanthemum crassum under salinity stress and the initiation of morphological changes

Chrysanthemum crassum is a decaploid species of Chrysanthemum with high stress tolerance that allows survival under salinity stress while maintaining a relatively ideal growth rate. We previously recorded morphological changes after salt treatment, such as the expansion of leaf cells. To explore the underlying salinity tolerance mechanisms, we used an Illumina platform and obtained three sequencing libraries from samples collected after 0 h, 12 h and 24 h of salt treatment. Following de novo assembly, 154,944 transcripts were generated, and 97,833 (63.14%) transcripts were annotated, including 55 Gene Ontology (GO) terms and 128 Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. The expression profile of C. crassum was globally altered after salt treatment. We selected functional genes and pathways that may contribute to salinity tolerance and identified some factors involved in the salinity tolerance strategies of C. crassum, such as signal transduction, transcription factors and plant hormone regulation, enhancement of energy metabolism, functional proteins and osmolyte synthesis, reactive oxygen species (ROS) scavenging, photosystem protection and recovery, and cell wall protein modifications. Forty-six genes were selected for quantitative real-time polymerase chain reaction detection, and their expression patterns were shown to be consistent with the changes in their transcript abundance determined by RNA sequencing.


Introduction
Salt can diffuse from underground soil, wind and rain can carry salt from the ocean, and human overuse of chemical fertilizers also causes salty soils [1,2]. These factors have increased the amount of land affected by salt. According to a study by the Food and Agriculture Organization (FAO) in 2011, more than 800 million hectares of land were estimated to be affected by salt around the world.

Plant materials and growth conditions
Plants were obtained from the Chrysanthemum Germplasm Resource Preservation Center of Nanjing Agricultural University, China. Shoot cuttings of C. crassum were rooted and grown in matrix (vermiculite: perlite: nutritive soil = 1:1:1). Rooted seedlings at the 6-7 leaf stage were selected and then transplanted into plastic containers (volume 20 L) filled with diluted (1:2) Hoagland nutrient solution, with aeration for 24 h•d -1 . The plants were maintained in a greenhouse with a normal photoperiod and sunlight, an average temperature of 25˚C and relative humidity of 70%; all solutions were renewed every two days. After one week, one set of plants grown in Hoagland solution alone was kept as a control (CK), and salt treatment (S) was performed by supplementing the nutrient solution with 120 mM NaCl.

Recording of morphological changes after salt treatment and paraffin sectioning
At 0 d, plants that showed similar sizes and shapes were selected from the salt treatment and control groups and photographed. Photographs of these plants were subsequently taken at 5 d, 10 d, and 20 d. At the same time points, the fourth unfolded leaves of the plants were traced on paper without injuring them.
The fourth unfolded leaves numbered from the tops of both salt-treated and control plants were picked (at 0 d, 5 d, 10 d, and 20 d), cut into pieces of 0.3 Ã 0.3 cm, then immediately fixed in FAA solution. The samples were subsequently dehydrated through an alcohol series, infiltrated with xylene, and embedded in paraffin wax. Transverse sections were cut to a thickness of 25 μm and stained with 0.1% toluidine blue solution (w/v). The sections were finally observed under an Olympus BX41 microscope [16].

Determination of Na + and K + levels
The unfolded leaves numbered from the top of both salt-treated and control plants were collected (approximately 3 g, at 0 d, 5 d, 10 d, 20 d, and 30 d). The samples were then oven-dried at 70˚C for 48 h, after which 50 mg DW was digested in 35% HNO 3 . Na + and K + was resuspended in 10 ml of HCl (0.1 N) and the solutions were filtered. Ions were quantified via flame atomic absorption spectrophotometry (VARIAN spectra-300) [17]. Three replicates of each sample were examined.
bp. The fragments were used for the synthesis of first-strand cDNA employing random hexamer primers. Second-strand cDNA was synthesized with the SuperScript Double-Stranded cDNA Synthesis kit (Invitrogen, Camarillo, CA). The cDNA was then purified using a Qia-Quick PCR Extraction Kit (Qiagen, Hilden, Germany) and resolved with EB buffer for end reparation and poly(A) addition. The products were ligated with each another using sequencing adapters, and after agarose gel electrophoresis, a suitable size range of fragments was selected for PCR amplification. The resulting library was sequenced using the Illumina HiSeqTM 2000 system [18].

Data filtering and de novo assembly
The image data output from the sequencing system was transformed into raw reads and stored in FASTQ format. These data were filtered to remove raw reads that included adapter sequences, had more than 5% unknown nucleotides or that were of low quality (the percentage of reads with quality value ! 10 was more than 20%). Transcriptome de novo assembly was carried out with Trinity [19]. The resulting Trinity sequences were considered unigenes. After sequence splicing and redundancy removal, gene family clustering was performed, and the unigenes were divided into two classes. In one cluster (with the prefix "CL" followed by the cluster id), there were several unigenes showing similarity of greater than 70%. The others were singletons (which have the prefix "Unigene"). BLASTX [20] alignment between each unigene sequence and those registered in the Nr (non-redundant protein database, NCBI), Nt (non-redundant nucleotide database, NCBI), Swiss-Prot, and GO (gene ontology, http://www. geneontology.org/) databases was performed, and the best alignments for inferring the directionality of each unigene were determined. For the best alignments in which the outcome from the various databases conflicted with one another, the priority order applied was Nr, followed by Swiss-Prot. The software tool ESTScan [21] was used to assign directionality.

Gene annotation and analysis
Functional annotation was assigned using the protein (Nr and Swiss-Prot) and GO databases. BLASTX was employed to identify related sequences in the protein databases based on an Evalue of less than 10 −5 . The annotations acquired from Nr were processed by using the Blas-t2GO program [22] to obtain the relevant GO terms, which were then analyzed with WEGO software [23] to assign a GO functional classification and to illustrate the distribution of gene functions.

Prediction of unigene coding regions
Unigenes were aligned to protein databases using BLASTX with the following order of priority: Nr, Swiss-Prot, and KEGG. The proteins with highest ranks in the BLAST results were selected to determine the coding region sequences of the unigenes, and the sequences of the coding regions were then translated into amino sequences using the standard codon table. Thus, both the nucleotide sequences (5' -> 3') and amino acid sequences of the Unigene coding region were acquired. Unigenes that could not be aligned to any database were scanned by ESTScan to determine the direction of the nucleotide sequence (5' -> 3') and the amino acid sequence of the predicted coding region [21].

Gene validation and expression analysis
Samples were prepared, and total RNA was extracted using the method indicated above. A total of 16 unigenes that responded to salt stress were chosen for validation. Three independent biological replicates of each sample were used in the analysis. A set of gene-specific primer pairs was designed using Primer3 software [24]. Reverse transcription was performed with M-MLV reverse transcriptase (TaKaRa). The relative expression of these genes was determined through qRT-PCR analysis using a SYBR 1 Green reaction kit (TaKaRa), with the elongation factor 1-alpha gene as a reference. Then, a Mastercycler ep realplex device (Eppendorf, Hamburg, Germany) was used to run the qPCR assays. The transcription data were calculated using the −ΔΔCt method [25].

Results
Morphological changes of C. crassum and Na + accumulation in leaves under salinity stress C. crassum showed remarkable tolerance under salinity stress, as the salt-treated plants grew almost as fast as the controls (Fig 1A), with only the leaves near the roots wilting, revealing the damage due to salinity stress. As shown in Fig 1A, the leaves treated with salt grew differently from the controls; therefore, we chose salt-treated and control leaves of the same size and followed their growth by tracing the outlines of the leaves without injuring them at 0 d, 5 d, 10 d, and 20 d. The shapes of the leaves are shown in S1 Fig. It was apparent that the leaves that were treated with salt were larger than those of the controls, and the analysis of paraffinic slices of leaves ( Fig 1B) showed that the leaves that were salt-treated were thicker than the leaves of the controls. This incrassation was mainly caused by cell expansion and not cell replication, as the number of cell layers was not increased. The morphological changes that occur in C. crassum under salinity stress have rarely been studied; however, they may be relevant to the remarkable salinity tolerance of C. crassum.
The accumulation of Na + in the leaves of salt-treated C. crassum was slow, requiring days to reach a peak, and a balance was maintained around the peak as time passed (Fig 2A). The K + content of leaves also increased under salinity stress. These behaviors are beneficial to the survival and growth maintenance of C. crassum under salinity stress. The transcriptomes of the three samples helped us to determine the implications of these behaviors and identify which strategies C. crassum employs under salinity stress.

De novo assembly and quantitative assessment of RNA sequences
To investigate the global gene expression profiles of the C. crassum transcriptomes at different stages of salinity stress, we constructed three cDNA libraries using leaves from hydroponic cutting seedlings treated with 120 mM NaCl for 0 h, 12 h, and 24 h.
Sequencing was carried out on an Illumina platform; a total of 29,831,115,240 nt were generated, and 346,363,148 raw reads were obtained. After removing low-quality regions, adapters, and contamination, we obtained a total of 331,456,836 clean reads with Q20 > 98.34% and a GC percentage between 43.38% and 43.75% (Table 1).
In the assembly results, 154,944 unigenes were detected; the total length of the unigenes was 167,656,676 nt; the average length was 1,082 nt; and the N50 was 1,661 nt. We generated 160,301, 164,718, and 174,748 contigs for the samples obtained at 0 h, 12 h, and 24 h after salt treatment, with average lengths of 328, 321, and 318 nt, respectively, and average N50 lengths of 603, 556, and 546 (Table 1).

Sequence annotation and functional categories of annotated sequences
For the annotation of functional information in the assembled unigenes, we used several applications, including information on protein sequence similarities, GO and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. All of the sequences were successively subjected to BLAST searches in the NCBI Nr database, the Nt database, and the Swiss-Prot protein database, with E-values < 1e-5. Using the best hits provided by BLAST, the most reasonable functions were assigned to the sequences. After functional annotation analysis, we obtained 93,296, 68,322, 65,127, 60,906, and 55,754 unigenes annotated to the NR, NT, Swiss-Prot, KEGG, and GO databases, respectively, totaling 97,833 annotated unigenes.  Regarding the functional categories of the 154,944 unigenes, a total of 55,754 unigenes were assigned at least one GO term, which are mainly divided into three categories: biological processes, molecular functions, and cellular components. Among the 22 GO terms corresponding to biological processes, the terms, "cellular processes", "metabolic processes", "single-organisms", "response to stimulus", and "biological regulation" were significantly overrepresented, while "cell", "cell part", and "organelle" were significantly overrepresented among the 17 GO groups of cellular components, and "catalytic activity" and "binding" were significantly overrepresented among the 16 GO groups of molecular functions (Fig 3).
Additionally, in the prediction of protein coding regions, the number of CDSs that were mapped to the protein database was 92,001, and the number of predicted CDSs was 7,012. The total number of CDSs was 99,013 (S3 Fig).

Determination and functional categorization of differentially expressed genes (DEGs)
To determine the genes showing altered expressions under salinity stress conditions, DEGs were defined as those with an FDR 0.001, a log 2 ratio ! 1, and an RPKM > 2 in the Note: The Q20 percentage is the proportion of nucleotides with a quality value > 20. The N percentage is the proportion of unknown nucleotides in clean reads. The GC percentage is the proportion of guanidine and cytosine nucleotides among total nucleotides. N50 is 50% of the assembled bases that were incorporated into sequences with a length of N50 or longer. Therefore, there were many genes that responded to salinity stress in a very early stage, whose levels gradually returned to normal conditions. All of the DEGs were categorized into 55 functional groups. The GO classification maps of the two comparisons were similar to each other. In these maps, among biological processes, the dominant groups were "cellular process", "metabolic process", "single-organism process", and "response to stimulus". In the cellular component category, the dominant groups were "cell", "cell part", "organelle", and "membrane". In the molecular function category, the dominant groups were "catalytic activity" and "binding" (S5 Fig).

Verification of differential gene expression through qRT-PCR
To validate the results obtained from the Illumina sequencing data, 46 relevant unigenes were selected for quantitative real-time PCR (qRT-PCR) analysis with samples that were treated with 120 mM NaCl for 0 h, 12 h, and 24 h. The expression trends of the unigenes from the qRT-PCR and RNA sequencing analyses were consistent (S6 Fig). Ten of these Unigenes are shown in the text (Fig 5a), and their transcript abundances at 0 h, 1 h, 6 h, 12 h, 24 h and 48 h were also measured by qRT-PCR, demonstrating that these unigenes respond strongly to salt (Fig 5b). These results demonstrated that the transcriptomic profiling data accurately reflected the responses of C. crassum to salt stress.

Global patterns of transcription in response to salt treatment
To improve the salinity tolerance of plants, researchers worldwide have performed myriad studies to elucidate how plants respond to salinity stress. However, salinity tolerance mechanisms are too complex for researchers to have achieved a full understanding as yet. Despite being limited by technologies and knowledge, studies on elite salinity-tolerant germplasms are helpful for improving the salinity tolerance of other germplasms. For example, over-expression of the C. crassum plasma membrane Na + /H + antiporter gene CcSOS1 improved the salinity tolerance of chrysanthemum 'Jinba' [12]. Therefore, functional genes from C. crassum could be used to improve the salinity tolerance of other species of chrysanthemums.
There have been many studies on the salinity tolerance of model plants, such as Arabidopsis and Nicotiana tabacum [26], and on some valuable crops with known genomes, such as Oryza sativa and Zea mays [27]. These studies have provided instructional theories and established components related to plant salinity tolerance, such as transcription factors, plant hormones, photosynthesis, osmolyte accumulation, ion elimination and toxicity alleviation. In cases where whole-genome sequencing data are lacking, NGS technology provides a powerful tool for transcriptome analysis of plants, and the de novo assembly of transcript sequences offers a rapid approach for determining expressed gene catalogs of non-model plants [28]. In this study, we performed sequencing analysis of C. crassum in response to salinity stress using the Illumina platform and de novo transcriptome assembly. Many DEGs were identified after salt treatment (Fig 4A). To better understand the function of the DEGs in leaves of C. crassum under salinity stress, the DEGs were categorized into 52 functional groups based on GO classification terms (S5 Fig). The ten most-enriched, dominant GO terms were "cell", "cell part", "metabolic process", "cellular process", "catalytic activity", "organelle", "binding", "single-organism process", "response to stimulus", and "membrane". These categories are also dominant in the transcriptomes of Halogeton glomeratus [29] and Reaumuria soongorica [30]. This result shows a similarity to other salt-resistant plants and emphasizes the remarkable salinity tolerance of C. crassum.
The DEGs determined in the two comparisons were also mapped to 128 pathways with KEGG annotations. By analyzing the obtained maps, we could not only identify which metabolic functions responded to salinity stress but also determined the particular genes in those metabolic pathways that were up-or downregulated. This provided us with conventional references for determining the functional genes and metabolic pathways involved in responses to salinity stress. For example, key enzymes involved in energy metabolism pathways (such as "gluconeogenesis (ko00010)", "citrate cycle (ko00020)", and "oxidative phosphorylation (ko00190)"); amino acid synthesis pathways (such as "alanine, aspartate and glutamate metabolism (ko00250)", and "valine, leucine and isoleucine biosynthesis (ko00290)"); and DNA-, RNA-and protein-related pathways (such as "DNA replication (ko03030)", "RNA transport (ko03013)", "ribosome (ko03013)" and "protein export (ko03030)") were upregulated; the genes encoding these enzymes were regarded as functional genes that respond to salinity. Using these tools, we selected functional genes and drew conclusions about strategies employed by C. crassum under salinity stress.

Earliest-induced genes under salinity stress, such as transcription factors and signal molecules
Genes encoding transcription factors are some of the first genes to be induced under stress conditions [8]. Transcription factors are considered to have the greatest effect on crop salinity tolerance because these key transcription factors may regulate the induction or repression of a range of salinity tolerance genes. In the transcriptomes obtained in the present study, many reported transcription factors responding to salinity stress were found, such as bHLHs (Unigene40838, Unigene37206, and Unigene40839), DREB1A (CL6489.Contig1 and CL6489.Contig2), AP2/ ERF (Unigene38328 and CL12828.Contig2), NAC (CL4758.Contig1, CL14612.Contig1, and Unigene60684), AREB (Unigene41165), WRKY (Unigene40047 and CL15510.Contig1), GTL (CL10054.Contig1) and MYB2 (CL1046.Contig3 and CL1046.Contig1). All of these genes were upregulated at 12 h and returned to a low level at 24 h, and these genes have been reported to respond to stresses and regulate the growth conditions of plants under stress [31][32][33].
Additionally, some signal transduction molecules were found to be upregulated, such as CaMs (Unigene69398, Unigene70126, and Unigene70474), YWHAs (Unigene60260 and Uni-gene8553), and MAPKs (CL7742.Contig3, Unigene71806, and Unigene46385). These molecules have been reported to play a regulatory role in signal transduction [34]. These earliestinduced genes under salinity stress trigger a series of reactions related to stress endurance.

Genes related to plant hormone transduction pathways were differentially altered under salinity stress
In the obtained transcriptomes, genes involved in hormone transduction pathways were differentially expressed. The four key enzymes in the abscisic acid transduction pathway (PYL  [35] and key enzymes in the salicylic acid transduction pathway (TGA [Unigene21774], and PR-1 [Uni-gene33280]) [36]. These three well-studied hormones respond to stresses, can trigger many reactions, and may react with growth-related hormones, resulting in morphological changes in plants [37].
Some of the identified genes in hormone signal transduction pathways were considered to present functions regulating the morphological changes observed in C. crassum. For example, key enzymes in the auxin transduction pathway (TIR1 [Unigene46121], AUX [CL1780.Con-tig1], ARF [Unigene53204], and SAUR [CL16292.Contig2]) [38] were upregulated. Auxin plays an important role in coping with salinity and could help to maintain an ideal growth rate [39]. In the cytokinin transduction pathway [40], A-ARR (CL9399.Contig1) was downregulated, which could affect cell division [41], and leaf cells did not divide under salinity stress (Fig 1B). In the brassinosteroid transduction pathway [42], TCH4 (Unigene38906, Uni-gene13917, and Unigene54097) was upregulated, which could contribute to cell elongation [43], while CYCD3 (CL15069.Contig1, and Unigene61981), a cyclin that may contribute to cell division [44], was downregulated. These findings suggested that these three hormones played an important role in leaf incrassation and succulence in C. crassum.
Genes that contribute to overcoming osmotic stress caused by salinity in C. crassum Osmotic stress must accompany salinity stress [8]. In the transcriptomes of C. crassum, many water retention-related genes were identified. For example, aquaporins can regulate water transport and may have an important function (i.e., the detection of osmotic and turgor pressure gradients) [45]. There are many types of aquaporins in the plasma membrane and vacuolar membrane of plants, such as PIPs, TIPs, and NIPs. These aquaporins present distinct cell type-and tissue-specific expression patterns, and some of them respond to stress. However, their expression profiles vary under salinity stress, and some posttranscriptional mechanisms regulate aquaporin trafficking to the plasma membrane [46]. In our data, PIPs (Unigene72479, CL11177.Contig14, and CL392.Contig3) and NIPs (Unigene43528, CL6774.Contig4, and Uni-gene54541) were upregulated after salt treatment, while TIPs (CL1377.Contig1 and CL6342. Contig4) were downregulated after salt treatment. This means that different aquaporins showed different expression profiles under salinity stress in C. crassum and contribute to the regulation of water transport.
In addition, LEAs (Unigene52031, CL16208.Contig1, Unigene69354, Unigene70057, Uni-gene70137, and Unigene69772) and DHNs (CL3765.Contig11, Unigene44626, and Uni-gene33157) were significantly upregulated and were widely expressed after salt treatment. These LEAs and DHNs have been reported to be abundant proteins in plants that are highly hydrophilic and provide protection to plants under osmotic stress [47].
To improve osmotic pressure, some osmolytes must be synthesized to endure osmotic stress [8]. In our transcriptomes, a number of osmolyte genes were upregulated after salt treatment, such as TPS (Unigene42563, and Unigene14695), TPP (CL14450.Contig3, and Unigene45055), and GMDS [CL9344.Contig1]) [49] were also upregulated. Amino acids and monosaccharide are types of osmolytes [50]. We also found that some enzymes that participate in the synthesis of ascorbate (including GME [CL9645.Contig2], VTC4 [Unigene72075], and D-threo-aldose 1-dehydrogenase [Unigene71897]) [51] were upregulated after salt treatment. Ascorbate is a soluble carbohydrate that is abundant in leaves, and its accumulation would protect photosynthesis and improve the osmotic pressure of cells. Moreover, ascorbate is a cofactor for a large number of key enzymes and can influence both cell wall synthesis and hormone transduction [52]. Through the expression of the above genes, the plants were prevented from lacking water, and the ability of water to be taken up and retained was improved.
Genes that contribute to the exclusion and compartmentalization of Na + To prevent ion toxicity caused by salt, Na + must be kept from accumulating to excessive levels in protoplasts. The results in Fig 2 show that there must be some metabolic pathways that control the exclusion of Na + and the passive assimilation of K + . In our transcriptomes, F-type H +transporting ATPases (Unigene65029, Unigene69416, Unigene69674, CL11916.Contig1 and Unigene69648) and V-type H + -transporting ATPases (Unigene69565, Unigene70189, Uni-gene70388, Unigene70530, and Unigene71306) were significantly upregulated after salt treatment. These H + -transporting ATPases establish an electrochemical H + gradient across membranes to energize the membranes and participate in the synthesis of ATP [53,54], thereby providing energy and electric potential for Na + exclusion. Some ion transporters that have been frequently reported were identified in our transcriptomes, such as HKTs (CL15304. Contig2, and CL15304.Contig1), SOS1 (Unigene38015), NHXs (Unigene12739, CL15766.Con-tig3, and Unigene12656), and AVPs (Unigene61687, Unigene61689, and Unigene61692). They were upregulated after salt treatment and shown to contribute to the reduction of the Na + concentration in the cytoplasm [6].
Genes contributing to the supply of energy and materials to the reactions triggered by salinity Salinity causes many types of reactions in cells, and these reactions require energy and materials to function properly. Searching our transcriptomes revealed many energy metabolismrelated genes, such as AAC (CL17351.Contig1, Unigene69746, and Unigene68723), ALDH (Unigene68789, Unigene69865, Unigene69866, and CL4198.Contig4), GTP (Unigene69854, Unigene56731, and Unigene68499), GAPDH (CL2444.Contig1, CL370.Contig1, CL2444.Con-tig2,and Unigene68703), and ADHs (CL14096.Contig2, Unigene71447, and Unigene69264), that were significantly upregulated after salt treatment and are helpful in supplying energy to other metabolic pathways to endure salinity stress. This is a common reaction of plants under stress [56,57].
The above genes were significantly upregulated at 12 h and returned to a low level at 24 h. Thus, mechanisms that scavenge genes that have performed their functions and have no need to remain active are expected to function at this stage. We found that some proteasome subunits (Unigene71130, Unigene70519, Unigene71062, Unigene70385, and Unigene71486) and some ubiquitin-conjugating enzymes (Unigene69417, Unigene30464, and CL15637.Contig1) were upregulated after salt treatment, which work as "cleaners" to digest unwanted mRNA and proteins [60] and return the expression of certain genes to normal levels (Fig 4 and S4  Fig).

Behavior of the photosynthetic systems under salinity stress
The growth of plants must be supported by photosynthetic systems. Generally, in plants with less salinity tolerance, there is more damage to the photosynthetic systems under salinity stress [3]. After salt treatment, various parts of the photosynthetic systems of C. crassum showed different expression profiles. We found that components such as RPS4e (Unigene1713 and Unigene41867) and LHCA4 (CL12376.Contig1 and CL13565.Contig1) of the "LHC (light-harvesting chlorophyll protein complex)", which is responsible for trapping and transporting light energy to "photosystem I" and "photosystem II" [61], were downregulated, whereas components of "photosystem I" and "photosystem II" (such as psbA  [62], were downregulated. These changes may have been caused by the lack of water under osmotic stress. However, components of "cytochrome b6f complex", "photosynthetic electron transport", and "F-type ATPase" (such as petA [Unigene26782], FAD [CL14326.Contig2 and CL14326.Contig1], and F-type H + -transporting ATPase subunit gamma [Unigene43756]), which are responsible for mediating electron transport between PSII and PSI and converting the redox energy into part of the proton gradient used for ATP formation [62], were upregulated, as were key enzymes of the "carbon fixation in photosynthetic organisms" pathway (such as FBA [Unigene71038, CL13604.Contig2, Unigene69763, and Unigene35247], PK [Unigene71285 and Uni-gene62599], and RPIA [Unigene69514]) [63]. Thus, because of the difficulty in water uptake under salt shock, some components responsible for oxidizing H 2 O in the photosystem of C. crassum were downregulated, whereas other components were stimulated by the stress. At 24 h, all of these genes showed a tendency to recover to control levels. The previously mentioned ROS-scavenging enzymes provide protection to photosynthetic systems [55] and aid in the recovery of the photosynthetic systems of C. crassum. The recovery of the photosynthetic systems provides nutrition for the survival and growth of C. crassum under salinity stress.

Explanation of the morphological changes in C. crassum based on transcriptome analysis
The growth condition under salinity stress is considered a main norm for judging the salinity tolerance of a plant, with a better growth condition indicating better salinity tolerance [39]. The growth rate of C. crassum (Fig 1A) after salt treatment showed no obvious decrease, indicating that this species exhibits remarkable salinity tolerance. Incrassation and succulence of salt-treated leaves that are often observed in halophytes such as Salicornia europaea [64] and Populus euphratica [10]. The leaves of C3 plants such as Arabidopsis have also been reported become slightly succulent after salt treatment [39]. As shown in Fig 1B, the leaves of C. crassum became obviously succulent after salt treatment, which mainly due to the expansion of each cell. It is common for plants to show some morphological changes under salinity stress, which will improve the adaptation of plants to stress. The morphological changes displayed in C. crassum are rarely observed in chrysanthemums, and it would be very novel and interesting to elucidate how these changes were initiated. Here, we discuss how the morphological changes observed in C. crassum are initiated based on transcriptome analysis. The genes contributing to morphological changes are listed in S5 Table. As previously discussed, these morphological changes are mainly regulated by hormones such as auxins, cytokinins, and brassinosteroids (genes such as TIR1, AUX, ARF, A-ARR, TCH4, CYCD3 and SAUR are hypothesized to be responsible). Additionally, the appropriate functioning of photosynthetic systems provides nutritional preconditions for the ideal growth of C. crassum (functional genes are SOD, MDHAR, CAT, LHCs, photosystem proteins, FBA, and RPIA). additionally, stress-induced hormones trigger certain components of metabolism, such as energy metabolism and transporters, to disrupt growth by depleting nutritional resources [65].
Plant cell expansion is often characterized as being the product of the opposing forces of intracellular turgor pressure (growth-promoting) and the turgor resistance of the cell wall (growth-inhibiting) [66]. As previously mentioned, we found that osmolytes accumulated to high levels and water absorbability and retention were improved in C. crassum. In this part, the mentioned genes are TIPs, NIPs, PIPs, TPS, TPP, mt1D, MIP, BADH, CHDH and P5CS. Therefore, as time went by, the "intracellular turgor pressure" became stronger.
We next consider the cell wall to examine "turgor-resisting forces". The plant cell wall, and especially the cell wall proteins, undergoes adjustments under stress, as reported previously [67,68]. The cell wall proteins exhibit functions involved in the morphogenesis and signal transduction of plants, and changes in the cell wall are a sign of cell differentiation [11]. In our transcriptomes, many types of cell wall proteins were found to be upregulated after salt treatment, such as HRGPs (Unigene69380, CL1283.Contig2, Unigene68945, and Unigene70662), which have an expansion-like domain and show a strong correlation with final cell length [69]. These expansions cause the cellulose microfibrils to slide apart and are considered to be the primary determinants of wall elongation [70], together with a group of enzymes known as XETs (CL10820.Contig2, CL10820.Contig3, and Unigene33390) [71], which were also upregulated after salt treatment. PRPs (Unigene69599, Unigene70083, Uni-gene70589, and Unigene70081) were also upregulated after salt treatment, which are extension proteins involved in the response to stress [72]. The adjustment of cell wall proteins of C. crassum under salinity stress showed a signal interaction and served as a marker for cell morphological changes. Furthermore, the expression of many extension proteins weakened "turgor-resisting forces", making it easier for mesophyll cells to expand. Additionally, the cell wall proteins called AGPs (Unigene69390, Unigene68964, Unigene70408, and Uni-gene72883) were also upregulated after salt treatment; these proteins do not have structural function, but instead act as cell positional markers or as messengers in cell-cell interactions [73] Other upregulated genes included GRPs (Unigene68520, Unigene39763, and CL6740. Contig11), which play important roles in the development of vascular tissues and response to stresses [74].
With the expansion of cells, the membranes of the cells must also expand. In this context, the LTPs (Unigene70157, CL1283.Contig3, CL16785.Contig3, CL16785.Contig4, CL14640. Contig2, Unigene41619, and CL9481.Contig2) were also significantly upregulated and highly expressed after salt treatment, accounting for as much as 4% of the total soluble proteins of higher plants [75]. LTPs are thought to participate in membrane biogenesis, regulation of the intracellular fatty acid pools, cutin synthesis and responses to various stresses [76]. The abundance of this type of protein would benefit C. crassum in maintaining membrane stabilization and facilitating fatty acid-related pathways that would be helpful in stress tolerance.

Conclusions
Using the Illumina platform, we surveyed three transcriptomes of C. crassum leaves. We aimed to list the functional genes and pathways that contributed to salinity tolerance. Combined with the morphological changes observed under salinity stress, the results of the present study allowed us to infer and list the strategies adopted by C. crassum under salinity stress (Fig 6). Briefly, salinity stress first caused osmotic stress, which led to a lack of water. To maintain the water status of cells and cell functions effectively in C. crassum, large amounts of highly hydrophilic proteins, such as LEAs and DHNs were expressed, along with a reduction of water waste and regulation of the expression of aquaporins to sense turgor pressure and regulate water transportation. C. crassum also synthesized osmolytes to improve osmotic pressure. Osmotic stress was overcome by these strategies. Salinity also caused toxicity and perturbed the normal function of cell reactions. This toxicity mainly consisted of the accumulation of Na + and the production of ROS. Accordingly, C. crassum expressed ion transporters to exclude and compartmentalize Na + and expressed ROS-scavenging enzymes to remove ROS. The damage due to salinity was controlled by these strategies. The applied salinity stress triggered some transcription factors and plant hormones (especially some stress-induced hormones), including abscisic acid, jasmonic acid, and salicylic acid. C. crassum improved the activity of energy metabolism and transporters to supply energy and materials to other reactions to endure salinity stress. Additionally, C. crassum activated cleaners, such as proteasomes and ubiquitin, to remove unwanted mRNA and proteins. As time passed, the plants exhibited morphological changes to adapt to salinity stress, and these changes were initiated during an early stage of the stress treatment. The salinity treatment initiated changes in some transcription factors and plant hormones, including growth-related hormones such as auxin, cytokinins, and brassinosteroids. C. crassum expressed many types of cell wall proteins to expand the cell wall, and osmolyte accumulation and water storage increased turgor pressure. Therefore, the mesophyll cells of C. crassum expanded over time under salinity stress. To maintain normal growth, the photosynthesis systems of C. crassum were protected and recovered well under salinity stress. Through these strategies, C. crassum adapted to salinity and grew well.