Transcriptome and Small RNA Sequencing Reveals the Basis of Response to Salinity, Alkalinity and Hypertonia in Quinoa (Chenopodium quinoa Willd.)

Quinoa (Chenopodium quinoa Willd.) is a dicotyledonous cereal that is rich in nutrients. This important crop has been shown to have significant tolerance to abiotic stresses such as salinization and drought. Understanding the underlying mechanism of stress response in quinoa would be a significant advantage for breeding crops with stress tolerance. Here, we treated the low-altitude quinoa cultivar CM499 with either NaCl (200 mM), Na2CO3/NaHCO3 (100 mM, pH 9.0) or PEG6000 (10%) to induce salinity, alkalinity and hypertonia, respectively, and analyzed the subsequent expression of genes and small RNAs via high-throughput sequencing. A list of known/novel genes were identified in quinoa, and the ones responding to different stresses were selected. The known/novel quinoa miRNAs were also identified, and the target genes of the stress response ones were predicted. Both the differently expressed genes and the targets of differently expressed miRNAs were found to be enriched for reactive oxygen species homeostasis, hormone signaling, cell wall synthesis, transcription factors and some other factors. Furthermore, we detected changes in reactive oxygen species accumulation, hormone (auxin and ethylene) responses and hemicellulose synthesis in quinoa seedlings treated with stresses, indicating their important roles in the response to saline, alkaline or hyperosmotic stresses in quinoa. Thus, our work provides useful information for understanding the mechanism of abiotic stress responses in quinoa, which would provide clues for improving breeding for quinoa and other crops.


Introduction
Quinoa (Chenopodium quinoa Willd.), a dicotyledonous chenopodiaceae crop, is called the "Mother of Grain", as its grains contain balanced nutrients, including protein with all essential amino acids, resistant starch, lipids, mineral elements, vitamins, dietary fiber and many functional components such as saponin and polyphenol. Quinoa originates from the Andes, which have poor soil quality and a harsh climate, and thus possesses a higher tolerance to abiotic stresses, such as drought, cold, UV irradiance and salinization, than most crops grown on plains [1,2]. As a facultative halophyte, some quinoa species can tolerate high salinity levels up to about 400 mM NaCl (as found in sea water) [3]. As a drought-tolerant crop, quinoa can grow and produce seed grains in the semi-desert conditions of Chile with less than 200 mm of annual rainfall [4]. Thus, quinoa is not only a great cereal for optimizing the utilization of saline-alkali farmland, but also a good model for improving stress tolerance breeding in other crops.
Recently, salinity, alkalinity and drought have hampered agricultural production all around the world. They sometimes occur separately, but more often in combination [5]. Drought stress harms plants through increasing osmotic pressure. Saline stress additionally causes ion toxicity, arising from high levels of Na + and Cl − [6]. In addition to these toxicities, alkali stress also induces toxic effects due to high pH and CO 3 2− levels [7]. Thus, plants respond to these three types of stresses through some shared and some specific mechanisms. Plants usually respond to drought through closing the stomata, wilting and other behaviors, and increasing cell osmotic pressure through accumulating carbohydrates [8]. When suffering from saline stress, plants also remove or store harmful ions using ion pumps [9]. While in a high-pH soil environment, proton pumps in plant root cells and synthesis of organic acids play important roles [10,11]. In addition, many abiotic stresses disturb the normal function of mitochondria and chloroplasts and produce reactive oxygen species (ROS), which can be purified via the oxidation cleaning system or neutralized via reducing substances [12]. Furthermore, many plant hormones, such as ABA (abscisic acid), IAA (indole-3-acetic acid) and so on, are involved in the regulation of these stress responses [13]. As important post-transcriptional regulators, miRNA and other non-coding RNA have also been found to participate in the regulation of these stress responses [14].
Compared with other plants, quinoa displays a higher tolerance to many abiotic stresses. When water is limited, quinoa grows deeper roots, rapidly closes its stomata and produces ROS scavengers, antioxidants (e.g., ornithine and raffinose) and osmoprotectants (e.g., soluble sugars and proline) [4]. Expression analysis revealed that ABA signaling, heatshock proteins (HSPs) and the biosynthesis of flavonoids are involved in drought response in quinoa [4,15]. When suffering from saline stress, quinoa maintains K + and excretes Na + through various strategies, including Na + sequestration in leaf vacuoles via Na + /H + exchanger 1 (NHX1) [16], Na + elimination in roots via the salt overly sensitive 1 (SOS1) pathway [17] and K + retention via high-active H-ATPase [18]. Similar to drought stress, salt stress in quinoa also leads to stomatal closure, decreased photosynthetic efficiency and ROS accumulation [3], which can be relieved through increasing pore density and up-regulating the activity of ROS scavengers such as superoxide dismutase (SOD), peroxidase (POD), ascorbate peroxidase (APX), catalase (CAT) and antioxidants such as rutin [19][20][21]. In addition, expression analyses have indicated that ABA signaling and betaine metabolism may play important roles in the salt stress response in quinoa [22,23]. On the other hand, although alkalization is widespread and often combined with salinization, there is limited information on the alkaline stress response in quinoa. Moreover, thus far, no comparative analysis of responses to salt, alkali and drought in quinoa have been performed.
Here, using high-throughput sequencing, we aimed to analyze the similarities and differences of the molecular responses to salt, alkali and drought at the transcriptional level and the small-RNA-mediated post-transcriptional level in quinoa. Our work will provide a better understanding of the response to alkali, as well as the specific and common mechanisms of saline, alkaline and hyperosmotic stress responses in quinoa, which can provide guidance for the abiotic stress tolerance breeding of quinoa and other crops.

