Jimai 22 is a high-yielding variety selected by the Crop Research Institute of Shandong Academy of Agricultural Sciences, containing wax in wild leaf epidermis. Combined with the experimental experience of our research group for many years, 0.6% (v/w) EMS was used to mutagenize Jimai 22 in this study, and then waxless and waxy were selected from the separated high-generation population. Compared with waxy, waxless showed less wax in the whole plant (Figures 1A,B). In October 2020, the mutants were planted in the Hancun Experimental Base of Wheat Research Institute of Shanxi Agricultural University (36°N, 111°E), with three replicates. The length was 2 m, the row spacing was 0.3 m, and 30 seeds per row were sown evenly. Wax content and transcriptome analysis were performed on flag leaves at the heading stage. Three biological repeats were performed for both waxless and waxy, and two samples with good repeatability were selected in transcriptome analysis.

In order to analyze the effect of leaf epidermis wax on the water loss rate of wheat leaves, in the same period as transcriptome sequencing, the flag leaves at the heading stage of waxy and waxless mutants were fully soaked and dehydrated in the dark at room temperature for 10 h (stomata were completely closed). The leaf water loss rate was calculated by 0.001 mg analysis balance (AUW320, Japan) every hour. Each sample was repeated three times.

Liposome was detected by liquid chromatogramphy-tandem mass spectrometry (LC-MS). The leaves were immersed in chloroform for 30 s to dissolve the epidermis wax and dried with nitrogen. Lipids were extracted according to the MTBE method (Pizarro et al., 2013). Briefly, samples were spiked with a suitable amount of internal lipid standards and then homogenized with 200 µl water and 240 µl methanol. After that, 800 µl of MTBE was added and the mixture was ultrasound-sonicated 20 min at 4°C followed by keeping still for 30 min at room temperature. The solution was centrifuged at 14,000 g for 15 min at 10°C and the upper organic solvent layer was obtained and dried under nitrogen. Reverse-phase chromatography was selected for LC separation using the CSH C18 column (1.7 µm, 2.1 × 100 mm, Waters). The lipid extracts were re-dissolved in 200 µl 90% isopropanol/acetonitrile, centrifuged at 14,000 g for 15 min; finally, 3 µl of the sample was injected. Solvent A was acetonitrile-water (6:4, v/v) with 0.1% formic acid and 0.1 Mm ammonium formate and solvent B was acetonitrile-isopropanol (1:9, v/v) with 0.1% formic acid and 0.1 Mm ammonium formate. The initial mobile phase was 30% solvent B at a flow rate of 300 μl/min. It was held for 2 min and then linearly increased to 100% solvent B in 23 min, followed by equilibrating at 5% solvent B for 10 min. Mass spectra were acquired by Q-Exactive Plus in positive and negative modes, respectively. ESI parameters were optimized and preset for all measurements as follows: source temperature was set at 300°C, capillary temp at 350°C, the ion spray voltage at 3000 V, S-Lens RF level at 50% and the scan range of the instruments at m/z 200–1800.

The clean leaves of the two mutants were taken at the heading stage, wrapped in tin foil, frozen in liquid nitrogen quickly, and then frozen in the refrigerator at −80°C. The samples were sent to Beijing Baimaike Biotechnology Co., Ltd., for transcriptome sequencing. Total RNA was extracted using the Trizol method, and the concentration and purity of RNA were measured by NanoDrop 2000 (Thermo Fisher Science, Wilmington, DE). RNA integrity was evaluated using the RNA Nano 6000 analysis kit of Agilent Biological Analyzer 2,100 system (Agilent Technologies, CA, United States ). A total amount of 1 μg RNA per sample was used as input material for the RNA sample preparations. Sequencing libraries were generated using NEBNext Ultra TM RNA Library Prep Kit for Illumina (NEB, United States ) following the manufacturer’s recommendations and index codes were added to attribute sequences to each sample. In order to select cDNA fragments of preferentially 240 bp in length, the library fragments were purified with the AMPure XP system (Beckman Coulter, Beverly, United States ). Then, 3 μl USER Enzyme (NEB, United States ) was used with size-selected, adaptor-ligated cDNA at 37°C for 15 min followed by 5 min at 95°C before PCR. PCR was performed with Phusion High-Fidelity DNA polymerase, Universal PCR primers, and Index (X) Primer. At last, PCR products were purified (AMPure XP system) and library quality was assessed on the Agilent Bioanalyzer 2,100 system.

