Abstract
Plant roots perceive heat stress (HS) and adapt their architecture accordingly, which in turn influence the yield in crops. Investigating their heterogeneity and cell type-specific response to HS is essential for improving crop resilience. Here, we generate single-cell transcriptional landscape of maize (Zea mays) roots in response to HS. We characterize 15 cell clusters corresponding to 9 major cell types and identify cortex as the main root cell type responsive to HS with the most differentially expressed genes and its trajectory being preferentially affected upon HS. We find that cortex size strongly correlated with heat tolerance that is experimentally validated by using inbred lines and genetic mutation analysis of one candidate gene in maize, providing potential HS tolerance indicator and targets for crop improvement. Moreover, interspecies comparison reveals conserved root cell types and core markers in response to HS in plants, which are experimentally validated. These results provide a universal atlas for unraveling the transcriptional programs that specify and maintain the cell identity of maize roots in response to HS at a cell type-specific level.
Similar content being viewed by others
Introduction
Global warming and climate change increased the frequency of associated abiotic stress in major food crops worldwide. Among abiotic stresses, heat is currently considered the most challenging limitation to global agricultural productivity1. Selection for good root traits is an effective strategy to improve yield under abiotic stress2. Heat stress causes changes in root system architecture, leading to reduced primary root length, fewer lateral roots, less root branching, and diminished water and nutrient uptake3. It has been reported that heat stress shrinks the root meristem and blocks phloem unloading by inducing callose buildup at the plasmodesmata between sieve elements and phloem pole pericycle4. In addition, heat stress may inhibit the formation of lateral root primordia (LRP) in maize5. Cortical tissue with fewer but large cell files is beneficial for anatomy that confers tolerance under drought stress6,7. The combined stress of drought and heat reduces maize primary root length more significantly than either stress alone8. Previous studies have revealed that cells tightly regulate gene expression by altering the transcriptional capacity under stress. Heat stress interferes with protein folding, ion channel activity, cell membrane integrity, and enzyme functions. All of these signals could conceivably contribute to root growth arrest during heat stress9. However, knowledge about the transcriptional regulatory network underlying the heat stress response at the cellular level in maize has been largely unexplored in previous studies because of technical bottlenecks10.
Plant roots provide anchorage, uptake, storage, and transport of minerals and water by serving as an interface between plants and the environment. Quantifying gene expression in specific cell types is crucial for dissecting the genetic regulatory networks that modulate root development in response to biotic and abiotic stimuli, as gene expression is both a product of and a contributor to these regulatory networks11. Maize plays a significant role in global food, feed, and fuel production due to its domestication, global diversification, and C4 photosynthetic pathway, with its unique root system contributing to its adaptability in diverse environments12. Since roots show high developmental plasticity and often adapt to the environment, investigating their anatomical and architectural phenotypes is essential for improving crop resilience under heat stress. Roots are radially symmetrical organs composed of three fundamental tissue types: the epidermis on the outside, the ground tissue at the middle, and a core of vascular elements plus a pericycle that lies in a central cylinder known as the stele13. The ground tissue is further divided into two different cell types, the endodermis and cortex, which are arranged as concentric layers around the stele13. Variations in ground tissue patterning, particularly the number of cortex cell layers, are observed across different plant species and might be one of the key factors contributing to root morphological diversity7. This diversity allows plants to adapt to biotic and abiotic stress. For example, the maize cortex plays a role in drought and flood tolerance and hosts mycorrhizal symbiosis that reduces stress and improves nutrient uptake14,15,16. Mutant studies have identified some cell type- or domain-specific regulators that orchestrate root architecture or stress response by influencing distinct developmental zones within the root13,17,18. Haem Oxygenase-1 (ZmHO-1) and Giberellic Acid-Stimulated Like-1 (ZmGSL-1) have been found to be involved in lateral root development, and promoter-associated histone acetylation may regulate the expression of these genes, which are linked to heat-induced inhibition of lateral root primordium formation in maize5. Similarly, heat shock protein 16.9 (ZmHSP16.9), a cytosolic class I small heat shock protein (HSP) in maize, enhances heat and oxidative stress tolerance in transgenic tobacco by improving seed germination, root length, and antioxidant enzyme activities compared to wild-type plants19. However, these insights are limited by genetic redundancy and pleiotropy, thus, a high-resolution expression atlas of specific cell types and domains is needed to gain further insights into the gene networks that control root development. In addition, most stress response mechanisms were elaborated at the average levels of cell mixtures of the whole plant, rather than the cellular level. However, the expression of stress response genes tends to be highly variable from cell to cell20, which can result in different phenotypes at the single-cell level and varying stress responses21,22. The underlying mechanisms controlling functional differentiation at the single-cell level in response to abiotic stress are largely unknown.
Single-cell RNA sequencing (scRNA-seq) can generate high-resolution cell type-specific expression signatures, revealing the developmental trajectories of new cell types and cell lineages. Quantifying gene expression in individual cells or cell types is critical to understanding the complex gene regulatory networks that control root development. In recent years, rapid developments in scRNA-seq techniques have allowed the record of robust transcriptional signatures of cell types in a wide spectrum of eukaryotic organisms23,24,25,26,27,28. The technical advances have also made it possible to analyze the transcriptome of plant roots at single-cell resolution in dicot Arabidopsis, monocots rice, and maize, demonstrating a high rate of heterogeneity of roots and the expression signatures for the major cell type13,25,29,30,31,32,33,34. Recent genome-wide analyses of transcription with scRNA-seq data also provide a new picture of the molecular basis of gene expression regulation for different cell types in plant response to stress35. Therefore, further understanding of maize root responses to temperature at the cell type will be critical for developing strategies to improve maize root stress resistance.
In this study, we applied scRNA-seq analysis of maize root cells to examine cell type-specific responses to heat stress and constructed a single-cell transcriptional atlas comprising more than 35,000 valid cells. In total, 9 major cell types were classified, and a group of cell-type-specific marker genes in these heterogeneous cell populations were identified. We construct continuous developmental trajectories of columella and cortex cell lineages that were derived by differentiation of a subset of initiative cells with a continuous pseudotime series along their own developmental trajectory. Through the co-expression network analysis, we identified 15 regulatory modules that were activated or repressed in nine cell types for HS response. In particular, comparative analyses of single-cell expression data among maize, Arabidopsis, and rice allow us to examine the conserved features of root cell-type transcriptomes and identify key orthologous genes in response to HS. Moreover, we identified cortex cell type-specific genes and confirmed their function to regulate the apical cortex of maize taproot that are tightly correlated with heat tolerance. Our results profile a single-cell level map of maize root tip in response to HS, offering insights into the cellular mechanisms underlying heat tolerance.
Results
Construction of a single-cell transcriptome atlas in maize root under HS
To investigate the single-cell level transcriptional response to HS, we performed a well-established scRNA-seq on root tips from B73 seedlings exposed to HS (Fig. 1A). 4 mm root tips from 4-day-old plants were finely chopped to generate protoplasts for transcriptome analysis using the 10x Genomics platform (Supplementary Fig. 1). HS samples were prepared by exposing seedlings to 42 °C for 2 h, maintaining the same growth time as the control group. To ensure reproducibility, two replicates were included for both control and HS samples. Highly correlated results were observed between two replicates in both groups (CK-1 vs. CK-2, HS-1 vs. HS-2) and across groups (CK vs HS) (Supplementary Fig. 2A–C), indicating the reproducibility of our experimental procedure. Data were filtered at the cell and gene level, resulting in 7741, 9001 and 9202, 9159 root cells from two control samples and two HS samples, respectively. Median genes per cell ranged from 2796 to 3492 in four samples, and median unique molecular identifiers per cell varied from 5894 to 8137 (Supplementary Table 1), representing a high gene coverage for further analysis, which was comparable with previous studies13,25. To assess the robustness of scRNA-seq and the effect of protoplasting, we conducted conventional bulked RNA-seq on root tips without protoplasting. A strong correlation (R = 0.77 ~ 0.81, p < 2.2 × 10−16) was observed between the combined scRNA-seq data with bulk RNA-seq from independently isolated protoplasted root tips, suggesting the validity of scRNA-seq and the efficiency of maize root protoplast extraction (Supplementary Fig. 2D–G). Although the protoplasting process is generally considered to have minimal impact on gene expression, we performed bulk RNA-seq on protoplasted cells and identified genes that were significantly affected by this process. These genes were excluded from our scRNA-seq analysis to minimize the influence of protoplasting-induced artifacts34 (Supplementary Data 1). Subsequently, the two scRNA-seq expression matrices were merged using the Seurat integration method for cell type identification and subsequent analysis (Supplementary Fig. 3).
A Overview of the maize scRNA-seq and cluster annotation workflow: Protoplasts were isolated from 4 mm root tips of maize under control and heat shock (42 °C for 2 h), respectively. ScRNA-seq libraries were generated using the 10x Genomics platform, followed by high-throughput sequencing. B Visualization of nine broad populations and their spatial distribution in maize root tips using UMAP. C Expression patterns of specific marker genes for nine cell types. Dot diameter represents the proportion of cluster cells expressing a given gene, while color indicates the mean expression across cells in these cell clusters. D RNA in situ hybridization of cell-type marker genes for the nine putative clusters and response to HS. Scale bars, 100 µm. The right panel of each gene indicates UMAP visualization of the expression of cell-type marker genes in the maize root tip. At least three repeats, each consisting of individual plants were performed for each gene.
Following linear dimensional reduction, the Uniform Manifold Approximation and Projection (UMAP) algorithm and t-Distributed Stochastic Neighborhood Embedding (t-SNE) tools36,37 were employed to visualize the data, revealing the grouping of cells into 15 clusters (Fig. 1B and Supplementary Fig. 3). The maize cortex marker gene Zm00001d017508 and xylem marker gene Zm00001d035689, previously identified through in situ hybridization experiments13,18,38, were also detected in our single-cell experiments (Supplementary Data 2). Due to the limited availability of validated maize root cell-type-specific genes, we identified potential marker genes with high and specific expression in individual clusters and also used available information from recent relevant references13,17,18,38 to select cell type-specific marker genes for annotating the clusters (Fig. 1C, D and Supplementary Data 2)13. Through in situ hybridization, we found Zm00001d032822 is specific to epidermal cells, Zm00001d012081 (PLT29) to cortex cells, Zm00001d050168 to the endodermis, Zm00001d021192 (umc2686b) to the stele, Zm00001d032672 to the xylem, Zm00001d037032 to the phloem, Zm00001d005472 to the pericycle, Zm00001d004089 (PRP18) to the columella, and Zm00001d040390 to the initials (Fig. 1C, D and Supplementary Fig. 4). The five validated marker genes, including Zm00001d032822, Zm00001d012081, Zm00001d050168, Zm00001d004089, and Zm00001d040390, were identified in the same cell types as in our study (Supplementary Data 3) in previous studies13,18,38. In addition, a total of 939 marker genes were identified, of which 331 genes were similarly mapped in Ortiz-Ramírez et al.’s study13, 55 in Cao et al.’s study18, 4 in Marand et al.’s study17, and 33 in Li et al.’s study38 (Supplementary Data 2 and 3). Our markers, combined with previously known ones13,17,18, enabled us to annotate 15 clusters corresponding to 9 cell types in the maize root tips (Fig. 1B). To gain insights into each cell type’s functions, we conducted Gene Ontology (GO) analysis (Supplementary Fig. 5A). For instance, the columella cluster featured genes associated with serine-type endopeptidase inhibitor activity, which has an important role in plant defense against pests and phytopathogenic microorganisms39,40. Our analyses show that distinct maize root cell types might contribute uniquely to root development and stress response. The initial cluster, representing metabolically active cells, is enriched in genes related to microtubule binding, motor activity, and DNA protection, possibly indicating readiness for rapid growth and stress adaptation41. The epidermis cluster, enriched in genes linked to transcription and protein folding, might be crucial for responding to environmental signals42. Cortex cells show enrichment in proton transporter and oxidative stress response genes, possibly indicating their role in metabolic adjustments under stress (Supplementary Fig. 5A and Supplementary Data 4). In summary, our scRNA-seq dataset revealed 15 distinct clusters, corresponding to nine major cell types, validated by marker genes and in situ hybridization utilizing.
While nine cell types were identified in root tips from control and HS samples (Supplementary Fig. 6), their relative proportions were different but not statistically significant for HS (Supplementary Fig. 5B, 5C, and Supplementary Table 2). In addition, we compared the relative representation of root cell types with estimates from microscopy studies (Supplementary Fig. 5D). Regardless of the annotation method used, the proportions of 4 cell types (stele, columella, xylem, and phloem) measured by microscopy matched our scRNA-seq data, while 3 cell types showed slight overrepresentation (epidermis), or underrepresentation (pericycle and cortex), possibly reflecting a bias in the microscopy analysis or in the protoplast preparation of the whole root tissue (Supplementary Fig. 5D). In addition, the expression levels of some cell type-specific genes were altered upon HS. The specific marker genes were significantly up-regulated after HS, including Zm00001d032822 (specifically expressed in epidermal cells), Zm00001d040390 (initials cells), Zm00001d032672 (xylem cells), and Zm00001d037032 (phloem cells). While the specific marker genes of other clusters were significantly down-regulated after HS, such as PLT29 (cortex cells), Zm00001d050168 (endodermis cells), and Zm00001d005472 (pericycle cells) (Supplementary Fig. 4), indicating potential differential stress responses among cell types.
Cell type-specific responses to HS in distinct root cell types
After HS exposure, there are still nine same cell types identified in maize roots (Supplementary Fig. 5B). We aimed to identify differentially expressed genes in each cell type and understand conserved and cell type-specific variations of the responses among cell types. Among the nine cell types, cortex clusters exhibited the highest number of differentially expressed genes (290 DEGs) upon HS (Fig. 2A, Supplementary Fig. 7, and Supplementary Data 5), highlighting the essential role of this cell type in HS response, followed by the pericycle (229), endodermis (206) and xylem (162). We also performed pseudobulk differential expression analysis using DESeq243 and edgeR44, applying consistent filtering criteria, and found that the cortex exhibited the highest number of DEGs (Supplementary Fig. 7A, and Supplementary Data 6).
A The upset plot depicts the count of up-regulated DEGs in response to HS for each cluster (bars on the left panel) and the intersections of DEGs between clusters (bars on the top panel). B Heatmap shows the correlation of DEGs among nine cell types from root tips response to HS. Here, DEGs of each cell type after HS were used to perform correlation analysis. C Volcano plot shows the DEGs of the cortex after HS. D The GO terms of up-regulated DEGs in the cortex under HS. Hypergeometric tests in ClusterProfiler software were used to test statistical significance. E, F Visualization of the ZNRD15 expression in cortex cells under HS by UMAP (E) and violin plot (F) (ScRNA-seq, ****represents p < 0.0001). G, H Visualization of the HSF17 expression in cortex cells under HS by UMAP (G) and violin plot (H) (ScRNA-seq, ****represents p < 0.0001). I, J Visualization of the MAX1b expression in cortex cells under HS by UMAP (I) and violin plot (J) (ScRNA-seq, ****represents p < 0.0001). K–M Fluorescence in situ hybridization shows the changes of ZNRD15 (K), HSF17 (L), and MAX1b (M) genes after heat shock in maize root tips. Bars = 20 μm. The white arrow indicates a positive signal. Statistical significance was tested using a two-tailed t test in (F, H, and J). At least three samples were analyzed in (K, L, and M) with similar results for each experiment.
Gene expression changes upon HS exhibited both similarity and heterogeneity across cell types. The Person correlation coefficients (PCCs) for DEGs among the epidermis, cortex, stele, and columella ranged from 0.71 to 0.81, indicating a high degree of similarity in their gene expression responses to HS (Fig. 2B). Initials and epidermis also exhibited a highly relevant response, with a PCC of 0.71 (Fig. 2B). However, DEGs in phloem, endodermis, and pericycle showed low correlation with other clusters (0.1–0.18), indicating distinct heat response patterns (Fig. 2B). After HS, all nine cell types exhibited high expression of numerous HSP-related genes, which were enriched in GO terms like peptide metabolic process and response to stress (Fig. 2C, D). This suggests the involvement of cortex cells in responding to environmental signals. Some HSP-related genes exhibited elevated expression across all cell types during HS (Supplementary Figs. 7 and 8), while others were highly expressed in cortical cells (Fig. 2E–M). Genes like zinc-ribbon domain protein 15 (ZNRD15, Zm00001d010183), heat shock transcription factors 17 (HSF17, Zm00001d033987), and cytochrome P450 24 (ZmMAX1b, Zm00001d039697) in cortex cells significantly increased after HS, confirmed by in situ experiments (Fig. 2E–M). HSP genes in all clusters exhibited elevated expression during HS, such as HSP19 (Zm00001d048073), HSP24 (Zm00001d007271), and HSP20 (Zm00001d025508), being prominently induced (Supplementary Fig. 8), validated by RT-qPCR experiments (Supplementary Fig. 9). In summary, scRNA-seq of maize root tips unveiled distinct cell-type-specific responses to HS, with certain cell types exhibiting shared stress-response pathways, while others displayed unique response features.
Two subtypes of columella cells inmaize roots
We identified 15 cellular clusters and 9 major cell types, consistent with published maize root scRNA-seq studies13,17,18. There exists heterogeneity among cells even in the same cell type. Interestingly, we found two subtypes of columella cells (Columella 1 and Columella 2) with distinct gene expression profiles (Fig. 3A and Supplementary Fig. 10A), which were not reported in previous studies17,18,38. Columella cells, crucial for graviperception and abiotic stress response, prompted an investigation into their distinct developmental trajectories45,46. Pseudotime analysis using Monocle 226 revealed that columella1 cells primarily belong to branch2 and pre-branch, whereas branch1 cells predominantly consist of columella 2 cells (Fig. 3A, B). The distinct distribution across the reconstructed developmental trajectory of the two columella subtypes was also confirmed by Monocle3 (Supplementary Fig. 10). To further characterize the two columella subtypes, branched point expression analysis identified 1997 DEGs across the pseudotime branch, representing transcriptional differences between Columella 1 and Columella 2 and categorizing these genes into three clusters indicative of transcriptional rewiring during columella development (Fig. 3C and Supplementary Data 7). Cluster1 genes, predominantly expressed in branch2, were associated with responses to various stimuli (cell redox homeostasis, cysteine-type endopeptidase, and protein kinase binding) (Fig. 3C). Cluster2 genes, preferentially expressed in branch1, were related to protein metabolism (protein polymerization, electron transport chain, and unfolded protein binding) (Fig. 3C). Cluster3 genes, expressed in columella pre-branch cells, were associated to glucose metabolism (cellular amino acid biosynthetic process, glycolytic process, and pyruvate kinase activity) (Fig. 3C). Notably, VQDC (VQ domain-containing protein, Zm00001d018844) exhibited significant enriched expression at the collumella2 cluster, while VPP1 (vacuolar proton pump homolog1, Zm00001d037492) showed prominent expression in the columella pre-branching stage (Fig. 3D). Despite columella cells’ known role in responding to low-temperature stress45, pseudotime analysis revealed the developmental trajectory of the two subtypes remained unchanged after HS (Fig. 3E). Collectively, pseudotime analysis provided insights into the dynamic process of columella development and essential gene expression rewiring during cell state transitions.
A Columella cells were reclustered to two putative subclusters. B Pseudotime analysis of two subtypes of columella cells, cells from the two subtypes were unequally distributed across the reconstructed trajectory. C Heatmap showing the differentially expressed genes in two branches during the differentiation trajectory between two subclusters of columella cells. Gene functions revealed by GO enrichment analysis for genes in the three clusters of expression patterns are shown at right. D Illustrative instances include an early (pre-branch), an earlier (columella1), and a late (columella2) expressed gene specific to columella cells. The gene expression in each cell is overlaid onto the UMAP cluster and trajectory, with deeper red colors signifying higher gene expression. E Distribution difference of cells along the pseudotime trajectory of columella cells under HS.
Cortex cells are essential in response to heat shock
The analysis of HS responses across root cell types revealed both conserved and cell type-specific gene changes. Intriguingly, cortical cells exhibited the highest number of DEGs after HS, emphasizing their crucial role in the stress response (Fig. 2A). Maize contains 8–15 layers of parenchymatous cortex tissue47, which can be divided into multiple cell subtypes with distinct functions13. Therefore, we further delved into the cortex cells’ HS response, uncovering two subtypes with distinct gene expression features: cortex fate1 (cortex1) and cortex fate2 (cortex2) (Fig. 4A and Supplementary Fig. 10C). To further characterize the biological significance of the two cortex subtypes, we performed pseudotime analysis with both Monocle2 and Monocle 3 (Fig. 4B and Supplementary Fig. 10D, E), and predicted their developmental potential by CytoTRACE (Fig. 4C). We found that the developmental trajectory started from cortex2 cells and transitioned to cortex1 cells (Fig. 4B, C). GO enrichment analysis shows that the 2000 DEGs of cortex1 and cortex2 fell into two functional clusters, highlighting the inherent transcriptional differences between these two subtypes (Fig. 4D and Supplementary Data 8). Cluster 1, with a cortex2 preference, enriched in GO terms related to DNA-templated transcription, cellular response to DNA damage, and protein folding (Fig. 4D). These GO terms are closely related to plant responses to environmental stimuli and physical injury48, suggesting that cortex2 cells play a crucial role in response to external stress. In contrast, cluster 2, cortex1-preferred, enriched in vesicle-mediated transport, proton transmembrane transport, GTPase activity, and clathrin adapter complex (Fig. 4D). These GO terms are closely related to the physiological and developmental functions of plants, suggesting that cortex1 cells are mainly involved in the developmental functions of maize. Notably, CCDP1 (Cortical cell-delineating protein, Zm00001d047114) that is a cortex-specific gene exhibited significant expression at the cortex1 cluster, while CSU632b (Zm00001d022421) was prominently expressed at the cortex2 cluster (Fig. 4E).
A UMAP visualization of putative subclusters from cortex cells. B Pseudotime trajectory of two subclusters within cortex cells by Monocle 2. The color of the dots indicates the subtypes. C Combined application of CytoTRACE and Monocle 2 to cortex cell differentiation. Automatic selection of the correct root (developmental trajectory from the left to the right) by CytoTRACE. D Heatmap showing the expression of the most important genes in two branches during the differentiation trajectory between two subclusters of cortex cells. Gene functions revealed by GO enrichment analysis for genes preferentially expressed in the two clusters are shown at right. E Examples of a cortex1 and cortex2 expressed cortex-cell–specific gene. F Developmental potential of cells along the pseudotime trajectory of cortex cells under HS by CytoTRACE. G Distribution difference of cortex cells along the pseudotime trajectory of two HS-related genes.
Interestingly, the two cortex subtypes respond differentially to HS, with individual cortex cells projected into both ends of the pseudotime backbone, signifying two distinct final states after HS (Fig. 4F and Supplementary Fig. 10). The gene expression pattern of cortex cells after HS more closely resembled that of cortex2 cells (Fig. 4C, F), which was consistent with the observation that the highly expressed gene in cortex2 are enriched in GO terms related to stress response (cellular response to DNA damage stimulus and protein folding) (Fig. 4D). Post-HS, cortex cells showed significant alterations in gene expression patterns across the reconstructed trajectory, exemplified by Zm00001d001951 (encoding a PHD finger-like domain-containing protein 5 A, also named PHF5A) in cortex2 (HS branch) and Zm00001d004340 (encoding amino acid/auxin permease 12, also named Aaap12) in cortex1 (control branch) (Fig. 4G and Supplementary Fig. 11A, B). Notably, two cortex subcluster marker genes also showed completely opposite distribution patterns after HS, especially the CCDP1 gene (Supplementary Fig. 11C, D). The CCDP1 gene encodes an extensin-like protein with its possible role in cell expansion49. Moreover, the expressions of CCDP1 and Aaap12, which are specifically expressed in cortex1, are decreased after HS, while the expressions of CSU632b and PHF5A, which are specifically expressed in cortex2, are increased after HS (Supplementary Fig. 11). This is consistent with cortex2 being a cell subcluster that responds to external stimuli and cortex1 being a developmental cell subcluster (Fig. 4D). Collectively, these data suggest there are two subtypes of cortex cells, with differential expression features, and play distinct roles under HS. This insight into cortex-specific gene dynamics affirms their vital role in HS response, especially cortex2 cells, and sheds light on the gene expression reconfiguration essential for post-HS cell state transitions.
Co-expression modules of DEGs in heat response at single-cell level
To reveal dynamic patterns and the single-cell transcriptional regulatory networks for HS response in maize root, a gene co-expression network analysis of genes across all nine cell types was performed using Hotspot software50. The resultant co-expression network was composed of 15 modules (M1 to M15), containing 167 to 1745 genes each (Fig. 5A, Supplementary Figs. 12, 13, and Supplementary Data 9). Module trait correlation revealed twelve modules (M1–7, M9, M11, and M13–15) exhibited significant correlations (r values: 0.2–0.9) with maize root tip cell types under both HS or normal (CK) conditions (Fig. 5A). Eleven modules (M1–7, 9, and 13–15) were associated with CK, while ten modules (M1-5, 7, 9, 11, 13, and 15) were linked to HS (Fig. 5A and Supplementary Data 9). Notably, modules 6 and 14 were specifically tied to CK and module 11 to HS. Module 11’s gene network included HS-related genes in the stele, such as Zm00001d014463 (B4FTK0), tassel seed2 (TS2), dehydrin15 (DHN15), and FCS-like zinc finger35 (FLZ35) (Fig. 5B). These genes, associated with osmotic stress, jasmonic acid biosynthesis, cold stress tolerance, and SnRK1 complex regulation, likely form a network responding to external high-temperature stimulation in the stele51,52. Modules 1 and 2 in HS-responsive modules contained 47 HS-related genes (Supplementary Data 9), with module 2 linked to the cortex under HS (Fig. 5A). This module featured ribosomal proteins (RPS5, RPS3, and RPS4) potentially interacting with ZmMAX1b, known for maize resistance to witchweed (Fig. 5B)53. ZmMAX1b’s involvement in strigolactone biosynthesis alterations suggests its role in the cortex’s thermal response. This gene network analysis illuminates the power of identifying HS resistance genes through co-expression regulatory networks using single-cell transcriptome data.
A 15 co-expression modules associated with each cell types of maize root tips under HS. Red and blue represent high and low correlation coefficient values (Pearson), respectively. The p-value is typically calculated using Moran’s I statistic combined with permutation testing to assess spatial autocorrelation significance. *represents p < 0.05; **represents p < 0.01; ***represents p < 0.001. B 5 gene association modules harboring key DEGs in response to HS. Each node represents a DEG under HS, and DEGs exhibiting a high relationship are connected with edges.
Conservation of HS responses at Single-Cell Level among Different Plant Roots
For a comprehensive understanding of plant root biology, we compared scRNA-seq profiles of maize, rice, and Arabidopsis roots under normal growth conditions, using previously published datasets25,54. Utilizing orthologous genes, we identified common cell types (cortex, stele, endodermis, and epidermis) in maize and rice (Fig. 6A), with higher PCCs observed for cortex cells (Fig. 6A, C). In maize and Arabidopsis, at least 7 common cell types were identified (Fig. 6B). PCC for cortex cells remained higher when comparing maize and Arabidopsis (Fig. 6B, D). Pairwise comparisons revealed a substantial percentage of overlapping orthologous genes for four cell types among Arabidopsis, rice, and maize (Fig. 6E). Homologous genes between maize and rice were enriched in GO terms related to hydrolase activity, hydrolyzing O-glycosyl compounds, vesicle-mediated transport, and transporter activity, suggesting potentially common functions in root nutrient transport for both maize and rice (Supplementary Fig. 14). However, homologous genes between maize and Arabidopsis were enriched in GO terms related to GTPase activity, methyltransferase activity, kinase activity, and vesicle-mediate transport (Fig. 6F). In particular, the gene expression patterns of core marker genes such as ARF-like protein 1b (ARL1b, Zm00001d000379), phytol kinase 2 (also called VTE5, Zm00001d001896) (Supplementary Fig. 15) might have potential overlap and homogeneity in across species.
A, B UMAP visualization of distinct cell types in root tips between maize and rice (A) or Arabidopsis (B). C, D Heatmap showing the correlation among different cell types in root tips between maize and rice (C) or Arabidopsis (D). The r value represents correlation coefficients. E Sankey diagram illustrating that maize shares a high degree of similarity in cortex and endodermis with rice and Arabidopsis. The number of overlapped orthologous genes for each cell type is indicated on the top of each cluster. F GO enrichment analysis for orthologous genes conserved between maize and Arabidopsis. G Volcano plot showing the common DEGs of cortex in maize and Arabidopsis under HS. H, I Venn diagram of upregulated genes (H) and downregulated (I) genes in cortex cluster between maize and Arabidopsis under HS. J GO enrichment results for co-upregulated cortex genes between maize and Arabidopsis. K, L Violin plot illustrating the HS response of two homologous genes, Cypbc7 and HSA32, in the cortex of maize and Arabidopsis (ScRNA-seq, ****represents p < 0.0001). Statistical significance was tested using a two-tailed t test. M GO enrichments for 1288 genes specifically upregulated in maize under heat shock. Hypergeometric tests in ClusterProfiler software were used to test the statistical significance of GO terms in (F, J, and M).
Comparing maize and Arabidopsis root responses to HS at the single-cell level33, 75 upregulated homologous genes, including HSP genes (HSP19, HSP26, HSP27, HSP29), cytochrome b-c1 complex subunit 7 (Cypbc7, Zm00001d005114), Zm00001d014486 (HSA32) and ZNRD15, were identified in both cortical cells (Fig. 6G, H and Supplementary Fig. 16). These 75 genes were associated with proteostasis regulation under HS, i.e., protein folding, response to stress, and RNA helicase activity (Fig. 6J). In addition, 83 homologous genes were downregulated in both cortex types, related to GTPase activity and NAD binding (Fig. 6I and Supplementary Fig. 16J). Interestingly, maize cortex exhibited 1288 specifically upregulated genes in response to HS, involved in methyltransferase activity, iron-sulfur cluster binding, RNA processing, pseudo uridine synthesis, and RNA modification (Fig. 6H, M), contributing potentially to maize’s higher HS tolerance. In summary, the cortex exhibits a certain overlap of heat-related gene expression across plant species in response to HS, with dynamic changes in gene expression potentially contributing to heat-related processes.
Maize enlarge their cortex to promote heat tolerance
The maize cortex plays a pivotal role in abiotic stress7,15,55,56. In a previous study, it was demonstrated that larger root cortical cell size caused by various environmental factors enhances drought tolerance in maize55. Prompted by the cortex-specific gene expression rewiring, we asked whether cortex size variation caused by genetic factors could affect heat tolerance. We examined cortex cell size in certain maize lines, including Mo17, B73, M35, ZJ618, W22, B104, B73, and TF11, and identified a strong association between larger cortex and increased heat tolerance in these associations (Fig. 7A–E and Supplementary Fig. 17). Among 40 randomly selected inbred lines (Supplementary Table 3), root cortex size negatively correlated (R = − 0.70, p = 4.26 × 10−7) with heat-induced damage, emphasizing cortex’s essential role in heat tolerance (Fig.7F).
A–C The primary root phenotypic (A) and apical microstructure (B) of Mo17 and B73 of maize after 4 days germination and quantitative of cortical cell area (C) in apical sections of Mo17 and B73 maize lines. The red box points to amplify the phenotypic of the root apical. Statistical significance was tested using a two-tailed t test (n = three biological replicates; mean ± SD; ***represents p < 0.001). D The phenotypic of Mo17 and B73 of maize seedlings upon HS (Bars = 10 cm). E The relative injury rate of the maize seedlings after HS for Mo17 and B73 (n = three biological replicates; mean ± SD; **represents p < 0.01, ***represents p < 0.001). F The correlation of root cortex size and plant relative injured rate of inbred lines after HS. G Sequence alignment (partial) of the ZmMAX1b gene in B73 wild-type (WT) and zmmax1b-1 mutant plants. The blue arrow marks the mutation site in the zmmax1b-1 mutant. H–J The primary root apical phenotypic (H) and microstructure (I) of B73 and zmmax1b of maize after 4 days germination and quantitative of cortical cell area (J) in apical sections in zmmax1b mutants (n = three biological replicates; mean ± SD; * represents p < 0.05). Statistical significance was tested using a two-tailed t test. The red box points to the amplified phenotypic of the root apical (Bars = 1 cm). K The mRNA relative expression of ZmMAX1b gene in B73 and zmmax1b−1 mutants (n = three biological replicates; mean ± SD; ***represents p < 0.001). Statistical significance was tested using a two-tailed t test. L, M The phenotypic (L) and relative injury rate (M) of B73 and zmmax1b-1 seedlings under HS (Bars = 10 cm) (n = three biological replicates; mean ± SD; **represents p < 0.01; ***represents p < 0.001). N A proposed model for cortex’s functions in response to HS. One-way ANOVA and Bonferroni multiple comparison tests were used to test for statistical significance in (E, M). Source data and exact p-values are provided in the Source Data file.
Recognizing the significance of cortex structure, we hypothesized that altering cortex size could impact maize’s heat tolerance. ZmMAX1b, a gene specifically upregulated in the cortex under HS, exemplifies this relationship (Fig. 2I, J, M). Two EMS mutants of ZmMAX1b were utilized to validate the impact of cortex-specific gene alterations on maize response to HS. The EMS mutant of zmmax1b-1 (EMS4-1f1aff) harbored an early termination site due to the replacement of Guanine deoxynucleotide with Adenine deoxynucleotide at position 786 (Fig. 7G), while the zmmax1b-2 (EMS4-045ac9) with a G to A mutant got the acceptor splice site mutation (Supplementary Fig. 18A). The F1 generation plants were produced by backcrossing the zmmax1b-1 or zmmax1b -2 mutants with the B73. In the F2 plants, derived from self-pollinating the F1 generation, the phenotypic segregation ratio was approximately 3:1 (high temperature insensitive: sensitive = 69:21 for zmmax1b-1 and 45:14 for zmmax1b-2) (Supplementary Fig. 18B). This suggests that the high-temperature-sensitive phenotype in zmmax1b results from a single recessive mutation in B73. In the root tips of 4-day-old seedlings, both zmmax1b-1 and zmmax1b-2 plants exhibited the lower expression of the ZmMAX1b gene, thinner root tips, a smaller cortical area, correspondingly lower heat tolerance, and a higher HS damage ratio compared to the wild-type B73 plants (Fig. 7H–M and Supplementary Fig. 18C–H).
The ZmMAX1b gene contains a P450 domain (Supplementary Fig. 19A) and primarily catalyzes the conversion of zealactol into zealactone53. In maize, mutations in this gene lead to a significant reduction in zealactone and its derivative strigolactones, thereby enhancing resistance to witchweed53. However, its effects on downstream gene expression have not yet been reported. To elucidate the regulatory mechanism of ZmMAX1b in response to HS, we separated the root stele and the outer part of the roots as described previously56, and the outer part of root tips was used for RNA-sequencing (Supplementary Fig. 19B). Heatmap analysis demonstrated clear distinctions between the zmmax1b and wild-type (WT) groups, with high correlation observed among the three biological replicates (Supplementary Fig. 19C). Utilizing the RNA-seq data, we identified 1428 upregulated and 411 downregulated differentially expressed genes (DEGs) (|log2 fold change | ) (Supplementary Fig. 19D and Supplementary Data 10). Gene Ontology (GO) analysis indicated that the DEGs were predominantly associated with processes such as translation, regulation of gene expression, flower development, and developmental processes (Supplementary Fig. 19E). Notably, HSP genes, such as HSFTF28, HSP70-12, HSP7, HSP70-5, HSFTF30, HSP60, and HSP17.4, were identified in 1839 DEGs of zmmax1b mutants (Supplementary Data 10). Remarkably, approximately 21.7% (63/290, p = 1.62 × 10–14) of the DEGs in the cortex cell types under HS overlapped with DEGs in zmmax1b mutants (Supplementary Fig. 19F). These 63 DEGs included a substantial number of genes encoding ribosomal proteins (Supplementary Fig. 19G), a finding validated by quantitative polymerase chain reaction (qPCR) (Supplementary Fig. 19H). This outcome aligns with the gene network of ZmMAX1b associated with ribosomal proteins (Fig. 5B). In addition, the zmmax1b mutant induced the downregulation of three key stress response genes, namely proteolipid membrane potential regulator5 (pmpm5, also called ZmPMP3)57, hsp12 58, abscisic acid stress ripening1 (aasr1) 59. Overall, we demonstrated that the mutant of ZmMAX1b reduced cortex size, decreasing heat tolerance emphasizing the cortex size’s crucial role. The single-cell transcriptome atlas aids in identifying cell fate regulators for resilient crops under abiotic stress.
Discussion
A comprehensive understanding of stress response across distinct cell types and developmental trajectories is essential for uncovering the determinants of cellular identity and the sources of regulatory variation. Single-cell sequencing allows us to distinguish the specific cell from different types, which is essential for deciphering the molecular mechanisms underlying stress responses in maize root systems, as specific expression patterns required for stress response may be induced in the certain cell type. It is worth noting that recent scRNA-seq studies on root or shoot cell atlases have made great progress in cell-type identification and cell trajectories construction13,17,38, similar to a high-resolution picture under one condition. However, it is unresolved to distinguish cells from different types and correlate the statuses of individual cells with certain functions upon stress. Here, we reported a comprehensive single-cell atlas of maize root tips in response to HS and identified nine major cell types by using known and additionally defined cell markers in these heterogeneous cell populations. HS affects gene expression in a cell type-specific manner, in which the cortex was the cluster with the most number of differentially changed genes. We demonstrate that the cortex was mainly affected by HS compared to other cell types. Notably, we found that cortex size was tightly correlated with relative injured rate and heat tolerance, which was further experimentally confirmed using variant inbred lines and mutation maize plants. These data resources will contribute to a more comprehensive understanding of the heat stress at early developmental stages and breeding of heat-resistant maize varieties.
Although scRNA-seq has been utilized in investigating gene pathways involving abiotic stress in plants33,35,38,60, there are few studies on HS in crop plants. One possible explanation for this is that despite protoplast-based plant single-cell sequencing has been well established in the study of plant developmental17,33,54, this preparation is often difficult for complex tissues after stress treatment. Here, we profiled the high-quality single-cell transcriptomes on over 35,000 cells (~ 17,000 in each condition), compared to ~ 3121 cells (~ 1500 in each condition) in the previous study in Arabidopsis33, which enabled us to comprehensively analyze the response of root to HS with plenty number of cells in each cluster. In addition, the number of cell clusters and cell types observed in our study were consistent with published studies of maize root scRNA-seq13,18. Notably, we identified that cortex and columella both have two subcellular types, and these two cell types exhibit distinct cellular trajectories under HS, which has not been described previously (Figs. 3 and 4). Furthermore, among 939 identified cell markers, 331 genes matched those mapped in the previous study by Ortiz-Ramírez et al.13, 55 in the study by Cao et al.18, 4 in Marand et al.’s study17, and 33 in Li et al.’s study38. Furthermore, five validated marker genes, Zm00001d032822, Zm00001d012081, Zm00001d050168, Zm00001d004089, and Zm00001d040390, were found in the same cell types in our study as in previous studies13,18,38. (Supplementary Data 3).
By identifying plant cell types that display dramatic responses to abiotic stresses such as heat, drought, and nutrient starvation, ultimately, we may genetically manipulate relevant cell types to generate stress-tolerant crops without pleiotropically affecting plant fitness and yield. The single-cell atlas of Chinese cabbage, Arabidopsis, and rice showed that abiotic stresses, such as heat stress and salt stress, often affect gene expression in a cell-type-specific manner33,61,62. We would like to note that our scRNA-seq data also reveal predominantly expressed genes in cell type-specific patterns in maize roots under stress treatment. We found that the patterns of cortex, xylem and stele in response to HS were similar, while the patterns of initials, epidermis, columella1 and columella2 were similar (Fig. 2B). Under HS, some genes were generally higher in all clusters of HS samples, while some genes were found to be specifically high in cortical cell populations (Fig. 2A). During a rapid heat response, the synthesis and accumulation of specific proteins are ensured, and these proteins are designated as heat shock proteins (HSPs)63. We found that most of the heat-shocked cells indeed showed gene expression changes, especially some typically canonical heat-shock genes (Supplementary Figs. 8 and 9). Given that research on abiotic resistance has now advanced to the cellular level, we have discovered some cell-specific heat-corresponding marker genes that could be used to identify plant heat responses in different cells to stimulate future heat stress studies.
Interestingly, the highest percentage of heat-shock-response specific genes was found in cortex clusters (Fig. 2A), which is consistent with the finding that the cortex is highly associated with heat tolerance in maize (Fig. 4). We also found the change percentage of specific genes of pericycle, columella2, endodermis, and columella1 were also higher after HS compared to other cell types (Supplementary Data 5). Since these cell types are the outermost cell layer of the root, they may be more directly exposed and respond more quickly to thermal shock. This is consistent with that the heat stress response is fully integrated as a system-wide gene network coordinated across a variety of cells33. Because environmental factors tend to change faster than internal factors, the expression of genes that respond to external stress in the outermost cell layer may diverge faster than that of inner cells64. scRNA-seq has been developed to record the diversity of cell types, lineage relationships, and gene expression patterns in plant roots of Arabidopsis, rice, and maize, demonstrating the high heterogeneity of these tissues and unprecedented expression signatures of major cell types. Previous reports indicate that gene expression in different cell types is not conserved between monocot and dicot65. However, one recent study showed that data from Arabidopsis are sufficient to annotate the cell types of Medicago, indicating that gene expression in the roots is conserved enough to allow cross-species annotation66. To explore the degree of regulatory conservation in root development, interspecies analysis of root scRNA-seq data was conducted from maize and rice25, and the UMAP visualization and Pearson’s correlation coefficient analysis revealed that clusters of cortex cells were highly correlated (Fig. 6A, C). In addition, we compared single-cell datasets of root tip cells between maize and Arabidopsis54 and found that clusters of cortex and columella cells were highly correlated (Fig. 6B, D). Therefore, our results revealed that cell-type clusters of roots are conserved among monocotyledonous and dicotyledonous plants, especially cortex, epidermis, and stele were discovered as the conserved cell clusters in C3 and C4 plants. Notably, through interspecies comparison of conserved cell types, we identified 75 orthologous as core marker genes specifically expressed in cortex cells, such as VTE5, and ARL1b (Supplementary Fig. 15). These results indicated the functional similarity of conserved orthologous genes among three plant species, and the conserved core genes may contribute to the conserved regulation network and similarities in functions of root cell types among species. Most importantly, these core genes could be used as biomarkers that predict maize root development at single-cell resolution.
Although previous studies have revealed that cells tightly regulate gene expression by altering the transcriptional capacity under stress18,33,38,55, the relationship between specific cell types and stress response in plants, especially in crops, remains unclear due to technical bottlenecks. In our work, we discovered the link between the cortex and heat stress tolerance in maize roots at early developmental stages that have not been described previously. First, our single-cell data show that the cortex is the group with the most numbers of up-regulated genes in all clusters upon heat shock (Fig. 2A). The GO analysis showed that these altered genes are mainly enriched in response to stress (Fig. 2D), suggesting the essential role of this cell type in response to HS. Second, the interspecies comparison between our data and Arabidopsis data under HS33 reveals that cells of the same type tended to cluster together (Fig. 6B). Notably, we identified 75 orthologous genes involved in stress response with increasing expression levels in the cortex cells of Arabidopsis and maize after heat treatment (Fig. 6H, J), suggesting the conserved role of apical cortex in response to HS in monocots and dicots. Third, we found that cortex size was negatively correlated with the relative injury rate in maize inbred lines treated by HS. Fourth, the mutation of the cortex-specific gene ZmMAX1b led to a reduction in cortex size and diminished heat tolerance, confirming the positive correlation between cortex area and maize heat tolerance (Fig. 7). These results greatly expand understanding of how heat affects plant development at the cell-type and genome-wide levels (Fig. 7N).
In conclusion, taking advantage of the single-cell extraction method, we successfully established the single-cell level transcriptomic atlas of maize roots in response to heat stress, suggesting spatiotemporal symbiotic perception and early response during development. We finely dissected heat stress-responsive genes, cell trajectories, and core markers at single-cell resolution, leading to an improved understanding of maize roots in response to stress. The obtained information could serve as a valuable resource for comprehensively investigating the molecular mechanism of cell development and differentiation in response to heat stress in crop plants, and also facilitate the breed better crop plants with improved high-temperature tolerance and markedly improve cropping systems.
Methods
Maize growth conditions and protoplast isolation
Maize seeds were wrapped in germination paper sheets, placed in water in the dark until germination, and subsequently cultured for 4 days (28 °C with 10-h-light/14-h dark cycles, with the light intensity of 300 lux). For the heat shock, maize seedlings of germination for 4 days were transferred from 28 °C to 42 °C for 2 h (42 °C with 10-h-light/14-h dark cycles, with the light intensity of 300 lux) as previously described67,68, and the roots were harvested immediately after. One hundred 4 mm tips of primary roots were pre-treated following the methods outlined by the previous study69. Subsequently, the tips were cut into small pieces of approximately 1 mm and incubated in a digestion solution (2% cellulase RS, 0.75% macerozyme R10, 1% hemicellulose, 0.5% pectolase Y-23, 10 mM MES (2-(N-morpholino) ethanesulfonic acid), 0.6 M mannitol, 1 mM CaCl2, 0.1% BSA, and 0.04% mercapto-ethanol) at room temperature in a 28 °C incubator at 60–70 rpm for 2.5 h to release root cell protoplasts. The protoplasts were passed through an 80-μm nylon mesh, followed by centrifugation at 130 × g for 5 minutes. They were then washed three times with 8% mannitol at room temperature. After discarding the supernatant, a minimal amount of 8% mannitol solution was added to resuspend the protoplasts. Trypan blue staining was conducted to verify the quality (< 10% stained) of the protoplasts. The concentration of the protoplasts was adjusted to 1000–1200 cells/ml before being processed with the 10x Genomics Single Cell Protocol (CG00052, RevC).
ScRNA-seq library preparation and sequencing
The isolated protoplasts underwent processing using the 10x Chromium 3’ Single Cell Platform following the manufacturer’s instructions (10x Genomics). Specifically, a batch of 25,000 cells was loaded onto a single-cell chip for library construction using the Chromium Single Cell 3’ v2 Reagent Kit. Protoplasts were partitioned into single-cell GEMs (gel beads in emulsion), followed by cell lysis and barcoded reverse transcription of RNA within the droplets. Quantification of the DNA library was performed using an Agilent 2100 Bioanalyzer, and sequencing was carried out using an Illumina NovaSeq6000 sequencer.
Cell clustering and identification of marker genes
The scRNA-seq data were aligned to the Maize B73 RefGen_v4 genome assembly (http://ftp.ensemblgenomes.org/pub/plants/release-50/fasta/zea_mays/dna/Zea_mays.B73_RefGen_v4.dna_sm.toplevel.fa.gz). The counting process was conducted using the Cell Ranger pipelines (version 6.1.1, 10x Genomics). Subsequently, the raw count matrix data were imported into R utilizing the Seurat (v4.3.0) package for further data analysis70. We applied the following criteria for filtering the low-quality cells and genes: (1) each gene must be expressed in at least 3 cells (Seurat parameter: min.cells = 3); (2) each cell must have a minimum of 200 genes expressed and a maximum of 9000 genes expressed (min.features = 200, nfeature_RNA < 9000); and (3) cells with mitochondrial gene expression levels exceeding 10% were excluded. The DoubletFinder(v.2.0.3)71 tool was utilized to eliminate doublets from the dataset. Normalized transcript abundance for each gene was calculated by dividing counts by the total counts per cell, multiplying by a scaling factor (10,000) and log transforming the result using log1p (NormalizeData function in Seurat). The variable genes detecting and data scaling were respectively processed using the FindVariableFeatures function (“vst” method, 2000 features) ScaleData function. Principal component analysis (PCA) dimensionality reduction was performed on the top 2000 highly variable genes. The first 30 principal components (PCs) were selected based on the PCA elbow plot and used for clustering with a resolution parameter of 0.6. Clusters were visualized and explored using t-SNE36 and UMAP70. Marker genes (cluster-enriched genes) were identified using the FindAllMarkers function in the Seurat package. Differential expression analysis was conducted using the Wilcoxon rank-sum test. Cluster-enriched genes were identified with parameters min. log2fc = 0.25. The final cluster annotation was validated through in situ hybridization experiments.
RaceID, employed as an independent clustering method for analyzing maize scRNA-seq data, presents enhancements in detecting cell types with a small cell number. It also provides utilities for batch effect removal and optional imputation of gene expression. RaceID was executed with the following parameters: mintotal = 200, minexpr = 5, minnumber = 5, metric = “pearson”, clustnr = 3071.
Bulk RNA-seq and analysis
Bulk RNA-seq was conducted on RNA extracted independently from 4-day-old root tips (4 mm) of primary roots and protoplasts. The samples, including primary roots and independently isolated protoplasts, were obtained from 4-day-old root tips (4 mm in length) of primary roots of maize subjected to both control (28 °C) and heat treatment (42 °C for 2 h). The RNA-seq procedure was carried out as previously described72. In brief, total RNA extraction was performed using the RNeasy Plus Mini Kit (Qiagen). The quantity and quality of RNA were assessed using the Qubit NA Assay Kit (Life Technologies), gel electrophoresis, and the Agilent Bioanalyzer 2100 system. The mRNA purification was achieved using oligo(dT)-attached magnetic beads. The library was constructed with the Nextera XT DNA Library Preparation Kit (Illumina), and the amplified library was subjected to sequencing on the NovaSeq6000 platform. Sequence reads were trimmed using fastp (version 0.23.0)73 and aligned to the Maize B73_v4 reference genomes, with HISAT2 (version 2.0.4)74. Gene expression values were calculated on uniquely mapped reads using HTSeq (version 0.11.2)75, and DEseq2 (version 1.4.5) was used to calculate differential expression (fold change > 2 and P < 0.05)43.
Toluidine blue staining
The 4-day-old root tips (4 mm) of B73 primary roots were excised, fixed in a formaldehyde acetic acid solution consisting of 50% ethanol, 5% acetic acid, and 3.7% formaldehyde in deionized water, and subsequently dehydrated through an ethanol series and Van-Clear (Van-Clear used as an alternative to xylene is environmental protection transparent agent and high purity organic solvent of reagent grade with solubility in ethanol and paraffin, but not in water.) (Huntz, Wuhan, China). The dehydrated tissue samples were then embedded in paraffin, cut into 10 µm-thick paraffin sections, and mounted on microscope slides (Fisher Scientific). Sections from both control and heat-treated primary roots were stained with 0.01% toluidine blue. Following deparaffinization and hydration to distilled water, sections were stained in toluidine blue solution for 90 seconds, washed in distilled water, dehydrated quickly through an ethanol series, and covered with a coverslip using a resinous mounting medium. After staining, the samples were sealed with coverslips and observed under a light microscope (Leica DM6 B).
RNA in situ hybridization
The 4 mm root tips from 4-day-old plants were carefully sectioned for RNA in situ hybridization experiments. A cDNA fragment (100–200 bp) was amplified with gene-specific primers and cloned into the pGEM-T Easy vector under 24 °C for 1 h. The plasmid was used as a template for PCR with T7 and SP6 primers. The PCR products were purified and used as templates for in vitro transcription using the SP6/T7 transcription with a mix of transcription buffer, RNase inhibitor, DIG RNA labeled mixtures and T7/SP6 RNA Polymeras. According to the orientation of the cDNA insertion, either T7 or SP6 RNA polymerase was used to synthesize antisense and sense RNA probes. Subsequently, the hybridization and immunological detection procedures followed established protocols76,77, and microscopy was performed in bright-field mode using Leica DM6 B. The primer details are provided in Supplementary Data 11.
RNA fluorescence in situ hybridization (FISH)
For FISH, sample preparation was similar to chromogenic ISH. After blocking, the sections were incubated with anti-digoxigenin-POD antibody (1:200, Roche, 11633716001) for digoxigenin-labeled probes. The hybridization was kept at room temperature for 3 h, and the signals were then detected using TSA Plus Cy5 Fluorescence System (Perkin Elmer). Images were taken with Zeiss LSM 700 confocal microscopy. The Cy5 emission uses a long-pass LP 640 filter, and the wavelength range for detection was 640–700 nm78.
Pseudotime analysis
Pseudotime analysis for cell differentiation and cell fate determination was conducted using the Monocle2 (v.2.22.0) and CytoTRACE R packages26,79. To investigate the developmental trajectory of specific cell types, we analyzed a subset of raw data containing target clusters. The process involved calculating the variance in gene expression across cells using the dispersion Table function. Variable genes were selected based on their average expression levels to define the developmental trajectory. Next, dimensionality reduction to two dimensions was performed using the ‘DDRTree’ method (set max_components = 2). The state transition of single cells in the lower-dimensional space was described using the orderCells function. The trajectory was visualized using the plot_cell_trajectory function in Monocle 2. To specify the beginning of the trajectory, the orderCells function was run again with the root_state argument. Genes dynamically expressed along the pseudotime were clustered and visualized using the plot_pseudotime_heatmap function. To identify genes contributing to the branching of the developmental trajectory, a branch point was selected. The differentialGeneTest function was employed to analyze pseudotime-dependent or branch-dependent genes. Genes that were significantly branch-dependent were visualized using the plot_genes_branched_heatmap function.
Identification of cell type-specific HS genes
The normalized scRNA-seq data from both control and heat-treatment conditions were consolidated into a single object. Variable genes were selected using the FindVariableFeatures function with default parameters. The FindIntegrationAnchors function, with default parameters (dims = 1:30), was employed to identify integration anchors across all samples. The IntegrateData function was then executed on the anchor set with default arguments. Subsequently, ScaleData and RunPCA were applied to the integrated assay to compute 30 principal components (PCs). Clustering analysis was conducted on the integrated assay with a resolution of 0.6, and UMAP was utilized for cluster visualization. The FindMarkers function in Seurat was employed to analyze differentially expressed genes with |avg_log2FC | > = 0.75 and p_val_adj < = 0.01. The intersection size of DEGs across cell types was visualized using the upset function in the UpSetR R package (v. 1.4.0)80.
To explore the transcriptional regulatory network for HS response at the single-cell level, we conducted gene co-expression analysis using Hotspot50. In brief, the data were initially normalized based on the total unique molecular identifiers (UMIs) number for each gene bin. Subsequently, a K-Nearest Neighbors (KNN) graph of genes was constructed using the ‘create_knn_graph’ function with the following parameters: n_neighbors = 8. Genes exhibiting significant spatial autocorrelation (FDR < 0.05) were retained for further analysis. Module identification was carried out using the ‘create_modules’ function with specific parameters: a minimum gene threshold = 15 and fdr_threshold = 0.05. For the cell type enrichment analysis of each module, we calculated the cell type with high module scores (the module score calculation method is referenced in Hotspot50).
Comparative analysis of root scRNA-seq data from maize and rice or Arabidopsis
ScRNA-seq data for rice roots were acquired from Liu et al. 25, while Arabidopsis root scRNA-seq data originated from studies conducted by Shahan R et al. 54 and Jean-Baptiste K et al. 33 Homology relationships between maize and rice or Arabidopsis were evaluated using orthofinder81. The one-to-one orthologous gene pairs were employed for interspecies analysis. Variable genes were chosen through the SCtransform function with default parameters. The FindIntegrationAnchors function command, with normalization.method = “SCT,” was executed to discover integration anchors across all samples. Subsequently, the IntegrateData function was applied to the anchor set with normalization.method = “SCT.” RunPCA was then conducted on the integrated assay to calculate 20 PCs. UMAP dimensionality reductions were performed, and an SNN graph was constructed using dimensions 1:20 as input features and default PCA reduction. Average expression values for each cluster were computed using the AverageExpression function in Seurat. Spearman’s correlation coefficients were calculated to compare maize with rice or Arabidopsis. The FindMarkers function in Seurat was utilized to identify differentially expressed genes across homologous cell types between maize and Arabidopsis.
Functional analysis
Cluster-enriched genes with statistical significance were applied to clusterProfiler (v.4.2.2)82. The R package was used for GO enrichment analysis, and GO terms were visualized using ggplot2 (v.3.4.2).
Identification of zmmax1b mutants
The EMS mutant zmmax1b-1 (Mut_Sample EMS4-1f1aff) contains a mutation that alters the tryptophan codon TGG at amino acid position 262 to an early stop codon (TGA). In addition, the EMS mutant zmmax1b-2 (Mut_Sample EMS4-045ac9) has a mutation where the base G at position 2659 is changed to A, resulting in a mutation of the splice site acceptor. These mutants were obtained from the Maize EMS-induced Mutant Database (http://maizeems.qlnu.edu.cn/)83. Subsequent to sequencing for mutation confirmation (Supplemental Data 11), the EMS mutant zmmax1b homozygous mutants were purified through three generations of backcrossing and subsequently propagated in Hebei and Hainan province, China. Identification of the homozygous mutant zmmax1b and wild-type (WT) plants was performed. These experiments were followed by detailed phenotyping, including anatomical analysis of the mutant root tips, measurement of root cortex size, assessment of the heat shock phenotype, and evaluation of injury rate after HS.
Reporting summary
Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.
Data availability
The raw data of Illumina sequencing, including single-cell RNA-seq and bulk RNA-seq, have been deposited in the National Genomics Data Center (NGDC) database with the accession number CRA015077. The other data generated in this study are provided in the Supplementary Information/Source Data file. Source data are provided with this paper.
References
Wahid, A., Gelani, S., Ashraf, M. & Foolad, M. R. Heat tolerance in plants: An overview. Environ. Exp. Bot. 61, 199–223 (2007).
Tuberosa, R. et al. Identification of QTLs for root characteristics in maize grown in hydroponics and analysis of their overlap with QTLs for grain yield in the field at two water regimes. Plant Mol. Biol. 48, 697–712 (2002).
de Lima, C. F. F., Kleine-Vehn, J., De Smet, I. & Feraru, E. Getting to the root of belowground high temperature responses in plants. J. Exp. Bot. 72, 7404–7413 (2021).
Liu, J., Liu, Y., Wang, S., Cui, Y. Q. & Yan, D. W. Heat stress reduces root meristem size via induction of plasmodesmal callose accumulation inhibiting phloem unloading in Arabidopsis. Int. J. Mol. Sci. 23, 2063 (2022).
Zhang, H. et al. The role of promoter-sssociated histone acetylation of haem oxygenase-1 (HO-1) and giberellic acid-stimulated like-1 (GSL-1) henes in heat-induced lateral root primordium inhibition in maize. Front. Plant Sci. 9, 1520 (2018).
Colombi, T., Herrmann, A. M., Vallenback, P. & Keller, T. Cortical cell diameter is key to energy costs of root growth in wheat. Plant Physiol. 180, 2049–2060 (2019).
Chimungu, J. G., Brown, K. M. & Lynch, J. P. Reduced root cortical cell file number improves drought tolerance in maize. Plant Physiol. 166, 1943–1955 (2014).
Yang, X. Y. et al. Primary root response to combined drought and heat stress is regulated via salicylic acid metabolism in maize. BMC Plant Biol. 22, 417 (2022).
Karlova, R., Boer, D., Hayes, S. & Testerink, C. Root plasticity under abiotic stress. Plant Physiol. 187, 1057–1070 (2021).
He, J. et al. Genome-wide transcript and small RNA profiling reveals transcriptomic responses to heat stress. Plant Physiol. 181, 609–629 (2019).
Petricka, J. J., Winter, C. M. & Benfey, P. N. Control of Arabidopsis root development. Annu. Rev. Plant Biol. 63, 563–590 (2012).
Rogers, E. D. & Benfey, P. N. Regulation of plant root system architecture: implications for crop advancement. Curr. Opin. Biotech. 32, 93–98 (2015).
Ortiz-Ramírez, C. et al. Ground tissue circuitry regulates organ complexity in maize and Setaria. Science 374, 1247–1252 (2021).
Gunawardena, A. H. L. A. N., Pearce, D. M., Jackson, M. B., Hawes, C. R. & Evans, D. E. Characterisation of programmed cell death during aerenchyma formation induced by ethylene or hypoxia in roots of maize (Zea mays L.). Planta 212, 205–214 (2001).
Zhu, J. M., Brown, K. M. & Lynch, J. P. Root cortical aerenchyma improves the drought tolerance of maize (Zea mays L.). Plant Cell Environ. 33, 740–749 (2010).
Liu, F. et al. Identification and functional characterization of a maize phosphate transporter induced by mycorrhiza formation. Plant Cell Physiol. 59, 1683–1694 (2018).
Marand, A. P., Chen, Z. L., Gallavotti, A. & Schmitz, R. J. A cis-regulatory atlas in maize at single-cell resolution. Cell 184, 3041–3055 (2021).
Cao, Y. Y. et al. Single-cell RNA sequencing profiles reveal cell type-specific transcriptional regulation networks conditioning fungal invasion in maize roots. Plant Biotechnol. J. 21, 1839–1859 (2023).
Sun, L. P. et al. ZmHSP16.9, a cytosolic class I small heat shock protein in maize (Zea mays), confers heat tolerance in transgenic tobacco. Plant Cell Rep. 31, 1473–1484 (2012).
Locke, J. C. W., Young, J. W., Fontes, M., Jiménez, M. J. H. & Elowitz, M. B. Stochastic pulse regulation in bacterial stress response. Science 334, 366–369 (2011).
Sánchez-Romero, M. A. & Casadesús, J. Contribution of phenotypic heterogeneity to adaptive antibiotic resistance. Proc. Natl. Acad. Sci. USA 111, 355–360 (2014).
El Meouche, I., Siu, Y. & Dunlop, M. J. Stochastic expression of a multiple antibiotic resistance activator confers transient resistance in single cells. Sci. Rep. 6, 19538 (2016).
Tang, F. C., Lao, K. Q. & Surani, M. A. Development and applications of single-cell transcriptome analysis. Nat. Methods 8, S6–S11 (2011).
Macosko, E. Z. et al. Highly parallel genome-wide expression profiling of individual cells using nanoliter droplets. Cell 161, 1202–1214 (2015).
Liu, Q. et al. Transcriptional landscape of rice roots at the single-cell resolution. Mol. Plant 14, 384–394 (2021).
Trapnell, C. et al. The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells. Nat. Biotechnol. 32, 381–386 (2014).
Zeisel, A. et al. Cell types in the mouse cortex and hippocampus revealed by single-cell RNA-seq. Science 347, 1138–1142 (2015).
Cao, J. Y. et al. Comprehensive single-cell transcriptional profiling of a multicellular organism. Science 357, 661–667 (2017).
Efroni, I. & Birnbaum, K. D. The potential of single-cell profiling in plants. Genome Biol. 17, 65 (2016).
Shulse, C. N. et al. High-throughput single-cell transcriptome profiling of plant cell types. Cell Rep. 27, 2241–2247 (2019).
Zhang, T. Q., Xu, Z. G., Shang, G. D. & Wang, J. W. A single-cell RNA sequencing profiles the developmental landscape of root. Mol. Plant 12, 648–660 (2019).
Rich-Griffin, C. et al. Single-cell transcriptomics: A high-resolution avenue for plant functional genomics. Trends Plant Sci. 25, 186–197 (2020).
Jean-Baptiste, K. et al. Dynamics of gene expression in single root cells of Arabidopsis thaliana. Plant cell 31, 993–1011 (2019).
Denyer, T. et al. Spatiotemporal developmental trajectories in the root revealed using high-throughput single-cell RNA sequencing. Dev. Cell 48, 840–852 (2019).
Li, P. et al. Transcriptional landscape of cotton roots in response to salt stress at single-cell resolution. Plant Commun. 5, 100740 (2023).
van der Maaten, L. & Hinton, G. Visualizing data using t-SNE. J. Mach. Learn Res. 9, 2579–2605 (2008).
Becht, E. et al. Dimensionality reduction for visualizing single-cell data using UMAP. Nat. Biotechnol. 37, 38–44 (2019).
Li, X. H. et al. Single-cell RNA sequencing reveals the landscape of maize root tips and assists in identification of cell type-specific nitrate-response genes. Crop J. 10, 1589–1600 (2022).
Quilis, J. et al. A potato carboxypeptidase inhibitor gene provides pathogen resistance in transgenic rice. Plant Biotechnol. J. 5, 537–553 (2007).
Lingling, L., Jianjun, L., Ming, S., Liyun, L. & Bihao, C. Study on transformation of cowpea trypsin inhibitor gene into cauliflower (Brassica oleracea L. var. botrytis). Afr. J. Biotechnol. 4, 45–49 (2005).
Pesquet, E., Korolev, A. V., Calder, G. & Lloyd, C. W. The microtubule-associated protein AtMAP70-5 regulates secondary wall patterning in Arabidopsis wood cells. Curr. Biol. 20, 744–749 (2010).
Sun, J. L., Li, J. Y., Wang, M. J., Song, Z. T. & Liu, J. X. Protein Quality Control in Plant Organelles: Current Progress and Future Perspectives. Mol. Plant 14, 95–114 (2021).
Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550 (2014).
Robinson, M. D., McCarthy, D. J. & Smyth, G. K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26, 139–140 (2010).
Hong, J. H. et al. A sacrifice-for-survival mechanism protects root stem. Cell 170, 102–113 (2017).
Moore, R. & Evans, M. L. How roots perceive and respond to gravity. Am. J. Bot. 73, 574–587 (1986).
Hochholdinger, F., Woll, K., Sauer, M. & Dembinsky, D. Genetic dissection of root formation in maize (Zea mays) reveals root‐type specific developmental programmes. Ann. Bot. 93, 359–368 (2004).
Liu, Y. et al. Unfolded protein response in balancing plant growth and stress tolerance. Front. Plant Sci. 13, 1019414 (2022).
Dembinsky, D. et al. Transcriptomic and proteomic analyses of pericycle cells of the maize primary root. Plant Physiol. 145, 575–588 (2007).
DeTomaso, D. & Yosef, N. Hotspot identifies informative gene modules across modalities of single-cell genomics. Cell Syst. 12, 446–456 (2021).
Acosta, I. F. et al. tasselseed1 Is a Lipoxygenase affecting jasmonic acid signaling in sex determination of maize. Science 323, 262–265 (2009).
Alamin, M. et al. Rice SDSFL1 plays a critical role in the regulation of plant structure through the control of different phytohormones and altered cell structure. J. Plant Physiol. 231, 110–123 (2018).
Li, C. et al. Maize resistance to witchweed through changes in strigolactone biosynthesis. Science 379, 94–99 (2023).
Shahan, R. et al. A single-cell Arabidopsis root atlas reveals developmental trajectories in wild-type and cell identity mutants. Dev. Cell 57, 543–560 (2022).
Chimungu, J. G., Brown, K. M. & Lynch, J. P. Large root cortical cell size improves drought tolerance in maize. Plant Physiol. 166, 2166–2178 (2014).
Saleem, M. et al. Specification of cortical parenchyma and stele of maize primary roots by asymmetric levels of Auxin, Cytokinin, and Cytokinin-regulated proteins. Plant Physiol. 152, 4–18 (2010).
Luo, M. J. et al. Molecular dissection of maize seedling salt tolerance using a genome-wide association analysis method. Plant Biotechnol. J. 19, 1937–1951 (2021).
Qi, H. H. et al. Genome-wide identification and characterization of heat shock protein 20 genes in maize. Life Basel 12, 1397 (2022).
Jeanneau, M. et al. Improvement of drought tolerance in maize: towards the functional validation of the Zm-Asr1 gene and increase of water use efficiency by over-expressing C4-PEPC. Biochimie 84, 1127–1135 (2002).
Feng, D. et al. Chromatin accessibility illuminates single-cell regulatory dynamics of rice root tips. Bmc Biol. 20, 274 (2022).
Wang, Y., Huan, Q., Li, K. & Qian, W. F. Single-cell transcriptome atlas of the leaf and root of rice seedlings. J. Genet. Genom. 48, 881–898 (2021).
Sun, X. X. et al. Single-cell transcriptome reveals dominant subgenome expression and transcriptional response to heat stress in Chinese cabbage. Genome Biol. 23, 262 (2022).
Tan, W., Meng, Q. W., Brestic, M., Olsovska, K. & Yang, X. H. Photosynthesis is improved by exogenous calcium in heat-stressed tobacco plants. J. Plant Physiol. 168, 2063–2071 (2011).
Galsurker, O. et al. Differential response to heat stress in outer and inner onion bulb scales. J. Exp. Bot. 69, 4047–4064 (2018).
Zhang, T. Q., Chen, Y., Liu, Y., Lin, W. H. & Wang, J. W. Single-cell transcriptome atlas and chromatin accessibility landscape reveal differentiation trajectories in the rice root. Nat. Commun. 12, 2053 (2021).
Liu, Z. J. et al. Single-nucleus transcriptomes reveal spatiotemporal symbiotic perception and early response in Medicago. Nat. Plants 9, 1734–1748 (2023).
Zhao, Y. L. et al. Physiological and proteomic analyses reveal cAMP-regulated key factors in maize root tolerance to heat stress. Food Energy Secur. 10, e309 (2021).
Si, W. A. et al. Ectopic overexpression of maize heat stress transcription factor ZmHsf05 confers drought tolerance in transgenic rice. Genes 12, 1568 (2021).
Ortiz-Ramírez, C., Arevalo, E. D., Xu, X., Jackson, D. P. & Birnbaum, K. D. An efficient cell sorting protocol for maize protoplasts. Curr. Protoc. Plant Biol. 3, e20072 (2018).
Butler, A., Hoffman, P., Smibert, P., Papalexi, E. & Satija, R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. Nat. Biotechnol. 36, 411–420 (2018).
McGinnis, C. S., Murrow, L. M. & Gartner, Z. J. DoubletFinder: Doublet detection in single-cell RNA sequencing data using artificial nearest neighbors. Cell Syst. 8, 329–337 (2019).
Liang, Z. et al. DNA N6-Adenine methylation in arabidopsis thaliana. Dev. Cell 45, 406–416 (2018).
Chen, S. F., Zhou, Y. Q., Chen, Y. R. & Gu, J. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 34, 884–890 (2018).
Kim, D., Paggi, J. M., Park, C., Bennett, C. & Salzberg, S. L. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat. Biotechnol. 37, 907–915 (2019).
Anders, S., Pyl, P. T. & Huber, W. HTSeq-a Python framework to work with high-throughput sequencing data. Bioinformatics 31, 166–169 (2015).
Zhao, Y., Hu, Y. F., Dai, M. Q., Huang, L. M. & Zhou, D. X. The WUSCHEL-Related homeobox gene is required to activate shoot-borne crown root development in rice. Plant Cell 21, 736–748 (2009).
Liang, Z., Brown, R. C., Fletcher, J. C. & Opsahl-Sorteberg, H. G. Calpain-mediated positional information directs cell wall orientation to sustain plant stem cell activity, growth and development. Plant Cell Physiol. 56, 1855–1866 (2015).
Yang, W. B. et al. Visualization of protein coding, long noncoding, and nuclear RNAs by fluorescence in situ hybridization in sections of shoot apical meristems and developing flowers. Plant Physiol. 182, 147–158 (2020).
Gulati, G. S. et al. Single-cell transcriptional diversity is a hallmark of developmental potential. Science 367, 405–411 (2020).
Lex, A., Gehlenborg, N., Strobelt, H., Vuillemot, R. & Pfister, H. UpSet: Visualization of intersecting sets. IEEE Trans. Vis. Comput. Graph. 20, 1983–1992 (2014).
Emms, D. M. & Kelly, S. OrthoFinder: phylogenetic orthology inference for comparative genomics. Genome Biol. 20, 238 (2019).
Yu, G. C., Wang, L. G., Han, Y. Y. & He, Q. Y. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics 16, 284–287 (2012).
Lu, X. D. et al. Gene-indexed mutations in maize. Mol. Plant 11, 496–504 (2018).
Acknowledgements
This work was supported by funding from Biological Breeding-National Science and Technology Major Project (2023ZD04071) (L.L.), National Natural Science Foundation of China (32261143757) and Sustainable Development International Cooperation Program from Bill & Melinda Gates Foundation (2022YFAG1002) (to L.P), Sci-Tech Innovation 2030 Agenda (2022ZD0115703) (to L.P. and L.L.), Biological Breeding-National Science and Technology Major Project (2023ZD0407304 and 2023ZD04076) (to L.P.) and Fundamental Research Funds for Central Non-Profit of Chinese Academy of Agricultural Sciences and the Agricultural Science and Technology Innovation Program of CAAS (CAAS-CSCB-202403 and Y2024PT07).
Author information
Authors and Affiliations
Contributions
L.P. and L.L. conceived and designed the study. T.W., F.W., F.X., J.Y., and Y.W. performed the experiments. T.W., F.W., S.D., K.W., D.F., W.G., H.W., S.L., X.G., L.L. and L.P. analyzed data. T.W., F.W., L.L., and L.P. wrote the manuscript.
Corresponding authors
Ethics declarations
Competing interests
The authors declare no competing interests.
Peer review
Peer review information
Nature Communications thanks Peng Yu and Josh Cuperus for their contribution to the peer review of this work. A peer review file is available.
Additional information
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Supplementary information
Source data
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.
About this article
Cite this article
Wang, T., Wang, F., Deng, S. et al. Single-cell transcriptomes reveal spatiotemporal heat stress response in maize roots. Nat Commun 16, 177 (2025). https://doi.org/10.1038/s41467-024-55485-3
Received:
Accepted:
Published:
DOI: https://doi.org/10.1038/s41467-024-55485-3
This article is cited by
-
Recent progress in single-cell transcriptomic studies in plants
Plant Biotechnology Reports (2025)
