Introduction

Maize (Zea mays) represents one of the most important cereal crops for global human and animal consumption, with a production of over 1.1 billion tons in 2018 (FAOSTAT 2020). It is considered as a staple food especially in Africa, Central and South America (Ranum et al. 2014). Because of its agronomic significance, maize has been extensively investigated. Its genome as well as transcriptome are well annotated (Schnable et al. 2009; Sekhon et al. 2011; Jiao et al. 2017), rendering it a well-suited and convenient model organism for both aboveground and belowground investigations of cereal crop genetics, physiology, development, and interactions with the environment (Strable and Scanlon 2009; Monaco et al. 2013; Hochholdinger et al. 2018).

The maize root system consists of embryonal primary and seminal roots, as well as post-embryonal brace and crown roots with mostly thin but numerous lateral roots (Hochholdinger et al. 2004). Root hairs, specialized epidermal cells that enhance the root surface, develop on all these root types, and are particularly important for plant growth in nutrient-poor soils (Bates and Lynch 2000, 2001; Gahoonia et al. 2001; Klamer et al. 2019). The maize root hair elongation mutant roothairless3 (rth3) shows unimpaired growth compared to wildtype (WT) (Wen and Schnable 1994), but exhibits significantly lower grain yields in field trials (Hochholdinger et al. 2008). While it is known that the defective rth3 gene encodes for a COBRA-like protein involved in cell wall expansion and biosynthesis (Hochholdinger et al. 2008), it remains unclear whether and how the rth3 defect affects the expression of other genes in maize roots.

Crop health and yield are inherently linked to soil properties (Jin et al. 2017; Correa et al. 2019), and in recent years, research focus has been increasingly allocated to the belowground portion of the plants. In this context, complex physical, chemical and biotic interactions take place between the root system and the rhizosphere, the soil volume influenced by roots (Hiltner 1904). In order to dissect these delicate interactions and their drivers, numerous studies have investigated the role of individual components of the soil-plant-biome system, such as nutrient (Finzi et al. 2015; Daly et al. 2016) and water distribution (Carminati et al. 2010), plant genotype (Chugh et al. 2013; Gomes et al. 2018), soil properties (Correa et al. 2019) and biotic interactions (Berendsen et al. 2012; Mendes et al. 2013).

Soil texture is an important determinant for physico-chemical properties such as soil porosity, compaction, and water and nutrient fluxes (Jones 1983; Rich and Watt 2013; Rogers et al. 2016), that are associated with changes in root system architecture, rhizodeposition, and microbial communities in the rhizosphere (Berg and Smalla 2009; Shahzad and Amtmann 2017; Helliwell et al. 2019; López-Coria et al. 2019; Liu et al. 2020). Under equilibrium conditions, soil water content increases with depth due to the interplay of gravity and matric forces. The steepness of this gradient relates to the particle size distribution and can be derived from the water retention curve. Sand exhibits much steeper water retention curves compared to loam (Vetterlein et al. 2020b). In addition to water content, depth gradients in root architecture (Koebernick et al. 2014) can be expected and related to this, as well as soil biota compositional gradients, such as colonization by microbes (Engelhardt et al. 2018).

In past studies, soil system components have been often viewed in isolation, and the complex interactions of rhizosphere processes that result in modulated feedback loops between the plant, the soil and soil biota are frequently overlooked (York et al. 2016; Vetterlein et al. 2020a). In this context, maize transcriptomics studies have primarily focused on how plants respond to nutrient or water scarcity (Bi et al. 2014; Yu et al. 2018), but have so far not investigated how in general different soil properties such as soil texture and depth gradients affect maize root transcriptomes.

For this purpose, the research priority program SPP2089 of the German Research Foundation (“Spatiotemporal Organization – A Key to Rhizosphere Functions”) recently established standardized experimental platforms for laboratory and field studies (Vetterlein et al. 2020b). The aim of the experiments presented in this study was to compare maize gene expression patterns in a multi-factorial design with RNA sequencing, featuring soil texture, root hair genotype and sampling depth. For this, we used the substrates loam and sand, maize wildtype and rth3 mutant with strongly reduced root hair elongation, and sampling at three depths in a SPP2089 soil column experiment (Online Resource 1, Supplementary Fig. S1).

Our research was directed by three working hypotheses, the first being related to soil texture. Soil texture influences various physico-chemical parameters, root architecture, and microbial composition. Plants respond to changes in physico-chemical parameters, such as differences in soil moisture and nutrient distribution, for example by regulating expression levels of cell wall-modifying enzymes and adjusting root growth (Houston et al. 2016). These abiotic cues e.g. modify transcript abundances of xyloglucan endotransglycosylases, expansins, pectin modifiers, arabinogalactan protein and cellulose synthase gene families (Harb et al. 2010; Humbert et al. 2013; Wang et al. 2016a; Zhang et al. 2018). Changes in water availability and biotic interactions may induce similar changes in gene expression (Coolen et al. 2016), as e.g. genes for glycosyl hydrolases which target a range of polysaccharides and oligosaccharides of plant-interacting microbes, such as beta-glucanases and chitinases, are induced upon both drought and biotic interactions (Houston et al. 2016). Therefore, we anticipated that i) soil texture differences between loam and sand substrates predominantly affect gene expression related to plant cell wall and cytoskeleton, water channels, as well as plant immunity.