Primers were designed by Primer Premier Version 5.0 (Premier Biosoft International, Palo Alto, CA) and combined by Sangon (www.sangon.com) (Supplementary Table S1). LA-Taq enzyme from TaKaRa (www.takara.com.cn) was used for PCR amplification. PCR was performed in total volumes of 15 μl, including 3 pmol of each primer, 120 μM of each dNTP, 80 ng template DNA or cDNA, 0.75 unit La-Taq, and 7.5 μl of 2 × buffer (TaKaRa Biotechnology (Dalian) Co., Ltd., Product Code: DRR20AG). PCR was performed as follows: 95°C for 4 min, followed by 35 cycles of 95°C for 30 s, annealing (55–62°C) for 30 s, extension at 72°C (30 s–3 min), and 72°C for 30 s, with a final extension of 72°C for 10 min. Annealing temperatures and extension times depended upon individual primer sets and the length of expected PCR products.

Quantitative real-time PCR was performed using SYBR® Premix Ex Taq™ II (Takara) according to the manufacturer’s instructions on a 7,300 Real-time PCR System (Applied Biosystems), where the relative expression of each gene was calculated according to the 2^−ΔΔCT method (Zheng et al., 2020). The glyceraldehyde-3-phosphate dehydrogenase gene was used as an endogenous reference for real-time PCR, and all analyses were performed with three technical and three biological replications.

“Lipid Search” is a search engine for the identification of lipid species based on MS/MS math. Lipid Search contains more than 30 lipid classes and more than 1,500,000 fragment ions in the database. Both mass tolerances for precursor and fragment were set to 5 ppm. Under the positive and negative ion modes, the OPLS-DA model was constructed by SAMIC14.1 software, and the prediction rate Q² of the model was obtained by 7-fold cross-validation. Q² > 0.5 was taken as the reliable standard of the model; the multiple of variation analysis (Fold Change Analysis, FC) > 4 or < 0.25 and the variable importance for the projection (VIP) > 1, p value < 0.05, were used as the screening criteria to compare the overall differential expression multiples of lipid ions in leaves of waxless and waxy. Lipid differences were analyzed with SPSS19.0 and Microsoft Excel 2019, plotted with Origin 2018. The transcriptome data were further processed by the online platform of BMK Cloud (www.Biocloud.net). FPKM (Fragments Per Kilobase of exon model per Million mapped fragments) method was used to standardize the gene expression level. DESeq R package was used for differential analysis, and Fold Change ≥2 and FDR < 0.01 were used as screening criteria to determine the differentially expressed genes (DEGs) between waxless and waxy mutants.

The leaf water loss rate of the two mutants was faster in the darkroom temperature environment. The wax content had a significant effect on the leaf dehydration rate. After natural dehydration for 1 h, the water loss rate of low-wax mutant leaves reached 24.4% and that of multi-wax mutant leaves reached 17.4%. The leaves of the two mutants were dehydrated rapidly within 1–5 h, and the water loss rate of waxless was significantly higher than that of waxy. After dehydration for 5 h, the water loss rate of waxless leaves reached 63.5%, while that of waxy was only 44.3%. There was a significant difference between the two mutants. Although the dehydration rate of the two mutants slowed down after 5 h of dehydration, the water loss rate of low-wax mutants was still significantly higher than that of wax-rich mutants at 10 h of dehydration (Figure 2). This shows that the mutant waxy has a slow dehydration rate and strong water retention capacity and can maintain higher leaf water content, indicating that leaf wax content plays an important role in maintaining leaf water content.