Response of Quinoa Strain CM499 to Saline, Alkaline and Hyperosmotic Stresses
As the seedling stage is one of the most sensitive periods to salt and drought in quinoa [4], we transferred 3-day-old seedlings of CM499, a lowland quinoa cultivar suitable for growth in East China ( Figure S1), into 1/4 hoagland solution alone (pH 6.5) or with varying concentrations of either NaCl (100 mM and 200 mM), NaHCO 3 /Na 2 CO 3 (100 mM, pH 8.0 and pH 9.0) or PEG6000 (5% and 10%) to induce saline, alkaline and hyperosmotic stresses, respectively ( Figure 1). Five days later, both saline and alkaline stresses led to the abnormal growth of seedlings. In 100 mM NaCl, a slight growth-promoting phenotype was observed, suggesting a certain level of salt tolerance in CM499, while at concentrations of 200 mM NaCl, both shoots and roots became significantly shorter and the fresh weight was decreased (Figure 1a,d,g). Similarly, the growth of the roots was inhibited and the fresh weight was decreased after treatment with alkaline solution, starting at pH 8.0 (Figure 1b,e,h). Effects on shoot and root length were observed at pH 9.0. Seedlings treated with hyperosmotic stress also exhibited phenotypic changes, including a shorter shoot and a lower plant fresh weight, especially under the treatment of 10% PEG. But, noteworthily, the length of roots was not affected so significantly as the ones under saline or alkaline stresses (Figure 1c,f,i). Thus, we inferred that all three stresses might lead to similar physiological responses in seedlings, such as decreased shoot length and fresh weight. But, the growth inhibition of roots might be noticeable under alkaline stress, while the suppression of shoots was obvious under hyperosmotic treatment.
(100 mM, pH 8.0 and pH 9.0) or PEG6000 (5% and 10%) to induce saline, alkaline and hyperosmotic stresses, respectively ( Figure 1). Five days later, both saline and alkaline stresses led to the abnormal growth of seedlings. In 100 mM NaCl, a slight growthpromoting phenotype was observed, suggesting a certain level of salt tolerance in CM499, while at concentrations of 200 mM NaCl, both shoots and roots became significantly shorter and the fresh weight was decreased (Figure 1a,d,g). Similarly, the growth of the roots was inhibited and the fresh weight was decreased after treatment with alkaline solution, starting at pH8.0 (Figure 1b,e,h). Effects on shoot and root length were observed at pH 9.0. Seedlings treated with hyperosmotic stress also exhibited phenotypic changes, including a shorter shoot and a lower plant fresh weight, especially under the treatment of 10% PEG. But, noteworthily, the length of roots was not affected so significantly as the ones under saline or alkaline stresses (Figure 1c,f,i). Thus, we inferred that all three stresses might lead to similar physiological responses in seedlings, such as decreased shoot length and fresh weight. But, the growth inhibition of roots might be noticeable under alkaline stress, while the suppression of shoots was obvious under hyperosmotic treatment. Physiological response of quinoa strain CM499 to saline, alkaline and hyperosmotic stresses. (a-c) Growth phenotype of CM499 seedlings under saline, alkaline and hyperosmotic stresses, respectively. Scale bars are 1 cm. (d-f) Statistical analysis of shoot length and root length under saline, alkaline and hyperosmotic stresses, respectively. (g-i) Statistical analysis of fresh weight per plant under saline, alkaline and hyperosmotic stresses, respectively. The statistical difference between stress treatment groups and control group was analyzed using Student's T-Test. * means p < 0.05 and ** means p < 0.01.

Transcriptome Sequencing of CM499 Seedlings
As the differential phenotype might be the result of early response to treatment of stresses, in this work, we focused on the early transcription response of quinoa seedlings in stress treatment. To compare the molecular responses to saline, alkaline and hyperosmotic stresses in quinoa, we extracted total RNA from 3-day-old CM499 Figure 1. Physiological response of quinoa strain CM499 to saline, alkaline and hyperosmotic stresses. (a-c) Growth phenotype of CM499 seedlings under saline, alkaline and hyperosmotic stresses, respectively. Scale bars are 1 cm. (d-f) Statistical analysis of shoot length and root length under saline, alkaline and hyperosmotic stresses, respectively. (g-i) Statistical analysis of fresh weight per plant under saline, alkaline and hyperosmotic stresses, respectively. The statistical difference between stress treatment groups and control group was analyzed using Student's T-Test. * means p < 0.05 and ** means p < 0.01.

