Single-cell transcriptomics reveal how root tissues adapt to soil stress

Plant materials and growth conditions

The rice line used in this study is Xkitaake, a Kitaake line transformed with the XA21 gene driven by the maize (Zea mays) ubiquitin promoter^15,34. To ensure that the presence of XA21 does not influence the gene-expression trends observed in our scRNA-seq analysis, we also included the non-transgenic Kitaake line for comparison. Seeds were dehulled, sterilized with 50% bleach for 30 min and rinsed five times with sterilized water. Rice seedlings were then inserted into Yoshida’s media solidified with 0.15% gellan gum (Gelzan, Caisson)³⁵, with the embryo positioned facing upwards. Rice seeds were kept at 30 °C in the dark for 2 or 3 days until they germinated.

For gel-based growth conditions, germinated rice seedlings were then transferred to a Percival growth chamber set to 28 °C, and constant light (2,000 LUX) for 2–3 days before harvesting.

To establish the soil-based growth conditions, Wedowee sandy loam soils from Johnston County, NC, USA (15% clay, 75% sand, 15 g kg⁻¹ organic C, 1 g kg⁻¹ total N, cation exchange capacity (CEC) 6.4 meg per 100 g, base salt, 83%, P = 199 g m⁻³, K = 78 g m⁻³, Ca = 804 g m⁻³, pH 5.8) was used. Soils were air dried, crushed and then passed through a sieve with a 2-mm mesh size. To allow packing of soils to certain bulk densities, the soils were lightly sprayed with sterilized water and mixed thoroughly. Non-compacted soil condition was packed up to 1.2 g cm⁻³ (1.2 bulk density (BD)), Compacted soil was pressed to make 1.6 g cm⁻³ (1.6 BD). Soils were packed in three-dimensionally printed mesocosms at bulk density of 1.2 or 1.6 BD and then saturated with sterilized water¹⁴. In the soils (both compacted and non-compacted) used in our experiments, excess water was drained through gravitational pull to mimic the near-field capacity conditions. Germinated rice seedlings (maximum of four seedlings per mesocosm) with equivalent length radicles (roughly 0.5 cm) were placed below the soil surface (10 mm), and then grown in a Percival growth chamber set to 28 °C and constant light (2000 LUX) for 2–3 days before harvesting.

Bulk RNA-seq profiling of rice root sections from meristem, elongation and maturation zones

Sections (root tip to end of lateral root cap, meristem; end of lateral root cap to the start of root hair elongation, elongation; 1 mm each beyond the start of root hair elongation, maturation 1 and maturation 2) were collected into 10 μl of RNA-later (Ambion) in the lid of a 1.5-ml tube. Samples were frozen in liquid nitrogen and stored at −80 °C, and then processed by grinding with a blue homogenization pestle. RNA was isolated using the Zymo MagBead RNA Isolation kit according to the manufacturer’s protocol (Zymo). RNA was used as input into the Lexogen QuantSeq 3′ FWD RNA-Seq library preparation procedure according to the manufacturer’s protocol, using the unique molecular identifier (UMI) PCR add-on kit. Libraries were indexed and pooled on an Illumina NextSeq. Reads were aligned to Michigan State University Rice genome v.7 with the STAR aligner³⁶, deduplicated using UMI-Tools³⁷ and counted with HTSeq-Count.

scRNA-seq profiling of rice root protoplasts using the 10X Genomics Chromium system

For rice seedling harvesting, gel-grown rice seedlings were directly pulled out from the growth media and root tips were cut in the enzyme solution within the optimal osmotic conditions. For soil-grown rice seedlings, the three-dimensionally printed mesocosms were opened and rinsed with gentle water flow. The seedlings were exposed and further rinsed with a gentle water flow to remove attached soil particles. Gentle brushing with a small paint brush was also carried out to remove the remaining soil particles. The root tips were then cut in the enzyme solution with the optimal osmotic conditions.