Positive and negative ion patterns of electrospray ionization (ESI) were used in this study. By UPLC analysis, 31 lipid subclasses and 1,367 lipid molecules were identified in the leaf epidermis of the two mutants (Figure 3A). The lipids with a high number of lipid species are triacylglycerol (TG), ceramide (Cer), diacylglycerol (DG), wax ester (WE), (O-acyl)-1-hydroxy fatty acid (OAHFA), monogalactosyl diacylglycerol (MGDG), phosphatidylethanolamine (PE), phosphatidylglycerol (PG), cholesterol (ChE), cardiolipin (CL), phosphatidylinositol (PI), and phosphatidylserine (PS). Among them, TG, Cer, DG, WE, and OAHFA account for 83% of the total number of molecules (Figure 3B). Moreover, it was found that the main lipid components of wheat leaf epidermis were TG, Cer, DG, WE, and OAHFA.

To test the repeatability of samples, we performed OPLS-DA analysis on all samples (Figure 4B). Each group of samples gathered closely and located in the middle of each group, indicating the liposome analysis showed high reproducibility. By comparing the contents of lipid subclasses between low-wax mutants and multi-wax mutants at the heading stage, significant differences in lipids were screened out. As shown in Figure 4A, 11 kinds of lipid substances with significant differences were found between the two mutants (VIP >1, p < 0.05). Among them, there are six lipid subclasses with high lipid content in waxless, namely, Cer, ZyE, PE, MGDG, CL, and digalactosyl diacylglycerol (DGDG). On the other hand, in waxy, five kinds of lipid classes were found to be DG, ChE, OAHFA, WE, and monoacylglycerol (MG) (VIP >1, p < 0.05). We infer that the more abundant and different lipids in the waxy are the main components of the wax in the epidermis of wheat compared with the waxless. That is the main lipid components of wheat leaf epidermis at the heading stage are DG, ChE, OAHFA, WE, and MG.

We analyzed the differences of 1,367 detected lipid ions (FC > 4 or <0.25, VIP >1, p < 0.05). Compared with the mutant waxy, the mutant waxless downregulated 177 lipid ions and upregulated 246 lipid ions. The results were shown in the form of a volcanic map (Figure 4C). 177 lipid ions related to wax metabolism in the epidermis of wheat leaves were screened (FC < 0.25, VIP >1, p < 0.05), belonging to eight lipid subclasses, including WE (C19-C50), DG (C27-C53), TG (C29-C73), MG (C31-C35), OAHFA (C31-C52), Cer (C29-C48), ChE (C₄₆H₈₆O₂N₁), and CL (C₆₅H₁₂₀O₁₇P₂) (Supplementary Table S2).

The main differential lipid ions in DG are DG(45:2)+H, DG(28:0e)+H and DG (37:4e)+H; ChE(19:0)+NH₄ in ChE; OAHFA(16:0_31:1)-H, OAHFA(30:1)-H, and OAHFA(31:0)-H in OAHFA; WE(24:0_14:1)+H, WE(28:1_22:0)+H, WE(28:1_18:0)+H, WE(26:0_21:3)+H, and WE(12:0_19:3)+NH₄ in WE; MG(32:1)+NH₄, MG(28:1)+H, and MG(32:1)+NH₄ in MG (Table 1). These may be the main lipid ions that affect the gloss phenotype of wheat leaf epidermis.

In addition to the content of lipids, the carbon chain length of lipids is also an important factor that cannot be ignored. We added the lipid ions with the same carbon chain length and counted the lipid ions with diverse carbon chain lengths under different lipid subclasses to further analyze the main discrepant lipids in wheat epidermis wax (Figure 5). Compared with waxy mutants, the lipid ions significantly decreased in waxless mutants were WE (C29), WE (C31), WE (C44), WE (C46), DG (C47), DG (C45), DG (C43), DG (C37), DG (C36), DG (C31), OAHFA (C51), OAHFA (C50), OAHFA (C48), OAHFA (C31), MG (C35), MG (C31), and TG (C64).

These lipid ions affect the synthesis, transcription, and transport of wax through the change of carbon chain length, leading to the smooth green phenotype of wheat waxy deletion mutants. We found that WE in wheat leaf epidermis is mainly concentrated in C44-dominated wax ester molecules and also widely distributed in C46, C31, and C29 wax ester molecules but less in other chain lengths. DG in the wax mixture of the wheat epidermis is mainly concentrated in C47, C45, C37, and C31 molecules. In OAHFA, C48 plays a leading role, accompanied by a large number of C50-C51 and C31 molecules. C35 and C31 play significant roles in MG. On the other hand, in TG, the length of the carbon chain varies from C29 to C64, in which C48 plays a major role.