In order to complement the defective root hairs in the rth3 mutant, we expected that ii) compensatory gene expression might occur, reflected in biological processes related to epidermal differentiation, cell wall components and nutrient uptake. The presence/absence of root hairs has a soil-type specific effect on the structure of the rhizosphere’s microbial community (Robertson-Albertyn et al. 2017), therefore we also anticipated that the root hair-related maize gene expression patterns differ between soil textures, and include genes expressed during biotic interactions.

Lastly, we evaluated the relative influence of the depth of the experimental system. Here, we presumed that gradients of water (more in the bottom of the soil column), root age (younger in the bottom of the soil column), and soil biota form during the experiment. We expected that these gradients are more profound in the coarse-textured sand than the fine-textured loam. Correspondingly, their formation should be reflected by maize gene expression patterns. Aquaporins, membrane channels which guide water uptake and osmotic homeostasis of plants (Patel and Mishra 2021), are not only highly regulated by water availability in the soil (Bárzana et al. 2014), but also by phytohormones, pH, reactive oxygen species, and interactions with microbes (Javot and Maurel 2002; Mahdieh and Mostajeran 2009; Bárzana et al. 2014). Aging-related processes in plant roots are associated with changes in developmental and physiological processes relevant for growth, differentiation and stress tolerance; e.g. increased expression of genes for cell division and growth in young roots, cell wall formation in developing roots, and oxidative stress in older roots (Stelpflug et al. 2016; Hill et al. 2016; Liu et al. 2019). Furthermore, plant responses to changes in soil biota include differential expression of microbe-associated molecular pattern receptors and immune response-related genes (Hacquard et al. 2017). The third hypothesis therefore states that iii) the uppermost depth of the column is enriched in transcripts related to stress, the middle part in transcripts related to cell wall formation, and the lowest part in transcripts related to cell growth.

Materials and methods

The experiment was conducted within the framework of the research program “SPP2089 – Spatiotemporal Organization” following the experimental protocol for a standardized soil column experiment. For detailed information on the experimental setup and choice of variables, we refer to Vetterlein et al. 2020b. We applied a multi-factorial design with two substrates of differing soil textures (loam/sand) and two root genotypes (wildtype, WT/ roothairless3, rth3) at three depths (D1: 4.5–6.1 cm, D2: 9.0–10.6 cm, D3: 13.5–15.1 cm) (Online Resource 1, Supplementary Fig. S1).

Plant material and soil column preparation

In this study, we investigated two substrates with different soil textures. The substrate loam was derived from the 0–50 cm depth of a haplic Phaeozem near Schladebach, Germany (51°18’ N; 12°6’ E). The substrate sand was obtained by mixing 16.7 % w/w loam with quartz sand (Quarzwerke Weferlingen, Germany). To establish similar conditions in terms of nutrient availability, loam and sand were fertilized differently according to Vetterlein et al. 2020b. The substrates were sieved to ≤ 1 mm and filled into acrylic columns (3.5 mm radius, 5 mm wall thickness, 250 mm height). Filled columns were compacted in the same manner by stamping on a flat surface, resulting in a bulk density of 1.26 and 1.47 g cm− 3 for loam and sand, respectively.

To investigate the effect of root hair formation, we compared the Zea mays B73 wildtype (WT) to the root hair mutant roothairless3 (rth3) which exhibits defective root hair elongation (Hochholdinger et al. 2008). The rth3 mutant which is used in this study has been crossbred to the B73 line for > 8 generations and is therefore highly homozygous (Hochholdinger et al. 2008; Vetterlein et al. 2020b). Seeds were provided by the Hochholdinger lab (Universität Bonn, Germany). The surface-sterilized seeds were sown 1 cm below the soil surface and planted soil columns were irrigated to maintain 22 % (loam) / 18 % (sand) volumetric water content throughout the growth period. The planted columns were cultivated for 22 days in a climate chamber under standardized conditions (350 µmol m2 s− 1 photosynthetically active radiation, 65 % humidity, day/night cycle of 12 h/12 h at a temperature of 22 °C/18°C). The cumulative irrigation volume per plant, shoot height and leaf number were determined at harvest (Online Resource 1, Supplementary Fig. S1).

Sampling of shoot, roots and rhizosphere

Sampling took place within a time frame of approximately 2–6 h after onset of light in the climate chamber, in order to minimize potential variation due to circadian gene expression. Biological replicates were chosen from each treatment referring to identical sampling time points.

Root samples were collected from soil cross sections (1.6 cm height) at three depths (D1: 4.5–6.1 cm; D2: 9.0–10.6 cm; D3: 13.5–15.1 cm, distance from the soil top). Roots were washed for 3 × 5 s in sterile 0.3 % NaCl solution to remove adherent soil, shock-frozen in liquid nitrogen, and stored at − 80 °C until RNA extraction.

RNA extraction and sequencing

The roots were homogenized to a fine powder with mortar and pestle in liquid nitrogen. Approximately 50 mg was used for RNA extraction with the NucleoSpin RNA plant kit (Macherey-Nagel, Germany). Quantity and quality of the obtained total RNA was assessed with Nanodrop (ThermoFisher Scientific, US) and on an RNA 6000 Nano Chip on Bioanalyzer 2100 (Agilent, US). All RNA samples passed quality control with RIN values of > 8.8. Library preparation with DNBSEQ and sequencing for three biological replicates was conducted by the Beijing Genomics Institute (BGIseq-500 platform, single-end, 50 bp). RNA sequencing data were deposited to the NCBI Short Read Archive under the BioProject PRJNA666515.