For gel-based scRNA-seq protoplasting sample, roughly 1-cm root tips were harvested from 15–40 roots, chopped with sharp razor for 1 min and then placed into a 35-mm petri dish containing a 70-μm cell strainer and 4.5 ml enzyme solution (4% [w/v] cellulase (ONOZUKA R-10, GoldBio), 1.5% Macerozyme R10 (GoldBio), 0.24% Pectolyase (Sigma P3026), 0.4 M mannitol, 20 mM MES (pH 5.7), 20 mM KCl, 10 mM CaCl₂, 0.1% bovine serum albumin and 0.000194% (v/v) 2-mercaptoethanol). The digestion was incubated on an 85-rpm shaker at 28 °C for 2.5 h with extra pipette mixing every 30 min. The resulting cell solution was filtered twice through 40-μm cell strainers, transferred into a Falcon Round-Bottom Polystyrene Test Tubes and then centrifuged for 5 min at 500g in a swinging bucket rotor. The pellet was washed with 2 ml of washing solution (0.4 M mannitol, 20 mM MES (pH 5.7), 20 mM KCl, 5 mM CaCl₂, 0.1% bovine serum albumin and 0.000194% (v/v) 2-mercaptoethanol), and centrifuged again at 500g for 3 min. The washing step was repeated for one more time and the pellet resuspended in the washing solution (normally 50–80 µl) without CaCl₂ at a concentration of roughly 2,000 cells per µl. Cell concentration was counted using a C-chips disposable hemocytometer (Fuchs Rosenthal, 20 µl, VWR, catalogue no. 22-600-102).

For soil-based scRNA-seq protoplasting, the procedure mirrors that of gel-based RNA-seq protoplasting, with modifications to chopping time (reduced to 45 s) and digestion time (extended from 2.5 to 3 h). These adjustments aim to enhance protoplast yield without introducing excessive debris. Despite careful washing of soil from root tips, a significant number of epidermal cells were probably removed, potentially altering the proportions of trichoblast and atrichoblast cells under different growth conditions. We conducted root trichoblast cell-specific reporter image analysis in gel, non-compacted and compacted soil conditions and we did not see a difference in the number of cells expressing the proCSLD1:VENUS-N7 reporter³⁸.

For chromium-based droplet production, we loaded 16,000 (32,000) cells, with the aim to capture 10,000 (20,000 for High Throughput version) cells per sample with the 10X Genomics Chromium 3′ Gene expression v.3 (for sc_7), v.3.1 (for sc_108, sc_109, sc_115, sc_116, sc_192, sc_193, sc_194, sc_195 and sc_196) or v.3.1 High Throughput (for sc_199, sc_200, sc_201 and sc_202, sc_303, sc_304, sc_305 and sc_306) kits.

scRNA-seq data preprocessing

Raw sequencing reads underwent demultiplexing from Illumina BCL files to generate FASTQ files for each sample using CellRanger mkfastq (v.3.1.0, 10X Genomics). Subsequently, reads were aligned to the Oryza sativa genome BSgenome object (BSgenome.Osativa.MSU.MSU7) along with the MSU7 gene annotation file. This alignment was carried out using the scKB script within the COPILOT preprocessing pipeline⁹, which integrates kallisto³⁹ and bustools^40,41. Quality filtering of cells was performed with the R package COPILOT⁹. COPILOT uses a non-arbitrary scheme to eliminate empty droplets and low-quality cells, using a 5% mitochondrial expression threshold as the criterion for searching the initial cut-off defining low-quality cells (parameter mt.threshold set to 5). A single iteration of COPILOT filtering (parameter filtering.ratio set to 1) was applied, effectively segregating high-quality cells from the background, as indicated by barcode rank plots. To address issues related to doublets and outliers, the resulting high-quality cells underwent further filtering, removing the top 1% of cells based on UMI counts. Putative doublets were identified and removed using DoubletFinder⁴² with the estimated doublet rate from the 10X Genomics Chromium Single Cell 3′ Reagent Kit user guide.

Normalization, annotation and integration of scRNA-seq datasets

