Front. Plant Sci.Frontiers in Plant ScienceFront. Plant Sci.1664-462XFrontiers Media S.A.10.3389/fpls.2025.1540179Plant ScienceOriginal ResearchCombined multi-omics and physiological approaches to elucidate drought-response mechanisms of durum wheatArriagadaOsvin1MenesesClaudio2PedreschiRomina34Núñez-LilloGerardo3MaureiraCarlos5RevecoSamantha5VillarroelValentina5SteinfortÚrsula5AlbornozFrancisco5Cabas-LühmannPatricia5SilvaManuela6MatusIván7SchwemberAndrés R.5*1Centro de Estudios en Alimentos Procesados (CEAP), Talca, Chile2Facultad de Ciencias Biológicas, Pontificia Universidad Católica de Chile, Santiago, Chile3Escuela de Agronomía, Facultad de Ciencias Agronómicas y de los Alimentos, Pontificia Universidad Católica de Valparaíso, Quillota, Chile4Millennium Institute Center for Genome Regulation, Santiago, Chile5Departamento de Ciencias Vegetales, Facultad de Agronomía y Sistemas Naturales, Pontificia Universidad Católica de Chile, Santiago, Chile6Linking Landscape, Environment, Agriculture and Food Research Center (LEAF), Associate Laboratory TERRA, Instituto Superior de Agronomia, Universidade de Lisboa, Lisboa, Portugal7Instituto de Investigaciones Agropecuarias (INIA)-Quilamapu, Programa de Mejoramiento Genético de trigo, Chillán, Chile
Edited by: Ali M. Missaoui, University of Georgia, United States
Reviewed by: Muluken Birara Enyew, Washington State University, United States
Moez Amri, Mohammed VI Polytechnic University, Morocco
*Correspondence: Andrés R. Schwember, aschwember@uc.cl
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
Introduction
Durum wheat is the most important cereal in the Mediterranean regions, where drought negatively affects grain yield. Therefore, our objective was to perform a multi-omics and integration analysis in conjunction with physiological trials to improve our understanding of drought tolerance mechanisms of durum wheat.
Methods
Genome-wide association study (GWAS) for yield components was performed on a panel of 225 elite durum wheat genotypes evaluated in eight sites under irrigated and rainfed conditions. Based on physiological parameters (net photosynthesis, intracellular CO2 content, transpiration and stomatal conductance) and grain yield, contrasting genotypes (susceptible and tolerant) to drought were identified. A transcriptomic (RNA-seq), metabolomic and integration analyses were performed to identify genes and metabolites associated with tolerance in durum wheat.
Results
Nine marker-trait associations were detected across 8 environments, and they were grouped into three QTL clusters (QTL_2A_TGW/GPS.1, QTL_2A_TGW/GPS.2, and QTL_2B_TGW/GPS.1), explaining between 5.15% and 14.29% of the phenotypic variation. One drought tolerant (QUC 3678-2016) and one susceptible (BRESCIA) genotype were identified based on physiological parameters. RNA-seq analysis showed that the genes regulated were mainly enriched in processes such as response to salicylic acid, plant organ senescence, synthesis of secondary metabolites, and immune response. Metabolic analysis showed that drought increased the contents of amino acids, sugars, and organic acids. The integration analysis identified 30 genes and six metabolites in the root and 30 genes and 10 metabolites in leaves as the primary variables in the drought-tolerant genotype, in which L-Proline was an important metabolite that allowed differentiating those two contrasting genotypes. A WRKY transcription factor was also positioned on the stable QTL QTN_2A_TGW/GPS.1 associated with the GENE-1342_238 SNP marker.
Discussion
These results open an opportunity to use new biomarkers in durum wheat breeding programs to develop resilient and high-yielding cultivars and ensure food security under water deficit conditions.
Triticum durum L.prolineWRKYQTLGWASRNA-seqAgencia Nacional de Investigación y Desarrollo10.13039/501100020884section-in-acceptancePlant BreedingIntroduction
Durum wheat (Triticum turgidum L. ssp. durum; 2n = 4 × = 28, AABB) is an essential crop in the world with a cultivated area of 17 million ha and a production of 33.6 million t in 2020 (Royo et al., 2021). This crop is the most important cereal in the Mediterranean region since it is deeply connected with the history and culinary traditions of this region (Martínez-Moreno et al., 2020). However, only half of the worldwide production is produced in this area (Zaïm et al., 2024). Among the main durum wheat-producing countries in the world are Spain, France, Italy, and Greece in southern Europe; Morocco, Algeria, and Tunisia in northern Africa; Turkey and Syria in southwest Asia; Canada, the United States, and Mexico in North America (Martínez-Moreno et al., 2020); and Argentina and Chile in South America (Colasuonno et al., 2019). Durum wheat plays a key role in the human diet because it is primarily used for making pasta and other semolina-based products, such as frike, couscous, bourghul, and unleavened bread, which are widely consumed in many countries of the world (Sharma et al., 2019), providing more than 20% of the nutritional demands for human consumption (Shewry and Hey, 2015).
Durum wheat is commonly grown in arid and semi-arid regions under rainfed conditions, where the water scarcity combined with high temperatures during grain-filling period significantly affects key plant metabolic pathways, causing grain yield losses of up to 50% to farmers (Dettori et al., 2017; Soriano et al., 2021). Additionally, global climate models predict a continued reduction in precipitation in the main durum wheat production areas (Schlaepfer et al., 2017; Daramola and Xu, 2022), which will significantly affect the production. Therefore, it is essential to study and understand the mechanisms that regulate drought tolerance in durum wheat to develop high-yielding varieties that are drought tolerance.
Plants have developed different physiological, biochemical, and molecular mechanisms to escape, tolerate, or avoid drought stress (Bashir et al., 2021; Dodd and Kudoyarova, 2021). Drought causes a reduction in leaf water potential (Ψ). The plant responds by closing its stomata, reducing the stomatal conductance, which affects transpiration and net photosynthesis rates and, ultimately, the grain yield (Gupta et al., 2020; Akhiyarova et al., 2023). The phytohormone abscisic acid (ABA) plays a key role in the regulation of stomatal closure under drought conditions (Hsu et al., 2021), where two signaling pathways in response to drought have been described, which are classified as ABA-dependent and ABA-independent (Yoshida et al., 2014). From a biochemical standpoint, drought signals induce stress-protective compounds with antioxidant activity, osmoprotectants, and secondary metabolites (Bashir et al., 2021). Finally, the drought response strategies are genetically controlled by numerous genes involved in osmolyte, redox, and hormonal metabolism (Singh PK. et al., 2019), such as dehydration-responsive element-binding factors (DREBs), basic leucine zipper (bZIP) proteins, MYBs, MYC, WRKY, NAC, sucrose non-fermenting1-related protein kinases (SnRKs) and SOS2-like protein kinases (CIPKs), among others (Joshi et al., 2016; Singh PK. et al., 2019; Bashir et al., 2021).
Genome-wide association studies (GWAS) have identified thousands of loci associated with complex traits in durum wheat (Sukumaran et al., 2018; Zaïm et al., 2024). However, the functional interpretation of GWAS results remains challenging due to large candidate regions and linkage disequilibrium (Knoch et al., 2024). Implementing and integrating different high-throughput omics techniques, including genomics, transcriptomics, proteomics, and metabolomics, have produced very valuable results in discovering specific genes, proteins, metabolites, and signaling processes in different tissues that are involved in the stress response in wheat (Shah et al., 2018; Goel et al., 2020; Da Ros et al., 2023; Le Roux et al., 2024). In addition, robust gene regulatory networks associated with drought stress have been constructed using multi-omics analysis in cereals (Baldoni et al., 2021), including wheat (Chen et al., 2025), which have had a major impact on accelerating crop improvement in the era of omics. However, multi-omics studies are still very limited for wheat. In fact, to our knowledge, multi-omics studies have only been performed in common wheat, opening the possibility of exploring and identify new markers, genes and metabolites that are associated with the response to drought in durum wheat. Therefore, this study aimed to perform physiological analyses and integrate different omics (genomics, transcriptomics, and metabolomics) approaches for a better understanding of the interaction between the different tolerance mechanisms involved in response to drought stress. These results are novel and valuable for developing resilient durum wheat cultivars, which are necessary to ensure production sustainability and global food security.
Materials and methods
All the experimental setup of this work is summarized in Supplementary Figure S1.
Genome-wide association approachPlant material and field management
A panel of 225 durum wheat genotypes comes from Mexico (International Maize and Wheat Improvement Center [CIMMYT]), Chile (Chile’s National Agricultural Research Institute [INIA], and Pontifical Catholic University of Chile [PUC]), Spain (Institute of Agrifood Research and Technology [IRTA]), and Tunisia (National Institute of Agronomic Research [INRAT]) breeding programs were used as plant material. The complete list of the germplasm used in this study is available in Supplementary Table S1. The field trials were established at Santa Rosa (36° 32′ S, 71° 55′
W; 220 m.a.s.l.) and Los Tilos (33° 70′ S, 70° 70′ W; 530 m.a.s.l.) INIA experimental stations in southern and central Chile, respectively. Considering that durum wheat is grown in Chile from the Valparaíso (central) to Biobío (south) Regions, the sites evaluated in this study represent contrasting production environments. The trials were established during 2021-22 (hereinafter referred to as 22) and 2022-23 (hereinafter referred to as 23) growing seasons under rainfed (r) and irrigated (i) conditions, giving a total of eight environments. The monthly precipitation and temperature data for each site and growing season are shown in Supplementary Table S2. For the irrigation treatments, furrow irrigation was used at Santa Rosa (SR) and Los Tilos (LT) locations; one irrigation with 50 mm of water at the end of tillering (Z22) and at the flag leaf stage (Z37) was applied in the two growing seasons in each location. The experimental design was α-lattice with two replications and 15 incomplete blocks per replicate, each containing 15 genotypes. Each entry was planted in plots of five rows of 2 m length and 0.2 m distance between rows. The sowing rate was 22 g m2, and sowing dates were 13 and 18 of August in 2021 and 2022, respectively, at Santa Rosa; 30 of July and 28 of July in 2021 and 2022, respectively, at Los Tilos. The harvest dates were 17 and 19 January 2022 and 2023, respectively, at Santa Rosa, and 10 and 20 January 2022 and 2023, respectively, at Los Tilos.
Plots were fertilized with a mix of 700 kg ha−1 with the following element components: 18% P2O5, 4% N, 15% of K2O, 6% of S, 5% of Ca, 0.2% of B and 0.3% of Zn. Fertilizers were incorporated with a cultivator before sowing. At Z20, and Z30, 111 kg ha−1 of N was applied. Weeds were controlled with the application of flufenacet + flurtamone + diflufenican (96 g a.i.; Bacara Forte 360 SC, Bayer) as pre-emergence controls and a further application of MCPA (525 g a.i.; Anasac) + metsulfuron-metil (5 g a.i.; Anasac) as post-emergent. To keep the genotypes free of leaf disease as rust (Puccinia striiformis and Puccinia triticina), two applications were made of the foliar fungicide Miravis Aeon, Syngenta (75 mL ha−1; 100.1 g/L azoxystrobin, 125 g/L propiconazole, 75.1 g/L pydiflumetofen). These applications were made before symptoms appeared, to avoid any interference of these diseases in the development of the plants.
Analysis of phenotypic data
Five yield-related traits such as spikes per square meter (SPM), grains per spike (GPS), thousand grain weight (TGW), grain yield (GY), and harvest index (HI) were evaluated. For each yield-related trait, the best linear unbiased predictors (BLUP) were predicted for each site, combined by treatment, and across environments using the META-R software version 6.0 (Alvarado et al., 2020), using the restriction maximum likelihood method (REML). The normal distribution was tested with the Shapiro–Wilk test using the shapiro.test function in R software. Most of the traits in the different data sets followed a normal or approximately normal distribution (Supplementary Table S3). The raw data and BLUP values are provided in the supplementary Supplementary Table S4 and Supplementary Table S5, respectively. The linear mixed models used in META-R are implemented in the lme4 package (Bates et al., 2014) and adjusted as follows:
Yijk=u+Repi+Blockj(Repi)+Genk+ϵijk
Where Yijk is the trait of interest, u is the overall mean effect, Repi is the effect of the ith replicate, Blockj(Repi) is the effect of the jth incomplete block within the ith replicate, Genk is the effect of the kth genotype and ϵijk is the effect of the error associated with the ith replication, jth incomplete block, and kth genotype. The distribution was assumed standard with mean zero and variance σϵ2. All effects, except the overall mean, were treated as random. Regarding the analysis between environments by treatment or between all environments, the terms associated with the environment and the environment x genotype interaction were added to the model (Alvarado et al., 2020). Additionally, phenotypic and genetic correlations among the traits combined across all environments were performed in META-R software version 6.0 (Alvarado et al., 2020).
Analysis of genotypic data
The total genomic DNA of the 225 genotypes was extracted from fresh leaves using Plant/Fungi DNA Isolation Kit (Norgen Biotek, Canada), according to the manufacturer’s instructions. DNA quality and concentration were measured using agarose gel electrophoresis (1%) and Qubit 4.0 (Thermo Fisher Scientific, Inc.). SNP genotyping was performed by Trait Genetics GmbH (http://www.traitgenetics.com/) using the 90K SNP array developed by Wang et al. (2014), of which a total of 51,423 SNPs were distributed on the A and B genome of durum wheat. Monomorphic SNP and SNPs with minor allele frequency (MAF) of less than 5% and those SNPs and accessions containing 10% or more of missing values were removed using TASSEL 5.2 software (Bradbury et al., 2007). Additionally, the genotypes with 10% or more heterozygotes were also eliminated.
Population structure and kinship analysis
SNPs were pruned for linkage disequilibrium (LD) to remove linked loci at the r2 > 0.2 using PriorityPruner software (http://prioritypruner.sourceforge.net/) to infer the population structure (Q) and kinship (K) matrices. The number of subgroups among the 225 genotypes was inferred based on a Bayesian model approach implemented in Structure v2.3.4 (Pritchard et al., 2000). An admixture ancestry model with correlated allele frequencies and no prior information of population origin was used. The putative number of subpopulations (K) ranged from 1 to 6. Ten replicates were performed for each K with a burn-in period of 100,000 steps followed by 1,000,000 MCMC iterations. The optimal K value was determined in Structure Harvester (https://taylor0.biology.ucla.edu/structureHarvester/) using the ad hoc statistic ΔK (Evanno et al., 2005). To confirm the results, a discriminant analysis of the principal components (DAPC) was implemented in ade4 (Chessel et al., 2004) and adegenet (Jombart, 2008) and visualized with ggplot2 (Wickham, 2010) packages in R software. The K matrix was estimated using the ‘Scald_IBS method’ in the TASSEL 5.2 software (Bradbury et al., 2007).
Association analysis
The marker-trait association analysis was carried out in TASSEL 5 (Bradbury et al., 2007) with the following mixed linear model (MLM) proposed by Yu et al. (2006):
y=Sa+Qv+Zu+ϵ
Where y is the vector of BLUPs of each trait at each environment, combined for the two water treatments (irrigated and rainfed), and across all environment, a is the vector of SNP effect (fixed), v is the vector of population structure (Q; fixed), u is the vector of kinship effects (K; random), and ϵ is the vector of residual effects. S, Q, and Z are incidence matrices relating y to av and u, respectively. A threshold of p< 0.001 [−log10(p) > 3] was used to indicate the significant SNP-trait association (MTA). Manhattan plots were drawn to visualize significant markers using SRplot (Tang et al., 2023), and quantile-quantile (Q-Q) plots to show important p-value distributions. Furthermore, for comparison purposes of the confidence of the MTAs detected using MLM, the general linear model (GLM) and the fixed and random circulating probability unitization model (FarmCPU) were performed on the data sets combined by treatment using the R package rMVP (Yin et al., 2021). MTAs detected for the same trait in at least two environments were reported as stable MTAs and MTAs detected for the same trait in all environments were considered constitutive. Moreover, MTAs were clustered into a unique QTL using the genome LD decay value, where 8.71 Mbp and 7.37 Mbp were used for genome A and B, respectively. The MTAs not in LD were considered as independent QTL.
Physiological studyYield tolerance index and experimental design
The yield tolerance index (YTI) was calculated for the 225 genotypes established at Santa Rosa
during the 2021–22 growing season. This index combines the yield performance of a genotype
under drought with its potential yield under irrigated conditions (Ober et al., 2004). The three most tolerant (highest YTI; QUC 3418-2020, QUC 3592-2020, and QUC 3678-2016) and susceptible (lowest YTI; QUC 3319-2020, QUC 3572-2018, and BRESCIA) genotypes were coded from G1 to G6and selected for the evaluation of physiological parameters (Supplementary Table S6).
The experiment was performed in a greenhouse under a natural photoperiod with control of temperature and relative humidity at the Pontifical Catholic University of Chile (33°29’ 46.117” S and 70°36’ 27.454” W) using the protocol outlined in Vargas et al. (2016). The experiment consisted of two factors associated with drought treatments: well-watered (WW, 100% of container capacity) and water deficit (WD, 50% of container capacity). The experimental layout was a split plot with five replications, arranging drought treatments as the main factor (2) and the genotypes (6) as the subfactor. The experimental unit was a pot containing 20 seeds arranged randomly for each block. The seeds were sown in rectangular pots (18 x 48 x 14 cm) with a substrate mixture based on mature compost and sand (7:3 v/v). A dose of 10 g of slow-release fertilizer Basacote® plus 6 M (16% N, 8% P2O5, 12% K2O) was applied both at sowing and at the tillering stage, as described by Vargas et al. (2016).
The initial irrigation was the same for all pots. Then, the pots were weighed, and their volumetric water content (VWC) was monitored daily using a Procheck reader with a GS2 METER Group, Inc. USA FDR (Frequency Domain Reflectometry) probe. Watering to replenish the VWC was performed when the volumetric moisture content of the WW 100% treatment decreased to 75% and that of the DW 50% treatment decreased to 15% of the VWC in the pots. The total amount of water used was calculated as the difference between the initial and final weight of the container. The WD treatment began from the full stage of plant tillering (stage 32 of the scale of Zadoks) until the end of the crop.
Gas exchanges and leaf pigments
The net photosynthesis rate (Pn; μmol CO2 m-2 s-1), intracellular CO2 content (Ci; μmol CO2 mol air-1), transpiration (T; mol H2O m-2 s-1), and stomatal conductance (gs; mol H2O m-2 s-1) were measured using a portable Infrared Gas Analyzer (IRGA) LI-6400/LI-6400XT across treatment between 8:00 to 13:00 on sunny days. All measurements were taken using the fully expanded flag leaf from the appearance of the first awns to the fully emerged spike. The variables CO2 (400 ppm), temperature (20°C), and relative humidity (~20%) were monitored. Moreover, the greenhouse radiation at the start of the measurements was close to 200 umol m-2 s-1, therefore, this radiation was established as a reference. Leaf pigments (chlorophyll content [CC] and nitrogen balance index [NBI]) were non-destructively measured in three flag leaves per pot using a handheld Dualex® Scientific instrument (Force A DX16641, Paris, France). These parameters were measured at three phenological stages: at the end of stem elongation (corresponding to stage 49 of the Zadoks scale, Z49, Zadoks et al., 1974), at the beginning of flowering (stage 61 of the Zadoks scale, Z61), and at the beginning of grain filling (stage 73 of the Zadoks scale, Z73). The leaf water potential (LWP; Scholander et al., 1965) and the relative water content (RWC; Barrs and Weatherley, 1962) were also calculated.
Yield components and data analysis
The different plant tissues were harvested separately per pot and dried at 70°C for two days, and weights were recorded until they remained constant. The following yield components were measured: number of grains, number of empty and productive spikes, number of tillers, number of spikes, 1000-grains weight, pot yield, harvest index, and total dry biomass. A general linear model (ANOVA) analyzed physiological parameters and yield component data. For comparisons of means, the LSD-Fisher test was used with a significant level of 0.05. Additionally, the t-student test was performed to compare treatments. All analyses were performed using the R program (R Core Team, 2023).
Multi-omics approach
Two contrasting genotypes for drought tolerance were identified based on their physiological response and yield components under water stress conditions using data from the 2021–2022 season. Therefore, a new experiment using the drought-tolerant (DT) and drought-susceptible (DS) genotypes was carried out under greenhouse conditions as described above in the “yield tolerance index and experimental design” section to perform the transcriptomic and metabolomic analyses.
RNA extraction and transcriptomic analysis
For the RNA extraction, flag leaf and root tissues were collected at three sampling times: 0, 14, and 30 days after anthesis (61, 75, and 85 on the Zadoks scale, respectively) from the two contrasting genotypes under well-watered (WW) and water-deficient (WD) conditions. Total RNA was extracted from each sample using the NucleoSpin RNA Plant Kit (Macherey-Nagel, Düren, Germany), following the manufacturer’s instructions. RNA quality was evaluated using the Qubit™ RNA IQ Assay Kit with a Qubit™ 4 Fluorometer (Thermo Scientific) and RNA quantity was assessed with a NanoDrop™ Lite spectrophotometer (Thermo Scientific). RNA samples were sequenced using Illumina NovaSeq6000 with an average of 20 million reads per sample (Macrogen, Inc. Seoul). The raw reads were subjected to quality filtering using the FastaQC (v0.11.9) and MultiQC (v1.12) software. Low-quality sequences and adapters were trimmed with the TrimGalore v.0.6.5 (https://github.com/FelixKrueger/TrimGalore), considering those with a Q value > 30 and a sequence size greater than 50 bp. Subsequently, contamination was removed using the Kraken2 software, and a library was built based on high-quality sequences. Filtered reads were aligned against the Triticum aestivum reference transcriptome IWGSC RefSeq v1.0 (https://www.wheatgenome.org/) employing the STAR v2.07 software. The differential expression of genes (DEG) was proceeded by DESeq2 (Anders and Huber, 2010). Genes were differentially regulated with a Log2-fold change ≤ - 1 or ≥ 1 and false discovery rate (FDR)< 0.05. The functional annotations of DEGs, the gene sets of Kyoto Encyclopedia of Genes and Genome (KEGG) pathways, and Gene ontology (GO) terms based on biological processes (BP) were employed in gene set enrichment analysis.
Polar metabolite profiling
Leaf root and seed tissues were collected at three sampling times: 0, 14, and 30 days after anthesis (61, 75, and 85 on the Zadoks scale, respectively) in the two contrasting genotypes under WW and WD conditions. The extraction and derivatization of polar metabolites were determined according to the protocol of Hatoum et al. (2014), slightly modified by Fuentealba et al. (2017). Briefly, 100 mg of frozen tissue powder was mixed with 500 μL of cold methanol, and 20 μL of 2910 ng L−1 phenyl β-D-glucopyranoside, and the sample was incubated at 70°C for 15 min with shaking (Labnet International Inc., Edison, NJ, USA). After centrifugation at 17–000 g for 20 min, 100 μL of supernatant was dried under a stream of nitrogen gas. The derivatization consisted of methoximation and trimethylsilylation reactions. GC-MS methods were performed, one for more concentrated compounds, such as sugars, and the other for less concentrated compounds, such as organic and amino acids. For both methods, the injector and interface temperatures were 220°C and 280°C, respectively, and 1 μL of sample was injected. Helium was used as carrier gas with a constant flow of 1 mL min-1. Mass spectra in the 50–600 m/z range were recorded at a scanning speed of 2.66 scan cycles per second. The MS ion source and quadrupole temperatures were 230°C and 150°C, respectively. The method for more abundant compounds included an injection with a split ratio of 1:150, and the oven temperature was programmed to start at 120°C (for 1 min), increased to 300°C at 10°C min−1, and then held for 6 min. The method for less abundant compounds used a splitless injection mode, and the oven temperature was set to start at 50°C (for 1 min), increased to 310°C at 10°C min−1, and then held for 13 min. The chromatographic peaks were deconvolved and identified by comparing retention times and mass spectra to a home-built library of commercial standards and the NIST14 library using Mass Hunter Quantitative software (Agilent Technologies, Santa Clara, CA, USA). To obtain a relative response of each compound, the peak area data were corrected using the peak area of phenyl β-D-glucopyranoside (as internal standard), the sample dry weight, and a quality control (QC) sample representative of all samples.
Partial Least Squares Discriminant Analysis (PLS-DA) was performed on the normalized dataset using the MetaboAnalyst 5.0 software (https://www.metaboanalyst.ca/MetaboAnalyst/). PLS-DA analysis was used with metabolites as predictor variables and the drought and genotype treatments as response variables. Variables were mean-centered and weighted by the standard deviation to assign an equal variance. Variable Important in Projection (VIP) scores were used to filter PLS-DA results to determine essential features.
Transcriptomics and metabolomics integration
Normalized expression/abundance values of transcripts and metabolites were scaled considering the mean-centered and divided by the standard deviation of each variable, which was used to perform multi-omics integration analysis. The R package mixOmics v6.20.0 (Rohart et al., 2017) was used to perform the multiblock PLS-DA, considering the average representation space between all datasets with the block.plsda function described for the DIABLO software (Singh A. et al., 2019). Correlations between omics datasets (transcriptomics and metabolomics) were calculated and plotted using the plotDiablo function of the same package. Finally, the tune.block.splsda function was used to identify the minimum number of variables that explain the sample dispersion of each variate, and correlations between selected variables were plotted using the circosPlot function considering a correlation cutoff = 0.8.
ResultsPhenotypic variation for yield-related traits
A panel of 225 durum wheat genotypes were grown under rainfed and irrigated conditions in two
locations (Los Tilos and Santa Rosa, central and south of Chile, respectively) for two consecutive years (the 2022 and 2023 seasons), and five yield-related traits were measured, such as spikes per square meter (SPM), grains per spike (GPS), thousand grain weight (TGW), grain yield (GY), and harvest index (HI). The combined analysis of variance (ANOVA) across all sites showed significant effects (p< 0.05) of genotypes (G), environments (E), and their interactions (G×E) for all traits, except for grains per spike (GPS) and spike per square meter (SPM), in which the G×E factor was not significant (Supplementary Table S7). The E effect explained most of the variation for grain yield (GY; 79.8%) and 1000-grains weight (TGW; 63.2%), harvest index (HI; 55.8%), and GPS (54.7%), whereas the G×E effect explained the most significant variation for SPM (45.5%). Moreover, the G×E interaction contributed more to the total variability than the G effect for all traits except for TGW, which is in accordance with the value of its broad-sense heritability value of H2 = 0.882. Under rainfed conditions, the E effect explained over 52% of the variation for all yield-related traits, indicating that the environment has a key impact on durum wheat yield. Phenotypic variation was found for all traits across the 8 environments (Supplementary Figure S2). The highest average GY was recorded in the SR_2022 site (9.46 t ha-1) under
irrigated conditions, whereas SR_2023 had the highest average GY (4.3 t ha-1) under rainfed conditions. On average, a 36% reduction in yield was observed at the two locations during the two growing seasons due to drought (Supplementary Table S8). The grand mean, variance components, and heritability for each trait are shown in Table 1. The phenotypic variances were higher than the genotypic variances in all the traits in the combined analysis across all environments. However, the genotypic variances were higher than their corresponding environmental variances for GPS (7.95) and TGW (10.02), indicating that these traits have greater genetic variability within the traits evaluated in the durum wheat panel. In addition, the GPS and TGW were the ones that presented the highest heritability with values of H2 = 0.825 and H2 = 0.882, respectively. These results show that GPS and TGW have greater potential for genetic gains due to selection in breeding programs using this durum wheat panel. On the other hand, the harvest index (HI) was the trait that presented the lowest genetic variation (0.000) and heritability (H2 = 0.416), being a trait with no potential to be used as a selection criterion to obtain genetic gains.
Grand mean, variance components, and heritability for grain yield (GY), spikes per square meter (SPM), harvest index (HI), grains per spike (GPS), and thousand grain weight (TGW) under irrigated, rainfed, and across all environments.
The genetic and phenotypic relationship among traits under irrigated, rainfed, and across
environments showed significant (p< 0.05) correlations among most of the traits (Supplementary Table S9). Under well-water (WW) conditions, the strongest genetic positive associations were obtained between GPS and HI (r = 0.999; p< 0.05) and GPS and GY (r = 0.611; p< 0.05), while the strongest negative correlation was found between TGW and HI (r = -0.771; p< 0.05). The highest positive phenotypic correlation was found between GPS and GY (r = 0.515; p< 0.05), while the negative one was between TGW and GPS (r = -0.542; p< 0.05). Under water deficit (WD) conditions, the highest positive genetic correlation was between HI and GY (r = 0.924; p< 0.05), while the highest negative correlation was between TGW and GPS (r = -0.771; p< 0.05). The highest positive and negative phenotypic correlations were found between HI and GY (r = 0.445; p< 0.05) and TGW and GPS (r = -0.553; p< 0.05), respectively. Overall, GPS recorded one of the highest and most consistent genetic and phenotypic correlations with GY in both water regimes and across all environments.
Markers distribution and population structure
12,498 SNP markers distributed across the 14 chromosomes, primarily positioned in the regions close to the telomeres, were retained after filtering criteria (Supplementary Figure S3). Out of the total markers, 5,969 and 6,529 SNPs were located on the genomes A and B,
respectively. On genome A, the number of markers varied from 631 to 1,061 on chromosomes 4A and 7A, respectively, with an average of 852 SNPs per chromosome. In the B genome, the number of markers ranged from 644 to 1,134 on chromosomes 4B and 2B, respectively, with an average of 932 SNPs per chromosome. The physical chromosome length ranged from 593.5 to 829.5 Mpb on chromosomes 1A and 3B, respectively. The highest density of markers was found on chromosome 1B with ~1.57 SNP/Mpb (Supplementary Table S10).
Based on the structure and DAPC analyses, using 5,258 SNPs retained after filtering by LD (r2< 0.2), the 225 durum wheat lines were divided into two groups. Population structure analysis indicated K = 2 as the most probable number of subpopulations. Based on the membership coefficient (Qi > 0.8), 18 lines (7 from INIA [Chile’s National Agricultural Research Institute] and 11 from CIMMYT [International Maize and Wheat Improvement Center]) were grouped into group 1, while 194 lines (87 from INIA, 98 from CIMMYT, 3 from PUC, 2 from INRAT, and 4 IRTA) were classified into group 2. This result indicates that the durum wheat lines were not grouped according to their genebank origin. Conversely, five sub-populations were inferred according to DAPC analysis by the origin of the durum wheat lines (Supplementary Figure S4). The first and second principal components (PC) explained 58.3% and 20.5% of the total variation, respectively, in which the PC1 allows the differentiation of lines from the INIA and CIMMYT gene banks. Finally, the results from both Structure and DAPC analyses showed an admixture between the genotypes of the different germplasm banks, mainly those of INIA and CIMMYT, common in durum wheat lines developed in breeding programs. The subsequent association analysis was performed using the two principal components.
Linkage disequilibrium and genome-wide association
The LD decay was similar in both genomes. In genome A, the LD decayed at values less than r2< 0.2 at an approximate distance of 8.71 Mpb; in genome B, the LD decayed decay at 7.37 Mbp (Supplementary Figure S5).
Based on the individual analysis by environment, a total of 772 significant SNPs
(-log10P > 3) distributed across the 14 chromosomes of the durum wheat genome were associated with the five yield-related traits through the eight environments using the MLM model (Supplementary Table S11). Marker-trait associations (MTAs) were 331 and 441 under rainfed and irrigated conditions, respectively. The chromosomes with the highest number of MTAs were 2A (351) and 2B (102), whereas the lowest number of MTAs were 4B (2), 4A (18), and 5A (18) (Figure 1). In this context, 300, 206, 115, 88, and 63 MTAs were identified for TGW, GPS, HI, SPM, and
GY, respectively. The phenotypic variation explained (PVE) for each MTA ranged from 5% (BS00000297_51 on chromosome 2A) to 24.3% (Excalibur_c12446_155 on chromosome 3A), both for HI. A total of 57 and 33 MTAs were stable in more than one environment under irrigated and/or rainfed conditions, respectively. Within the stable MTAs, eight were detected exclusively under rainfed conditions, which were clustered into two genomic regions associated with TGW on chromosomes 2A (QTN_2A.TGW.1_R) and 3A (QTN_3A_TGW.1_R), explaining between 5.11% and 11.07% of the phenotypic variance. Under irrigated conditions, 32 stable MTAs were detected, clustered into twelve genomic regions on chromosomes 1A, 2A, 2B, and 5A. All genomic regions were associated with TGW and GPS, explaining between 5.15% and 10.41% of the phenotypic variation, except a quantitative trait loci (QTL) (QTN_2B_GY.1_I) located on chromosome 2B at 56.78 Mpb that explained between a 5.80% and 7.42% of the GY variation (Supplementary Table S12). Finally, nine constitutive MTAs were identified across all environments (Table 2), which were grouped into three QTL associated with TGW and GPS located on chromosomes 2A (QTN_2A_TGW/GPS.1; QTN_2A_TGW/GPS.2) and 2B (QTN_2B_TGW/GPS.1), explaining between 5.15% and 14.29% of the variation.
Manhattan plot showing the SNPs significantly associated with the yield-related traits evaluated across eight environments. SR, Santa Rosa; LT, Los Tilos; 2022, 2021/22 growing season; 2023, 2022/23 growing season; Rainfed and Irrigated with 50 mm at the end of tillering and the flag leaf stage.
Nine constitutive MTAs were identified across all eight environments, grouped into three drought-tolerant QTLs.
Marker
Trait
Chr
Pos (Mpb)
PVE (%)
Name QTL
GENE-1342_238
TGW;GPS
2A
146.959587
5.15 - 13.6
QTN_2A_TGW/GPS.1
RAC875_rep_c74200_1415
TGW;GPS
2A
146.962374
5.15 - 13.6
wsnp_CAP11_c1787_968410
TGW;GPS
2A
146.962374
5.15 - 13.6
GENE-1319_146
TGW;GPS
2A
151.25509
5.62 - 14.29
TA001351-1193
TGW;GPS
2A
151.258908
5.62 - 14.29
BobWhite_c2532_966
TGW;GPS
2A
154.415628
5.62 - 14.29
Tdurum_contig27887_55
TGW;GPS
2A
154.415628
5.62 - 14.29
wsnp_BG607088A_Ta_1_1
TGW;GPS
2A
156.762886
5.39 - 13.21
QTN_2A_TGW/GPS.2
GENE-1181_694
TGW;GPS
2B
209.088615
5.62 - 14.29
QTN_2B_TGW/GPS.1
Chr, Chromosome; Pos, position; PVE, percentage of the phenotypic variation explained.
Based on the combined analysis by water treatment, the GLM model, as measured by the Q-Q plot, did not show a pronounced deviation from the expected distribution of P values in the tail area, indicating poor control of false positives and negatives. In contrast, while the FarmCPU model better controlled false positives and negatives, the MLM model used in this study also provided adequate control (Supplementary Figure S6). Therefore, our MTAs detected using MLM are reliable. A total of 152 significant MTAs
(-log10P > 3) were found for the five yield-related traits under irrigated conditions, with two for GY, 13 for SPM, 24 for HI, 54 for GSP, and 55 for TGW (Supplementary Table S13). The percentage of phenotypic variation explained by each MTA varied from 5.23%
(Kukri_rep_c108875_85 marker on chromosome 2B) to 16.49% (Tdurum_contig27887_55 marker on chromosome 2A), both markers associated with TGW and GPS. Under rainfed conditions (Supplementary Table S14), a total of 107 MTAs were detected for GY (10), SPM (2), HI (7), GPS (48), and TGW (40).
The marker BS00066873_51 with HI located on chromosome 2B explained the lowest percentage of phenotypic variance, while the Tdurum_contig27887_55 (2B) marker associated with TGW and GPS explained the highest percentage (15.69%). Finally, under the combined analysis of both conditions (Supplementary Table S15), a total of 133 MTAs were detected for GY (2), SPM (2), HI (5), GPS (64), and TGW (60). The percentage of phenotypic variation explained by each MTA varied from 5.13% (Excalibur_c4847_631 marker on chromosome 2A) to 18.03% (Tdurum_contig27887_55 marker on chromosome 2A), both markers associated with TGW and GPS. Interestingly, the Tdurum_contig27887_55 marker located at 156,366,746 bp on chromosome 2B, which was associated with GPS and TGW, it was consistently the one that explained the greatest phenotypic variation under irrigated, rainfed, and in the combined analysis of both conditions.
Physiological effects of drought
The summary of the physiological parameters evaluated for the six genotypes can be visualized in Table 3. According to the results, the G3 (BRESCIA) and G4 (QUC 3678-2016) genotypes were the most and least affected, respectively, due to water deficit (WD) treatment. Therefore, these genotypes are appropriate candidates for studying mechanisms associated with susceptibility (BRESCIA) and tolerance (QUC 3678-2016) in the durum wheat panel.
Effects of water treatment on different physiological parameters evaluated in six contrasting drought tolerant/susceptible durum wheat genotypes.
Category
Parameter
Measure
G1
G2
G3
G4
G5
G6
Gas exchange
Pn
Fully-expanded leaf
*
T
Fully-expanded leaf
*
*
*
Gs
Fully-expanded leaf
*
*
*
Ci
Fully-expanded leaf
*
*
Leaf pigment
CC
End of stem elongation
*
*
Beginning of flowering
*
*
Beginning of grain filling
*
*
*
NBI
End of stem elongation
*
Beginning of flowering
*
*
Beginning of grain filling
*
*
Water status
LWP
End of stem elongation
*
Beginning of grain filling
*
*
RWC
End of stem elongation
Beginning of grain filling
Yield
N° grains
Harvest
*
*
*
*
TGW
Harvest
*
*
Yield
Harvest
*
*
Number of tillers
Harvest
*
*
*
*
*
*Significant differences between treatments according to the Tukey test (p<0.05).
Regarding the gas exchange parameters studied, the genotype and the genotype x water treatment interaction were significant only for stomatal conductance (gs) according to the ANOVA (Supplementary Table S16), in which the genotypes that showed the highest gs under water-deficit (WD) conditions were G1 (0.26 mol H2O m-2 s-1) and G2 (0.24 mol H2O m-2 s-1), while the lowest values of gs were found in G3 (0.03 mol H2O m-2 s-1) and G5 (0.07 mol H2O m-2 s-1). The genotypes that most and least reduced gs due to WD treatment were G3 (75%) and G6 (7.14%), respectively. The G4 genotype had a gs reduction of 66.66%, and, surprisingly, the G1 genotype increased gs by 136.36%. The highest net photosynthetic rate (Pn) was observed in G4 (Pn = 6.3 μmol CO2 m-2 s-1) and G2 (Pn = 8.73 μmol CO2 m-2 s-1) for well-watered (WW) and WD conditions, respectively. The G3 and G4 genotypes reduced Pn due to WD treatment by 32.14% and 37.14%, respectively. According to the t-student test performed, significant differences (p< 0.05) were found between treatments for Ci and T in G3 and G4 genotypes. The G3 genotype was the most affected with respect to the intracellular CO2 concentration (Ci) and transpiration (T) with a reduction of 23.8% and 69%, respectively. While the G4 genotype had a reduction of 14.29% and 60.94% for Ci and T, respectively.
At the end of the stem elongation stage (Z 49), nitrogen balance index (NBI) did not show
significant differences between genotypes and water treatments according to the ANOVA (Supplementary Table S17). The G3 genotype was the most affected with a reduction of 25.38%, while the G4 genotype had a 4% reduction. Chlorophyll content (CC) was significantly reduced due to WD treatment, with genotype G6 standing out with a 26% reduction, followed by G3 with 16%. The G4 genotype was little affected due to WD treatment and had a reduction of 6.32% in CC. During flowering (Z 61, Zadoks et al., 1974), no differences were found between genotypes and between water treatments. However, NBI and CC increased due to WD treatment in the G3 and G4 genotypes. The G3 had an increase of 3% and 5% for NBI and CC, respectively. While the G4 genotype increased by 4.5% and 6% for NBI and CC, respectively. At the beginning of grain filling (Z 73), only 4 genotypes (G1, G2, G3, and G4) were evaluated due to technical limitations. Significant differences were found for NBI and CC only for the factor treatment. The G3 genotype experienced a more significant reduction for NBI (39.76%) and CC (42.69%) due to WD treatment. It is worth noting that G4 had the lowest reduction between treatments associated with NBI (11%) and CC (18%).
At the end of stem elongation (Z 49) and during beginning of flowering or grain filling (Z 60),
no significant differences were observed between genotypes for leaf water potential (LWP). Compared to the control, the G3 genotype significantly reduced its LWP by 57.62% and 53.96% at stages Z 49 and Z 60, respectively. The G4 genotype had a 1.5% reduction at stage Z 49 and an 11.59% increase at stage Z60 for LWP (Supplementary Table S18). For relative water content (RWC), no significant differences were identified between
genotypes and treatments at stages Z 49 and Z 60. However, the genotype x treatment interaction was substantial at stage Z 49. The G3 genotype had a reduction of 1% and 19% at Z 49 and Z 60 stages, respectively. While the G4 genotype had a reduction of 1% in both stages of development. Finally, for the yield-related traits, there were significant differences between genotypes and between treatments, in which the G1, G5, and G3 genotypes were the most affected under WD conditions. In contrast, the G6 and G4 genotypes were the least affected. The grain yield was reduced by 67.5% and 38.09% for G3 and G4 genotypes, respectively (Supplementary Table S19).
Overall gene expression profile under drought stress
Gene transcription from leaf and root tissues at 0-, 14-, and 30- days post anthesis (DPA) under WW and WD conditions were compared to investigate the molecular responses of the two contrasting genotypes of durum wheat to drought previously reported. An average of 27,278,116 reads were generated, of which 26,294,562 were retained after the quality and trimming filters, meaning 96.3% sequencing integrity. In total, 85% of sequences were aligned to the genome of Triticum turgidum Svevo v1.0, within which 81.43% had an alignment to a single locus while 4% showed an alignment to multiple loci. Overall, there was a higher amount of differentially expressed genes (DEGs) in the root compared to the leaf for comparisons between genotypes and treatments at the three developmental stages studied. 11,176 DEGs were identified in the leaf, of which 4,819 were up- and 6,357 down-regulated. In root, 12,229 DEGs were identified, of which 5,489 were up- and 6,740 down-regulated. In the leaf, a more significant number of DEGs were determined by the comparison between genotypes (9,173) than between treatments (2,003). In root, the number of DEGs identified by comparing genotypes and treatments was 6,949 and 5,280, respectively. It was confirmed through principal component analysis (PCA) that the transcriptome data of the leaf could separate the G3 (drought-susceptible, DS) and G4 (drought-tolerant, DT) genotypes by PCA1 (39% of variance) and the drought treatment by PC2 accounting for 22% of variance. At root, there was no clear separation between genotypes and treatments, where PC1 and PC2 explained 36% and 15% of the variance, respectively (Supplementary Figure S7).
The summary of the total number of up-and down-regulated genes identified by the comparison of WD/WW in the two contrasting genotypes and DT/DS under different water conditions at the three developmental stages in the leaf and roots is shown in Figure 2. A large proportion of drought-responsive genes were induced in leaf at 0-DPA in both genotypes after drought treatment, and the number of drought-responsive genes in the DT genotype (1,174 DEGs) was much higher than in the DS genotype (84 DEGs) after drought stress exposure. 151 and 1852 DEGs were identified by comparing WD/WW, with 92 and 569 up-regulated genes and 59 and 1283 down-regulated genes for DS and DT genotypes, respectively (Figure 2A). Moreover, four DEGs were constitutive across 0 to 30-DPA in the DT genotype (Figure 2B), which were TRITD2Bv1G146530 (encodes a plastid-localized arogenate dehydratase involved in phenylalanine biosynthesis), TRITD4Bv1G179760 (encodes a chloroplast RNA-binding protein), TRITD5Av1G000200 (encodes a grain softness protein), and TRITD5Bv1G186950 (Serine protease inhibitor (SERPIN) family protein). Interestingly, the TRITD5Bv1G186950 and TRITD5Av1G000200 genes were down-regulated (Log2FC< -20) at 0-DPA, and then at 14- and 30-DPA were highly up-regulated (Log2FC > 20). Additionally, 583 DEGs were compared to the two genotypes under WD conditions.
A total number of up-and down-regulated genes identified by the comparison of well-water (WW) and water deficit (WD) conditions in the drought-tolerant (DT) and drought-susceptible (DS) genotypes and the comparison DT/DS genotypes under WW/WD conditions at 0-, 14-, and 30-days post anthesis (DPA) in leaf (A, B) and root (C, D).
Different expression patterns across the three developmental stages studied were observed in the root. In the DS genotype, a more significant amount of DEGs were observed at 0-DPA (1,096) and 14-DPA (907), whereas a greater amount of DEGs were observed at 0-DPA (2,041) and 30-DPA (857) in the DT genotype, which were much higher compared to the DS genotype. 2192 and 3088 DEGs were identified by comparing WD/WW, with 841 and 1429 up-regulated genes, 1351 and 1659 down-regulated genes for DS and DT genotypes, respectively (Figure 2C). 11 DEGs were constitutive in the DT genotypes (Figure 2D). The total number of DEGs under WD was higher than that under WW conditions. Finally, 569 and 550 DEGs were constitutive in the three developmental stages under WW and WD conditions, respectively.
Gene classification by gene ontology
Since the reference genome of T. turgidum does not have a robust functional annotation, the functional classification was performed from an alignment against the genome of A. thaliana. Gene ontology (GO) analysis was based on classifying DEGs into various biological processes by overrepresenting these genes within the genome. GO analysis showed genotype-, time-, and tissue-specific results, revealing that some important biological processes occur differently in the two contrasting genotypes evaluated under water stress. However, some biological processes were common, highlighting cell death, response to salicylic acid, and biosynthesis of secondary metabolites.
In leaf tissue, the DS genotype only showed enrichment at 14 DPA, whereas the DT genotype exhibited enrichment at 0 DPA and 30 DPA (Supplementary Figure S8). Additionally, genes associated with the DS genotype were classified into fewer biological processes than the DT genotype, with most genes related to gibberellin response and ribosome biogenesis. In contrast, the functional annotation of genes associated with the DT genotype involved a wide range of biological processes, highlighting cell death, salicylic acid response, immune response, and secondary metabolite biosynthesis. Comparing the two genotypes under WW conditions, the overrepresentation of genes associated with senescence, monocarboxylic acid biosynthesis, aromatic compound catabolism, and cell death was found. Under WD conditions (WD_DT/DS), the enriched processes were senescence, hypoxia, cell death, and salicylic acid response at 0-DPA. In later developmental stages, the persistence of genes associated with senescence and citric acid response was observed, as well as the appearance of genes related to secondary metabolite biosynthesis and metabolic processes associated with hydroxyl group-containing compounds. The DS genotype showed gene enrichment across all three times studied in root tissue. Many genes associated with various biological processes were identified, including salicylic acid response, signaling regulation, secondary metabolite biosynthesis, senescence, response to external stimuli, ABA signaling pathway, and cellular communication regulation (Figure 3A). In the DT genotype, gene enrichment was observed in processes unrelated to drought stress response. In early developmental stages, genes were primarily associated with DNA replication, chromosomal organization, microtubule movement, ribosome biogenesis, cell wall biogenesis, monocarboxylic acid biosynthesis, and ABA signaling pathway. In contrast, in later developmental stages, genes were mainly classified into catabolic processes and senescence (Figure 3B). In both genotypes under WW conditions, gene enrichment was observed in flavonoid biosynthesis, secondary metabolite biosynthesis, and processes related to the plant immune system (Figure 3C). Under WD conditions, gene overrepresentation was associated with immune response, cell death, and salicylic acid response (Figure 3D).
Gene ontology (GO) enrichment analysis in the root. Comparison between well-water (WW) and water deficit (WD) treatments for drought-susceptible (DS; A) and drought-tolerant (DT; B) genotypes. Comparison between DS and DT genotypes under WW (C) and WD (D) conditions.
Metabolite profiling under drought stress
To investigate the metabolic response under drought in selected genotypes (G3 and G4), the relative abundance of polar metabolites in leaf, root, and seed tissues was determined at 0-, 14-, and 30-DPA by gas chromatography/mass spectrometry (GC–MS) (Figure 4). Overall, 43 metabolites were identified and quantified, 37 were identified as known compounds, and 6 metabolites were unknown. Of these metabolites, most were amino acids, sugars, and organic acids. Moreover, the highest number of known compounds were present in leaf (42), then in root (41) and seeds (38). A Partial Least Squares-discriminant Analysis (PLS-DA) showed that the metabolite profiles between DT and DS genotypes differed (Figure 4A), which suggests a distinction in the metabolite accumulation response due to genotypic differences. The score plots explained 29.1% (leaf), 45.5% (root), and 73.1% (seed) of the variability with the two first components. Variable importance of projection (VIP) scores were determined for displaying the 15 most important metabolites involved in the discrimination between genotypes and treatment (Figure 4B). Among the five most contrasting metabolites in the leaf are D-Mannitol, shikimic acid, erythrono-1,4-lactone, (Z), 2-oxoglutaric acid, and L-Serine. The root had mainly amino acids such as L-threonine, L-Serine, L-valine, L-glutamic acid, and shikimic acid. In contrast, a disaccharide unknown, D-Mannitol, shikimic acid, unknown sugar, and malic acid were the five most contrasting metabolites in seeds. Only four metabolites showed higher abundance in the DT genotype under WD conditions in the leaf (citric acid, 2-Keto-l-gluconic acid, glycol, and malic acid). Interestingly, vanylglycol was the only metabolite with a high abundance in roots, while no metabolite had the highest abundance in the DT genotype under WD conditions in seeds.
Partial least squares-discriminant analyses (PLS-DA) with the two first components (A) and variable’s importance (VIP scores) of secondary metabolite profiles (B) of leaf, root, and seed of the drought-tolerant (DT) and drought-susceptible (DS), genotypes under well-water (WW) and water deficit (WD) conditions at 0-, 14- and 30-DPA.
Further analysis of individual metabolites based on fold-change (FC) of metabolite abundance between the comparison of WD/WW in the two genotypes and DT/DS under the two water conditions at different developmental stages in leaf and roots is shown in Supplementary Figures S9, S10. Across the three time points in the leaf, the total number of differential abundance of metabolites (DAMs) was 30 and 35 in the DT and DS genotypes, respectively. Only two (asparagine and quininic acid) and one (L-Proline) metabolites were identified at all three times (constitutive) in DT and DS genotypes due to drought treatment, respectively. Further, 33 and 43 DAMs were determined by comparing DT/DS genotypes, with 24 and 25 up-concentrated and 9 and 18 down-concentrated metabolites for WW and WD conditions, respectively. Under WD conditions, a similar number of DAMs was observed at 0-DPA (13) and 14-DPA (12), while the number was higher (18) in the later stage (30-DPA). Moreover, five (D-mannitol, asparagine, aspartic acid, L-glutamic acid, malic acid) metabolites were constitutive across 0 to 30-DPA in the comparison between the two genotypes under WD conditions, and only citric acid was under WW conditions.
A higher number of DAMs was determined in the root than in the leaf. In this sense, 41 and 49 DAMs were determined across the three phenological stages in the DT and DS genotypes under WD conditions, respectively. Furthermore, more constitutive DAMs were identified in roots compared to leaves. Nine metabolites were constitutive in the DS genotype, which were amino acids (L-Proline, Asparagine, L-Glutamine, L-Threonine, L-Valine, and L-Phenylalanine), sugars (D-Fructose, and D-Glucose) and the polyol D-Mannitol. In the DT genotype, the metabolites constitutive were organic acids (Quininic acid, Aconitic acid, 2-Oxoglutaric acid, and Quininic acid), the amino acid L-Proline, and an unknown disaccharide. Comparing both contrasting genotypes, a similar number of metabolites were up-concentrated under WW (14) and WD (15) conditions.
Interestingly, more metabolites were down-concentrated under WD (41) than WW (8) conditions. Four metabolites were constitutive under WD conditions: organic acids (aconitic, aspartic, and glyceric acid) and L-valine. In contrast, the D-Mannitol and 5-O-Feruloylquinic acid were constitutive under WW conditions. These results indicate that the metabolic response of durum wheat under water stress conditions is modulated at a specific time, indicating the differential abundance of metabolites at different developmental stages. However, constitutive metabolites are up-and down-concentrated at different times (0- to 30-DPA), which can be used as biomarkers associated with drought stress in leaves and roots.
Integration of transcriptomic and metabolomic data
To integrate transcriptomic and metabolomic data and find the lowest number of variables that allow explaining the separation between the DT and the DS genotypes under WD conditions and both genotypes under WW conditions, a sPLS-DA was performed using each sampling time as replicas (Figure 5). The number of DEGs used for the analysis was obtained from the union of a total of DEGs identified that respond to drought in both genotypes and a total of DEGs identified due to a differential response between genotypes. Therefore, the total list of DEGs was 697 and 900 DEGs for leaf and root, respectively. A clear separation was observed between the genotypes under WW condition and the DT genotype under the WD condition with variate 1. In contrast, a separation was observed between the genotypes under WW condition and the DS genotype under the WD condition with variate 2 (Figure 5A). Considering a two-variate analysis, the tuning function was performed to determine the optimal number of variables to select for the leaf and root dataset. For variable 1, which separates the DT genotype, 60 (30 in root and 30 in leaf) transcripts and 16 (6 in root and 10 in leaf) metabolites that best explained the separation were selected (Figure 5B), and their expression/abundance is shown in Supplementary Table S20. Regarding variable 2, which separates the DS genotype, 70 (40 in root and 30 in leaf) transcripts and 28 (15 in root and 13 in leaf) metabolites were selected (Figure 5C), and their expression/abundance is shown in Supplementary Table S21. A circosPlot with the correlations between the selected transcripts and the metabolites for each variate in the DT and DS genotypes is shown in Figures 5B, C. Finally, L-Proline was an important metabolite in both tissues (leaf and root) that allowed differentiating the tolerant genotype from the susceptible one (Figure 5D), suggesting it as a good biomarker.
Leaf and root transcriptomic and metabolomic integration analysis of durum wheat considering two water treatments (WW and WD) at three developmental stages (0-, 14-, and 30-DPA) in susceptible (DS) and tolerant (DT) genotypes. (A) Average sparse partial least square discriminant Analysis (sPLS-DA) with selected variables using leaf metabolomic (ML), leaf transcriptomic (TL), root metabolomic (MR), and root transcriptomic (TR) as omic datasets for integration. Circosplots showing the root and leaf selected metabolites and transcripts for the DT genotype (variate 1; B) and DS genotype (variate 2; C). The expression values of each transcript are plotted on the outside part of the circosplots, considering the average value of all replicates for each condition (WW-mean, DS-WD-mean, and DT-WD-mean). Positive and negative correlations between omic data sets were represented by green and red lines in the center, respectively. In red, transcripts were located on chromosome 2A, where constitutive QTLs were identified. (D) L-Proline abundance in tolerant (DT) and susceptible (DS) genotypes in leaf and root under well-water and water deficit conditions.
To identify possible biomarkers associated with drought tolerance/susceptibility in durum wheat, the variables that exhibited a consistent differential expression/abundance across the three sampling times in each genotype were determined. Regarding the transcriptomic data, 15 and 24 genes had higher differential expression in leaf and roots, respectively. Interestingly, the TRITD6Av1G194890 gene that encodes a Galactoside 2-alpha-L-fucosyltransferase and the L-proline metabolite were significant in both tissues in the DT genotype. In addition, the metabolites such as malic acid, glyceric acid, Myo-inositol, and L−proline were appropriate candidates for leaf, and the sugars such as D-fructose and D-glucose and the amino acid L-proline were appropriate for roots. These results confirm that the L-proline is a very good candidate as a biomarker of tolerance to drought, independent of tissue and developmental stages. In addition, our study allowed us to determine new tissue-specific biomarkers in durum wheat. The catalase was an interesting gene in the DS genotype that responds to drought stress in both tissues. Therefore, it could be a proper biomarker of susceptibility to drought in durum wheat. The metabolites such as 2-oxoglutaric acid, D-mannitol, and an unknown sugar were adequate candidates for leaf. The amino acids such as L-phenylalanine, L-leucine, asparagine, L-isoleucine, L-valine, L-threonine, and b-sitosterol were suitable for roots.
DiscussionGenetic dissection of drought tolerance in durum wheat
In recent decades, climate change has increased drought frequency, duration, and extent during plant development, causing significant grain yield losses worldwide and affecting food security (Bailey-Serres et al., 2019). Uncertainties about the future effects of climate change require breeders to develop varieties that perform well in diverse environments, which is challenging since broad adaptation is hampered by high genotype-environment interactions for most yield-related traits (Zaïm et al., 2024). The ANOVA analysis across eight environments showed a larger contribution of the G×E interaction to the total variability compared to the G effect for all traits except for TGW, which agrees with their heritability as reported in previous studies (Brinton and Uauy, 2019; Arriagada et al., 2020). Thus, developing high-yielding varieties adapted to diverse environments requires the accumulation of loci with additive genetic effects on yield components. In this sense, numerous studies have been carried out to identify quantitative trait loci through QTL or associative mapping related to the main yield components under irrigated and rainfed conditions in durum wheat (Sukumaran et al., 2018; Maccaferri et al., 2019; Arriagada et al., 2022; Zaïm et al., 2024). Most common wheat loci have been localized on chromosomes 2B, 3A, 4A, 4B, 7A, and 7B (Gupta et al., 2017; Goel et al., 2020). In durum wheat, the most significant MTAs associated with yield components under irrigated, rainfed, and heat stress were identified on chromosomes 2A and 2B (Sukumaran et al., 2018). This agrees with our results, where the chromosomes with the highest number of MTAs were 2A (351), 2B (102), and 7B (53), and the traits most reported were TGW and GPS (Zaïm et al., 2024).
The QTLs identified across different environments (stable QTLs) can be considered robust candidates for improving drought tolerance and increasing grain yield (Soriano et al., 2021; Arriagada et al., 2022). This study identified nine constitutive MTAs across eight environments under rainfed and irrigated conditions. These MTAs were grouped into three QTLs associated with TGW and GPS, two QTLs located on chromosomes 2A (QTN_2A_TGW/GPS.1; QTN_2A_TGW/GPS.2) and one on chromosome 2B (QTN_2B_TGW/GPS.1), explaining between 5.15% and 14.29% of the variation. In a previous Meta-QTL (MQTL) analysis, Arriagada et al. (2022) identified three MQTL associated with yield components on chromosomes 2A (2) at 71.4 and 88.7 cM and 2B (1) at 28.9 cM associated with TGW and HI under rainfed and irrigated conditions. Moreover, Soriano et al. (2021) identified an MQTL on chromosome 2A at 50.8 cM related to normalized difference vegetation index (NDVI) and chlorophyll content (SPAD), traits that are associated with grain yield under drought stress. Sukumaran et al. (2018) also identified a stable QTL under irrigation and drought conditions on chromosome 2A at 66–70 cM associated with TGW. These results suggest that regions on chromosomes 2A and 2B are useful in durum wheat breeding programs, in which the goal is to develop varieties adapted to different environments with high grain yields.
Physiological response to drought stress in durum wheat
Reductions of 1% in net photosynthesis (Pn), 8% in intercellular carbon concentration (Ci), and 21% in transpiration (T) were observed in the water deficit (WD) treatment, with stomatal conductance (gs) 9% lower across genotypes. Zhao et al. (2020) reported a 20% decrease in Pn in bread wheat under moderate deficit irrigation, while Bogale et al. (2011) noted significant reductions in Pn and T in durum wheat under WD conditions, aligning with earlier findings (Condon et al., 2002; Allahverdiyev, 2015). Measurements in this study were taken two weeks after starting the WD treatment, and gas exchange and yield components response vary based on the frequency, intensity of water stress, and wheat growth stage and may not have been reflected to a full extent.
The G3 genotype exhibited significant susceptibility to water deficit (WD), with all gas exchange parameters declining and lower NBI and Chlorophyll Content (CC). In contrast, G4 genotype showed a 66% reduction in gs and a 61% decrease in T but maintained stable Ci and Pn and higher NBI and CC, linked to the chromosome 2 findings previously discussed (Sukumaran et al., 2018; Soriano et al., 2021; Zaïm et al., 2024). Pn reduction is primarily due to stomatal closure under water stress, linked to increased ABA in the xylem (Liu et al., 2015). Still, non-stomatal limitations also arise from lowered enzymatic efficiency or activation, including RuBP carboxylation and regeneration (Flexas and Medrano, 2002). Water stress can also damage the oxygen-evolving complex and PSII reaction centers (Sainz et al., 2010). Severe drought slows ATP renewal as proton concentration rises in the thylakoid lumen, leading to more PSII heat dissipation. Research indicates that Rubisco activity and RuBP reduction do not limit Pn until drought becomes severe (Boutraa et al., 2010). Also, gas exchange parameters may not change significantly until the volumetric water content of the soil decreases by more than two-thirds (Sinclair, 2005).
Transcriptomics response to drought in durum wheat
The stress signaling mechanism and its regulatory networks that control the expression of critical genes are the basis of plant tolerance to drought (Janiak et al., 2016). In general, it was observed that the response to drought stress in durum wheat acted in a tissue- and time-specific manner since neither the number of DEGs nor their functional annotation was standard between leaf and root in the different sampling times studied in this work. There was a higher amount of DEGs in the root compared to the leaf for comparisons between genotypes and treatments at 0-, 14- and 30- DPA. Moreover, most drought response genes were regulated at 0-DPA in both genotypes, where the number of DEGs in the DT genotype (1,174) was much higher than in the DS genotype (84). This result is expected since the plant responds to drought stress by sensing the water shortage in the soil through its roots (Janiak et al., 2016). Then, a signaling cascade is activated that sends chemical signals to the tissues of the aerial part, such as shoots and leaves (Kuromori et al., 2022), to initiate a response that leads to physiological, biochemical, and morphological changes that allow plants to tolerate stress conditions (Janiak et al., 2016).
Drought stress modulates molecular responses associated with the biosynthesis of antioxidants, photosystem components, ABA metabolism, biosynthesis of secondary metabolites, and metabolism of carbohydrates and tricarboxylic acids in durum wheat (Zuluaga et al., 2023). This study’s GO analysis showed genotype-, time-, and tissue-specific results, revealing that some important biological processes occur differently in the two contrasting genotypes [G3 (BRESCIA, drought-susceptible) and G4 (QUC 3678-2016, drought-tolerant)] evaluated under water stress. However, some biological processes were common, including cell death, response to salicylic acid, and biosynthesis of secondary metabolites. Under drought stress, the salicylic acid (SA) phytohormone regulates the antioxidant defense system, transpiration and photosynthesis rates, nutrient uptake, stomatal movement, nitrogen metabolism, cell elongation, proline metabolism, inhibiting lipid oxidation, and deactivating reactive oxygen species (ROS) (Khan et al., 2015; Khalvandi et al., 2021; González-Villagra et al., 2022). In wheat, the SA induced responses against drought by increasing antioxidative enzymes and osmolytes such as proline and total soluble sugars (Sharma et al., 2017; Khalvandi et al., 2021).
Additionally, plants with exogenous applications of SA have higher biomass, higher chlorophyll content, higher water content, higher Rubisco carboxylase activity, and higher superoxide dismutase (SOD) activity relative to untreated plants (Singh and Usha, 2003). On the other hand, plants synthesize secondary metabolites such as flavonoids, spermidine, and plant hormones to scavenge ROS and protect themselves from lipid peroxidation when subjected to adverse environmental conditions (Yadav et al., 2021). Direct regulation of genes in the biosynthetic pathway of secondary metabolites or indirect regulation through phytohormone regulation provides the plant with drought tolerance by avoiding oxidative stress (Jogawat et al., 2021). Therefore, phytohormones such as salicylic acid and secondary metabolites play a key role in the drought tolerance of durum wheat. Finally, a recent study indicates that post-anthesis drought stress greatly influences wheat grain development and programmed cell death in the endosperm. The expression of genes related to the glutathione-ascorbic acid cycle pathway (DNA repair) and the mitogen-activated protein kinase (MAPK) signaling cascade pathway (membrane functions) were key pathways involved in regulating the programmed cell death at two key stages of endosperm development, 7 days and 14 days post anthesis (Li et al., 2024), which agrees with our results.
At the gene level, four DEGs showed a constitutive response to drought stress from 0- to 30-DPA in leaves of the DT genotype. However, only two genes (TRITD2Bv1G146530 and TRITD4Bv1G179760) were up-regulated all three times sampled. The TRITD2Bv1G146530 gene located on chromosome 2B (434782763–434784101 bp) encodes a plastid-localized arogenate dehydratase (ADT) involved in phenylalanine biosynthesis. The amino acid phenylalanine serves as a building block for proteins and as a main precursor for the biosynthesis of phenylpropanoids that play a key role in plant growth, development, and the response to environmental signals (El-Azaz et al., 2020). In Arabidopsis, the arogenate dehydratase gene (ADT2) plays an essential role in seed development because an ADT2 deficiency causes seed abortion (El-Azaz et al., 2018), which in wheat can translate into a loss of grain yield. In common wheat, the arogenate dehydratase 5 (ADT5) expression was increased in leaf tissues of a sensitive and tolerant cultivar after 4 h drought stress. However, the expression level of ADT5 in the leaf was increased by eight-fold and decreased by two-fold in the tolerant and susceptible cultivars, respectively, after 8 h of drought stress (Cevher-Keskin et al., 2023). The TRITD4Bv1G179760 gene on chromosome 4B (607000032 – 607006–701 bp) encodes a chloroplast RNA-binding protein (RBP). The RBPs govern many aspects of RNA metabolism, including pre-mRNA processing, transport, stability/decay, and translation, and are also involved in regulating the response to abiotic stresses (Ambrosone et al., 2012; Muthusamy et al., 2021; Yan et al., 2022), including drought in wheat (Xu et al., 2014; Sun et al., 2020). Transgenic plants were generated in rice to investigate the potential use of RBPs in the development of drought-tolerant cultivars. The expression of the Arabidopsis glycine-rich RNA-binding protein (AtGRP2 or AtGRP7) caused an increase in the number of grains per panicle, therefore increasing grain yield under drought conditions (Yang et al., 2014). Our results show that the ADT and RBP genes were up-regulated at 0-, 14-, and 30 DPA in the DT genotype under WD conditions, indicating that these genes can be keys in the drought response in durum wheat, which has not been previously reported.
In root, 11 DEGs were constitutive in the DT genotypes under WD conditions (Supplementary Table S22). Among them, three genes were highly up-regulated (Log2FC > 30). The TRITD1Av1G135510 gene encodes a bifunctional inhibitor/lipid-transfer protein/seed storage 2S albumin superfamily protein related to the lipid metabolic process (Quinga et al., 2018). This gene was up-regulated under osmotic and drought stress in barley (Paul et al., 2023) and wheat (Pál et al., 2018), respectively. In addition, it is highly up-regulated in phellem cells, participating in suberization and/or secondary cell wall development in the roots of Arabidopsis (Leal et al., 2022). The TRITD2Bv1G203390 gene that encodes a proline-rich protein DC2.15. The proline-rich proteins (PRP) are stress-induced and are involved in cell wall modification and organization (Fan et al., 2014). In a transcriptomic analysis of DEGs in two alfalfa cultivars with different salt tolerance, the MS.gene011517 (14 kDa proline-rich protein DC2.15) had higher and consistent expression under salt stress (Bhattarai et al., 2021). In addition, PRP may regulate free cellular proline levels during drought stress to provide drought tolerance in tomatoes (Gujjar et al., 2018). Finally, the TRITD4Av1G258040 gene that encodes a Serine/threonine-protein kinase (S/T PK) was also highly up-regulated in the root of durum wheat. Protein kinases are essential enzymes that regulate many metabolic pathways through protein phosphorylation/dephosphorylation to regulate cellular functions such as signal transduction and response to abiotic stress (Mao et al., 2010). The S/T PK is important in maintaining root morphology during stress conditions (Muhammad et al., 2019). In wheat, several S/T PKs involved in drought tolerance have been characterized, such as TaSnRK2.7 (Zhang et al., 2011), TaSnRK2.4 (Mao et al., 2010), PK4 (Rahimi et al., 2021), RPK1 (Rahim et al., 2022). In durum wheat, the only reported and characterized gene that encodes for an S/T PK and is associated with drought tolerance is Td4IN2 (Rampino et al., 2017).
Metabolic response to drought in durum wheat
To explore the metabolic response of durum wheat grown under drought, GC–MS in leaf and root tissues determined polar metabolites relative abundances. Metabolites play a fundamental role in drought tolerance since they are considered signaling molecules (Obata and Fernie, 2012). Overall, 43 metabolites were identified and quantified, most of them were amino acids, sugars, and organic acids. This result agrees with other metabolic studies in wheat under drought stress conditions, where the most abundant identified metabolites are amino acids, organic acids, and sugars (Solomon and Labuschagne, 2005; Guo et al., 2018; Marček et al., 2019; Fadoul et al., 2021). Amino acids are essential plant metabolites for protein synthesis and cellular function (Ahanger et al., 2018). Free amino acids act as osmoprotectants and improve drought-induced nitrogen uptake through reassimilation and maintaining protein homeostasis (Živanović et al., 2020). In addition, they also function as scavengers of reactive oxygen species (ROS) generated in plants under drought conditions (Kumar et al., 2021). Organic acids are important in energy metabolism because they are involved in the oxidative and synthetic pathways of carbohydrates and fatty acids (Fadoul et al., 2021). Besides this, they are precursors of amino acids and increase the osmotic potential of the cell, which implies an increased water-holding capacity (Martinelli et al., 2013; Zeng et al., 2021). Finally, sugars play an important role in plant metabolism as sources of carbon and energy in cells, the maintenance of osmotic balance, as signaling molecules for photosynthesis, and detoxification of ROS by acting as metabolic signals in stress conditions (Li C. et al., 2016; Ahmad et al., 2020).
Five metabolites showed differential abundance across 0- and 30-DPA when genotypes under WD conditions in the leaf were compared. These metabolites were D-mannitol, asparagine, aspartic acid, L-glutamic acid, and malic acid. The amino acid asparagine is central to nitrogen storage and transport due to its relatively high nitrogen-to-carbon ratio and non-reactive nature (Curtis et al., 2018). Further, it accumulates in high concentrations during processes such as seed germination and in response to various abiotic stresses due to its potential role as an osmolyte (Carillo et al., 2005; Lea et al., 2007). In this study, asparagine levels were higher in the DS than in the DT genotype under drought stress, consistent with rice observations (Degenkolbe et al., 2013) and wheat (Yadav et al., 2019). The asparagine accumulation may be associated with a biochemical response to dehydration, likely senescence, that has negative consequences for yield (Yadav et al., 2019), in agreement with our results where senescence was an enriched process in the transcriptomic analysis comparing genotypes under WD conditions. Glutamic acid participates in different biological processes, playing an important role in the biosynthesis of chlorophyll and proline (Franzoni et al., 2021). It has been shown that the exogenous application of glutathione protected wheat seedlings from water deficit damage, increasing chlorophyll pigments and antioxidant enzyme activity and promoting the root system, shoot growth, and relative leaf water content (Suliman et al., 2024). Moreover, the overexpression of the cytosolic glutamine synthetases (GS) 1 and plastidic GS2 genes from wheat enhanced drought tolerance of both root and leaf tissues in tobacco (Yu et al., 2020). At 0-DPA, levels of L-glutamic acid were higher in the DS than in the DT genotype. However, at 14- and 30-DPA the levels were higher in DT than in the DS genotypes. Therefore, the higher proline contents in the DT genotype could indicate the greater tolerance to drought stress, induced by glutamic acid.
In root, four DAMs were constitutive in the comparison between contrasting genotypes under WD conditions: organic acids (aconitic acid, aspartic acid, and glyceric acid) and L-Valine. Aconitic acid is an intermediate in the glyoxylate and tricarboxylic acid (TCA) cycles produced during the conversion of citric acid to isocitrate (Min et al., 2018). In wheat, its positive regulation increases energy for the defense against oxidative stress (Liu et al., 2020). In this study, aconitic acid levels were higher in the DS genotype across the three developmental stages under WD conditions. However, the aconitic acid levels decreased significantly after 0-DPA in DS genotype with drought stress, which agrees with previous findings in common wheat (Guo et al., 2018). Aspartic acid is obtained from a transamination reaction between glutamate and oxaloacetate, which then produces the amino acids lysine, threonine, methionine and isoleucine, in a series of reactions known as the aspartic acid metabolic pathway (Alfosea-Simón et al., 2021). In common wheat, exogenous application of aspartic acid has alleviated salt stress-induced growth decline by enhancing antioxidants, compatible solutes, and reducing reactive oxygen species (Sadak et al., 2022). Finally, Weng et al. (2024) showed that eight metabolites (L-tryptophan, L-valine, L-leucine, glycerol 1-phosphate, L-threonine, gluconic acid lactone, malonic acid, D-malic acid) were increased and two (quinic acid, aspartic acid) were decreased constitutively in a DT genotype of common wheat under drought stress conditions, similar as it was found in this study.
Integration of transcriptomic and metabolomic data
Tolerance to drought involves activating a series of complex mechanisms by the plant. This begins with the detection of stress, which activates a signal transduction pathway that activates stress response genes that are translated into proteins and compounds whose main function is to prevent water loss in the cells. Therefore, the integration of multi-omics data from genomics, transcriptomics, proteomics and metabolomics analysis allows to obtain a comprehensive understanding of these complex biological processes and to identify the key molecular players involved in the drought response of the plant (Roychowdhury et al., 2023). In this sense, to identify possible key biomarkers associated with drought tolerance/susceptibility in durum wheat, the variables that showed consistent differential expression/abundance across 0- to 30-DPA in each genotype were determined.
Regarding the transcriptomic data of the DT genotype on leaf, the TRITD2Av1G062130 gene located on chromosome 2A (140,485,202-140,488,307 bp) that encodes a WRKY transcription factor (WRKY3) is positioned on the stable QTL “QTN_2A_TGW/GPS.1”, whose SNP marker in LD is GENE-1342_238. The WRKY proteins constitute a large family of transcription factors involved in the regulation of development, senescence and stress tolerance, including drought (Zhu et al., 2013; Wang et al., 2023). Several studies in wheat have shown that overexpression of WRKY transcription factors such as TaWRKY1, TaWRKY10, TaWRKY33, TaWRKY40, TaWRKY75 enhances drought tolerance (Yu et al., 2023). In root, four identified genes were located on chromosomes 2A and 2B, where stable QTLs were located. The TRITD2Av1G187060 (Amino acid permease), TRITD2Av1G178200 (Uridylate kinase), TRITD2Av1G265880 (Organic cation transporter protein), and TRITD2Bv1G172660 (DNA-directed RNA polymerase II, putative) genes, however, none of them was present within the regions of the detected stable QTLs. Interestingly, TRITD6Av1G194890 gene that encodes a Galactoside 2-alpha-L-fucosyltransferase was important in both tissues in the DT genotype, in agreement with studies showing that this enzyme participates in drought tolerance (Li S. et al., 2016; Moradi Tarnabi et al., 2020; Sui and Wang, 2020). Regarding metabolites in the DT genotype, the malic acid, glyceric acid, myo-Inositol, and L−proline were good candidates for leaf, and the sugars such as D-fructose and D-glucose and the amino acid L-proline were appropriate for roots. This result confirms that the L-proline is an outstanding candidate as biomarker for drought tolerance, independent of tissue and developmental stages. Proline accumulation functions as an electron sink mechanism that can reduce the amount of singlet oxygen present due to drought stress, contributing to cellular redox balance (Ullah et al., 2017). In the DS genotype, the catalase (CAT2) was an interesting gene that responds to drought stress in both tissues. Therefore, it could be a promising biomarker of susceptibility to drought in durum wheat. The CAT gene expression regulates growth, development, and response to environmental stress (Zhang et al., 2022). The CAT gene family encodes an antioxidant enzyme that modulates the hydrogen peroxide (H2O2) contents in cells by translating this toxic compound into water (H2O) and O2− to reduce ROS contents in cells, produced because of drought stress (Ghorbel et al., 2023). Finally, the metabolites such as 2-oxoglutaric acid, D-mannitol and an unknown sugar were good candidates for leaf, and the amino acids such as L-phenylalanine, L-leucine, asparagine, L-isoleucine, L-valine, and L-threonine, and the b-Sitosterol were appropriate for roots.
To our knowledge, this is the first study evaluating different omics approaches and their integration to decipher drought tolerance in durum wheat. The genomic approach identified nine stable QTLs under rainfed and irrigated conditions across eight environments. These QTLs were located on chromosomes 2A and 2B and associated with TGW and GPS traits, suggesting that these regions may be useful in durum wheat breeding programs, where the goal is to develop varieties adapted to different environments with a high yield. Through independent transcriptomic and metabolomic analysis, it was possible to identify genes and metabolites that show differential expression/abundance across the different developmental stages (0 to 30 DPA) in the drought-tolerant genotype of durum wheat. Finally, by integrating omics, a set of leaf- and root-specific genes and metabolites were identified, allowing differentiation between the tolerant and susceptible genotypes of durum wheat. This result permits the determination of new tissue-specific biomarkers for drought-stress responses in durum wheat through their characterization and validation in future studies. These results are valuable and novel in durum wheat breeding programs for developing resilient and high-yielding cultivars to ensure food security.
Data availability statement
The original contributions presented in the study are publicly available. This data can be found here:NCBI, PRJNA1261096.
The author(s) declare that financial support was received for the research and/or publication of this article. This work was supported by the Agencia Nacional de Investigación y Desasrrollo de Chile (ANID) grants Fondecyt Regular n° 1210092 (Andrés R. Schwember, Iván Matus, Carlos Maureira, Samantha Reveco and Valentina Villarroel) and Fondecyt Postdoctorado n°3200981 (Osvin Arriagada). This research was also partially funded by the project ANID – Millennium Science Initiative Program – ICN2021_044 (Romina Pedreschi and Claudio Meneses).
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declare that no Generative AI was used in the creation of this manuscript.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2025.1540179/full#supplementary-material
ReferencesAhangerM. A.GulF.AhmadP.AkramN. A. (2018). “Environmental stresses and metabolomics-deciphering the role of stress responsive metabolites,” in Plant Metabolites and Regulation under Environmental Stress. Eds. AhmadP.AhangerM. A.SinghV. P.TripathiD. K.AlamP.AlYemeniM. N. (Academic Press, Cambridge, MA, USA), 53–67.AhmadF.SinghA.KamalA. (2020). “Osmoprotective role of sugar in mitigating abiotic stress in plants,” in Protective chemical agents in the amelioration of plant abiotic stress: biochemical and molecular perspectives. Eds. RoychoudhuryA.TripathiD. K. (John Wiley & Sons), 53–70. doi: 10.1002/9781119552154.ch3AkhiyarovaG.VeselovD.IvanovR.SharipovaG.IvanovI.DoddI. C.. (2023). Root ABA accumulation delays lateral root emergence in osmotically stressed barley plants by decreasing root primordial IAA accumulation. Int. J. Plant Biol.14, 77–90. doi: 10.3390/ijpb14010007Alfosea-SimónM.Simón-GraoS.Zavala-GonzalezE. A.Cámara-ZapataJ. M.SimónI.Martínez-NicolásJ. J.. (2021). Physiological, nutritional and metabolomic responses of tomato plants after the foliar application of amino acids aspartic acid, glutamic acid and alanine. Front. Plant Sci.11. doi: 10.3389/fpls.2020.581234AllahverdiyevT. I. (2015). Effect of Drought Stress on Some Physiological Traits of Durum (Triticum durum Desf.) and Bread (Triticum aestivum L.) Wheat Genotypes. J. Stress Physiol. Biochem.11, 29–38.AlvaradoG.RodriguezF. M.PachecoA.BurgueñoJ.CrossaJ.VargasM.. (2020). META-R: A software to analyze data from multi-environment plant breeding trials. Crop J.8, 745–756. doi: 10.1016/j.cj.2020.03.010AmbrosoneA.CostaA.LeoneA.GrilloS. (2012). Beyond transcription: RNA-binding proteins as emerging regulators of plant response to environmental constraints. Plant Sci.182, 12–18. doi: 10.1016/j.plantsci.2011.02.004AndersS.HuberW. (2010). Differential expression analysis for sequence count data. Nat. Prec. doi: 10.1038/npre.2010.4282.1ArriagadaO.GadaletaA.MarcotuliI.MaccaferriM.CampanaM.RevecoS.. (2022). A comprehensive meta-QTL analysis for yield-related traits of durum wheat (Triticum turgidum L. var. durum) grown under different water regimes. Front. Plant Sci.13. doi: 10.3389/fpls.2022.984269ArriagadaO.MarcotuliI.GadaletaA.SchwemberA. R. (2020). Molecular mapping and genomics of grain yield in durum wheat: A review. Int. J. Mol. Sci.21, 7021. doi: 10.3390/ijms21197021Bailey-SerresJ.ParkerJ. E.AinsworthE. A.OldroydG. E. D.SchroederJ. I. (2019). Genetic strategies for improving crop yields. Nature575, 109–118. doi: 10.1038/s41586-019-1679-0BaldoniE.FrugisG.MartinelliF.BennyJ.PaffettiD.ButiM. (2021). A comparative transcriptomic meta-analysis revealed conserved key genes and regulatory networks involved in drought tolerance in cereal crops. Int. J. Mol. Sci.22, 13062. doi: 10.3390/ijms222313062BarrsH. D.WeatherleyP. E. (1962). A re-examination of the relative turgidity technique for estimating water deficits in leaves. Aust. J. Biol. Sci.15, 413–428. doi: 10.1071/BI9620413BashirS. S.HussainA.HussainS. J.WaniO. A.Zahid NabiS.DarN. A.. (2021). Plant drought stress tolerance: Understanding its physiological, biochemical and molecular mechanisms. Biotechnol. Biotec. Eq.35, 1912–1925. doi: 10.1080/13102818.2021.2020161BatesD.MächlerM.BolkerB.WalkerS. (2014). Fitting linear mixed-effects models using lme4. J. Stat. Softw.67, 1–48. doi: 10.18637/jss.v067.i01BhattaraiS.FuY. B.CoulmanB.TaninoK.KarunakaranC.BiligetuB. (2021). Transcriptomic analysis of differentially expressed genes in leaves and roots of two alfalfa (Medicago sativa L.) cultivars with different salt tolerance. BMC Plant Biol.21, 446. doi: 10.1186/s12870-021-03201-4BogaleA.TesfayeK.GeletoT. (2011). Morphological and physiological attributes associated to drought tolerance of Ethiopian durum wheat genotypes under water deficit condition. J. Biodivers. Environ. Sci.1, 22–36.BoutraaT.AkhkhaA.Al-ShoaibiA. A.AlhejeliA. M. (2010). Effect of water stress on growth and water use efficiency (WUE) of some wheat cultivars (Triticum durum) grown in Saudi Arabia. J. Taibah Univ. Sci.3, 39–48. doi: 10.1016/S1658-3655(12)60019-3BradburyP. J.ZhangZ.KroonD. E.CasstevensT. M.RamdossY.BucklerE. S. (2007). TASSEL: software for association mapping of complex traits in diverse samples. Bioinformatics23, 2633–2635. doi: 10.1093/bioinformatics/btm308BrintonJ.UauyC. (2019). A reductionist approach to dissecting grain weight and yield in wheat. J. Integr. Plant Biol.61, 337–358. doi: 10.1111/jipb.12741CarilloP.MastrolonardoG.NaccaF.FuggiA. (2005). Nitrate reductase in durum wheat seedlings as affected by nitrate nutrition and salinity. Funct. Plant Biol.32, 209–219. doi: 10.1071/FP04184Cevher-KeskinB.YıldızhanY.Sekmen-ÇetinelA. H.SezermanO. U.ÖzerB.FayetörbayR.. (2023). Deciphering drought response mechanisms: transcriptomic insights from drought-tolerant and drought-sensitive wheat (Triticum aestivum L.) cultivars. doi: 10.1101/2023.10.26.563210ChenB.LiuY.YangY.WangQ.LiS.LiF.. (2025). A system genetics analysis uncovers the regulatory variants controlling drought response in wheat. Plant Biotechnol. J. doi: 10.1111/pbi.14605ChesselD.DufourA. B.ThioulouseJ. (2004). The ade4 package-I-One-table methods. R news.4, 5–10.ColasuonnoP.MarcotuliI.BlancoA.MaccaferriM.CondorelliG. E.TuberosaR.. (2019). Carotenoid pigment content in durum wheat (Triticum turgidum L. var durum): An overview of quantitative trait loci and candidate genes. Front. Plant Sci.10. doi: 10.3389/fpls.2019.01347CondonA. G.RichardsR. A.RebetzkeG. J.FarquharG. D. (2002). Improving intrinsic water-use efficiency and crop yield. Crop Sci.42, 122–131. doi: 10.2135/cropsci2002.1220CurtisT. Y.BoV.TuckerA.HalfordN. G. (2018). Construction of a network describing asparagine metabolism in plants and its application to the identification of genes affecting asparagine metabolism in wheat under drought and nutritional stress. Food Energy Secur.7, e00126. doi: 10.1002/fes3.126DaramolaM. T.XuM. (2022). Recent changes in global dryland temperature and precipitation. Int. J. Climatol.42, 1267–1282. doi: 10.1002/joc.7301Da RosL.BollinaV.SoolanayakanahallyR.PahariS.ElferjaniR.KulkarniM.. (2023). Multi-omics atlas of combinatorial abiotic stress responses in wheat. Plant J.116, 1118–1135. doi: 10.1111/tpj.16332DegenkolbeT.DoP. T.KopkaJ.ZutherE.HinchaD. K.KöhlK. I. (2013). Identification of drought tolerance markers in a diverse population of rice cultivars by expression and metabolite profiling. PloS One8, e63637. doi: 10.1371/journal.pone.0063637DettoriM.CesaraccioC.DuceP. (2017). Simulation of climate change impacts on production and phenology of durum wheat in Mediterranean environments using CERES-Wheat model. Field Crop Res.206, 43–53. doi: 10.1016/j.fcr.2017.02.013DoddI. C.KudoyarovaG. (2021). “Water Uptake in Drying Soil: Balancing Supply and Demand,” in The Root Systems in Sustainable Agricultural Intensification. Eds. RengelZ.DjalovicI.. John Wiley & Sons Ltd, 245–270. doi: 10.1002/9781119525417.ch9El-AzazJ.CánovasF. M.ÁvilaC.de la TorreF. (2018). The arogenate dehydratase ADT2 is essential for seed development in Arabidopsis. Plant Cell Physiol.59, 2409–2420. doi: 10.1093/pcp/pcy200El-AzazJ.de la TorreF.PascualM. B.DebilleS.CanletF.HarvengtL.. (2020). Transcriptional analysis of arogenate dehydratase genes identifies a link between phenylalanine biosynthesis and lignin biosynthesis. J. Exp. Bot.71, 3080–3093. doi: 10.1093/jxb/eraa099EvannoG.RegnautS.GoudetJ. (2005). Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol. Ecol.14, 2611–2620. doi: 10.1111/j.1365-294X.2005.02553.xFadoulH. E.Martínez RivasF. J.NeumannK.BalazadehS.FernieA. R.AlseekhS. (2021). Comparative molecular and metabolic profiling of two contrasting wheat cultivars under drought stress. Int. J. Mol. Sci.22, 13287. doi: 10.3390/ijms222413287FanW.LouH. Q.GongY. L.LiuM. Y.WangZ. Q.YangJ. L.. (2014). Identification of early Al-responsive genes in rice bean (Vigna umbellata) roots provides new clues to molecular mechanisms of Al toxicity and tolerance. Plant Cell Environ.37, 1586–1597. doi: 10.1111/pce.12258FlexasJ.MedranoH. (2002). Drought-inhibition of photosynthesis in C3 plants: Stomatal and non-stomatal limitations revisited. Ann. Bot.89, 183–189. doi: 10.1093/aob/mcf027FranzoniG.CocettaG.FerranteA. (2021). Effect of glutamic acid foliar applications on lettuce under water stress. Physiol. Mol. Biol. Plants.27, 1059–1072. doi: 10.1007/s12298-021-00984-6FuentealbaC.HernandezI.SaaS.ToledoL.BurdilesP.ChirinosC.. (2017). Color and in vitro quality attributes of walnuts from different growing conditions correlate with key precursors of primary and secondary metabolism. Food Chem.232, 664–672. doi: 10.1016/j.foodchem.2017.04.029GhorbelM.ZribiI.BesbesM.BoualiN.BriniF. (2023). Catalase gene family in durum wheat: Genome-wide analysis and expression profiling in response to multiple abiotic stress conditions. Plants12, 2720. doi: 10.1016/j.foodchem.2017.04.029GoelS.SinghK.GrewalS.NathM. (2020). Impact of “omics” in improving drought tolerance in wheat. Crit. Rev. Plant Sci.39, 222–235. doi: 10.1080/07352689.2020.1778924González-VillagraJ.Reyes-DíazM. M.Tighe-NeiraR.Inostroza-BlancheteauC.EscobarA. L.BravoL. A. (2022). Salicylic acid improves antioxidant defense system and photosynthetic performance in Aristotelia Chilensis plants subjected to moderate drought stress. Plants11, 639. doi: 10.3390/plants11050639GujjarR. S.KarkuteS. G.RaiA.SinghM.SinghB. (2018). Proline-rich proteins may regulate free cellular proline levels during drought stress in tomato. Curr. Sci.114, 915–920. doi: 10.18520/cs/v114/i04/915-920GuoR.ShiL.JiaoY.LiM. X.ZhongX. L.GuF. X.. (2018). Metabolic responses to drought stress in the tissues of drought-tolerant and drought-sensitive wheat genotype seedlings. AoB Plants10, ply016. doi: 10.1093/aobpla/ply016GuptaP. K.BalyanH. S.GahlautV. (2017). QTL analysis for drought tolerance in wheat: present status and future possibilities. Agronomy7, 5. doi: 10.3390/agronomy7010005GuptaA.Rico-MedinaA.Caño-DelgadoA. I. (2020). The physiology of plant responses to drought. Science368, 266–269. doi: 10.1126/science.aaz7614HatoumD.AnnaratoneC.HertogM.L.A.T.M.GeeraerdA.NicolaiB. (2014). Targeted metabolomics study of ‘Braeburn’ apples during long-term storage. Postharv. Biol. Technol.96, 33–41. doi: 10.1016/j.postharvbio.2014.05.004HsuP. K.DubeauxG.TakahashiY.SchroederJ. I. (2021). Signaling mechanisms in abscisic acid-mediated stomatal closure. Plant J.105, 307–321. doi: 10.1111/tpj.15067JaniakA.KwaśniewskiM.SzarejkoI. (2016). Gene expression regulation in roots under drought. J. Exp. Bot.67, 1003–1014. doi: 10.1093/jxb/erv512JogawatA.YadavB.ChhayaLakraN.SinghA. K.NarayanO. P. (2021). Crosstalk between phytohormones and secondary metabolites in the drought stress tolerance of crop plants: a review. Physiol. Plant172, 1106–1132. doi: 10.1111/ppl.13328JombartT. (2008). adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics24, 1403–1405. doi: 10.1093/bioinformatics/btn129JoshiR.WaniS. H.SinghB.BohraA.DarZ. A.LoneA. A.. (2016). Transcription factors and plants response to drought stress: current understanding and future directions. Front. Plant Sci.7. doi: 10.3389/fpls.2016.01029KhalvandiM.SiosemardehA.RoohiE.KeramatiS. (2021). Salicylic acid alleviated the effect of drought stress on photosynthetic characteristics and leaf protein pattern in winter wheat. Heliyon7, e05908. doi: 10.1016/j.heliyon.2021.e05908KhanM. I. R.FatmaM.PerT. S.AnjumN. A.KhanN. A. (2015). Salicylic acid-induced abiotic stress tolerance and underlying mechanisms in plants. Front. Plant Sci.6. doi: 10.3389/fpls.2015.00462KnochD.MeyerR. C.HeuermannM. C.RieweD.PelekeF. F.SzymańskiJ.. (2024). Integrated multi-omics analyses and genome-wide association studies reveal prime candidate genes of metabolic and vegetative growth variation in canola. Plant J.117, 713–728. doi: 10.1111/tpj.16524KumarM.Kumar PatelM.KumarN.BajpaiA. B.SiddiqueK. H. (2021). Metabolomics and molecular approaches reveal drought stress tolerance in plants. Int. J. Mol. Sci.22, 9108. doi: 10.3390/ijms22179108KuromoriT.FujitaM.TakahashiF.Yamaguchi-ShinozakiK.ShinozakiK. (2022). Inter-tissue and inter-organ signaling in drought stress response and phenotyping of drought tolerance. Plant J.109, 342–358. doi: 10.1111/tpj.15619LeaP. J.SodekL.ParryM. A. J.ShewryP. R.HalfordN. G. (2007). Asparagine in plants. Ann. Appl. Biol.150, 1–26. doi: 10.1111/j.1744-7348.2006.00104.xLealA. R.BarrosP. M.ParizotB.SapetaH.VangheluweN.AndersenT. G.. (2022). Translational profile of developing phellem cells in Arabidopsis thaliana roots. Plant J.110, 899–915. doi: 10.1111/tpj.15691Le RouxM. S.KunertK. J.CullisC. A.BothaA. M. (2024). Unlocking wheat drought tolerance: the synergy of omics data and computational intelligence. Food Energy Secur.13, e70024. doi: 10.1002/fes3.70024LiC.BaoY.GuoW.LiC.LiC. (2024). Multi-omics analysis of response mechanism of programmed cell death in wheat endosperm to post-anthesis drought stress. Environ. Exp. Bot.218, 105614. doi: 10.1016/j.envexpbot.2023.105614LiC.NongQ.SolankiM. K.LiangQ.XieJ.LiuX.. (2016). Differential expression profiles and pathways of genes in sugarcane leaf at elongation stage in response to drought stress. Sci. Rep.6, 25698. doi: 10.1038/srep25698LiS. W.ShiR. F.LengY.ZhouY. (2016). Transcriptomic analysis reveals the gene expression profile that specifically responds to IBA during adventitious rooting in mung bean seedlings. BMC Genomics17, 43. doi: 10.1186/s12864-016-2372-4LiuQ.LiuM.WuS.XiaoB.WangX.SunB.. (2020). Metabolomics reveals antioxidant stress responses of wheat (Triticum aestivum L.) exposed to chlorinated organophosphate esters. J. Agric. Food Chem.68, 6520–6529. doi: 10.1021/acs.jafc.0c01397LiuH.SearleI. R.MatherD. E.AbleA. J.AbleJ. A. (2015). Morphological, physiological and yield responses of durum wheat to pre-anthesis water-deficit stress are genotype-dependent. Crop Pasture Sci.66, 1024–1038. doi: 10.1071/CP15013MaccaferriM.HarrisN. S.TwardziokS. O.PasamR. K.GundlachH.SpannaglM.. (2019). Durum wheat genome highlights past domestication signatures and future improvement targets. Nat. Genet.51, 885–895. doi: 10.1038/s41588-019-0381-3MaoX.ZhangH.TianS.ChangX.JingR. (2010). TaSnRK2. 4, an SNF1-type serine/threonine protein kinase of wheat (Triticum aestivum L.), confers enhanced multistress tolerance in Arabidopsis. J. Exp. Bot.61, 683–696. doi: 10.1093/jxb/erp331MarčekT.HamowK. A.VéghB.JandaT.DarkoE. (2019). Metabolic response to drought in six winter wheat genotypes. PloS One14, e0212411. doi: 10.1371/journal.pone.0212411MartinelliF.RemoriniD.SaiaS.MassaiR.TonuttiP. (2013). Metabolic profiling of ripe olive fruit in response to moderate water stress. Sci. Hortic.159, 52–58. doi: 10.1016/j.scienta.2013.04.039Martínez-MorenoF.SolísI.NogueroD.BlancoA.Özberkİ.NsarellahN.. (2020). Durum wheat in the Mediterranean Rim: Historical evolution and genetic resources. Genet. Resour. Crop Evol.67, 1415–1436. doi: 10.1007/s10722-020-00913-8MinK.ShowmanL.PereraA.AroraR. (2018). Salicylic acid-induced freezing tolerance in spinach (Spinacia oleracea L.) leaves explored through metabolite profiling. Environ. Exp. Bot.156, 214–227. doi: 10.1016/j.envexpbot.2018.09.011Moradi TarnabiZ.IranbakhshA.MehreganI.AhmadvandR. (2020). Impact of arbuscular mycorrhizal fungi (AMF) on gene expression of some cell wall and membrane elements of wheat (Triticum aestivum L.) under water deficit using transcriptome analysis. Physiol. Mol. Biol. Plants.26, 143–162. doi: 10.1007/s12298-019-00727-8MuhammadT.ZhangJ.MaY.LiY.ZhangF.ZhangY.. (2019). Overexpression of a mitogen-activated protein kinase SlMAPK3 positively regulates tomato tolerance to cadmium and drought stress. Molecules24, 556. doi: 10.3390/molecules24030556MuthusamyM.KimJ. H.KimJ. A.LeeS. I. (2021). Plant RNA binding proteins as critical modulators in drought, high salinity, heat, and cold stress responses: an updated overview. Int. J. Mol. Sci.22, 6731. doi: 10.3390/ijms22136731ObataT.FernieA. R. (2012). The use of metabolomics to dissect plant responses to abiotic stresses. Cell. Mol. Life Sci.69, 3225–3243. doi: 10.1007/s00018-012-1091-5OberE. S.ClarkC. J.Le BloaM.RoyalA.JaggardK. W.PidgeonJ. D. (2004). Assessing the genetic resources to improve drought tolerance in sugar beet: agronomic traits of diverse genotypes under droughted and irrigated conditions. Field Crops Res.90, 213–234. doi: 10.1016/j.fcr.2004.03.004PálM.MajláthI.NémethE.HamowK.Á.SzalaiG.RudnóyS.. (2018). The effects of putrescine are partly overlapping with osmotic stress processes in wheat. Plant Sci.268, 67–76. doi: 10.1016/j.plantsci.2017.12.011PaulM.TanskanenJ.JääskeläinenM.ChangW.DalalA.MoshelionM.. (2023). Drought and recovery in barley: key gene networks and retrotransposon response. Front. Plant Sci.14. doi: 10.3389/fpls.2023.1193284PritchardJ. K.StephensM.RosenbergN. A.DonnellyP. (2000). Association mapping in structured populations. Am. J. Hum. Genet.67, 170–181. doi: 10.1086/302959QuingaL. A. P.HeringerA. S.de Freitas FragaH. P.do Nascimento VieiraL.SilveiraV.SteinmacherD. A.. (2018). Insights into the conversion potential of Theobroma cacao L. somatic embryos using quantitative proteomic analysis. Sci. Hortic.229, 65–76. doi: 10.1016/j.scienta.2017.10.005RahimA. A.UzairM.RehmanN.RehmanO. U.ZahraN.KhanM. R. (2022). Genome-wide identification and characterization of receptor-like protein kinase 1 (RPK1) gene family in Triticum aestivum under drought stress. Front. Genet.13. doi: 10.3389/fgene.2022.912251RahimiY.IngvarssonP. K.BihamtaM. R.AlipourH.TaleeiA.Khoshnoodi Jabar AbadiS. (2021). Characterization of dynamic regulatory gene and protein networks in wheat roots upon perceiving water deficit through comparative transcriptomics survey. Front. Plant Sci.12. doi: 10.3389/fpls.2021.710867RampinoP.De PascaliM.De CaroliM.LuvisiA.De BellisL.PiroG.. (2017). Td4IN2: A drought-responsive durum wheat (Triticum durum Desf.) gene coding for a resistance like protein with serine/threonine protein kinase, nucleotide binding site and leucine rich domains. Plant Physiol. Bioch.120, 223–231. doi: 10.1016/j.plaphy.2017.10.010R Core Team (2023). R: A language and environment for statistical computing (R Foundation for Statistical Computing). Available at: https://www.R-project.org/ (Accessed July 14, 2024).RohartF.GautierB.SinghA.Lê CaoK. A. (2017). MixOmics: An R package for ‘omics feature selection and multiple data integration. PloS Comput. Biol.13, e1005752. doi: 10.1371/journal.pcbi.1005752RoychowdhuryR.DasS. P.GuptaA.PariharP.ChandrasekharK.SarkerU.. (2023). Multi-omics pipeline and omics-integration approach to decipher plant’s abiotic stress tolerance responses. Genes14, 1281. doi: 10.3390/genes14061281RoyoC.AmmarK.VillegasD.SorianoJ. M. (2021). Agronomic, physiological and genetic changes associated with evolution, migration and modern breeding in durum wheat. Front. Plant Sci.12. doi: 10.3389/fpls.2021.674470SadakM. S.SekaraA.Al-AshkarI.Habib-ur-RahmanM.SkalickyM.BresticM.. (2022). Exogenous aspartic acid alleviates salt stress-induced decline in growth by enhancing antioxidants and compatible solutes while reducing reactive oxygen species in wheat. Front. Plant Sci.13. doi: 10.3389/fpls.2022.987641SainzM.DiazP.MonzaJ.BorsaniO. (2010). Heat stress results in loss of chloroplast Cu/Zn superoxide dismutase and increased damage to photosystem II in combined drought-heat stressed Lotus japonicus. Physiol. Plant140, 46–56. doi: 10.1111/j.1399-3054.2010.01383.xSchlaepferD. R.BradfordJ. B.LauenrothW. K.MunsonS. M.TietjenB.HallS. A.. (2017). Climate change reduces extent of temperate drylands and intensifies drought in deep soils. Nat. Commun.8, 14196. doi: 10.1038/ncomms14196ScholanderP. F.BradstreetE. D.HemmingsenE. A.HammelH. T. (1965). Sap Pressure in Vascular Plants: Negative hydrostatic pressure can be measured in plants. Science148, 339–346. doi: 10.1126/science.148.3668.33ShahT.XuJ.ZouX.ChengY.NasirM.ZhangX. (2018). Omics approaches for engineering wheat production under abiotic stresses. Int. J. Mol. Sci.19, 2390. doi: 10.3390/ijms19082390SharmaM.GuptaS. K.MajumderB.MauryaV. K.DeebaF.AlamA.. (2017). Salicylic acid mediated growth, physiological and proteomic responses in two wheat varieties under drought stress. J. Proteom.163, 28–51. doi: 10.1016/j.jprot.2017.05.011SharmaJ. S.RunningK. L. D.XuS. S.ZhangQ.Peters HaugrudA. R.SharmaS.. (2019). Genetic analysis of threshability and other spike traits in the evolution of cultivated emmer to fully domesticated durum wheat. Mol. Genet. Genomics294, 757–771. doi: 10.1007/s00438-019-01544-0ShewryP. R.HeyS. J. (2015). The contribution of wheat to human diet and health. Food Energy Secur.4, 178–202. doi: 10.1002/fes3.64SinclairT. R. (2005). Theoretical Analysis of Soil and Plant Traits Influencing Daily Plant Water Flux on Drying Soils. Agron. J.97, 1148–1152. doi: 10.2134/agronj2004.0286SinghA.ShannonC. P.GautierB.RohartF.VacherM.TebbuttS. J.. (2019). DIABLO: an integrative approach for identifying key molecular drivers from multi-omics assays. Bioinformatics35, 3055–3062. doi: 10.1093/bioinformatics/bty1054SinghP. K.SrivastavaD.TiwariP.TiwariM.VermaG.ChakrabartyD. (2019). “Drought tolerance in plants: molecular mechanism and regulation of signaling molecules,” in Plant Signaling Molecules: Role and Regulation under Stressful Environments. Eds. KhanM. I. R.ReddyP. S.FerranteA.KhanN. A. (Woodhead Publishing, UK), 105–123. doi: 10.1016/B978-0-12-816451-8.00006-XSinghB.UshaK. (2003). Salicylic acid induced physiological and biochemical changes in wheat seedlings under water stress. Plant Growth Regul.39, 137–141. doi: 10.1023/A:1022556103536SolomonK. F.LabuschagneM. T. (2005). Differences in the level of D-glucose and sucrose among durum wheat (Triticum turgidum var durum L.) genotypes differing in their responses to drought stress. S. Afr. J. Plant Soil.22, 122–126. doi: 10.1080/02571862.2005.10634693SorianoJ. M.ColasuonnoP.MarcotuliI.GadaletaA. (2021). Meta-QTL analysis and identification of candidate genes for quality, abiotic and biotic stress in durum wheat. Sci. Rep.11, 11877. doi: 10.1038/s41598-021-91446-2SuiD.WangB. (2020). Transcriptome analysis reveals complex defensive mechanisms in salt-tolerant and salt-sensitive shrub willow genotypes under salinity stress. Int. J. Genom.2020, 6870157. doi: 10.1155/2020/6870157SukumaranS.ReynoldsM. P.SansaloniC. (2018). Genome-wide association analyses identify QTL hotspots for yield and component traits in durum wheat grown under yield potential, drought, and heat stress environments. Front. Plant Sci.9. doi: 10.3389/fpls.2018.00081SulimanM. S. E.ElradiS. B. M.ZhouG.MengT.ZhuG.XuY.. (2024). Exogenous glutathione protected wheat seedling from high temperature and water deficit damages. Sci. Rep.14, 5304. doi: 10.1038/s41598-023-47868-1SunJ.BieX. M.WangN.ZhangX. S.GaoX. Q. (2020). Genome-wide identification and expression analysis of YTH domain-containing RNA-binding protein family in common wheat. BMC Plant Biol.20, 351. doi: 10.1186/s12870-020-02505-1TangD.ChenM.HuangX.ZhangG.ZengL.ZhangG.. (2023). SRplot: A free online platform for data visualization and graphing. PloS One18, e0294236. doi: 10.1371/journal.pone.0294236UllahN.YüceM.Neslihan Öztürk GökçeZ.BudakH. (2017). Comparative metabolite profiling of drought stress in roots and leaves of seven Triticeae species. BMC Genomics18, 969. doi: 10.1186/s12864-017-4321-2VargasV. H.SchulthessA.RoyoC.MatusI.SchwemberA. R. (2016). Transcripts levels of Phytoene synthase 1 (Psy-1) are associated to semolina yellowness variation in durum wheat (Triticum turgidum L. ssp. durum). J. Cereal Sci.68, 155–163. doi: 10.1016/j.jcs.2016.01.011WangH.ChengX.YinD.ChenD.LuoC.LiuH.. (2023). Advances in the research on plant WRKY transcription factors responsive to external stresses. Curr. Issues Mol. Biol.45, 2861–2880. doi: 10.3390/cimb45040187WangS.WongD.ForrestK.AllenA.ChaoS.HuangB.. (2014). Characterization of polyploidy wheat genomic diversity using a high-density 90000 single nucleotide polymorphism array. Plant Biotechnol. J.12, 787–796. doi: 10.1111/pbi.12183WengY.MegaR.AbeF.TsujimotoH.OkamotoM. (2024). Metabolic profiles in drought-tolerant wheat with enhanced abscisic acid sensitivity. PloS One19, e0307393. doi: 10.1371/journal.pone.0307393WickhamH. (2010). A layered grammar of graphics. J. Comput. Graph. Stat.19, 3–28. doi: 10.1198/jcgs.2009.07098XuT.GuL.ChoiM. J.KimR. J.SuhM. C.KangH. (2014). Comparative functional analysis of wheat (Triticum aestivum) zinc finger-containing glycine-rich RNA-binding proteins in response to abiotic stresses. PloS One9, e96877. doi: 10.1371/journal.pone.0096877YadavA. K.CarrollA. J.EstavilloG. M.RebetzkeG. J.PogsonB. J. (2019). Wheat drought tolerance in the field is predicted by amino acid responses to glasshouse-imposed drought. J. Exp. Bot.70, 4931–4948. doi: 10.1093/jxb/erz224YadavB.JogawatA.RahmanM. S.NarayanO. P. (2021). Secondary metabolites in the drought stress tolerance of crop plants: A review. Gene Rep.23, 101040. doi: 10.1016/j.genrep.2021.101040YanY.GanJ.TaoY.OkitaT. W.TianL. (2022). RNA-binding proteins: the key modulator in stress granule formation and abiotic stress response. Front. Plant Sci.13. doi: 10.3389/fpls.2022.882596YangD. H.KwakK. J.KimM. K.ParkS. J.YangK. Y.KangH. (2014). Expression of Arabidopsis glycine-rich RNA-binding protein AtGRP2 or AtGRP7 improves grain yield of rice (Oryza sativa) under drought stress conditions. Plant Sci.214, 106–112. doi: 10.1016/j.plantsci.2013.10.006YinL.ZhangH.TangZ.XuJ.YinD.ZhangZ.. (2021). rMVP: A memory-efficient, visualization-enhanced, and parallel-accelerated tool for genome-wide association study. Genom. Proteom. Bioinform.19, 619–628. doi: 10.1016/j.gpb.2020.10.007YoshidaT.MogamiJ.Yamaguchi-ShinozakiK. (2014). ABA-dependent and ABA-independent signaling in response to osmotic stress in plants. Curr. Opin. Plant Biol.21, 133–139. doi: 10.1016/j.pbi.2014.07.009YuJ.PressoirG.BriggsW. H.Vroh BiI.YamasakiM.DoebleyJ. F.. (2006). A unified mixed-model method for association mapping that accounts for multiple levels of relatedness. Nat. Genet.38, 203–208. doi: 10.1038/ng1702YuY.SongT.WangY.ZhangM.LiN.YuM.. (2023). The wheat WRKY transcription factor TaWRKY1-2D confers drought resistance in transgenic Arabidopsis and wheat (Triticum aestivum L.). Int. J. Biol. Macromol.226, 1203–1217. doi: 10.1016/j.ijbiomac.2022.11.234YuH.ZhangY.ZhangZ.ZhangJ.WeiY.JiaX.. (2020). Towards identification of molecular mechanism in which the overexpression of wheat cytosolic and plastid glutamine synthetases in tobacco enhanced drought tolerance. Plant Physiol. Biochem.151, 608–620. doi: 10.1016/j.plaphy.2020.04.013ZadoksJ. C.ChangT. T.KonzakC. F. (1974). A decimal code for the growth stages of cereals. Weed Res.14, 415–421. doi: 10.1111/j.1365-3180.1974.tb01084.xZaïmM.Sanchez-GarciaM.BelkadiB.Filali-MaltoufA.Al AbdallatA.KehelZ.. (2024). Genomic regions of durum wheat involved in water productivity. J. Exp. Bot.75, 316–333. doi: 10.1093/jxb/erad357ZengT.FangB.HuangF.DaiL.TangZ.TianJ.. (2021). Mass spectrometry-based metabolomics investigation on two different indica rice grains (Oryza sativa L.) under cadmium stress. Food Chem.343, 128472. doi: 10.1016/j.foodchem.2020.128472ZhangH.MaoX.JingR.ChangX.XieH. (2011). Characterization of a common wheat (Triticum aestivum L.) TaSnRK2. 7 gene involved in abiotic stress responses. J. Exp. Bot.62, 975–988. doi: 10.1093/jxb/erq328ZhangY.ZhengL.YunL.JiL.LiG.JiM.. (2022). Catalase (CAT) gene family in wheat (Triticum aestivum L.): Evolution, expression pattern and function analysis. Int. J. Mol. Sci.23, 542. doi: 10.3390/ijms23010542ZhaoW.LiuL.ShenQ.YangJ.HanX.TianF.. (2020). Effects of water stress on photosynthesis, yield, and water use efficiency in winter wheat. Water12, 2127. doi: 10.3390/w12082127ZhuX.LiuS.MengC.QinL.KongL.XiaG. (2013). WRKY transcription factors in wheat and their induction by biotic and abiotic stress. Plant Mol. Biol. Rep.31, 1053–1067. doi: 10.1007/s11105-013-0565-4ŽivanovićB.Milić KomićS.TostiT.VidovićM.ProkićL.Veljović JovanovićS. (2020). Leaf soluble sugars and free amino acids as important components of abscisic acid—Mediated drought response in tomato. Plants9, 1147. doi: 10.3390/plants9091147ZuluagaD. L.BlancoE.ManginiG.SonnanteG.CurciP. L. (2023). A survey of the transcriptomic resources in durum wheat: Stress responses, data integration and exploitation. Plants12, 1267. doi: 10.3390/plants1206126