Read processing and differential gene expression analysis

The raw reads were processed as follows: Raw read quality was checked using ‘FastQC’ (Andrews 2010). ‘Trimmomatic’ removed adapter sequences and low quality bases/reads (Bolger et al. 2014), Clean reads were mapped with HISAT2 onto the maize B73_4.42 reference genome (Zerbino et al. 2018; Kim et al. 2015). Mapped reads were assigned to genes and quantified with ‘featureCounts’ (Liao et al. 2014). Gene descriptions were linked from the Ensembl database (Yates et al. 2020) and GO annotations from the maize-GAMER datasets (Wimalanathan et al. 2018). For down-stream analyses, we used the following packages within R (R Core Team 2017): ‘DESeq2’ (differential gene expression analysis), ‘GOseq’ (gene-set enrichment analysis), ‘pcaExplorer’ (principal component correlation), ‘variancePartitioning’ (contribution of variation to experimental factors), and ‘vegan’ (PERMANOVA). Genes with an adjusted (Benjamini-Hochberg correction for multiple testing) p-value < 0.05 were considered as differentially expressed.

qPCR arrays

Differential gene expression was verified on a 96.96 Dynamic Array chip (Fluidigm, US). 8 RNA samples were multiplexed with 19 genes. Primers used in this study are listed in Online Resource 2, Supplementary Table S1. The cDNA was synthesized from 750 ng RNA using the SuperScript IV First-Strand Synthesis System (Thermo Fisher Scientific, US) according to the manufacturer’s protocol. Target cDNAs were pre-amplified with the MyTaq DNA Polymerase (Bioline, UK) followed by an exonuclease I treatment (New England Biolabs, US). Gene expression values were normalized to three housekeeping genes (GAPDH, Actin1, and EF-1a) and analyzed with the ΔΔCT method (Schmittgen and Livak 2008).

Results

The soil column experiment was implemented successfully and plant growth parameters were assessed on the 22nd day after sowing. Plant growth in terms of leaf number, shoot height, and fresh weight was not affected by substrate, but rth3 height and water consumption were lower compared to WT (Online Resource 1, Supplementary Fig. 1).

RNA sequencing of maize roots

Transcriptomic profiles of maize roots with two substrates (loam/sand) and two root genotypes (wildtype/roothairless3) were investigated in three depths (D1: 4.5–6.1 cm, D2: 9.0–10.6 cm, D3: 13.5–15.1 cm). RNA-sequencing yielded on average 25,486,129 raw reads with a Q20 score above 96 %. After trimming adapter and low quality sequences, 99.1–99.8 % reads survived. Of these, 82.7–88.8 % mapped uniquely onto the B73 reference genome (Online Resource 3, Supplementary Table S2).

Comparative analysis of maize root transcriptomic profiles

Differentially expressed genes (DEG, adjusted p-value < 0.05) were determined by the ‘DESeq2’ package in R. Log2-fold changes (LFC) derived from RNAseq analysis were cross-compared for 19 differentially expressed genes to reverse transcription qPCR. The results were in accordance with each other (Online Resource 1, Supplementary Fig. S2). Differences in gene expression were compared pairwise for the main effect of the experimental variables – substrate (sand vs. loam), genotype (rth3 vs. wildtype) and column depth (D3 vs. D1, D2 vs. D1, and D3 vs. D2). Comparison between substrates provided the largest number of DEG (6236 up- and 5712 down-regulated), followed by that between depth D3 and D1 (5033 up- and 4526 down-regulated). In comparison, considerably fewer DEG were found between the rth3 mutant and the wildtype (657 up- and 485 down-regulated) (Online Resource 1, Supplementary Fig. S3; Online Resource 4, Supplementary Table S3).

Principal component analysis (PCA) of normalized gene expression levels revealed that samples were strongly separated by substrate (Fig. 1a). Comparing substrates, we observed a denser clustering of loam samples than of sand samples. In principal component correlation analysis, substrate contributed significantly to PC1, depth to both PC2 and PC3, and genotype to PC4 (Fig. 1b). Cumulatively, PC1 to PC4 explained 70.6 % of observed variance. To expand on the individual quantitative contributions of substrate, genotype, and column depth, a variance partitioning analysis was performed (Fig. 2). Column depth showed the highest median variance explained, followed by substrate, and the interaction between substrate and depth. Permutational analysis of variance (PERMANOVA, Table 1) confirmed significant effects of substrate (p < 0.001), depth (p < 0.001) and the interaction between substrate and depth (p < 0.001) on the gene expression profiles, but also the effects of genotype (p < 0.01) and interaction between substrate and genotype (p < 0.05).

Fig. 1
figure 1

Principal component analysis (PCA) of the RNAseq data for 36 maize root samples. a PCA plot of variance stabilizing transformed counts. Normalization was performed with the R package ‘DESeq2’. b Principal component correlation matrix. Correlation p-values for the experimental variables were determined with the R package ‘pcaExplorer’. Field labels represent Kruskal-Wallis rank sum statistic. Cumulatively, PC1 to PC4 contributed to 71% variance (data not shown)

Full size image
Fig. 2
figure 2