Downstream analyses were conducted using Seurat v.3.1.5. Individual processing and examination of samples were performed, followed by data normalization using SCTransform⁴³. As a standard step in scRNA-seq data processing, we identified protoplasting-induced genes using bulk RNA-seq (Supplementary Table 2). These genes were excluded from our analysis. Specifically, we conducted bulk RNA-seq comparisons between intact roots and digested roots to identify general protoplasting-induced genes. Furthermore, we compared roots digested for 2.5 versus 3 h to account for digestion time effects and further minimize their impact on gene-expression trends in the gel versus soil comparison (Extended Data Fig. 10).

All detected genes, excluding those associated with mitochondria, chloroplasts or affected by protoplasting (absolute log₂ fold change greater than or equal to 2), were retained for analysis (Supplementary Tables 2, 5, 8 and 10). PCA was executed by calculating 50 principal components using the RunPCA function (with approx = FALSE). Subsequently, uniform manifold approximation and projection (UMAP) nonlinear dimensionality reduction was computed by means of the RunUMAP function using all 50 principal components with default parameters.

These processing steps are detailed and documented in Jupyter notebooks (provided on GitHub at https://github.com/zhumy09/scRNA-seq-for-rice).

Data integration was carried out using Seurat v.3.1.5, following the Seurat reference-based integration pipeline^44,45. The sample with the highest median UMI/gene per cell and the highest number of detected genes was selected as the reference (sample name, tz2; Supplementary Data 1). The 12 WT replicates (tz2, tz1, sc_108, sc_109, sc_7, sc_115, sc_116, sc_192, sc_193, sc_194, sc_195, sc_196) were used to construct the WT atlas shown in Fig. 1, including two previously published samples (tz1, tz2; Supplementary Data 1). For the integrated object containing eight samples shown in Figs. 2 and 3, comprising gel-grown (sc_192, sc_193, sc_194, sc_195) and soil-grown samples (sc_199, sc_200, sc_201, sc_202), sample sc_201 was chosen as the reference. These processing steps are detailed and documented in Jupyter notebooks (provided on GitHub at https://github.com/zhumy09/scRNA-seq-for-rice).

The cell type annotation for both integrated objects was based on markers (Supplementary Table 3) that have been previously validated and show specific local expression on the atlas UMAP. Marker gene-expression z-scores were calculated depending on clustering. Clusters were defined using the Seurat FindClusters function by testing the modularity parameter, ranging from res = 2 (low) to res = 300 (high), until the reasonable cluster numbers were reached. Coarse and finely-resolved clusters were annotated by comparing average marker gene z-scores. Cells annotated with the same cell identity by both resolutions were considered confidently annotated, forming the consensus annotation. This combination effectively annotated rare cell types while capturing major cell types given that high resolution and low noise provided by low-resolution are balanced. New reference expression profiles for each cell type were built by averaging the expression values for cells in the consensus annotation. All cells were then re-annotated using the correlation-based approach, which calculates Pearson correlation coefficients between each cell and reference expression profiles for cell types, assigning each cell the cell type with the highest correlation coefficient.

To eliminate the potential occurrence of specific cell groups being filtered out during our COPILOT-based scRNA-seq data preprocessing, possibility as a result of induced cell stresses, we also conducted an examination of the cell type identities for the low-quality cells and found no enrichment of any particular cell type (Supplementary Data 9). This confirmed that we have inclusively incorporated high-quality cells representative of all major cell types in an unbiased manner.

For developmental-stage annotation, correlation annotation compared each cell from scRNA-seq to bulk data from morphologically defined sections (Supplementary Data 2) for both the 12-sample WT atlas and the 8-sample integrated object grown in gel versus soil.

Plotting gene-expression values on the UMAP projection

We examined the gene expression patterns by plotting the log-normalized, ‘corrected’ counts produced by the SCTransform function rather than the batch-corrected ‘integrated’ values. The UMAPs were generated with the ‘featureplot’ function in the Seurat package.

Jupyter notebooks illustrating the gene expression plotting process are available on GitHub at https://github.com/zhumy09/scRNA-seq-for-rice.

Pseudotime estimation and heatmaps of gene-expression trends

Rice root epidermal cells were extracted from the integrated Seurat objects (12 gel-grown Xkitaake). Pseudotime was then inferred on the SCT assay of the extracted epidermal cells using Monocle3 (ref. ⁴⁶). The learn_graph and order_cell functions in Monocle3 package were used to generate pseudotime metadata. Owing to the complexity of defining epidermal principal points, we opted to calculate pseudotime values separately for atrichoblast and trichoblast cells. Subsequently, these values were merged back into the pseudotime metadata. Furthermore, we manually delineated ten developmental groups. The construction of a UMAP representing the pseudotime trajectory and gene expression (SCT) was achieved using the ‘plot_cells’ and ‘plot_genes_in_pseudotime’ functions in the Monocle3 package. Differential expression analysis for genes was conducted using the ‘graph_test’ function within Monocle3. The modular expression trends of DEGs were visualized using the ComplexHeatmap package in R⁴⁷.

Jupyter notebooks illustrating the pseudotime analysis process are available on GitHub at https://github.com/zhumy09/scRNA-seq-for-rice.

Pseudobulk differential expression analysis

Pseudobulk methods, which aggregate cell-level counts for subpopulations of interest on a per-sample basis, have been identified as top performers for cross-condition comparisons in scRNA-seq^48,49. Hence, we used a pseudobulk approach implemented in muscat (multi-sample multi-group scRNA-seq analysis tools)⁴⁸.

Differential expression analysis was conducted for our non-compacted soil-based samples versus gel-based samples, as well as for our compacted soil-based samples versus non-compacted soil-based samples. Pseudobulk expression profiles for individual cell types in each sample were aggregated for these subpopulations by summing the raw counts (RNA assay) using the ‘aggregateData’ function. Subsequently, differential expression testing was performed using the edgeR method⁵⁰ incorporated in the ‘pbDS’ function. A gene was considered differentially expressed in a given subpopulation if the false discovery rate adjusted P value was less than or equal to 0.05, absolute fold change was greater than or equal to 1.5, and detection frequency was greater than or equal to 10% in any of the included conditions. GO enrichment analysis was carried out on the DEGs using the R package gprofiler2 (ref. ⁵¹). Visualizations were generated using Seurat⁴⁵, ComplexHeatmap⁴⁷ and ggplot2 (ref. ⁵²). The full tables containing gene-expression trends and GO term enrichment information for all detected genes and GO terms from the scRNA-seq data comparison across various growth conditions is available in Supplementary Data 10.

Jupyter notebooks illustrating the pseudobulk differential expression analysis process are available on GitHub at https://github.com/zhumy09/scRNA-seq-for-rice.

Spatial transcriptomic sample preparation

The spatial transcriptomic sample preparation followed the protocol provided by Resolve Biosciences, with minor adjustment. Root parts of rice seedlings were isolated and fixed in a paraformaldehyde (PFA)-Triton-X solution: 4% [w/v] PFA (Sigma, catalogue no. 158127) and 0.03% Triton-X (Fisher Sci, catalogue no. AC327371000) in 1× PBS solution. The fixation was conducted within a 20-ml glass scintillation vial (Fisher Sci, catalogue no. 03-340-25N). The vial, containing rice roots, was placed on ice under a vacuum chamber. Vacuum was applied to the rice roots for 10 min, and this was repeated four times. Subsequently, the rice roots were rinsed with 1× PBS and dehydrated with an ethanol gradient (15, 30, 50, 70, 80, 90 and 100%), each concentration for 1 h on ice. The roots were then kept in 100% ethanol overnight.

For clearing the roots, a mixture of ethanol and Histo-clear (VWR, catalogue no. 101412-878) was applied in the following concentrations: 100% ethanol, 75% ethanol + 25% Histo-clear, 50% ethanol + 50% Histo-clear, 25% ethanol + 75% Histo-clear and 2× 100% Histo-clear, each for 1 h. The Histo-clear was then aspirated, and the vial was filled halfway with a mixture of 100% Histo-clear and melted paraplast (Leica, catalogue no. 39601006). The roots were included overnight at precisely 60 °C. The top half of Histo-clear was later replaced with paraplast, following an embedding routine that involved exchanging the top half of the embedding solution twice a day for 2 or 3 days until the sample stayed at the bottom of the containers.

The embedded roots were then mounted into plastic tissue embedding moulds (VWR, catalogue no. 15160-339) with properly adjusted orientation using flamed forceps. Paraplast-embedded roots were cut into 10-µm sections. These root tissue sections were transferred to cover slips provided by Resolve Biosciences, and the cover slip was placed in a slide dryer at 42 °C overnight. To prevent detachment issues, the cover slip could be placed in a 60 °C incubator for 5–30 min before proceeding to the next step.

Tissue sections mounted were deparaffinized with Histo-clear (100% Histo-clear, 100% Histo-clear, 25% ethanol + 75% Histo-clear, 50% ethanol + 50% Histo-clear, 75% ethanol + 25% Histo-clear, 100% ethanol). This was followed by rehydration with an ethanol gradient (100, 90, 80, 70, 50, 30%). The tissue was then permeabilized with proteinase K (Invitrogen, catalogue no. 25530049) buffer: 10 µm ml⁻¹ Proteinase K, 100 mM Tris-HCl, 50 mM EDTA) and a 0.2% [w/v] glycine (Promega, catalogue no. H5073) solution. The tissue was also refixed with a 4% [w/v] PFA solution and acetylated with an acetylation solution: 0.1 M triethanolamine (Sigma, 90279), 0.5% [v:v] acetic anhydride (Sigma, catalogue no. 320102) and 0.4% [v:v] HCl in 1× PBS. Dehydration with an ethanol gradient (30, 50, 70, 80, 90, 100, 100%) followed.

Finally, SlowFade antifade Mountant (Invitrogen, catalogue no. S36967) was applied to the tissue, and the cover slip where the tissue sections were mounted was covered with another cover slip. A slide box was used to store the cover slips with root tissue, tightly sealed with parafilm and shipped with dry ice to Resolve Biosciences for messenger RNA (mRNA) detection and imaging, with Molecular Cartography technique.

In brief, preserved mRNA molecules were hybridized with specifically designed probes based on sequence complementarity. Each probe contained a long tail with many binding sites for various fluorescent dyes. These long tails facilitated several rounds of imaging of the same probe with different fluorescent colours, generating a unique barcode for each individual gene.

The probe–mRNA complexes were sequentially coloured, imaged and decoloured for several imaging rounds. Fluorescent signal images captured on the root tissue sections were processed to identify individual mRNA molecules. Detected mRNAs corresponding to the same gene were assigned a unified identity and false-coloured for clear visualization and presentation.

The raw data for the spatial transcriptomic data included in Figs. 1–3, Extended Data Figs. 2 and 5 can be found in Supplementary Data 4 (gel), Supplementary Data 6 (non-compacted soils) and Supplementary Data 8 (compacted soils).

Spatial transcriptomic data analysis

The Resolve Biosciences dataset comprises both stained root images and transcript detection profiles. Staining images using Calcofluor white to visualize cell boundaries were processed using the ImageJ app provided by Resolve Biosciences. The Molecular Cartography plugin facilitated the visualization of mRNA detection. Transcript information was stored in a .txt file, which could be loaded using the Molecular Cartography plugin. Specific genes with mRNA detection were selected, and each assigned unique colours and dot diameters. The resulting mRNA detection images were saved as screenshots. Subsequently, image brightness and contrast were adjusted using the auto-setting in ImageJ, for presentation.

It is notable that the detected mRNA levels in roots grown in compacted soils were considerably lower compared to those grown in gel and non-compacted soil conditions. We suspect this may be due to reduced fixation efficiency. Roots grown in compacted soils undergo radial expansion, enhanced barrier formation and increased mucilage secretion (data not shown), all of which probably hinder formaldehyde penetration into inner cell layers. As a result, mRNA preservation efficiency is diminished, particularly for markers in the stele tissue. Despite these challenges, we successfully identified roughly 20 robust cell type-specific markers under compacted soil conditions, as detailed in Supplementary Data 8.

Bulk RNA-seq of Xkitaake roots

Root tips (roughly 1 cm) from Xkitaake rice varieties were harvested from gel, non-compacted and compacted soil conditions, and flash-frozen in liquid nitrogen. For RNA isolation, root tips were ground to a fine powder in liquid nitrogen using a mortar and pestle, followed by the addition of 1 ml of RLT buffer to the powdered tissue. RNA was then isolated and purified using the RNeasy Mini Kit (Qiagen) according to the manufacturer’s protocol. Raw reads were processed by removing adapter sequences and filtering out low-quality nucleotides (base quality lower than 5). HISAT2 was used to align reads to the Oryza sativa (Japonica) genome, and gene-expression levels were quantified using the fragments per kilobase of transcript sequence per million mapped reads) method. Differential gene expression (log₂ fold change greater than or equal to 1.0) was analysed through read count normalization, model-dependent P value estimation (P_adj ≤ 0.05), and false discovery rate (FDR) adjustment.