Transcriptome Sequencing of CM499 Seedlings
As the differential phenotype might be the result of early response to treatment of stresses, in this work, we focused on the early transcription response of quinoa seedlings in stress treatment. To compare the molecular responses to saline, alkaline and hyperosmotic stresses in quinoa, we extracted total RNA from 3-day-old CM499 seedlings treated with either 200 mM NaCl, 100 mM NaHCO 3 /Na 2 CO 3 (pH 9.0), 10% PEG or normal condition (the CK group) for 6 h for high-throughput RNA sequencing (SRA Accession Number: PRJNA994889). There were 42.08-43.58 million clean reads detected in twelve samples (three replicates for the control and three stress treatment groups, respectively), of which 91.38-93.06% were mapped to quinoa genome and 79.63-84.71% to quinoa genes (Table S1). Among these genes, there were 49,138 known genes, of which 11,188 were detected from more than one transcript, and 5485 novel genes (Table S2). The lengths of the identified quinoa genes ranged from 201 bp to 23,672 bp, of which the transcripts were 195 bp to 121,276 bp (but the majority of transcript lengths were less than 10,000 bp). In addition, there were 37,938 genes, all belonging to known genes, with existing functional annotations.

Differential Expression Analysis of Quinoa Genes
Next, we identified the significantly differentially expressed genes (DEGs) whose expression changed significantly after stress treatment compared with the untreated group (the changes were over two-fold and both the p-values and the q-values should be under 0.05, Table S3). In total, 8624 genes had significantly altered expression across any of the treatments, among which 5165 genes were up-regulated and 3459 genes were downregulated ( Figure 2). Some of these genes were selected randomly for qRT-PCR and mostly displayed similar trends as the transcriptome, confirming the reliability of the data ( Figure S2). Further, we divided these DEGs into seven groups according to the stress they responded to. Over half of the identified DEGs were only found in salt stress, including 2859 genes with increased expression and 2633 genes with decreased expression. On the other hand, 495 DEGs (388 up and 107 down) were only found following PEG treatment and 261 DEGS (188 up and 73 down) were only found following alkaline treatment. There were 1527 DEGs that responded to two stresses, including 702 genes which were significantly altered in both saline and hyperosmotic stress, 768 genes found in both salt and alkaline stress and 57 in both alkaline and hyperosmotic stress. Additionally, 849 genes exhibited significantly changed expression when treated with all three stresses. Thus, the CM499 seedlings seemed to be more sensitive to 200 mM NaCl than 100 mM NaHCO 3 /Na 2 CO 3 (pH 9.0) or 10% PEG. PEG or normal condition (the CK group) for 6 h for high-throughput RNA sequencing (SRA Accession Number: PRJNA994889). There were 42.08-43.58 million clean reads detected in twelve samples (three replicates for the control and three stress treatment groups, respectively), of which 91.38%-93.06% were mapped to quinoa genome and 79.63%-84.71% to quinoa genes (Table S1). Among these genes, there were 49,138 known genes, of which 11,188 were detected from more than one transcript, and 5485 novel genes (Table S2). The lengths of the identified quinoa genes ranged from 201 bp to 23,672 bp, of which the transcripts were 195 bp to 121,276 bp (but the majority of transcript lengths were less than 10,000 bp). In addition, there were 37,938 genes, all belonging to known genes, with existing functional annotations.

Differential Expression Analysis of Quinoa Genes
Next, we identified the significantly differentially expressed genes (DEGs) whose expression changed significantly after stress treatment compared with the untreated group (the changes were over two-fold and both the p-values and the q-values should be under 0.05, Table S3). In total, 8624 genes had significantly altered expression across any of the treatments, among which 5165 genes were up-regulated and 3459 genes were down-regulated ( Figure 2). Some of these genes were selected randomly for qRT-PCR and mostly displayed similar trends as the transcriptome, confirming the reliability of the data ( Figure S2). Further, we divided these DEGs into seven groups according to the stress they responded to. Over half of the identified DEGs were only found in salt stress, including 2859 genes with increased expression and 2633 genes with decreased expression. On the other hand, 495 DEGs (388 up and 107 down) were only found following PEG treatment and 261 DEGS (188 up and 73 down) were only found following alkaline treatment. There were 1527 DEGs that responded to two stresses, including 702 genes which were significantly altered in both saline and hyperosmotic stress, 768 genes found in both salt and alkaline stress and 57 in both alkaline and hyperosmotic stress. Additionally, 849 genes exhibited significantly changed expression when treated with all three stresses. Thus, the CM499 seedlings seemed to be more sensitive to 200 mM NaCl than 100 mM NaHCO3/Na2CO3 (pH 9.0) or 10% PEG.