Variance partitioning analysis for variance stabilizing transformed counts. Violin plots describing relative explained variance by the experimental variables and their interactions and gene frequency. Variables are sorted sequentially by highest median of explained variance to lowest

Full size image
Table 1 Effect of substrate, depth and genotype on maize root gene expression profiles
Full size table

Substrate-related differential gene expression and GO terms in maize roots

In order to find relevant gene functions and pathways for large sets of differentially expressed genes, we conducted a gene set enrichment analysis and extracted individual genes from overrepresented gene ontology (GO) categories. The full dataset including p-values and log fold changes is given in Online Resources 4 and 5 (Supplementary Tables S3 and S4).

Higher expression in sand compared to loam was found for genes in functional categories relating to chitinase activity, water transport, secondary metabolism, and cell wall. In particular, we found chitinases, aquaporins, terpene synthases, expansins, protein EXORDIUM, formin-like and extensin-like proteins, as well as xyloglucan endotransferases (Online Resource 4, Supplementary Table S3; overview in Online Resource 9, Supplementary Table S8). These were enriched in several GO terms including chitinase activity, water channel activity, flavonoid and terpenoid biosynthetic process and plant type cell wall (Fig. 3; Online Resource 5, Supplementary Table S4). The reverse (lower expression in sand compared to loam) could be observed for genes relating to cellulose biosynthesis, such as cellulose synthases and fasciclin-like arabinogalactan proteins, albeit with a weak regulation.

Fig. 3
figure 3

GO term annotation and differentially expressed genes for substrate. a GO terms were determined with the ‘GOseq’ package and the MAIZE-GAMER annotation in R from differentially expressed genes (adjusted p-value (padj) <0.05, absolute log2 fold change (LFC)>1). +/– next to the GO terms indicates whether genes were up- or down-regulated. b Selected DEG from depicted functional categories and corresponding LFC. secmet: secondary metabolism

Full size image

Depth‐related differential gene expression is highly dependent on substrate

Next, we analyzed differentially expressed genes and corresponding GO terms between different sampling depths, using subsets of genes that were either up- or down-regulated by depth, and that were differentially expressed in sand only, loam only, or in both substrates. The full dataset of differentially expressed genes and GO terms with log2 fold changes and significance values for each substrate and depth contrast is given in Online Resource 7 and 8, Supplementary Tables S6 and S7. For better overview, discussed genes and GO terms are listed in Online Resource 9, Supplementary Table S8.

Comparison between the lowest and the uppermost depth (D3 vs. D1) identified 2262 genes that were differentially expressed in the samples of loam and sand (Fig. 4a, intersect of Venn diagram). Among these, we found in both substrates a strong up-regulation in D3 of genes relating to stress, defense, and cell proliferation (Online Resource 9, Supplementary Table 8). Among stress- and defense-related genes were for example: ethylene-response transcription factors, pathogenesis-related proteins, cytochrome P450, dehydration-responsive element binding proteins, lipoxygenases, and peroxidases. Cell growth and proliferation was indicated by the up-regulation of the following genes: cyclins, expansins, actin cross-linking protein, protein EXORDIUM, formin-like protein, tetraspanin7, yucca2, and xyloglucan endotransferases. These gene functions were supported by several overrepresented GO terms, such as respiratory burst involved in defense response, response to wounding, ethylene-activated signaling pathway, and cell proliferation.

Fig. 4
figure 4

Substrate-specific genes and GO terms for depth (a-c) and genotype (d).Venn diagrams show the numbers of differentially expressed genes (DEG) in either or both substrates. Heatmaps show mean normalized gene expression levels (variance-stabilized transformed counts, vst function of the DESeq2 R package). Red indicates higher expression, blue lower expression levels. Tables list relevant GO terms, where enrichment through up- or down-regulated genes is indicated by black triangles. For GO terms, DEG with a padj <0.05 were considered. A more detailed list can be found in Online Resource 7-9, Supplementary Tables 6-8

Full size image

Specifically in sand, the D3 layer in comparison to D1 showed a higher expression of peroxidases, as well as regulation of cytochrome c oxidases, supported by the enriched GO term hydrogen peroxide catabolic process (Fig. 4a). Genes related to water stress and transport were also regulated, with a lower expression of abscisic acid hydroxylase and stress ripening in D3 and a higher expression of aquaporins, with enriched GO terms response to water deprivation and water transport. The loam-specific D3 vs. D1 comparison indicated an elevated defense response, exemplified through genes with a higher expression in D3 such as gibberellin oxidases, 1-aminocyclopropane-1-carboxylate oxidase15, LRR family proteins, jasmonate-regulated gene21, and several ethylene-responsive transcription factors. GO terms relating to this were respiratory burst involved in defense response and regulation of plant-type hypersensitive response. Genes with a higher expression in D1 compared to D3 were mostly related to cell wall functions, and less strongly regulated. GO terms referred to plant type secondary cell wall biogenesis and xylan biosynthetic process, with genes such as cellulose synthases, extensin-like proteins, and fasciclin-like arabinogalactan proteins.