Cell area quantification

A 4% agarose gel was prepared and poured into a square petri dish, allowing it to cool for 2 min. Rice roots were then embedded in the gel for 45 min. Subsequently, the agarose block containing the root tips was radially sectioned with a razor. Transverse sections, each with a thickness of 500 µm at 0.7 cm from the root tips, were transferred to slides. Calcofluor white staining, at a concentration of 10 mg ml⁻¹, was applied to the transverse sections on slides for 1 min. After aspirating the staining solution, a drop of sterilized water was added on top of the sections. The root transverse sections were imaged using Zeiss 880 Confocal microscopy (excitation wavelength 405 nm, emission wavelength 410–585 nm). For data collection with the confocal microscopy, we used Zen 2009 v.6.0.0.303.

The acquired confocal images in CZI format were converted to TIF format and opened with MorphoGraphX⁵³. The images underwent the following processing steps: (1) Gaussian blur with x-sigma, y-sigma and z-sigma set to 1 µm. (2) Edge detect with a threshold of 100, multiplier of 2, adapt factor of 0.3 and fill value of 30,000. (3) Fill holes with the x and y radii both set to 10, threshold of 10,000, depth of 0 and fill value of 30,000. (4) Marching cube surface with a cube size of 5 µm and a threshold of 2,000. (5) Subdivide meshes and smooth meshes until the final vertices number was close to 700,000. (6) Project signal with minimum distance of 18, maximum distance of 22, minimum signal of 0 and maximum signal of 60,000.