Through the analysis of leaf transcriptional groups of low-wax mutants and multi-wax mutants at the heading stage, correlation coefficients among different biological repetitive samples ranged from 0.890 to 0.926 (Figure 6A). Compared with the mutant waxy, 6,840 DEGs were detected in the mutant waxless, of which 3,181 genes were upregulated and 3,659 genes were downregulated (Fold Change ≥2, FDR <0.01) (Figure 6B). To understand the biological significance between low-wax mutant and multi-wax mutant DEGs, these genes were enriched by Gene Ontology (GO) analysis. We found that the biological processes were mainly enriched in the metabolic process (28.89%), cellular process (24.73%), single-organism process (14.96%), biological regulation (8.35%), and response to stimulus (6.46%). These genes are mainly distributed in the membrane, membrane part, cell, component cell part, and organelle. The molecular functions mainly include binding and catalytic activity (Supplementary Figure S1).

Furthermore, we annotated the lipid function of 3,659 DEGs, screened 343 DEGs, and enriched them using the Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis (Figures 6C,D). It was found that these genes are mainly concentrated in glycerophospholipid metabolism, cutin, suberine and wax biosynthesis, glycerolipid metabolism, sphingolipid metabolism, terpenoid backbone biosynthesis, fatty acid degradation, and fatty acid elongation pathway, which are consistent with the main lipid components that control the glossy phenotype of wheat leaves in this study.

Combined with the determined lipid components, we marked the annotated genes in the corresponding pathway, constructed the metabolic model of the main lipid components in wheat epidermis according to the KEGG metabolic pathway (Figure 7), and screened 46 related genes (Table 2). Moreover, compared with the multi-wax mutant waxy, the genes annotated to the synthesis of wax ester (WE) in the low-wax mutant waxless were 3-ketoacyl-CoA synthase genes (KCS), very-long-chain enoyl-CoA reductase genes (TER), alcohol-forming fatty acyl-CoA reductase genes (FAR), wax-ester synthase genes (WSD1), aldehyde decarbonylase genes (CER1), and midchain alkane hydroxylase genes (MAH1). The downregulation of TER inhibits the conversion of long-chain 3-oxoacyl-CoA to long-chain acyl-CoA (n + 2), the downregulation of CER1 inhibits the transformation of upstream A long-chain aldehyde to A long-chain alkane, and the expression genes involved in KCS, FAR, WSD1, and MAH1 were both upregulated and downregulated.

In the process of OAHFA synthesis, the upregulated expression of acyl-coenzyme A thioesterase 1/2/4 genes (ACOT1_2_4) promotes the synthesis of long-chain fatty acid, while the expression level of acylglycerol lipase genes (MGLL) was both upregulated and downregulated during the conversion of monoacyl-glycerol to long-chain fatty acid. Downregulation of aldehyde dehydrogenase family 7 member A1 genes (ALDH7A1) expression inhibited the transition between A long-chain aldehyde and long-chain fatty acid, while aldehyde dehydrogenase genes (ALDH) upregulated the transition between them. During the transformation of long-chain fatty acid into omega-hydroxy fatty acid, the expression of long-chain fatty acid omega-monooxygenase genes (CYP704B1) was downregulated, while that of fatty acid omega-hydroxylase genes (CYP86) was upregulated in the mutant waxless. Glycerol-3-phosphate acyltransferase genes (GPAT) and ALDH7A1 were downregulated in the process of glyceride metabolism of DG, MG, and TG. The downregulated expression of GPAT inhibits the synthesis of 1-acyl-sn glycerol 3-phosphate in waxless, which in turn reduces the synthesis of MG. The downregulated expression of ALDH7A1 in the low-wax mutant decreases the synthesis of DG. The upregulated genes are ALDH, diacylglycerol diphosphate phosphatase genes (DPP1), and diacylglycerol kinase genes (dgkA), while the expressions of 1-acyl-sn-glycerol-3-phosphate acyltransferase genes (plsC), phospholipid: diacylglycerol acyltransferase genes (E2.3.1.158), MGLL, and alcohol dehydrogenase genes (AKR1AI) were both upregulated and downregulated. These genes work together to reduce the content of TG in the low-wax mutation waxless. In order to verify the reliability of our data for screening wax-related genes, we randomly selected 13 differentially expressed genes for real-time quantitative PCR detection, and the results of qRT-PCR of 13 genes were consistent with those of transcriptome data (Supplementary Figure S2), indicating that the screened wax-related genes are highly reliable, which is helpful for further functional verification and cloning of wax-related genes.