The Acquisition of Small RNA Sequence
Small RNAs, e.g., miRNAs, are important post-transcriptional regulators in both animals and plants but have not been explored much in quinoa. Using the same samples as used for transcriptome sequencing (three repeats for normal, saline, alkaline and hyperosmotic stresses, respectively), we analyzed small RNAs expressed in quinoa through sequencing the limited length components of the total RNA (18-30 nt). In total, 23.2-24.7 million raw reads were obtained from 12 samples, consisting of 21.9-24.0 million clean reads (SRA Accession Number: PRJNA994977). In the clean reads, the lengths ranged from 18 to 30 nt, with the majority being 24 nt (Table S5; Figure S4). Removing intergenic, intron, exon, unmapped sequences and other small RNA, such as tRNA and rRNA, 0.23-0.55 million reads mapping to miRNA were obtained (Table S5). After aligning the sequences in miRBase, there were 148 miRNAs matched to known miRNA families (Table S6), some of which are conserved across many species, e.g., miR166a-3p, miR160a-5p, miR396b-5p, while others were only identified in limited species, e.g., miR6173, miR166e-3p, miR156g ( Figure 4; Table S7). In addition, there were 87 novel miRNAs identified by potential stem-loop structures in the precursors (Table S6).

MiRNA Responses to Different Stresses in Quinoa
To understand the role of miRNAs in quinoa stress responses, differentially expressed miRNAs (DEMs) under saline, alkaline or hyperosmotic stresses were identified as follows: the absolute value of log 2 (stress treatment/CK) should be over 0.5-fold and both the p-values and the q-values should be under 0.05 (Table S8). There were, in total, 114 miRNAs with altered expression under any stress treatment, of which 46 were up-regulated and 68 were down-regulated ( Figure 5). Similar to the transcriptomic analysis, we also divided the DEMs into seven groups according to stress types. There were 21 (9 up and 12 down), 22 (8 up and 14 down) and 18 (9 up and 9 down) miRNAs that only responded to saline, alkaline or hyperosmotic stresses, respectively. Additionally, 20 (9 up and 11 down), 8 (3 up and 5 down) and 12 (3 up and 9 down) miRNAs displayed the same expression change trends when treated with saline and alkaline, saline and hyperosmotic or alkaline and hyperosmotic stresses, respectively. In addition, there were 13 miRNAs which responded to all three stresses similarly, including 5 up-regulated and 8 down-regulated ones. Noteworthy, some miRNAs were found to have altered expression under two stresses but in opposing directions (e.g., miR160b_1 was up-regulated under alkaline stress but down-regulated under hyperosmotic treatment), indicating the complex roles of miRNAs in responding to different stresses.    Table S6.  Table S6.
alkaline, saline and hyperosmotic or alkaline and hyperosmotic stresses, respectively. In addition, there were 13 miRNAs which responded to all three stresses similarly, including 5 up-regulated and 8 down-regulated ones. Noteworthy, some miRNAs were found to have altered expression under two stresses but in opposing directions (e.g., miR160b_1 was up-regulated under alkaline stress but down-regulated under hyperosmotic treatment), indicating the complex roles of miRNAs in responding to different stresses.