The resulting mesh files, representing the sample structure, were then manually segmented to identify individual cells. The mesh number in segmented cells facilitated the final quantification of cell areas.

Lignin and suberin imaging

Rice roots (WT (cv Nipponbare), mhz5, aba1 and aba2) grown for 3 days under ± compaction conditions were gently removed from the 3D-printed soil columns, cleaned using deionized water and a thin brush, and embedded in 4% melted agarose. The agarose blocks containing the roots were then positioned in a vibratome (Leica), cut into 100-µm thick primary root cross-sections (1–1.5 cm or 2–2.5 cm behind the root tip), and stored in 20% ethanol. For lignin staining, the cross-sections were incubated for 10 min in a 0.2% solution of basic fuchsin dissolved in ClearSee⁵⁴ and mixed 1:1 with aqueous calcofluor white to stain cell walls⁵⁵. The stained cross-sections were quickly rinsed with ClearSee and then washed for 1.5 h in fresh ClearSee, replacing the solution halfway through.

For suberin staining, the primary root cross-sections were stained for 10 min in a 0.01% fluorol yellow solution dissolved in pure ethanol, prepared from a 1% fluorol yellow solution dissolved in dimethylsulfoxide. The stained cross-sections were rinsed once with deionized water and incubated for 10 min in aqueous calcofluor white. Finally, the cross-sections were washed 2–3 times in 50% ethanol for 20 min.