Epidermal wax is a complex lipid mixture composed of VLCFAs and their derivatives, which plays a special role in plant resistance to drought, diseases, and insect pests (Bernard and Joubes, 2013). At present, the GC-MS method is mainly used to determine the wax of epidermis (Rolim et al., 2015) and primary alcohols, secondary alcohols, aldehydes, alkanes, ketones, esters, triterpenes, sterols, and flavonoids (Zhang et al., 2005; Tafolla-Arellano et al., 2018). Based on the limitations of the determination method, there may still be many epidermal wax-related lipid components that have not been detected, which limits the study of the epidermal wax biosynthesis pathway. Lipidomics as a research model based on high-throughput analysis can systematically analyze the changes in lipid composition and expression in organisms. Perez-Navarro et al. (2019) used LC-MS technology to determine grape lipids and found new information on the composition of free fatty acids such as glycerol, glycerol phospholipids, and triterpenes in grape skins and seeds. Broughton et al. (2018) used electrospray ionization tandem mass spectrometry (ESI-MS/MS) to analyze the acyl parts of wax esters (WE) and sterol esters (SE) in common and mutant sunflower oils with different fatty acid profiles and discovered the methylsterol components in sunflower oil sterol esters (SE) for the first time.

Bianchi et al. (1980) studied the epidermal wax of common wheat in Chinese spring and found that its main components are n-alkanes, esters, aldehydes, free alcohols, β-diketones and hydroxy-β-diketones. Racovita et al. (2016) detected new waxy substances such as 2-alkyl alcohol, benzyl alcohol, phenylethyl alcohol, and hydroxyphenylethanol in wheat flag leaves and peduncles. Lavergne et al. (2018) analyzed the epidermis of wheat leaves and stems and found that the waxy components of the epidermis were alkanes (C20-C42), fatty acids (C7-C34), ketones (C9-C35), and primary alcohols (C22-C33). In this study, UHPLC-MS/MS analysis technique was used for the first time to detect the leaf epidermis wax of low-wax mutant waxless and multi-wax mutant waxy of Jimai 22. We found that the main lipid components of leaf epidermis wax of Jimai 22 were WE (C19-C50), DG (C27-C53), MG (C31-C35), OAHFA (C31-C52), TG (C29-C73), Cer (C29-C48), ChE (C₄₆H₈₆O₂N₁), and CL (C₆₅H₁₂₀O₁₇P₂). Among them, DG, MG and OAHFA were detected in wheat leaf wax for the first time, which improved the composition of wheat wax components. It is clear that the liposome analysis of wheat epidermis wax by LC-MS can more systematically study the changes and functions of related lipid subclasses and lipid molecules in the process of epidermal wax metabolism, which is helpful to improve the pathway and mechanism of epidermal wax synthesis.

In this research, we found that the wax content was closely related to the water loss rate of leaves, which may be related to drought resistance (Zhang et al., 2013). Besides, DG, MG, and OAHFA were newly found in wheat epidermis wax. DG, as one of the main components of glycerol, plays a critical role in maintaining the stability of cell membrane (especially plasma membrane and chloroplast membrane) under high temperature and drought stress (Narayanan et al., 2016). It can maintain the fluidity of cell membrane by downregulating saturated DGs (4 and 10 double bonds) and u-regulating low unsaturated DGs (0–4 double bonds) (Navarro-Reig et al., 2019). MG, as a waxy substance, can form stable hydrated dispersions in water. Under salt stress, the ratio of MG unsaturated fatty acids to saturated fatty acids increased significantly in salt-tolerant varieties (Gogna et al., 2020). The MG in wheat epidermis wax found in this study is also mainly MG (32:1) + NH₄ and MG (28:1) + H with low saturation. OAHFA acts as a surfactant in human tear film lipids and plays a key role in stabilizing tear film (Schuett and Millar, 2013; Marshall et al., 2016). Interestingly, we found OAHFA in wheat leaf epidermis wax for the first time and C48 molecules play a leading role. Its effects on wheat stress resistance and related mechanism need to be further studied.