Functional Analysis of DEMs in Quinoa
MiRNAs usually bind to mRNA with a homologous sequence to induce cleavage or inhibit translation of the mRNA. Thus, we predicted the candidate targets of the DEMs to analyze the possible function of these miRNAs in quinoa stress responses. The potential targets of the identified quinoa miRNAs were predicted using TargetFinder and psRobot. In total, 1666 candidate targets (Table S9), belonging to 651 transcripts from 370 genes, were predicted using both two tools (score < 3.5 for TargetFinder, score < 2.5 for psRobot).
Some GO clusters appeared in multiple stress groups, but the miRNA-target pairs were rather distinguishable across different groups. For example, auxin-related GO clusters appeared in the "alkaline only", "hyperosmotic only" and "saline & alkaline & hyperosmotic" groups, which matched to quinoa genes annotated as ARF (auxin response factor) 17-like and ARF18-like ( Figure 6; Table S9). However, these ARF-like genes were putatively regulated by miR160b_1 and miR160a-5p (Table S9), the former of which was induced by alkaline treatment but inhibited by hyperosmotic stress, while the latter was up-regulated by all three stresses ( Figure 5; Table S8). As another example, the target genes corresponding to the "lipid binding" GO cluster were homeobox-leucine zipper proteins which matched to miR166m_2 in "saline only", to miR166 and miR166a in "alkaline only" and to miR166e, miR166e-3p and miR166m_2 in "alkaline & hyperosmotic" stresses ( Figures 6 and S5; Table S9). For the "sulfate adenylyltransferase activity" GO cluster, the targets were identified as ATP sulfurylase 1, which matched to miR395a_5 in "saline & alkaline" and miR395b_3 in "saline & alkaline & hyperosmotic" stresses ( Figures 6 and S5; Table S9). For "glycine biosynthetic process", the targets were serine hydroxymethyltransferase (mitochondrial) that matched to miR172d-5p_4 in "saline only" and serine hydroxymethyltransferase 4 to miR6300 in "saline & alkaline" stresses ( Figures 6 and S5; Table S9). These particular examples demonstrated the complexity of the functions of miRNA in saline, alkaline and hyperosmotic stress responses in quinoa.
In addition to the above examples, there were also some other interesting targets matched to the DEMs. Many ethylene-responsive transcription factors, such as AP2-like, RAP2-7-like, ERF113-like and WIN1-like were candidate targets of miR172a_3, miR172e-3p_1, miR172c_2 and miR395b_3, which belonged to the "saline only", "saline & alkaline", "alkaline & hyperosmotic" and "saline & alkaline & hyperosmotic" groups, respectively (Table S9). Some ubiquitin-related factors, such as "U-box domain-containing protein 33like" and "probable ubiquitin-conjugating enzyme E2", were candidate targets of miR1436, miR399f_3, miR399a_6 and miR399b, which were DEMs of the "saline only", "alkaline only", "saline & alkaline", "alkaline & hyperosmotic" and "saline & alkaline & hyperosmotic" groups. Additionally, many squamosa promoter-binding-like protein genes were found as the targets of miR156e_1 and miR156_2, which responded to saline and alkaline stresses, respectively. A calcium uptake protein gene was predicted as the target of miR166h-3p_1 and miR166u, which were up-regulated in "hyperosmotic only" but down-regulated in "saline & alkaline". There were also many epigenetic factors that were found as targets of DEMs, e.g., histone-lysine N-methyltransferase SUVR4 for miR172-3P_1 in "saline & alkaline" stress, while SWI/SNF-related matrix-associated actin-dependent regulator and histone deacetylase HDT1-like were found to be targets for miR169h-3p_1, which was differentially expressed in "saline & alkaline & hyperosmotic". In addition, kinesin-like proteins, which are related to microtubules, were targets of miR172d-5p in "saline only" and of miR396a-3p_4 in "alkaline only". Int. J. Mol. Sci. 2023, 24, x FOR PEER REVIEW 10 of 22 Figure 6. GO enrichment analysis of miRNA targets in quinoa. Figure 6. GO enrichment analysis of miRNA targets in quinoa.

ROS Accumulation in Quinoa Seedlings under Stress Treatment
In plants, ROS are usually produced by the disruption of oxidation-reduction reactions, such as photosynthesis and respiration, induced by abiotic stresses. When under stress treatment, ROS are not only a harmful byproduct that needs to be cleaned, but also an important signal to regulate downstream responses in plants [24]. As our GO analysis displayed, "oxidoreductase activity" was enriched in many DEG groups, including many ROS-related DEGs, such as respiratory burst oxidase homolog (RBOH), peroxidase and catalase (Table S3), which indicated a possible role of ROS homeostasis in the quinoa stress response. Thus, we visualized the levels of superoxide and hydrogen peroxide in 5-day-old quinoa seedlings following 1-day stress treatment via staining with NBT (nitrotetrazolium blue chloride) and DAB (diaminobenzidine) stain, respectively. For seedlings under normal cultivation, there was obvious accumulation of superoxide (stood by NBT stain) in the cotyledons, young euphylla, roots and especially the root tips (Figure 7a). When treated with saline stress, the NBT stain was lighter in the cotyledons and roots compared with the control group. The biggest difference was seen in seedlings treated with alkaline stress, in which the root tips were nearly without any stain signal. The seedlings under hyperosmotic stress exhibited a similar NBT stain as the control group, although their root tips showed a slightly lighter coloration. According to the results of DAB staining (Figure 7b), hydrogen peroxide was slightly accumulated in the cotyledons and roots of quinoa seedlings under normal condition and following hyperosmotic stress. However, after being treated with saline stress, the level of hydrogen peroxide increased in the root tips. And even a more greatly increased hydrogen peroxide level was observed in seedlings treated with alkaline stress, which exhibited an almost opposite trend to superoxide. Thus, we concluded that the ROS metabolism was involved in the responses of quinoa to salt and especially alkali. Notably, superoxide and hydrogen peroxide appeared to play different roles in this process, as there was a down-regulation of superoxide but an up-regulation of hydrogen peroxide.