For confocal imaging, primary root tip cross-sections stained for lignin or suberin were mounted in a drop of ClearSee or 50% glycerol, respectively, and positioned on a Leica SP5 inverted confocal microscope. The excitation (Ex) and emission (Em) settings used were as follows: basic fuchsin, 561 nm (Ex), 600–650 nm (Em); calcofluor white, 405 nm (Ex), 425–475 nm (Em); fluorol yellow, 488 nm (Ex) and 520–550 nm (Em). For both basic fuchsin or fluorol yellow with calcofluor white, a sequential scanning was configured with the corresponding settings mentioned above.

Lignin analysis from compacted root tips

The dried root tips (WT (cultivar Nipponbare) and mhz5 mutant root tips) were grinded into a fine powder using a microcentrifuge tube with two metal beads (3 mm) for 1 min 30 s at 20 Hz, and then solvent extracted with sequential extractions of water (1 ml, 30 min, 98 °C), ethanol (1 ml, 30 min, 76 °C), chloroform (1 ml, 30 min, 59 °C) and acetone (1 ml, 30 min, 54 °C). The extract-free samples were dried under vacuum (overnight, 50 °C) and considered as cell wall residue.

Acetyl bromide lignin was determined as previously described⁵⁶ with modifications. In brief, 1–2 mg of cell wall residue was incubated in 200 µl of acetyl bromide solution (25% acetyl bromide in glacial acetic acid) in a 2-ml Eppendorf for 3 h at 50 °C. After cooling the samples on ice, 360 µl of 2 M NaOH, 65 µl of 0.5 M hydroxylamine hydrochloride and 375 ml of glacial acetic acid were added. After centrifuging for 5 min at 14,000g, 50 µl of supernatant and 150 µl of acetic acid were added to wells of a 96-well ultraviolet transparent plate (Thermo Scientific). The absorption was measured at 280 nm with a microtitre plate reader (Microplate-reader SpectraMax 250, Sopachem), SoftMax Pro v.5 was used for collecting data and applying the extinction coefficient for grasses 17.75 g l⁻¹ cm⁻¹. Two technical replicates of each biological replicate were analysed.