The comparison of D2 vs. D1 (Fig. 4b) in sand revealed an up-regulation of peroxidases and cell wall-related genes such as xyloglucan endotransferases, extensin-like proteins, fasciclin-like arabinogalactans, and expansins. These related to enriched GO terms such as peroxidase activity and plant-type cell wall organization. Down-regulated genes in sand included multiple glutathione transferases and the enriched GO terms glutathione transferase activity and glutathione binding. In loam, we found a higher expression in D2 compared to D1 of genes referring to the enriched GO terms cell wall modification and plant-type cell wall organization (expansins, extensins, protein EXORDIUM). Several down-regulated genes encoded chitinases, enriching GO terms chitinase activity, defense response to fungus, and chitin binding.

For D3 vs. D2 (Fig. 4c) in sand peroxidases were up-regulated, with the corresponding GO terms regulation of hydrogen peroxide metabolic process and respiratory burst defense response. Negative regulation, i.e. a lower expression in D2, was observed for expansins, formin-like proteins and tetraspanins. For loam, we also found a higher expression of genes in D3 compared to D2 relating to oxidative metabolism, but referring mainly to cytochrome P450. Other stress-related genes were ethylene-responsive transcription factors and gibberellin oxidases. The enriched GO terms respiratory burst involved in defense response and regulation of plant-type hypersensitive response related to these genes. Lower expression in D3 compared to D2 was found for several growth-related genes, such as actin cross-linking protein, cellulose synthases, protein EXORDIUM and xyloglucan endotransferases.

Highly depth‐dependent genes relate to growth and defense/stress

To address the core gene expression difference according to column depth and unaltered by substrate, we compiled genes that were differentially expressed between the depths in both substrates (i.e. the intersections of Venn graphs 4a–c) to a ‘core depth’ dataset. ‘Depth’ explained 57–85 % of variance for these genes (Fig. 5; Online Resource 6, Supplementary Table S5). These DEG belonged to two categories: (a) growth-related genes (formin-like 18, formin homologue 4, actin cross-linking, xyloglucan endotransglucosylase/hydrolase 30, tetraspanin-7, yucca2), and (b) defense/stress-related genes (nematode resistance protein-like HSPRO2, Peroxidase 55, N-acetylglucosaminyltransferase III, dehydration-responsive element-binding 1D, dehydration-responsive element-binding 1 A, DREB-like, and WRKY transcription factors).

Fig. 5
figure 5

Core depth gene dataset. a 45 DEGs that are “shared” between substrates. i.e. that were differentially expressed in both substrates, were determined for each depth comparison. A ‘core depth’ dataset of 45 genes was generated from the intersection of the three subsets. b Heatmap of the 45 ‘core depth’ genes showing variance-stabilized normalized counts with hierarchical clustering. Color scale corresponds to normalized gene expression

Full size image

The impact of root hair formation on maize gene expression

The number of DEG between rth3 and WT was higher in sand than in loam (1352 and 856, respectively). The small overlap (160, only 7.8 %) indicates that the root hair effect is conveyed by different sets of genes in the two substrates (Fig. 4d, Online Resource 7/8, Supplementary Table S6, S7).

The rth3 vs. WT up-regulated genes in sand contained cell wall-related genes like expansins and cellulose synthase-like, and were enriched in GO terms anatomical structure formation involved in morphogenesis and cell wall organization. Down-regulated genes in sand were related to defense-related gene functions, with chitinases, chorismate mutase, disease resistance protein and LRR receptor-like kinases, and GO terms salicylic acid biosynthetic process and regulation of plant-type hypersensitive response.

In loam, the rth3 vs. WT the down-regulated genes related to cell wall and growth (brittle stalk, cellulose synthase, fasciclin-like arabinogalactans, tubulins), and were associated with GO terms plant-type cell wall biogenesis, and cellulose biosynthetic process. Up-regulated genes were few, among them aspartyl protease AED1, bZIP transcription factor 16, and zinc finger RING type family protein (Online Resource 8, Supplementary Table S7).

The genes that were differentially expressed between rth3 and WT in both substrates included up-regulated aspartyl protease AED1, chloride channel CLC-f, xyloglucan 6-xylosyltransferase 2, bZIP transcription factor 16, asparagine synthetase 1 and nicotianamine synthase 6, and down-regulated benzoxazinone synthesis 2, high affinity nitrate transporter 2.5, COBRA-like 7, NRT1/ PTR FAMILY 6.2 and potassium transporter 5. Cell wall or membrane-related genes were predominantly down-regulated in both substrates (Online Resource 8, Supplementary Table S7).

Discussion

Substrate is the main driver of gene expression profiles

Soil texture is defined as the proportion of sand, silt and clay particles. Higher proportion of sand particles leads to a larger proportion of macropores which drain at higher water potentials than meso- or micropores. On the other hand, fine textured clay soils with a substantial proportion of micropores have high water storage capacities, but due to strong matric forces, part of the water is held too tightly for plant uptake. Volumetric water content was chosen for the two substrates (22 % in loam, 18 % in sand) to supply equal amounts of plant available water with small gradients in water content across column height (Vetterlein et al. 2020b). Nevertheless, we observed higher gene expression levels for aquaporins such as PIP1-5, PIP2-6 and TIP1-1 (Zm00001d051872, Zm00001d019565, Zm00001d027652) for roots in sand, as well as a higher total water consumption per gram biomass for plants grown on sand substrate (Online Resource 1, Supplementary Fig. S1). Aquaporins are integral membrane proteins that facilitate passive water transport along a water potential gradient (Javot and Maurel 2002; Maurel et al. 2002). Water uptake is influenced by soil texture and associated root growth rate, root architecture, individual root morphology and aquaporin expression (Hukin et al. 2002), and a significant positive correlation occurs between the root water uptake rate and the length of thinner roots (Nosalewicz and Lipiec 2014). These parameters are reflected by the root hydraulic architecture, i.e. the ratio of radial and axial conductivity for individual root segments (Doussan 1998). This suggests that thicker roots might compensate potential increase in radial resistance based on anatomy by increased aquaporin expression. In line with this, the roots in sand, where the plasma membrane aquaporin genes were up-regulated, were thicker than the roots in loam (DV, Eva Lippold, Maxime Phalempin, Steffen Schlüter, submitted).