Auxin and Ethylene Were Involved in Quinoa Stress Responses
As described in Tables S3 and S9, many ethylene-responsive transcription factors as well as auxin-related genes (e.g., ARFs) or their miRNAs displayed altered expression when treated with saline, alkaline or hyperosmotic stresses in quinoa seedlings. To further analyze the relationship between these hormones and stress responses, a series of co-treatment assays were designed (Figure 8). When disturbing endogenous auxin signaling with NAA (naphthylacetic acid) or PEO-IAA (2-(1H-Indol-3-yl)-4-oxo-4phenyl-butyric acid, an IAA antagonist), the seedlings in the normal culture solution grew more slowly than those with no additives, with shorter shoot length, root length

Auxin and Ethylene Were Involved in Quinoa Stress Responses
As described in Tables S3 and S9, many ethylene-responsive transcription factors as well as auxin-related genes (e.g., ARFs) or their miRNAs displayed altered expression when treated with saline, alkaline or hyperosmotic stresses in quinoa seedlings. To further analyze the relationship between these hormones and stress responses, a series of co-treatment assays were designed (Figure 8). When disturbing endogenous auxin signaling with NAA (naphthylacetic acid) or PEO-IAA (2-(1H-Indol-3-yl)-4-oxo-4-phenyl-butyric acid, an IAA antagonist), the seedlings in the normal culture solution grew more slowly than those with no additives, with shorter shoot length, root length and lower fresh weight. With the addition of NAA, seedlings treated with any of the stresses displayed a significant decrease in shoot length or fresh weight compared with the NAA-treated no-stress group, exhibiting similar alteration trends as those with no additives. After treatment of PEO-IAA, unchanged or increased root length and fresh weight were observed compared to the no-stress, PEO-IAA-treated controls, especially for seedlings treated with saline stress, which was different from those with no additives. In addition, when disturbing endogenous ethylene signaling with ACC (1-aminocyclopropane-1-carboxylic acid, an ethylene precursor), the decrease in root length and fresh weight under hyperosmotic stress was partly alleviated. Thus, our results indicated potential roles of auxin and ethylene signaling pathways in the response to saline and alkaline stresses in quinoa.

Metabolism of Hemicellulose in Quinoa Stress Responses
According to the GO analysis of the DEGs, cell wall synthesis, which is associated with xyloglucosyl transferase activity (xyloglucan endotransglucosylase/hydrolase protein, XTH), was enriched when quinoa seedlings were treated with saline, alkaline and hyperosmotic stresses (Figure 3; Table S3). As the substrate for XTH, xyloglucan is the main type of hemicellulose, of which the content and fiber structure are closely related to cell growth. To clarify the effect of stresses on cell wall synthesis in quinoa, the levels of hemicellulose in seedlings with or without stresses were quantified. The absolute concentration of reducing sugars derived from hemicellulose in quinoa seedlings can be calculated from a standard curve created using D-xylose, and the

Metabolism of Hemicellulose in Quinoa Stress Responses
According to the GO analysis of the DEGs, cell wall synthesis, which is associated with xyloglucosyl transferase activity (xyloglucan endotransglucosylase/hydrolase protein, XTH), was enriched when quinoa seedlings were treated with saline, alkaline and hyperosmotic stresses (Figure 3; Table S3). As the substrate for XTH, xyloglucan is the main type of hemicellulose, of which the content and fiber structure are closely related to cell growth. To clarify the effect of stresses on cell wall synthesis in quinoa, the levels of hemi-cellulose in seedlings with or without stresses were quantified. The absolute concentration of reducing sugars derived from hemicellulose in quinoa seedlings can be calculated from a standard curve created using D-xylose, and the changed hemicellulose levels between seedlings under normal and stress conditions can be analyzed (Figure 9). Compared with the untreated CK group, the contents of hemicellulose were decreased significantly by all three stresses, with greater effects seen with the saline and alkaline stresses. Therefore, the metabolism of hemicellulose might play a role in the quinoa stress response; however, whether it influences stress tolerance is still unclear.

ROS Metabolism in Abiotic Stress Response of Quinoa
ROS is not only a harmful by-product during stress treatment but also an important molecular messenger in regulating stress response [12]. In our data, the GO cluster "oxidoreductase activity", which was associated with ROS homeostasis, appeared in all groups and was enriched in the "alkaline only", "hyperosmotic only" and "alkaline & hyperosmotic" groups ( Figures 3 and S3). As the function description exhibited, many RBOH members were induced in the "saline only", "hyperosmotic only" or "saline & alkaline & hyperosmotic" groups (Table S3), suggesting higher production of ROS [63]. Furthermore, many ROS scavengers, e.g., peroxidase and catalase, were also found in the DEG groups "saline only", "hyperosmotic only", "saline & alkaline" and "saline & alkaline &hyperosmotic" (Table S3). Thus, we inferred that ROS homeostasis, modulated by oxidoreductases such as RBOHs, peroxidase and catalase, played an important role in quinoa stress responses. In agreement with this hypothesis, NBT and DAB staining displayed changes in ROS accumulation in quinoa seedlings, especially those treated with alkaline and saline stresses (Figure 7). Interestingly, unlike some other crops [64], ROS levels did not increase in the root tips of stress-treated quinoa and even decreased under alkaline treatment, which might indicate the stronger activity of antioxidant enzymes (e.g., peroxidase) than ROS producers (e.g., RBOH). Further, lower superoxide levels appeared with a stronger H 2 O 2 accumulation, perhaps indicating a relatively lower activity of catalase than peroxidase. Thus, quinoa might have a special ROS homeostasis in abiotic stress responses.