Radial water loss assay

Rice seedlings (either WT or mhz5), grown for 3 days under ±compaction, were gently removed from the three-dimensionally printed soil columns. They were then delicately brushed with deionized water to remove soil particles, and the diameter of each seminal root was measured. The primary root of each seedling (4–6 seedlings were used for each replicate) was cut into a 3-cm segment, including the root tip. After gently blotting with paper towels, each segment was positioned inside a five-digit balance closed chamber (Automatic balance, Mettler Toledo) over a thin nylon mesh. The cut ends of the segments were sealed using vacuum grease (Dow Corning) before placing them in the balance.

After 1 min of equalization inside the chamber, the fresh weight was recorded and subsequently, the weight was recorded every 30 s for up to 25–30 min. A constant relative humidity was maintained by adding bags with silica gel, which maintained the relative humidity inside the chambers at 30–35%. The silica gel was replaced after every three replicates. The temperature and relative humidity were monitored using a digital logger. Following the measurements, the root segments were wrapped and preweighed in aluminium foil and placed inside a 65 °C oven for 48 h to obtain the dry mass. The dry mass was subtracted from the initial fresh mass to obtain the total water content of each replicate. Water loss at every time point was recorded to plot the cumulative water loss (percentage of total water content). The length and diameter of the roots were used to calculate the total lateral surface, and the water loss at each time point was divided by this value to obtain the radial water loss rates (µmol m⁻² s⁻¹).

Cell wall mechanical imaging in compacted soils (phonon imaging)