We also observed in sand extensive substrate-related expression patterns of genes related to cell wall expansion and cytoskeleton (e.g. expansins, exordium, xyloglucan endotransferases, cellulose synthases, tubulin gamma-2, fasciclin arabinogalactans, Fig. 3b), suggesting that root growth rate and root morphology may also be affected by the substrate. Indeed, plants in sand showed thicker roots, lower total root fresh biomass and total root length density than plants grown on loam (DV, Eva Lippold, Maxime Phalempin, Steffen Schlüter, submitted), with the rth3 mutant displaying the lowest root biomass and length density. Expansins loosen plant cell walls, and whereas the expression of some expansin genes is positively correlated with root elongation (Cheng et al. 2016), other expansins are induced during root swelling (Li et al. 2014). The exposure of maize roots to 200 mM NaCl leads to growth inhibition due to the reduction in root tip cell division activity (Li et al. 2014). An enlargement of the stele tissue and cortex cells contributed to root swelling in the elongation zone, accompanied with the up-regulation of expansin A1/A3/A5/B1/B2 genes in the swollen maize roots (Li et al. 2014). By contrast, another set of expansins, A2/A6/A11/B3/B4, were up-regulated in sand vs. loam, indicating that the formation of thicker roots in sand involves a different developmental program than root swelling by salt treatment. These results imply that root growth in sand leads to the induction of certain expansin genes, resulting in a more slowly root elongation, but enhanced cell wall expansion and increased root thickness.

In addition to regulation of water and growth-related genes, a large number of genes were enriched in GO terms that indicate plant immunity functions. Plants perceive beneficial and pathogenic microorganisms by sensing microbe-associated molecular patterns by pattern recognition receptors that initiate immune responses to control the microbial communities (Trdá et al. 2015). In sand, we found such up-regulated genes including L-type lectin-domain containing receptor kinases (Zm00001d019411, Zm00001d021897, and Zm00001d021896). These receptors recognize damage-associated molecular pattern in plants and start a defense response against various pathogenic microorganisms. Gene for ethylene-insensitive 3 transcription factor that initiates downstream transcriptional cascades for ethylene responses was also induced (Zm00001d022530), as well as pathogenesis related protein 5 PR5 (Zm00001d031158). Maize PR5 is induced during an incompatible interaction with a pathogen in maize, when plant resistance triggers a hypersensitive response to restrict progression of a pathogen (Morris et al. 1998). Additionally, induction of several chitinases (e.g. endochitinase A (Zm00001d003190), basic endochitinase B (Zm00001d027524), acidic endochitinase (Zm00001d018966), chitinase B1 (Zm00001d025753)) points to an elevated defense response against fungi. To test if defense responses are more effective in sand than loam, maize seedlings could now be cultivated in both substrates and confronted by root and leaf pathogens.

Plant interacts with the rhizosphere microbiomes by secondary metabolites, and their transcripts were in general induced on sand. These included genes for trans-cinnamate 4-monooxygenase, chalcone synthase 2, flavone 3’-O-methyltransferase 1 and benzoxazinone synthesis 14 (Zm00001d016471, Zm00001d052915, Zm00001d048087, Zm00001d004921). They catalyze steps in phenylpropanoid formation to produce more disease-related secondary metabolites, such as lignin, coumaric acid, and flavonoids. For instance, maize resistance to Fusarium verticillioides involves trans-cinnamate 4-monooxygenase and 4-coumaric acid coenzyme A ligase induction (Wang et al. 2016b). Maize also produces numerous terpenoid natural products which are not only involved in plant defense against insects and microbes (Schnee et al. 2002; Ding et al. 2020), but also influence rhizosphere microbiome assembly (Murphy et al. 2021). Most of the structural diversity among terpenes can be attributed to terpene synthases (TPSs). Up-regulated genes in sand included TPS1 (Zm00001d002351), TPS6 (Zm00001d024207), TPS 8 (Zm00001d029195), TPS11 (Zm00001d024210) and TPS12 (Zm00001d024208); and genes for cytochrome P450s ZmCYP71Z19 (Zm00001d014121) ZmCYP81A38 (Zm00001d034096), and ZmCYP81A39 (Zm00001d034097), which also participate in terpenoid biosynthesis. The maize TPS1 catalyzes the formation of farnesene, nerolidol, and farnesol, which are emitted after herbivore damage of leaves to attract herbivore parasitic insects and TPS1 transcript abundance is elevated after herbivory (Schnee et al. 2002). Maize TPS6 is induced systemically in roots during leaf herbivory, but also locally when infected with the leaf pathogen Ustilago maydis (Basse 2005). Together with the three P450s ZmCYP71Z19, ZmCYP81A38 and ZmCYP81A39, the terpene synthases TPS6, TPS12 and TPS11 contribute to combinatorial terpenoid antibiotic biosynthesis (Ding et al. 2020). Using rhizosphere samples of this very same experiment, Gebauer et al. 2021 reported that substrate is the main driver of the 1-aminocyclopropane-1-carboxylic acid (ACC) deaminase gene carrying part of the microbial community. Our data indicates that differences in the structure of rhizosphere microbiome may, at least in part, be explained by maize secondary metabolism and the expression of pathogenesis-related genes.