Cytoskeleton Dynamics Might Be Involved in Stress Responses in Quinoa
The cytoskeleton is a dynamic component of plant cells mainly composed of microtubules (MTs) and actin filaments (AFs), the former of which regulate the synthesis of cellulose to microfibrils for cell wall establishment and can effect stomatal movement in response to ABA signaling, while the latter play roles in delivering materials, expanding the cell wall and selectively enabling and/or inhibiting endocytosis and can regulate cell morphogenesis [65][66][67][68]. In our data, some GO clusters such as "microtubule motor activity", "microtubule-based movement", "microtubule-based process" and "microtubule binding" were enriched only under saline treatment (Figure 3). These clusters matched to 70 DEGs, mainly including kinesin-like proteins, tubulins and microtubule-associated proteins, surprisingly, all of which were down-regulated by saline stress (Table S3), perhaps indicating a functional inhibition of MTs by saline stress in quinoa. Kinesin is a motor protein which binds to MTs to transport intracellular materials and has been found as a regulator of drought stress response in Arabidopsis [69], while in our data the kinesin-like protein KIN-6 was predicted as the target of miR172d-5p, a DEM in the "saline only" group (Table S3). Some AF-related clusters, e.g., "actin filament binding", "actin filament-based process" and "actin cytoskeleton organization", also appeared in the DEG groups "saline & alkaline", "alkaline & hyperosmotic", "saline & hyperosmotic", "saline & alkaline & hyperosmotic" and mainly in the "saline only" group (Figures 3 and S3). AF-related DEGs included protein-tyrosine-phosphatase MKP1-like, actin-depolymerizing factor, fimbrin-like, SDA1 homolog, myosin-like, F-actin-capping protein and glucan 1,3-beta-glucosidase A (Table S3), some of which have been found as regulators in stress response [70]. Thus, in quinoa, cytoskeleton-related functions were possibly modified through stress treatment, especially under saline stress.

Cell Wall Organization Was Associated with Stress Response in Quinoa
The plant cell wall is a protective barrier composed of cellulose, hemicellulose, pectin, structural proteins and lignin, which can determine cell size and shape through controlling cell expansion, and is associated with the response to drought and salt in plants [71]. In our work, we found changed levels of hemicellulose in quinoa seedlings under stress treatment (Figure 9), which might be associated with the modulation of the cell wall. There were many GO clusters, such as "xyloglucosyl transferase activity", "xyloglucan metabolic process" and "cell wall biogenesis", enriched in the DEG groups "saline & alkaline", "saline & hyperosmotic", "alkaline & hyperosmotic" and "saline & alkaline & hyperosmotic" (Figures 3 and S3). These clusters mainly contained DEGs characterized as xyloglucan endotransglucosylase/hydrolases (XTHs), nearly all of which were suppressed by the stresses (Table S3). XTHs can split and reconnect xyloglucan to cellulose microfibrils, mediating the remodification of the cellulose/hemicellulose framework in the cell wall [72]. As a family with dozens of members, XTHs have complex effects on plant stress response. Over-expression of the hot pepper CaXTH3 enhances tolerance to salt and drought in Arabidopsis [73], while AtXTH30 negatively affects the salt tolerance of Arabidopsis [74]. In summary, in quinoa, the expression of many XTHs as well as the synthesis of hemicellulose were influenced by saline, alkaline and hyperosmotic stresses, which might ultimately influence the cell wall assembly.