Phonon microscopy is an optical elastography technique that uses the phenomenon of Brillouin scattering to probe mechanical information in biological specimens with subcellular resolution. Phonon microscopy photoacoustically stimulates GHz frequency coherent acoustic phonons that, as they propagate through the specimen, periodically modulate the local refractive index that induces resonant optical scattering of a probe laser⁵⁷. Through conservation of energy, the Brillouin scattered probe photons are frequency shifted by the phonon frequency (the so-called Brillouin frequency shift) and this can be detected either using a high-resolution spectrometer as with Brillouin microscopy⁵⁸, or interferometrically in the time domain⁵⁹.

Phonon microscopy is capable of measuring a specimen’s mechanical properties through the relationship between the measured Brillouin frequency shift (f_B) and the sound velocity (v):

for normal optical incidence where n is the refractive index and λ_probe is the optical probing wavelength. Provided n is known a priori, a measurement of the Brillouin frequency shift infers a measurement of the local sound velocity, which is determined by the elasticity of the specimen in the form of the longitudinal elastic modulus (\(M={\rho v}^{2}\)).

An absolute measurement of M requires knowledge about the mass density; however, refractive index and mass density of plant cells have been shown to vary substantially less than inter-specimen and inter-environmental variation in elasticity⁶⁰. In this work, we use the relative difference in Brillouin frequency shift (\(\Delta {f}_{{\rm{B}}}\)) between the cell wall and the water:ethanol filled cytoplasm as a proxy for the relative difference in cell wall elasticity in compacted and non-compacted conditions. It is worth noting that the longitudinal modulus should not be directly compared with the Young’s modulus, as the two describe elasticity at very different time and frequency scales (for example, Hz to kHz deformations compared with GHz); however, it has been shown that there is an empirical relationship between the two quantities⁶¹.

Sample preparation and signal processing for phonon microscopy

The harvested and cleaned root tips (1.5 cm) were embedded in 4% molten agarose within a three-dimensionally printed root tip cassette. Agarose blocks containing the root tips were sectioned transversely into 50-µm slices. These root cross-sections were fixed in 20% ethanol for phonon imaging experiments. A cross-section was laid flat onto a photoacoustic transducer (200-nm thick partially transparent metal:dielectric cavity on a 170-µm sapphire cover slip), covered in roughly 50–100 µl of water:ethanol medium and then topped with a glass cover slip. Residual medium was wicked away and the cover slip sandwich was sealed shut using varnish.

Once placed into the phonon microscope, a region of interest was selected (for example, the endodermis) and a 2D raster scan was performed. A phonon time-of-flight signal was detected at each spatial pixel position, and the relative Brillouin frequency shift (Δf_B) and the acoustic attenuation (α_B) were measured for each pixel using a fast Fourier transform and wavelet transform, respectively (Extended Data Fig. 7i,j and k,l, respectively). The spatial resolution of the technique will be determined by the optical diffraction limit (a function of optical wavelength and numerical aperture), and in this case was roughly 300 nm. This is greater than the expected thickness of the cell wall, and so the technique is probing the average elasticity of the sample volume weighted by the optical intensity distribution.

To isolate the endodermal region of interest, the Brillouin and attenuation maps were manually segmented based on positioning, morphology and size. From these segmented datasets, Δf_B versus α_B cluster maps were generated and then segmented using a two-component Gaussian mixture model. This grouped the data into two clusters that were labelled ‘background’ and ‘cell wall’. Intervals of roughly 70% confidence were determined within these clusters and mean Δf_B and α_B values were calculated. The distributions identified through the two-component Gaussian mixture model are in good agreement with the spatial positions of the cell walls and cytoplasm regions.

Using the above methodology, we report in Extended Data Fig. 7n that the relative Brillouin frequency shifts in compacted endodermal cell walls are statistically significantly greater than the equivalent cell walls grown in non-compacted conditions (P < 0.0001). Furthermore, the measurements extracted from the cytoplasm regions can be used as a control, and a Yuen’s t-test indicates that the two groups are not statistically significantly different (P > 0.05). These data indicate that the compacted cell walls have greater elasticity than those grown in non-compacted soils.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.