Depth gradient relates to gene expression patterns for growth and differentiation-related gene functions

In total, 45 genes were differentially expressed in D3 vs. D1, D3 vs. D2 and D2 vs., D1 in both loam and sand, and these genes were according to variance partitioning strongly associated with the depth of the column (Online Resource 8, Supplementary Table S7). These genes contributed mainly to gene functions relating to growth, differentiation and stress, supporting our third hypothesis. Several genes that showed an expression gradient from D1 (low) to D3 (high) encoded functions related to cell growth and differentiation, such as cell proliferation and differentiation related tetraspanin (Reimann et al. 2017), primary cell wall strength and extensibility modifying xyloglucan endotransglucosylase/hydrolase (Rose et al. 2002), actin and tubulin cytoskeleton and dimension of cell growth-affecting formin (Yi et al. 2005), and tryptophan-dependent auxin biosynthesis enzyme encoded by yucca 2 (Zhao et al. 2001). This is in accordance to the observation of increased abundance of young roots in the lower soil layers (DV, Eva Lippold, Maxime Phalempin, Steffen Schlüter, submitted). The second set of depth-dependent genes are related to stress genes such as peroxidase (Passardi et al. 2005), dehydration-response element binding (Lata and Prasad 2011), nematode-resistance protein-like HSPRO2 (Williamson 1999). These genes could be a response to the soil-physical factors that vary along the depth gradient (e.g. water distribution, oxidative conditions), as for example dehydration-responsive element binding type transcription factors have been described to respond to abiotic stresses such as cold, drought, and salt stress (Agarwal et al. 2006). On the other hand, it has been shown that the expression of these transcription factors is modulated by microbial presence leading to improved stress resistance (Singh et al. 2020). The analysis of rhizosphere microbial samples taken in the same experiment shows that community composition of ACC deaminase producing bacteria was strongly affected by the depth of the soil column (Gebauer et al. 2021). Such changes may have an impact on the growth and stress resistance of maize roots (e.g. Nadeem et al. 2020).

Root hair formation plays a subdued role under non-limiting conditions

According to our second hypothesis the root hair elongation mutant exhibited differentially expressed genes that are related to plant immunity, and these responses were strongly affected by substrate. Essentially, down-regulated genes in loam related to plant cell wall biogenesis with enriched GO terms plant-type cell wall biogenesis and cellulose biosynthetic process, including genes such as brittle stalk-like proteins (Zm00001d034049, Zm00001d013256), cellulose synthases (e.g. Zm00001d009795, Zm00001d037636, Zm00001d019149, Zm00001d034553), fasciclin-like arabinogalactan proteins (Zm00001d037540, Zm00001d043929, Zm00001d039010), protein exordium (Zm00001d018105), and tetraspanin-3 (Zm00001d035501). In contrast, in sand root hair formation affected cell wall functions were limited to expansins (up-regulated: e.g. Zm00001d038476, Zm00001d029899, Zm00001d032886, Zm00001d029907) and down-regulated defense-related genes, such as chitinases (Zm00001d032946, Zm00001d043988), chorismate mutase (Zm00001d015509), disease resistance protein RPM1 (Zm00001d035172), jasmonate-induced protein (Zm00001d041656), LRR protein kinases (e.g. Zm00001d041773, Zm00001d044656, Zm00001d045785, Zm00001d049286), pathogenesis-related protein 10 (Zm00001d028816), and regulatory protein NPR1 (Zm00001d012660).

This suggests that the formation of root hairs affects the structure or function of root associated microbiomes in a substrate-specific manner. This difference has been established in barley, where the reduction in root hair elongation leads to less diverse bacterial communities and to the enrichment of specific phyla (Robertson-Albertyn et al. 2017).

The effect of root hair formation on maize gene expression was much less pronounced than the effects of substrate and sampling depth displaying only small fold changes (Online Resource 1, Supplementary Fig. S3). This might be associated with the provision of sufficient irrigation and similar although marginal levels of nutrients through fertilization in the present study (Vetterlein et al. 2020b). Plant growth and phosphate uptake investigations of root hair formation mutants suggest that the absence of root hairs only has an impact on plant growth and nutrient uptake at dry or low mineral nutrient levels, especially low phosphorus (Bates and Lynch 2000; Klamer et al. 2019). Regarding nutrient uptake, gene functions relating to phosphate uptake were not regulated, but several genes involved in nitrogen transport were less expressed in the rth3 mutant, such as nitrate transport1, nitrate transport2, HA nitrate transporter 2.5 and HA nitrate transporter 2.7 (Zm00001d054057, Zm00001d054060, Zm00001d011679, Zm00001d044504). Although the regulation of the length and density of root hairs is not taken as a major player in maize nitrogen uptake (Lynch 2013), nitrate and ammonium transporters are expressed in the root hairs (Lauter et al. 1996; Wirén et al. 2000), and the lower expression of these transporters in in the rth3 mutant could indicate that they are predominantly localized to root hairs.