Plant Materials and Growing Conditions
Seeds of CM499 were germinated on wet filter paper at 22 • C for two days, then transferred into 1/4 hoagland for hydroponics. Three-day-old seedlings were transferred into 1/4 hoagland with no treatment or stress treatment for phenotypic observation, RNA extraction, or hemicellulose determination. The stress treatments included NaCl (100 and 200 mM), NaHCO 3 /Na 2 CO 3 (100 mM, pH 8.0 and pH 9.0) or PEG6000 (5% and 10%). For hormone treatment assays, the seedlings were transferred into 1/4 hoagland with no treatment, NaCl (200 mM), NaHCO 3 /Na 2 CO 3 (100 mM, pH 9.0) or PEG6000 (10% W/V), and each exposure condition was treated with either no additive, 5 µM NAA (CN7541,

Growth Phenotype Observation
For detection of quinoa's physiological response to the various stresses, five plants per treatment were cultured in stress solutions for 5 days. Subsequently, the shoot length, root length and fresh weight of each plant was measured. For hormone treatment assays, after treatment with both the hormone and stress for 5 days, the shoot length, root length and fresh weight of five plants for every group were measured.

Transcriptome and Small RNA Sequencing
Seedlings (20 plants per group) treated for 6 h under the various stresses (C1, C2, C3 for normal condition as the control groups; S1, S2, S3 for NaCl treatment as the saline groups; A1, A2, A3 for Na 2 CO 3 /NaHCO 3 treatment as the alkaline groups; P1, P2, P3 for PEG treatment as the hyperosmotic groups) were collected for total RNA extraction (TRIzol method, BGI, Shenzhen, China) and construction of the respective libraries for transcriptome (6GB per sample) or small RNA sequencing (20 MB per sample, PE150, Illumina HiSeq 2000 device, BGI, Shenzhen, China). Raw reads were initially edited with the Fastx-toolkit pipeline to remove adapter sequences, low quality reads, repetitive reads and, for small RNA, reads longer than 30 nucleotides or shorter than 18 nucleotides. The resulting reads were aligned against the quinoa expressed sequence tag database (https://ftp.ncbi.nlm.nih.gov/ genomes/all/GCF/001/683/475/GCF_001683475.1_ASM168347v1, accessed on 5 March 2021), the genomic sequence database and the Rfam database for removal of rRNA, tRNA, snRNA and snoRNA sequences and identification of transcripts. For small RNA sequencing, the reads were compared using BLAST analysis against the miRbase (www.mirbase.org, accessed on 6 March 2021) to identify known miRNAs. Potential novel miRNAs were identified from the remaining reads via miRA (https://github.com/mhuttner/miRA, accessed on 7 March 2021), which could screen miRNA precursors through exploring the characteristic hairpin structure [75]. Target prediction of the identified quinoa miRNAs was made via both TargetFinder and psRobot [76,77], with the default parameters (psRobot: -gl 17 -p 8 -gn 1, TargetFinder: -c 4). The sequence data of transcriptome and small RNA library was deposited in NCBI Sequence Read Archive (SRA, www.ncbi.nlm.nih.gov/sra, accessed on 15 July 2023), with accession numbers PRJNA994889 and PRJNA994977, respectively.

qRT-PCR Analysis
The total RNA from seedlings treated under various stresses (same as the samples for sequencing) for 6 h was extracted (FastPure Plant Total RNA Isolation Kit, RC411-C1, Vazyme, Nanjing, China) and used for reverse transcription (EasyScript One-Step gDNA Removal and cDNA Synthesis SuperMix, AE311-02, Transgen, Beijing, China). An amount of 0.5 g fresh seedlings of each sample was ground in liquid nitrogen, and we extracted the total RNA according to the kit instruction. For each sample, 1 µg total RNA was reverse transcripted to cDNAs following the description of the kit instruction. qRT-PCR was performed via the SuperReal PreMix Plus (SYBR Green) kit (FP205, TIANGEN, Beijing, China) and a real-time PCR detection system (480II, Roche, Switzerland). Quinoa CqTUB-6 (LOC110711758) was used as a reference gene. Abundances were normalized using the 2 −∆∆Ct method. The primers used were listed in Table S11, and synthesized by Beijing Tsingke Biotech Co., Ltd. (Beijing, China).

Hemicellulose Determination
Three-day-old seedlings (20 plants per group) treated under the various stresses (same as the samples for sequencing) for 3 days were collected and dried to constant weight in a bake oven at 65 • C. The samples were then ground and filtered using a sifter of 30-50 mesh. The powder was then processed using a Hemicellulose Content Detection kit (BC4445, Solarbio, Beijing, China), following manufacturer guidelines. The hemicellulose was hydrolyzed to reducing sugars, and then reacted with 3,5-dinitrosalicylic acid (DNS) to form a colored product, which could be detected at 540 nm via a spectrophotometer (UV-1800, SHIMADZU, Kyoto, Japan). The relative contents were calculated via a standard curve of various concentrations of xylose and standardized according to the control group.

Conclusions
Quinoa is a high-nutrient crop with broad stress tolerance. Our work analyzed its response to salt, alkali and hypertonic stresses at the transcriptome and small RNA levels, and we identified a series of stress-responsive genes and miRNAs. The function analysis of DEGs and DEMs revealed factors or pathways related to these stress responses, including plant hormones, TFs, ROS homeostasis, cytoskeleton dynamics and cell wall assembly. We further found that the signaling of ethylene and auxin, the accumulation of superoxide and hydrogen peroxide and the metabolism of hemicellulose were closely associated with the abiotic stress response of quinoa seedlings. In summary, quinoa perhaps shares many common mechanisms in response to saline, alkaline and hyperosmotic stresses, such as xyloglucan metabolism, and possesses some unique pathways such as ROS homeostasis under alkaline stress or cytoskeleton dynamics under saline stress. Our results would afford useful clues for understanding the mechanism of stress response in quinoa.