In addition, as an important oil substance, WE plays a significant role in plant resistance to drought because of its hydrolytic resistance (Ivarson et al., 2017; Wang et al., 2021). The WE wax found by predecessors in sunflower is mainly concentrated in C32–C48 molecules (Broughton et al., 2018). In our study, we found that the WE wax in the epidermis of wheat leaves is mainly distributed in the molecular range of C29–C46 chain length, and the wax ester molecule dominated by C44 plays a role, which is similar to the results of previous studies.

Previous studies have shown that the W1-W5 homologue of Arabidopsis CER protein in wheat increases the wax content of the epidermis by producing hydroxy-β-diketone and inducing the biosynthesis of β-diketone (Huang et al., 2017; Li et al., 2020). In this study, two CER1 genes, TraesCS1D03G0373900 and TraesCS1D03G0374000, were also detected in the low-wax mutants, which decreased the synthesis of A long-chain alkane by downregulating the expression of CER1. TaFAR1, TaFAR2, TaFAR3, TaFAR4, TaFAR5, and AtCER4 are homologous genes of CER4 in Arabidopsis. As a kind of alcohol-forming fatty acyl-CoA reductase, they induce the production of A long-chain primary alcohol and increase the wax content of wheat leaves (Wang W. et al., 2015; Wang et al., 2015b; Wang et al., 2015c; Wang et al., 2017). Besides, two FAR genes regulating primary alcohol synthesis, TraesCS4B03G0019500, and TraesCS7B03G1338900, were also found in the mutant waxless, which downregulated primary alcohol synthesis and promoted the smooth green appearance of the waxy wheat epidermis. In the research of low-wax mutants of wheat, Li et al. (2021) found that the expression levels of ACC1, LACS, KCS, and KCR were downregulated, which reduced the synthesis of VLCFAs. In our research, four KCS genes were also found to downregulate long-chain 3-oxoacyl-CoA synthesis in waxless. In addition, we found TER gene TraesCS6B03G1132900LC in wheat mutant waxless leaves for the first time, which inhibited the synthesis of long-chain acyl-CoA (n + 2) by downregulating its expression. This provides important reference information for enriching the molecular regulation mechanism of wheat epidermis wax synthesis.

Moreover, our study also found that the downregulation of two ALDH7A1 genes TraesCS5B03G0557800 and TraesCS5D03G0511400 inhibited the synthesis of long-chain fatty acid and DG, the downregulation of three CYP704B1 genes TraesCS3B03G0228000, TraesCS3B03G0938500, and TraesCS3D03G0167200 inhibited the synthesis of omega-hydroxy fatty acid, and the downregulation of GPAT genes TraesCS6A03G0083200 and TraesCS6B03G0120600 downregulated the expression of 1-acyl-sn glycerol 3-phosphate in DG, MG, and TG during glyceride metabolism, which inhibited the synthesis of 1-acyl-sn glycerol 3-phosphate in wax-free mutants. The downregulated expression of plsC gene TraesCS2B03G1160500 reduces the synthesis of 1,2-diacyl-sn glycerol 3-phosphate. Downregulation of phospholipid: diacylglycerol acyltransferase genes TraesCS7A03G0874000 expression inhibits the transformation from 1,2-diacyl-sn glycerol to TG. These genes are found for the first time to regulate waxy synthesis in the wheat epidermis. In-depth genomic comparison and functional identification of these genes are helpful in cloning wheat waxy functional genes and analyzing the regulatory mechanism of waxy metabolism. Our next research will focus on the functional verification of key genes to better explain the molecular regulation mechanism of wax in the leaf epidermis of wheat.