Analysis of Metabolomic Changes in Lettuce Leaves under Low Nitrogen and Phosphorus Deﬁciencies Stresses

: Nitrogen and phosphorus limitation affect the growth, development, and productivity of lettuce, which exert a marked influence on metabolites. To compare the influence of low-nitrogen and low-phosphorus stresses on various metabolites of lettuce leaves, experiments were performed under three conditions of treatment—low-nitrogen stress, low-phosphorus stress, and normal samples. Metabolomic analyses were conducted based on ultra-performance liquid chromatography-quadrupole time-of-flight mass spectrometry. Principle components analysis yielded distinctive clustering information among the holistic samples; fold change analysis, t -test and orthogonal partial least square discriminant analysis were used for the selection of metabolic biomarkers. Ten pathways were selected which were significantly enriched by metabolic biomarkers. Metabolic biomarkers were screened by fold change (FC) value, p -value and variable importance in the projection (VIP) value, low-nitrogen and low-phosphorus stresses caused an increase in 16 metabolites (FC > 2, p -value < 0.05, VIP > 1) and a decrease in 26 metabolites (FC < 0.5, p -value < 0.05, VIP > 1). Outside of these, our results showed that inositol, p-hydroxybenzoic acid, stachyose, dinoseb, and 7, 8-dihydroxycoumarin increase in low-nitrogen stress samples. Low-phosphorus stress caused accumulation of citrate, isocitrate, l-5-oxoproline, succinate, and histamine, which may be considered potential metabolic biomarkers. The metabolites could be used to monitor the nitrogen and phosphorus status of lettuce and to guide appropriate fertilization regimens. 95–65% (0.5–7 65–40% B (7–9 min), 40% B (9–10 min), B (10–11.1 (11.1–16 min). The injection temperature was 4 °C. The flow rate was 0.3 mL/min. The injection volume was 10 μ L. The column temperature was maintained at 25 °C. Metabolites were measured using an UPLC-QTOF MS platform. Experiments were performed with a triple time-of-flight 5600+ mass spectrometer (AB SCIEX, Framingham, MA, USA) in both electrospray ionization (ESI) positive and negative ion modes. The mass spectrometer was equipped with an electrospray ionization (ESI) source. The ESI source parameters in ionization mode were as follows:


Introduction
Mineral nutrients, such as nitrogen and phosphorus, frequently limit plant growth [1]. Nitrogen limitation directly influences the biosynthesis of compounds, such as proteins, vitamins, phytohormones, co-enzymes, chlorophyll, and nucleic acids [2,3]. Phosphorus is a component of many cellular molecules and plays an essential role in structural maintenance [4]. A key concern for growers is to be able to verify plant nutrition status early while nutrition deficiency can still be remedied [5]; thus, the early diagnosis of nutrition stress plays a key role in the regulation of plant nutrition. Studies have focused on topics such as the nutritional deficit of plants through the differences of surface macro characteristics, such as height, length, and color [6]; the differences of micro characteristics such as leaf surface roughness and texture [7]; the differences of the internal microstructure such as stomata, sponges and palisade tissues; and the differences of macromolecular compounds such as chlorophyll and lutein [8][9][10].
The diagnosis of cancer is based on tissue, cell, and molecule level. The detection of metabolic biomarkers at the molecular level can realize early warning [11], which is of great significance to improve the life span and quality of life. On research of plant physiology and nutrition has indicated that the reaction of crop nutrition stress is a complex mechanism, which involves macro characteristics, micro characteristics, internal organization, macromolecules, and small molecule compounds [12]. To resist nutrition stress, plants use their metabolic capacity to produce a large variety of small molecular compounds and reorganize their metabolic networks to maintain essential metabolism [13]. Plant metabolomics is now widespread and valuable biotechnology, which has been used to explore the resistance mechanisms of plants in adverse environments [14]. Metabolomics has largely focused on the study of all small molecular compounds and the related dynamic changes found in or produced by organisms [15,16]. In this regard, several recent studies have explored the change of small molecular compounds underlying the response to mineral nutrient stress. Sung et al. [17] conducted a detailed broad-scale identification of metabolic responses of tomato leaves and roots to nitrogen, phosphorus, or potassium deficiency, with deficiency of any of these three minerals affecting energy production and amino acid metabolism. Cevallos-Cevallos et al. [18] found significant differences in oxo-butanedioic acid, arabitol, and neo-inositol between zinc-deficient and healthy orange trees. The metabolomics of rice leaves under conditions of nitrogen and phosphorus deficiencies were analyzed by Shen [19] and Watanabe et al. [20], the stress-resistant substances and amino acid substitution products increased under nitrogen deficiency, with sinapate, benzoate, and glucuronate related to phosphorus deficiency.
Lettuce forms an important component of the human diet and is typically low in calories and fat, and high in protein, dietary fiber, iron, calcium, and phytochemicals. Lettuce metabolomic studies have focused only on associations with light quality and intensity conditions [21], temperature stress [22], saline stress [23], cultivation conditions [24], and lettuce varieties [25]. Few studies have focused on the metabolite variation of lettuce and the identification of biomarkers under mineral nutrient stress, particularly low-nitrogen and low-phosphorus stresses. In this study, the metabolites of lettuce leaves were analyzed by an ultra-performance liquid chromatography-quadrupole time-of-flight mass spectrometer (UPLC-QTOF MS) platform, insight into lettuce responses to low-nitrogen and low-phosphorus stresses were obtained by discovering new information about changes in metabolite abundance, aiming to provide a basis for further diagnosis of lettuce early nutrition stress and guide appropriate fertilization regimens.

Plant Material and Growth Conditions
Lettuce seed (Lactuca sativa, Italian) was obtained from Woshu Seeds Co. Ltd., Nanjing, China, and sown in sponge blocks, then the investigated plants were transplanted in Rockwool cubes (50 × 50 × 50 mm) during the period of five true leaves and grown in a micro plant factory (Hangzhou Shuolian instrument CO., LTD, Zhejiang, China) at Jiangsu University, China. This was a closed cultivation room with an air conditioner to control the temperature and artificial lighting. Sample cultivation is shown in Figure 1 4 3− in the nutrient solution. The EC and pH of the nutrient solution were adjusted to 1.6 dS m −1 and 6.0, respectively. Air temperature, relative humidity, and CO 2 concentration were maintained at 20 ± 5 • C, 60-80%, and 400 ± 10 µmol mol −1 , respectively [26]. After two weeks, eighteen plants having similar growth from each treatment were harvested between 10:00 and 12:00 to minimize diurnal effects on metabolite concentration. The leaves were rinsed briefly in deionized water and immediately frozen in liquid nitrogen, and stored at −80 • C until biochemical analysis.
Agriculture 2020, 10, x FOR PEER REVIEW 3 of 15 similar growth from each treatment were harvested between 10:00 and 12:00 to minimize diurnal effects on metabolite concentration. The leaves were rinsed briefly in deionized water and immediately frozen in liquid nitrogen, and stored at −80 °C until biochemical analysis.

Metabolite Extraction
Three groups of lettuce leaves were ground into liquid nitrogen, 1000 μL methanol/acetonitrile/water (2:2:1, V/V/V) was added to each sample, then the samples were ultrasonicated twice at low temperature for 30 min after vortexing, and incubated at −20 °C for 1 h. Acetonitrile (1499230-935) and methanol (144282) were purchased from Merck (Darmstadt, Germany). Subsequently, the samples were centrifuged at 16,626 g at 4 °C for 15 min. Supernatants were collected, lyophilized, and stored at −80 °C. Five quality control (QC) samples were prepared by pooling aliquots of all samples. The pretreatment of the QC samples was the same as that of the study samples. The QC samples were evenly inserted in each set of the analysis running sequence to monitor the stability of the large-scale analysis [27].
Metabolites were measured using an UPLC-QTOF MS platform. Experiments were performed with a triple time-of-flight 5600+ mass spectrometer (AB SCIEX, Framingham, MA, USA) in both electrospray ionization (ESI) positive and negative ion modes. The mass spectrometer was equipped with an electrospray ionization (ESI) source. The ESI source parameters in ionization mode were as follows: source temperature, 600 °C; ion source pressure, 60 psig; curtain pressure, 30 psig; ion spray voltage floating (ISVF), ±5.5 KV; TOF MS scan m/z range, 60-1200 Da; TOF MS scan accumulation time, 0.15 s/spectra; production scan m/z range, 25-1200 Da; production scan accumulation time, 0.03

Metabolite Extraction
Three groups of lettuce leaves were ground into liquid nitrogen, 1000 µL methanol/acetonitrile/water (2:2:1, v/v/v) was added to each sample, then the samples were ultrasonicated twice at low temperature for 30 min after vortexing, and incubated at −20 • C for 1 h. Acetonitrile (1499230-935) and methanol (144282) were purchased from Merck (Darmstadt, Germany). Subsequently, the samples were centrifuged at 16,626 g at 4 • C for 15 min. Supernatants were collected, lyophilized, and stored at −80 • C. Five quality control (QC) samples were prepared by pooling aliquots of all samples. The pretreatment of the QC samples was the same as that of the study samples. The QC samples were evenly inserted in each set of the analysis running sequence to monitor the stability of the large-scale analysis [27].

Data Processing, Statistical Analysis, and Metabolic Pathway Analysis
The raw data were converted to common data format (mzXML) files for peak alignment, retention time correction, and peak area extraction. The structure of metabolites was identified by mass number matching (<25 ppm) and secondary spectrum matching. The ion peaks with missing values >50% were deleted and normalized by Pareto-scaling. 13,046 and 13,028 metabolite ions were acquired separately in the negative ion mode of nitrogen and phosphorus, whereas 9761 and 9757 metabolite ions were acquired separately in positive ion mode for nitrogen and phosphorus.
Principal component analysis (PCA) was applied for unsupervised multivariate analysis using SIMCA-P 14.1 (Umetric, Umea, Sweden) software, which provided an overview of the similarities and differences among the samples. Screening metabolites was divided into two steps. Fold change analysis and t-test were used to search preliminarily for metabolic biomarkers, then fold change (FC) value and p-value were calculated, and volcano plots were generated. The processed data sets were imported to SIMCA-P 14.1 for orthogonal partial least squares-discriminant analysis (OPLS-DA). The variable importance in the projection (VIP) values were considered as differential variables for further screening metabolic biomarkers [28]. Metabolite structural identification was conducted based on UPLC-QTOF MS analysis and retention time, accurate molecular weight, and MS/MS data. The metabolomics pathway was analyzed by MetaboAnalyst (https://www.metaboanalyst.ca) based on the Kyoto Encyclopedia of Genes and Genomes (KEGG) metabolic pathway database.

Metabolic Biomarkers Screening
QC samples were visualized using total ion current, which indicated that the instrumental analysis of all samples was strong and the peak capacity was large, so the retention time was reproducible. The correlation between QC samples was greater than 0.9, which indicated the system was stable as shown in Figure 2. To observe the overall distribution of samples, the PCA was considered an unsupervised model to operate without any anthropogenic factors, which was conducive to understanding the holistic data and eliminating abnormal samples. The results showed that QC samples were aggregated, and stress and normal treatment groups of samples separated clearly. There were no outlier samples in the Hotelling's T2 region (95% confidence interval).
To examine the significant differences of metabolites between low-nitrogen stress, low-phosphorus stress, and normal samples, the univariate analysis with fold change analysis and t-test were used for screening metabolic biomarkers. The importance of metabolites in sample discrimination was visible in the volcano plots. The volcano plots were employed to recognize the downregulated (FC < 0.5) and upregulated (FC > 2) metabolites that showed significant differences (p-value < 0.05). As shown in Figure 3, red points were upregulated metabolites, green points were downregulated metabolites, the blue points were metabolites with no change. Low-nitrogen stress caused an increase in 4952 metabolites and a decrease in 6787 metabolites, low-phosphorus stress caused an increase in 4359 metabolites and a decrease in 7538 metabolites, the Venn diagrams are shown in Figure 4. According to structural features of metabolic biomarkers, 260 metabolites were identified for low-nitrogen stress (134 metabolites increased and 126 metabolites decreased). The low-phosphorus stress treatment caused an increase in 53 metabolites and a decrease in 116 metabolites.
The next step of the statistical analysis was to perform a supervised OPLS-DA to separate samples into two clusters and identify biomarkers between the normal and stress groups. To prevent the model from overfitting, the quality of the model was investigated by seven-fold interactive validation and a 200-response sequencing test. Validation parameters for the two OPLS-DA were R 2 Y = 0.948 and Q 2 Y = 0.885 for low-nitrogen stress, and for low-phosphorus stress R 2 Y = 0.921 and Q 2 Y = 0.892. 262 and 169 identified compounds of low-nitrogen and low-phosphorus stresses were further screened by Figure 5, in which the upper part was the V-plot figure, the lower part was the S-plot figure, and each box/point represented one compound. In the V-plot figure, the biomarkers were filtered by VIP > 1. In the S-plot figure, the closer the compound was to the lower left and the upper right corners the greater the contribution to each category. In sum, 55 and 76 compounds were identified as important variables that contribute to the low-nitrogen or low-phosphorus stress.       The next step of the statistical analysis was to perform a supervised OPLS-DA to separate samples into two clusters and identify biomarkers between the normal and stress groups. To prevent the model from overfitting, the quality of the model was investigated by seven-fold interactive validation and a 200-response sequencing test. Validation parameters for the two OPLS-DA were R 2 Y = 0.948 and Q 2 Y = 0.885 for low-nitrogen stress, and for low-phosphorus stress R 2 Y = 0.921 and Q 2 Y = 0.892. 262 and 169 identified compounds of low-nitrogen and low-phosphorus stresses were further screened by Figure 5, in which the upper part was the V-plot figure, the lower part was the S-plot figure, and each box/point represented one compound. In the V-plot figure, the biomarkers were filtered by VIP > 1. In the S-plot figure, the closer the compound was to the lower left and the upper right corners the greater the contribution to each category. In sum, 55 and 76 compounds were identified as important variables that contribute to the low-nitrogen or low-phosphorus stress.
(a)  The next step of the statistical analysis was to perform a supervised OPLS-DA to separate samples into two clusters and identify biomarkers between the normal and stress groups. To prevent the model from overfitting, the quality of the model was investigated by seven-fold interactive validation and a 200-response sequencing test. Validation parameters for the two OPLS-DA were R 2 Y = 0.948 and Q 2 Y = 0.885 for low-nitrogen stress, and for low-phosphorus stress R 2 Y = 0.921 and Q 2 Y = 0.892. 262 and 169 identified compounds of low-nitrogen and low-phosphorus stresses were further screened by Figure 5, in which the upper part was the V-plot figure, the lower part was the S-plot figure, and each box/point represented one compound. In the V-plot figure, the biomarkers were filtered by VIP > 1. In the S-plot figure, the closer the compound was to the lower left and the upper right corners the greater the contribution to each category. In sum, 55 and 76 compounds were identified as important variables that contribute to the low-nitrogen or low-phosphorus stress.

Cluster Analysis of Metabolic Biomarkers
To more intuitively and comprehensively display the expression patterns of the metabolites, six samples were randomly selected from each group. The metabolic biomarkers in the two groups were analyzed using heat maps, which showed that the low-nitrogen and low-phosphorus stresses responded differently compared with the normal treatment. The results showed the distribution of metabolic biomarkers in the two groups and demonstrated that metabolic biomarkers between the two groups had a clustering trend. As shown in Figure 6, low-nitrogen stress caused an increase in 21 metabolites and a decrease in 34 metabolites, 45 of which were negative ion mode and 10 were positive ion mode. In Figure 7, low-phosphorus stress caused an increase in 45 metabolites and a decrease in 31 metabolites, 54 metabolites in negative ion mode, 22 metabolites in positive ion mode.

Cluster Analysis of Metabolic Biomarkers
To more intuitively and comprehensively display the expression patterns of the metabolites, six samples were randomly selected from each group. The metabolic biomarkers in the two groups were analyzed using heat maps, which showed that the low-nitrogen and low-phosphorus stresses responded differently compared with the normal treatment. The results showed the distribution of metabolic biomarkers in the two groups and demonstrated that metabolic biomarkers between the two groups had a clustering trend. As shown in Figure 6, low-nitrogen stress caused an increase in 21 metabolites and a decrease in 34 metabolites, 45 of which were negative ion mode and 10 were positive ion mode. In Figure 7, low-phosphorus stress caused an increase in 45 metabolites and a decrease in 31 metabolites, 54 metabolites in negative ion mode, 22 metabolites in positive ion mode.

KEGG Metabolic Pathway Analysis
The KEGG enrichment results can be visualized by bubble diagrams, and the ten pathways with the smallest p-value were mapped. As shown in Figure 8, complex interactions between networks of metabolic pathways exist in many aspects of nitrogen and phosphorus metabolism, and these enriched pathways were ATP-binding cassette (ABC) transporters and Biosynthesis of plant secondary metabolites. ABC transporters manage the active transport of a wide range of molecules across biological membranes, including secondary metabolites, inorganic acids, lipids, and phytohormones [29]. In this study, ABC transporters were associated with phosphoric acid; Lglutamate; sucrose; xylitol; raffinose; maltotriose; and alpha, alpha-trehalose; among others. Secondary metabolites often act as defense molecules and protect plants in various adverse

KEGG Metabolic Pathway Analysis
The KEGG enrichment results can be visualized by bubble diagrams, and the ten pathways with the smallest p-value were mapped. As shown in Figure 8, complex interactions between networks of metabolic pathways exist in many aspects of nitrogen and phosphorus metabolism, and these enriched pathways were ATP-binding cassette (ABC) transporters and Biosynthesis of plant secondary metabolites. ABC transporters manage the active transport of a wide range of molecules across biological membranes, including secondary metabolites, inorganic acids, lipids, and phytohormones [29]. In this study, ABC transporters were associated with phosphoric acid; L-glutamate; sucrose; xylitol; raffinose; maltotriose; and alpha, alpha-trehalose; among others. Secondary metabolites often act as defense molecules and protect plants in various adverse conditions [30], the pathways included succinate; l-tryptophan; fumarate; l-asparagine; citrate; isocitrate; l-isoleucine; and cis-aconitate.

Discussion
Nitrogen and phosphorus deficiencies resulted in changes of morphology, color, and texture. For instance, nutritional deficiency caused the changes in foliar color due to the loss of chlorophyll [31,32], and changes in texture may occur because of changes in surface structure and yellowish appearance in the leaves [33]. There was some evidence that nitrogen and phosphorus may affect the endogenous levels of phytohormones, such as cytokinins, which are involved in the regulation of both cell division and cell elongation, known to be implicated in the regulation of plant morphogenesis [34,35]. Hence, there were two potential problems with the detection of nitrogen and phosphorus. First, plants had a similar performance in low-nitrogen and low-phosphorus stresses [36]. Second, changes in morphology, color, and texture are inconspicuous in the early stages of nutrient stress. When these characteristics show obvious differences, nutritional stress was severe [37].

Discussion
Nitrogen and phosphorus deficiencies resulted in changes of morphology, color, and texture. For instance, nutritional deficiency caused the changes in foliar color due to the loss of chlorophyll [31,32], and changes in texture may occur because of changes in surface structure and yellowish appearance in the leaves [33]. There was some evidence that nitrogen and phosphorus may affect the endogenous levels of phytohormones, such as cytokinins, which are involved in the regulation of both cell division and cell elongation, known to be implicated in the regulation of plant morphogenesis [34,35]. Hence, there were two potential problems with the detection of nitrogen and phosphorus. First, plants had a similar performance in low-nitrogen and low-phosphorus stresses [36]. Second, changes in morphology, color, and texture are inconspicuous in the early stages of nutrient stress. When these characteristics show obvious differences, nutritional stress was severe [37].
The change of internal small molecular compounds can occur earlier than changes in external characteristics [37]. Previous studies only focused on the changes of metabolites under single nutrient stress, such as nitrogen, phosphorus, potassium or zinc [17][18][19][20]38]. In fact, different nutritional stress can cause changes of the same metabolites. In this study, a large number of metabolic biomarkers were selected by FC value, p-value and VIP value. In these, low-nitrogen and low-phosphorus stresses caused an increase in 16   Since the content of downregulated metabolites was low in lettuce under low-nitrogen or low-phosphorus stress, the study focused on upregulated biomarkers. The 5 biomarkers of low-nitrogen stress, all in negative ESI mode, belonged to classes of organic compounds known as cyclohexanols, hydroxybenzoic acid derivatives, oligosaccharides, dinitrophenols, and phenylpropanoids. However, the biomarker panel of low-phosphorus stress included 29 metabolites, which is not beneficial for practical application. The top 5 biomarkers were further screened by VIP > 1.75, which were citrate, isocitrate, l-5-oxoproline, succinate, and histamine. Citrate and isocitrate belonged to the class of organic compounds known as tricarboxylic acids and derivatives. L-5-oxoproline, succinate, and histamine are included in alpha amino acids and derivatives, dicarboxylic acids and derivatives and 2-arylethylamines.

Conclusions
Our objective aimed to obtain insight into lettuce responses to low-nitrogen and low-phosphorus stresses by examining metabolites. We used LC-MS techniques to investigate the abundance and identities of metabolites. PCA, Fold change analysis, t-test and OPLS-DA were used to understand the holistic data and search for metabolic biomarkers. This study indicated that low-nitrogen stress caused an increase in 21 metabolites and a decrease in 34 metabolites; low-phosphorus stress caused an increase in 45 metabolites and a decrease in 31 metabolites. The most relevant pathways were ATP-binding cassette transporters and Biosynthesis of plant secondary metabolites. From our comparison and analysis, low-nitrogen stress caused accumulation of inositol, p-hydroxybenzoic acid, stachyose, dinoseb, and 7, 8-dihydroxycoumarin; low-phosphorus stress caused accumulation of citrate, isocitrate, l-5-oxoproline, succinate, and histamine. The metabolic biomarkers could be used to monitor the early nitrogen or phosphorus status of lettuce, presumably to guide appropriate fertilization regimens.

Conflicts of Interest:
The authors declare no conflict of interest.