Among 160 genes that were differentially expressed in maize genotype-related manner in both substrates, we mainly found down-regulated genes in the rth3 mutant that are involved in cell growth and differentiation. For instance, tetraspanins (Zm00001d028780) are intermembrane proteins that are involved in the control of proliferation and differentiation of plant cells (Reimann et al. 2017). In contrast, members of the COBRA-like gene family encode plant-specific glycosylphosphatidylinositol anchored proteins that play out in cell expansion and cell wall biosynthesis (Brady et al. 2007), and the down-regulation of maize COBRA-like 7 gene (Zm00001d028826) might be indicative of the rth3 mutation, since the rth3 gene also encodes a COBRA-like protein that belongs to a monocot-specific clade of the COBRA gene family (Hochholdinger et al. 2008). Finally, xyloglucan 6-xylosyltransferase 2 (Zm00001d028811) is involved in the biosynthesis the hemicellulose xyloglucan that is deposited during primary wall biosynthesis (Rui and Anderson 2016). In a drought stress experiment with the barley mutant rhl1.a, the root hair mutant also displayed a lower expression of cell wall-related genes including a COBRA-like 4 gene that belong to the same regulatory network that is generally induced in the onset of water stress (Kwasniewski et al. 2016). While we did not expose plants to drought, on the physiological scale, rth3 plants did consume significantly more water in relation to their biomass than their WT counterparts on the same soil (Online Resource 1, Supplementary Fig. S1). This reduced water use efficiency might be indicative of an impaired adaption to water stress.

Depth and root hair effects are highly dependent on substrate type

Analysis of maize gene expression patterns between the three depths D1 to D3 showed that the difference between D3 and D1 was the strongest. This effect can be related to root age distributions (DV, Eva Lippold, Maxime Phalempin, Steffen Schlüter, submitted). Sand root age distributions exhibit a clear depth gradient, with small, middle and large proportions of young roots for D1, D2, and D3; but in loam, depth layers D2 and D3 show relatively similar shares of young roots after 22 days cultivation compared to a low proportion of young roots in the D1 layer. In accordance to this more pronounced depth gradient in sand, the gene expression differences between D3 and D1 as well as between D2 and D1 were larger in sand than loam.

Moreover, DEG functions differed depending on the substrate: In sand substrate, deeper layers D2 and D3 displayed a higher level of gene expression compared to the top layer for oxidative metabolism, such as peroxidases (e.g. Zm00001d003707, Zm00001d040581, Zm00001d022354 in D3, Zm00001d022458, Zm00001d030199, Zm00001d014467 in D2). Of the three depths, D3 layer displayed highest expression of genes related to water availability and uptake, aquaporins (e.g. Zm00001d020383, Zm00001d051174, Zm00001d051872) and abiotic stress regulated abscisic acid stress ripening (Zm00001d025401, Zm00001d035409). Such changes in gene expression which could be interpreted as a response to differing water availability, might in fact be due to decreased gas diffusion with depth associated with an increase in water content. Under low oxygen levels, ethylene production activates oxidases (Yamauchi et al. 2017). Alternatively, the abundance of pathogenic microorganisms could be higher in D3 than D1, since defense responses to bacteria and fungi involve the accumulation of reactive oxygen species (Torres and Dangl 2005).

In loam, D3 compared to D1 defense and stress-related genes did not include peroxidases, but phytohormone-mediated defense response, including LRR family proteins (e.g. Zm00001d021836, Zm00001d051313, Zm00001d038153), jasmonate-regulated 21 (Zm00001d012456), and ethylene-responsive transcription factors (e.g. Zm00001d018191, Zm00001d027925). LRR-containing proteins play a vital role in the perception and defense against microbial and fungal pathogens (Marone et al. 2013; Trdá et al. 2015). Ethylene-responsive transcription factors are heavily implicated in response to both abiotic and biotic stresses (Park et al. 2001; Mizoi et al. 2012). This substrate-dependent response in consolidation with the above-mentioned discussion of differential gene expression in sand compared to loam regarding genes related to biotic defense, suggests the differing presence and abundance of soil microbes in the two substrates (Gebauer et al. 2021).

Perspectives

Soil texture, root hair formation, and sampling depth affected several activities that modulate rhizosphere microbiome. For instance, differentially expressed genes for microbiome perception, stress, phenylpropanoid and terpenoid metabolism were detected. The common experimental setup now allows cross-relating to other datasets, arising from the activity of the priority program SPP2089. While our sampling approach entailed a mixture of several root types, different root types of maize such as primary, seminal, and crown roots, have already been characterized by distinct physiological properties (Ahmed et al. 2018), transcriptomes (Tai et al. 2016) and root associated microorganisms (Yu et al. 2018). The experimental platform of the SPP2089 aims to unify these diverse study parameters and additionally offers the implementation of non-invasive methods such as X-ray computed tomography (Ganther et al. 2020) to associate root types with their spatial distribution. Gradients in the expression of plant defense response genes (Hill et al. 2016; Stelpflug et al. 2016) and the radial transport pathways (Duan et al. 2018) are not only present on a root system scale, but also along the length of individual roots which may relate to properties of the rhizosphere microbiome (Rüger et al. 2021). Therefore, future gene expression studies could target the zonation of the maize root system and its relation to the specialization of the rhizosphere microbiome.