Measurement of micro-scale soil deformation around roots using 4D synchrotron tomography and image correlation 

S. D. Keyes1, L. Cooper1, S. Duncan1, N. Koebernick1, D. M. McKay Fletcher1, C. P. Scotson1, A. V. van Veelen1 , I. Sinclair2, T. Roose1,*

1 – Bioengineering Sciences Research Group, Engineering Sciences Academic Unit, Faculty of Engineering and Environment, University of Southampton
2 – Materials Engineering Group, Engineering Sciences Academic Unit, Faculty of Engineering and Environment, University of Southampton. 

* Corresponding author: Professor Tiina Roose, Bioengineering Sciences Research Group, Engineering Sciences Academic Unit, Faculty of Engineering and Environment, University of Southampton, University Road, SO17 1BJ Southampton, UK. T.Roose@soton.ac.uk

Abstract
This study applied time lapse (4D) synchrotron X-ray computed tomography (SRXCT) to observe micro-scale interactions between plant roots and soil.   Functionally contrasting maize root tips were repeatedly imaged during ingress into soil columns of varying water content and compaction.  This yielded sequences of 3D densiometric data, representing time –resolved geometric soil and root configurations at the micron-scale.  These data were used as inputs for two full-field kinematic quantification methods, which enabled the analysis of 3D soil deformation around elongating roots.  Discrete object tracking was used to track rigid mineral grains, whilst continuum Digital Volume Correlation (DVC) was used to track grey-level patterns within local sub-volumes.  These techniques both allowed full-field soil displacements to be quantified at an intra-rhizosphere spatial sampling scale of <300 µm.  Significant differences in deformation mechanisms were identified around different phenotypes under different soil conditions.  A uniquely strong contrast was observed between intact and de-capped roots grown in dry, compacted soil.  This provides evidence that functional traits of the root cap significantly reduce the amount of soil disturbance per unit of root elongation, with this effect being particularly significant in drier soil.  
1.	Background
More efficient breeding and cultivation of crops requires improved understanding of how plants interact with their environments, with the interplay of roots and soils being particularly important [1,2].  All agricultural activity is linked directly to the functional efficiency of the plant root system [3], and in soils the co-dependent factors of mechanical impedance and water stress are strongly implicated in limiting root growth [3,4].  However, most understanding of how these factors influence root-soil interactions has been derived from studies in which the roots and soil were not directly observed.  The explicit mechanisms controlling soil deformation around roots are poorly understood for numerous reasons, including soil heterogeneity [5], the opacity of soils [6], and the difficulty of obtaining dynamic data [7].  As a result, little is known about the direct relationships between root traits, soil physical parameters, and the explicit mechanisms via which roots ingress soils.  Modelling approaches to plant and soil interactions usually overlook rhizosphere-scale heterogeneity by treating the soil as a granular continuum [3].    
There is thus a need for methods to quantify soil displacement at the single root scale, leading to mechanistic understanding of how phenotypic effects and soil parameters such as strength and water content influence root ingress [1].  Recent advances in imaging and image correlation have shown significant promise in this regard [8,9].  Vollsnes et al. used 2D image correlation (DIC) to show that contrasting conditions of soil compaction and root mucilage production led to quantifiable differences in local sand compaction around the root tips of maize seedlings [9].  Keyes et al. demonstrated a combination of 4D X-ray Computed Tomography (XCT) and digital volume correlation (DVC), allowing the 3D displacement fields around a growing maize seedling to be quantified in a field-derived soil [8].  However, neither was a strictly rhizosphere-scale approach, since both methods averaged displacements over sub-regions whose length scale was in the range of 2-3 mm.  Furthermore, pertinent questions exist as to the applicability of what are fundamentally continuum approaches at the grain scale, where continuum assumptions break down [10].    
Recent geotechnical research indicates the utility of applying discrete object tracking approaches to such questions.  These approaches to deformation analysis of granular media quantify the kinematics of distinct grains, thereby negating the need for continuum assumptions.  In sands, discrete approaches have quantitatively corroborated prior observations that important mechanical phenomena such as shear strain localisation are in part controlled by grain-scale phenomena [10,11]. 
Recent synchrotron studies at the rhizosphere scale have allowed visualisation of the soil structure surrounding roots in unprecedented detail [12–14].  At the micron-scale, the field-derived soils used in these studies can be described using a conceptual description as shown in Figure 1.  According to this description, rigid grains of homogenous, hard mineral composition (the ‘primary’ phase) are distributed within a compliant, porous matrix comprising the clay, silt and organic soil fractions (the ‘textural’ phase).  To understand how soil-root mechanical interactions are influenced by the features evident at this scale, it is necessary to conduct kinematic analyses utilising dynamic sequences of high-resolution SRXCT data.  These methods may improve modelling of root-soil interactions, and allow targeted screening of phenotypes for mechanisms of growth in different soil conditions.  Since root cap removal is known to increase penetration resistance in soil by up to 70% [15], the breeding of plants better adapted to soils with deleterious impedance to root growth is also likely to benefit from understanding how root cap traits such as mucilage production control root ingress performance.
In this study, discrete and continuum approaches to 3D kinematic quantification were applied to study soil deformation around the root tips of maize (Zea mays) seedlings.  In vivo, time-resolved (4D) SRXCT data were collected showing maize root tips growing in situ, using three soil treatments chosen to provide contrasts in soil compaction and/or water content.  Two different maize cultivars were chosen to contrast in root tip morphology and mucilage production.  The data permitted qualitative description of root soil interactions, quantification of grain-scale mechanics, and comparison of continuum and discrete approaches for a range of conditions.  

2.	Materials and Methods
2.1.	Growth chambers
A growth assay used in a previous kinematic study by Keyes et al. [8] was miniaturised (Figure 2), using 1 ml syringe barrels as soil columns, which are known to be well suited to the constraints of SRXCT imaging [12–14].  To mechanically stabilize roots for imaging, a plaster-of-Paris fixing protocol was adopted following Bengough and McKenzie [16], which  fixes the radicle in space, proximal to the seed, so that pressure generated in the meristematic zone is fully exerted against the soil.  A 3D printed root guide facilitated this step, and introduced the root tip into the soil chamber.  The basal 5.85 mm of each syringe was plugged with a rubber stopper to prevent water egress.  Between the top of this stopper and the base of the root guide, a 17.3 mm working length of chamber was able to be filled with growth media, giving a working volume of 0.3 ml.

2.2.	Growth media
The soil material was a sand-textured Eutric Cambisol [12,14,17] collected from the surface horizon of an experimental plot at Abergwyngregyn, North Wales, (53014’N, 4001’W), sieved to <5 mm, air dried for 2 days at 23±1°C then sieved between upper and lower bounds of 1180 μm and 600 μm.  This soil comprised 52.7% sand, 32.8% silt and 14.5% clay.  Three soil treatments were prepared: a loose texture with high water content, and two mechanically consolidated textures (with a bulk density of ~4/3 that of the unconsolidated treatment) prepared with contrasting (high and low) water contents.
To produce samples for the wet loose treatment, dry soil was introduced to each chamber, being lightly consolidated via tapping [18] to a volume of 0.3 ml.  The soil was hydrated to a volumetric water content of ~25% via addition 0.075 ml of Milli-Q water.
To produce samples for the wet compacted treatment, dry soil was introduced as before, to a volume of 0.2 ml, following which 0.0375 ml of Milli-Q water was added via pipette.  Another 0.2 ml of soil was introduced, followed by another 0.0375 ml of water.  The resulting 0.4 ml working volume was then compacted to a volume of 0.3 ml using a plastic piston, producing a volumetric water content of ~25%.
To produce samples for the dry compacted treatment, a ~3 mm depth of dry soil was first layered between two moistened filter papers (Whatman No. 2) of 185 mm diameter, which were sealed between sheets of polythene film, and allowed to equilibrate for 5 h at 21±10C.  The soil layer was gently homogenised by re-scattering at 0.5 h intervals, producing a medium with a gravimetrically determined bulk volumetric water content of 11.2%.  This soil was introduced to the working volume, consolidating lightly via tapping to a volume of 0.4 ml, following which it was compacted to a volume of 0.3 ml via a plastic piston, producing a volumetric water content of ~15%.
Following preparation, all samples were plugged with rubber bungs, and left at the imaging temperature of 21±10C for 6 h to allow mechanical settling and redistribution of soil water, being estimated via non-dimensionalisation of Richard’s equation [19].
2.3.	Growth conditions
Two maize cultivars were selected to provide a functional contrast in mucilage production and root cap morphology: a wild type (Zea mays cv. Earlybird), and a root-cap defective mutant (Zea mays cv. AGT1) [20].  AGT1 has a defective root cap junction, allowing the root cap to be removed.  These functionally distinct root types are henceforth referred to as decapped and intact.  
Seeds were sterilised in batches of 10, using 50 ml of sodium hypochorite solution diluted to 5% b.v. in Milli-Q water to which 3 drops of surfactant (Triton 100X) had been added.    Samples were agitated on rollers in 50 ml laboratory tubes for 15 minutes, following which they were rinsed 5 times and soaked for 24 h in Milli-Q water.  Seeds were then transferred to 120 mm culture dishes between moist filter papers (Whatman No. 2), foil-wrapped and germinated for 3 – 5 days at 21±10C.  In order to promote gravitropic elongation and root straightness, seedlings were transferred to sheets of filter paper sandwiched between sheets of foam rubber.  These were oriented vertically with root tips pointing vertically downwards.  Seedlings with optimal radicle length (~ 12 mm) were periodically taken from this assembly and transferred to a transparent growth-observation box.  Elongation was monitored for ~2 h immediately prior to imaging, with decapped specimens having root caps excised under a dissecting microscope [9] immediately before transfer to the observation box.
To conduct SRXCT imaging of each sample, a seedling was introduced to a root guide, such that ~3 mm of root tip was emerging from the basal aperture (Figure 2).  A small bead (~0.3 ml) of Plaster of Paris, mixed to 1:1 by mass with water, was applied via pipette to fix the seedling to the upper face of the guide, and the seed surrounded with moist fine sand housed in a plastic collar.  The root guide was then gently attached to the prepared soil chamber, such that the root tip penetrated 2 – 3 mm into the working soil volume (Figure 2).      
2.4.	Imaging
Tomographic data were acquired at the TOMCAT beamline of the Swiss Light Source (Villigen, Switzerland), under a 19 keV monochromatic beam condition.  To produce each tomogram, 750 2D projections were continuously acquired at a 75 ms exposure, with 30 dark- and 50 flat-field images acquired for radiograph correction.  Each full acquisition took ~60 s in total.   The attenuated beam was scintillated by a LuAG:Ce (Cerium activated Lutetium Aluminum Garnet) film of 20 µm thickness, and resultant images digitised by a PCO Edge 5.5 camera, following 4x optical magnification.  A sample-to-detector distance of 50 mm provided an intermediate degree of propagation phase contrast.  A tomographic acquisition was automatically triggered every 6 minutes, with a sequence of 8 tomograms collected for each sample over a total time of 48 minutes.
An implementation of the Paganin single-distance phase retrieval algorithm [21] was applied to all projections, which were then reconstructed to 16-bit volumes using a high-power computing cluster and an implementation [22] of the Gridrec algorithm [23].  This produced a time-series of 8 SRXCT volumes per sample, with a temporal spacing of 6 minutes and a voxel edge length of 1.6 μm.  Three biological replicates of each soil/cultivar combination were imaged in this manner.  Collecting sequences of 8 scans for each of 2 root types, three soil conditions and three replicates produced a total of 144 SRXCT volumes for analysis (representing ~3.5 TB of 16bit volume files).           
2.5.	Image analysis
2.5.1.	Pre-filtering
To quantify the mechanics of root-soil interaction, it was necessary to extract representations of soil and root morphology from the images.  Due to the computationally intensive nature of the kinematic analyses and the large volume of data, all volumes were first subsampled in FIJI [24]  by a factor of 2 along each axis, and downscaled to an 8-bit range, producing volumes of 1280 x 1280 x 960 voxels, at a voxel edge length of 3.2 μm.  
2.5.2.	Root morphology
3D root extraction
Due to an overlap in the grey-level distributions of the root and soil voxels, the morphology of the single root in each reconstructed volume was achieved using a manual approach [8].  Using the Segmentation Editor in FIJI, root cross-sections were defined at a vertical spacing 0f 25 voxels, with this sampling frequency increased slightly in the immediate vicinity of the root tip.  Automated 3D interpolation between cross-sections generated binary label fields representing the root morphology.  Root tip coordinates were exported during this step, using the point and measure tools in FIJI.  
Root distance transform
From each root morphology, a Euclidean distance transform [25] was computed, allowing all kinematic data to be related to distance from the root surface.  For every voxel in the soil domain, the algorithm returns the Euclidean distance () to the closest point on the root surface.
Root visualisation
For qualitative interpretation of the SRXCT data, a vertical sectioning plane was defined in each image volume, coincident with both the root tip coordinates and the coordinates of the root centroid in the uppermost slice of the volume.  
Soil morphology
The high-throughput approach necessitated a computationally inexpensive segmentation procedure. The Paganin filtering applied to the data provided sufficient contrast to allow the soil phases to be segmented via global thresholding between histogram peaks representing pore-space, mixed phase and primary mineral grains.  Prior to segmentation, a 3D median filter with a 3x3 kernel was used to smooth image noise and enhance peak separation in the histogram.  The details of this thresholding approach are detailed in the Supplementary Materials (S1).  Example slices following three-phase segmentation are shown in Figure 3.
Touching grains were separated using a 3D watershed algorithm in FIJI [26].  The soil segmentations were exported as three-phase 8-bit RAW volumes with pores, textural phase and primary grains represented by discrete grey-levels (Figure 1).
Spatiotemporal analysis of soil volume fractions
The spatial distributions of the volume fractions of primary, mixed and pore phases around the roots were quantified using an annulus approach (Koebernick et al 2017 in press).  The Euclidean distance transform of each root was used to subsample the three-phase soil segmentations, gridding the soil domain into 160 um steps with respect to both distance from the root and depth down the soil column (~3x larger than the diameter of the smallest grains).  The volume fraction of each phase at depth  from the top of the volume and distance from root  is denoted .  The change in volume of each phase at a given location is then given by: 
	,	(1)	
where and  are the volume fractions at location () at t = 0 min and t = 42 min respectively.
2.6.	Soil deformation kinematics
This study applied two distinct methods of kinematic analysis, discrete grain tracking and digital volume correlation, identified schematically in Figures 1 and 4.
2.6.1.	Discrete grain tracking
To quantify the kinematics of rigid mineral grain motion, a discrete tracking approach was implemented in MATLAB 2016b (Natick, MA), following Ando et al. [10].  This method tracks 3D objects using the assumption that they undergo rigid body motion only, with no appreciable strain within the population that is being tracked (i.e. the grains are treated as simple fiducials).  
To track grains between the time steps of each sample, each discrete grain must be uniquely identified in its initial (reference) and post-deformation (deformed) states.  Grain displacement is then defined as the vector linking the centroid coordinates between these states. Following Ando et al. [10], grains were matched on the basis of a morphometric description.  Discrete grain morphometry was quantified for the separated grain data using the BoneJ plugin in FIJI [24,27].  The parameters used were volume; surface area; 1st, 2nd and 3rd moments of inertia; major radius and minor radius, which were computed on a grain-wise basis for all grains in all 144 SRXCT volumes.  Due to memory limitation in the kinematic processing step, very small grains with a volume of  mm3 were filtered out.  Across all scans, this yielded a mean grain count of 1330 (SD=67), with a mean grain volume of  (SD= mm), and mean equivalent-sphere diameter of 0.22 mm (SD=0.05 mm).  A CSV file of morphometric data was exported for each step, with temporally sequential data pairs forming the basic input to the grain tracking algorithm.  
Grain matching
The tracking algorithm was implemented in MATLAB 2016b (Natick, MA).  For each time step, each of the grains in the reference set was paired to its most probable match in the deformed set via maximisation of morphological similarity.  
First, for eachof the 7 morphometric measures (), the morphometric difference was computed discretely for every permutation of reference and deformed grain pairs, producing a set of difference matrices () each of size ().  For a given matrix , each element is given by:

	,	(2)	
where  and  are the values of morphometric parameter  for the ith deformed grain and the jth reference grain respectively.

Each difference matrix was then converted to a normalised similarity matrix (), by dividing through by the maximum difference:
	,	(3)	
where  is the infinity norm of  (i.e. the magnitude of the maximum element value).  Each element () thus takes a value in the range, with 1 representing a perfect match between a deformed grain  and reference grain   for parameter ().
Next a correlation table was computed from the full set of similarity matrices, given by:

	,	(4)	
where  is the figure of merit for the morphometric similarity of reference grain  and deformed grain .  The indices of the element of greatest magnitude (i.e.  denote the pair of grains for which the overall morphometric match is maximised.
Next, to filter out erroneously large displacements, entries in  were masked if they would result in displacements greater than the root tip displacement ().  An () matrix of grain displacements () was defined, where each element  is the centroid-to-centroid distance between deformed grain  and reference grain .  To mask the correlation table, the condition was set over all  and  that if  then .
Grains were then matched.  In most cases,  as a consequence of image noise and segmentation errors [10].  Thus, to avoid double-matching of grains, the matching procedure always matches the larger set of grains against the smaller set.  For example, in the case of the  direction (), for each grain , the row  was searched for the maximum value.  The reference grain  corresponding to this element was taken to be the best morphometric match to grain .
The measurements exported for kinematic analysis were the centroid coordinates of each grain in the reference state, the displacement vector () between reference and deformed centroid coordinates, the grain volume, and the angle between root tip and grain displacement vectors (), given by:  

	,	(5)	
where  is the root tip elongation vector for the given step.  The displacement vector field showing all grain displacements was also exported in VTK format for 3D visualisation (Figure 4).
2.6.2.	Digital Volume Correlation
Alongside the discrete grain tracking approach, Digital Volume Correlation (DVC) was also applied to all SRXCT data.  DVC is a continuum, volume-averaged approach which tracks 3D grey-level distributions in local sub-volumes.  Following Keyes et al. [8], the unfiltered data of all 126 sequential scan pairs were processed with a fast-Fourier DVC algorithm in DaVis 2.3.2 (LaVision, Goettingen, Germany), using cubic sub-volumes in a regular array, with an overlap of 75%.  A multi-pass procedure was used to improve accuracy.
DVC calibration
A noise study was first conducted to determine the optimal sub-volume size, using the procedure outlined in [8].  A suitable sub-volume edge-length was determined to be 288 μm. In the interests of brevity, further details of this procedure are found in the Supplementary Materials (S2).
To verify that the sub-volume size was suitable for tracking the largest possible displacements between reference and deformed states, a calibration volume pair was defined by translating a reference volume by [0, 0, 80] voxels to produce a deformed volume.  This produced a uniaxial offset to the reference volume of 256 µm, equal to the maximal stepwise root elongation magnitude observed in the study.  The deviation of each sub-volume displacement vector from the displacement and orientation of the input translation vector was then computed and plotted.  The grain tracking approach was applied to the same data for comparison, and the results shown in the Supplementary Materials (S2).
Computing strain data
To process the exported DVC data and compute the strain tensor for each sub-volume, the MATLAB code used in [8] was applied, being modified to compute the average shear strain.  The strain tensor at each sub-volume location was decomposed into the dilation tensor (i.e the compressive components) and the deviator tensor, which characterises the deformation without the influence of the dilation.  The mean shear strain of each sub-volume was computed by taking the Euclidean norm of the deviator tensor.  Due to the variability in per-step root elongation across time, and within sets of biological replicates, the mean shear strain () of each sub-volume was normalised to the root tip elongation for the given step, giving a mean shear intensity () for each sub-volume, according to:
   
	=,	(6)	
where  is the mean shear strain, and  is the root tip displacement vector for a given step.  
DVC results from sub-volumes for which the correlation coefficient was in the range c<0.8 were masked out to reduce error [8,28].  The metrics exported for kinematic analysis were the centroid coordinates of each sub-volume in the reference state, the displacement vector between reference and deformed sub-volume centroids, the angle between root tip and sub-volume displacement vectors, and the mean shear strain for each sub-volume.  The 3D displacement vector field was also exported in VTK format for visualisation.
2.7.	Result visualisation
3D displacement vector fields from grain tracking and DVC approaches were visualised in Paraview 5.1.2 (Figure 4) [8,29]. 
Due to elongation rate heterogeneity, all displacement measurements (for grains and sub-volumes) were first normalised to the root elongation for the corresponding step, giving a displacement intensity (), computed according to:
	=,	(7)	
where  is the displacement vector of a given grain or sub-volume and  is the root tip growth vector at the given step.  Thus, a grain or sub-volume which was displaced by the same magnitude as the root tip growth returned a displacement intensity of .
Four key metrics were defined to describe and compare the kinematics of grains and sub-volumes: the displacement intensity (), the displacement intensities of the vertical and perpendicular displacement components ( and ), and the angle between the root growth direction and the grain or sub-volume displacement ().  The Euclidean distance transform of the root defined the distance of each grain or sub-volume centroid from the root surface (), permitting individual grain and sub-volume kinematics to be plotted with respect to distance from the root (Figure 4).
All results were plotted with respect to (), averaging over 0.3 mm steps (Koebernick 2017 in press), to a maximum of 3 mm, which represented the furthest extent of the imaged domain.  For each of the six experimental conditions (decapped dry compact, decapped wet compact, decapped wet loose, intact dry compact, intact wet compact, intact wet loose) the mean of all time-steps over all three replicates was taken in each case.
3.	Results
3.1.	Qualitative observations
The 4D SRXCT data reveal the dynamic geometrical interaction of rhizosphere soil and single roots at a magnification an order of magnitude higher than previous 4D imaging studies [8].  Examples of these time-lapse data are seen in Supplementary Video V1.  Discrete primary mineral grains were clearly distinguishable from the textural phase, and complex grain-scale interactions between the root and the soil were observable (Figure 5).   In the interests of brevity, further discussion of these observations, including the impingement of primary grains on the root epidermis, are given in Supplementary Materials (S3).
3.2.	Kinematic quantification
Welch’s t-tests, accounting for unequal variances and sample sizes, were used to compare results.  Comparisons were considered significant when p < 0.05.  Analyses were conducted using MATLAB 2016b (Natick, MA).  Only significant results are discussed, with other observations detailed in the Supplementary Materials.
3.2.1.	Root elongation rates
Figure 6 shows the average root elongation per time-step over the full SRXCT sequence for each of the 18 samples considered.  A generally decreasing trend was observed over time in all conditions.  Elongation rates of different specimens were highly variable, even within the same condition, and only one comparison was significant; with the decapped dry compact condition showing lower elongation per step than the intact wet compact condition over all time-steps except 12 mins.  Individual elongation data are detailed in the Supplementary Materials (S4).
3.2.1.	Comparison of DVC and grain tracking results
Statistical tests conducted between results for ,  ,  and  between the two sets of kinematic data, for each conditions at all distances, showed only one statistically significant difference in 240 comparisons.
3.2.2.	Soil densification
The results of the spatially-resolved soil densification analysis revealed no significant differences.  The results are reported and briefly discussed in the Supplementary Materials (S5).
3.2.3.	Discrete grain kinematics
Figure 7 shows the averaged grain kinematic behaviour for each experimental condition, plotted with respect to distance from the root ().  For each condition, the kinematic parameters at each distance were compared with the values immediately proximate to the root surface.  Only differences significant at p<0.05 are reported.   
The overall displacement intensity () decreased over distance in the decapped dry compact and decapped wet compact cases.  In the decapped dry compact case, this was observed for distances from the root in the range , and in the decapped wet compact case for .  The vertical component of displacement intensity () also exhibited a decreasing trend with distance from the root in the decapped dry compact and decapped wet compact cases.  In the decapped dry compact case, this was observed for , and in the decapped wet compact case for .  The intact wet compact case showed a decrease for .  The lateral component of grain displacement intensity () showed no decrease with distance in any condition.    The angle of grain displacement relative to the root growth direction () increased with distance in the intact wet compact and decapped wet compact cases, for which  was elevated relative to the root-adjacent value at .
Comparisons between conditions also revealed significant differences.  The overall displacement intensity () was higher for the decapped dry compact case compared to the intact dry compact case over of the majority of the range of  ().  The vertical component of displacement intensity () was higher for the decapped dry compact case compared to the intact dry compact case for .  The lateral component of grain displacement intensity () was higher for the decapped dry compact case compared to intact dry compact case for , and was elevated in the decapped wet compact case relative to the intact wet compact case for .    The angle between root and grain displacement () showed no significant differences between conditions.

3.2.4.	Sub-volume kinematics
Figure 8 shows the averaged sub-volume kinematic behaviour for each experimental condition, plotted with respect to distance from the root ().  For each condition, the kinematic parameters at each distance were compared with the values immediately proximate to the root surface.  Only differences significant at p<0.05 are reported
The overall displacement intensity () decreased over distance in the decapped dry compact, decapped wet compact, intact wet compact, and intact wet loose cases.  In the decapped dry compact case, this was observed for , in the decapped wet compact case for , in the intact wet compact case for , and in the intact wet loose case for .  The vertical component of displacement intensity () exhibited a decreasing trend in the decapped dry compact, decapped wet compact, intact wet compact, and intact wet loose cases.  In the decapped dry compact and decapped wet compact cases, this was observed for .  In the intact wet compact case, a decrease was observed for .  In the intact wet loose case, decreased decrease was observed for .  The lateral component of grain displacement intensity () only showed a decrease with distance in the intact wet loose case, observed for.   The angle of sub-volume displacement relative to the root growth direction () showed an increasing trend in the intact wet compact case, for which  was elevated at .
Comparisons between conditions also revealed differences.  As in the grain results, the overall displacement intensity () was higher for the decapped dry compact case compared to the intact dry compact case over of the majority of the range of  ().  As in the grain results, the vertical component of displacement intensity () was higher for the decapped dry compact case compared to the intact dry compact case, observed for for .  As in the grain results, the lateral component of grain displacement intensity () was higher for the decapped dry compact case compared to the intact dry compact case, observed for .  The angle between root and grain displacement () was elevated in the intact dry compact case compared to the decapped dry compact case for  .
The shear results showed large variance, and are detailed in the Supplementary Materials (S6).
4.	Discussion
4.1.	General
This study investigated the mechanical dynamics of soil and root interactions at the rhizosphere scale.  SRXCT is the only experimental technique that can provide such data at these fine spatiotemporal scales (~ 1 µm, and <<10 scans/h).  Triple-replicated samples were imaged through time, with the sample panel designed to include contrasts in root tip traits (tip shape and presence of mucilage), soil wetting and soil compaction.  The data reveal soil and root structures in detail similar to those observed in other SRXCT studies [12–14], but also show both soil deformation and root tip elongation over time thus allowing the link between root and soil kinematics to be observed and quantified.  Reconfiguration of grains, textural phase, pore gas and pore fluid were evident under the influence of root ingress (Figure 5).  
4.2.	Soil kinematic quantification
The particular strength of the 4D SRXCT method is that it provides suitable data to parameterise kinematic analyses.  In this study, 3D primary mineral grain geometry was extracted, separated, and classified using a discrete grain analysis.  This approach has previously been applied to SRXCT data of geotechnical sands [10], but not to soils and/or root-related deformation.  A DVC approach, previously applied to root-soil interactions [8], was applied in parallel, parameterised by unsegmented grey-level data, and using a spatial scale (i.e. sub-volume edge-length) an order of magnitude finer than the previous study.  The results indicated virtually no statistically significant differences between the discrete tracking and DVC results (<0.5% of all comparisons were significantly different), indicating that despite the methods tracking different aspects of the soil geometry (i.e. the grains vs the textural phase) the fundamental characteristics of soil deformation were preserved across both methods at this scale.  Both methods revealed significant differences both within and between different conditions of root type and soil condition, at substantially sub-millimetre scales of relevance.   Across all samples, the matching efficiency was 89% (SD = 4.75%), which is slightly lower than measured in Ando et al. (92-94%), and can be accounted for by the increased complexity of the segmentation problem due to the multiple phases present.
4.2.1.	Localisation of soil deformation mechanisms
In wet and dry compacted soil conditions, decapped roots induced a vertical displacement intensity () in the soil immediately proximal to the root () that was elevated compared to the values at distances from the root of  and greater (Figures 7 and 8).   The physical interpretation of this result is that due to frictional effects, soil immediately next to the root tends to be translocated in the dominant direction of root growth, with this effect localised to the <1 mm closest to the root surface.  The sub-millimetre scale is thus relevant for understanding at least some fundamental mechanics of root/soil interaction.  The DVC results also returned significant surface-localised vertical displacements at the root surface in the intact wet consolidated and intact wet loose conditions, but these were not significant in the grain tracking results at the 95% confidence interval.  Some differences of this nature are to be expected between the methods, given that the discrete grain method neglects the textural phase, whereas the DVC method predominantly tracks this phase.
In the direction lateral to the dominant root growth direction, the displacement intensity () quantified by the grain tracking approach did not vary significantly over distance from the root in any condition.  The DVC results only indicated a reduction for the intact wet loose condition for distances of .  Thus, the movement of soil in the dominant direction of growth (i.e. associated with friction at the boundary of the elongating root) appears to be a more highly localised mechanism than the displacement of soil in the lateral direction (i.e. that associated with radial compaction of soil during root ingress).
4.2.2.	Comparison of mechanisms across conditions 
Comparison of all conditions over the range of distance () from the root surface revealed the dominant difference to be between intact and decapped roots in the dry compacted soil condition (Figures 7 and 8). In the dry compact soil, the vertical displacement intensity () was significantly elevated for the decapped condition relative to the intact condition across virtually the entire range of  in both grain tracking and DVC results.  Lateral displacement intensity () was also significantly elevated for the decapped condition relative to the intact condition in the grain results, but only immediately proximal to the root ().  These observations indicate that in drier soil, decapped roots induced greater soil deformation per unit of root growth than intact roots, with soil displacements in the dominant direction of root growth being particularly significant.  This correlates with 2D DIC results from sand-grown maize roots [9].  The putative explanation of prior studies appear to be borne out by these results; namely that the lack of mucilage, border cells and pointed cap geometry necessitate a decapped root to induce more soil deformation per unit of root growth than an intact root [9,15].  Furthermore, the current study indicates a particularly strong influence of soil water content on the significance of the differences between intact and decapped specimens.  Because soil strength generally increases as soil dries [30], this corroborates the hypothesis that the presence of root cap is particularly important to penetrate high strength soils.   The shear strain results, discussed in Supplementary Materials (S6), lend more tentative support to this conclusion.   
Substantial heterogeneity was observed across the replicates in the study, evidenced by the relatively small number of significant comparisons between conditions.  Though the individual replicates for each condition were prepared in an identical manner, the specific configuration of grains and pores at the rhizosphere scale was observed to differ from sample to sample.  A particularly striking instance is seen in Figure 5c, which shows a dry, compacted sample, in which a large macropore was present immediately beneath the root tip in the initial condition.  This permitted elongation with little impedance by the soil, and may explain why this sample exhibited a high elongation rate compared to virtually all other samples across all conditions (Supplementary Figure S4).  As previously shown [31,32], heterogeneity at the macro-pore scale is particularly important for understanding root elongation.  This complexity and variability of dynamics at the rhizosphere scale is also illustrated in Figure 9, which shows orthogonal projections of grain displacements for a series of different samples.  4D SRXCT is the only technique for quantifying rhizosphere dynamics of root-soil mechanical interaction.  As with any technique, limitations exist.  The scale of samples is small (i.e < 10 mm) to achieve sufficient spatial resolution, and to allow sufficient transmission of X-rays [33].  Though the rigid walls of the root chambers possibly exert a constraint on soil displacement and root elongation, it was nonetheless possible to quantify significant deformation mechanisms with characteristic length scales of less than a tenth that of the chamber diameter.  It was not possible to determine if the observed decrease in elongation rate over time was influenced by factors beyond soil strength.  X-ray doses were minimised as far as was practical, and are unlikely to produce significant effects over the very short timescales considered [34].  The soil kinematics were normalised to root elongation to allow the influence of variable extension rates to be neglected in the analyses.  
The small samples required to conduct grain scale analyses entail greater heterogeneity in soil structure than at coarser scales.  These limitations are difficult to overcome without larger sample sizes or more idealised growth media.  Soil deformation has been recently been quantified at much coarser spatiotemporal scales using benchtop X-ray CT and DVC [8].  However, the poorer resolution and image noise precluded a grain tracking approach, and rhizosphere-scale dynamics could not be quantified.  This merits a complementary approach; using SRXCT data to probe the micro-scale fundamentals of root-soil interactions, and coarser-resolution XCT [19,35,36] to study deformation at whole root and root system scales.  One fruitful avenue for future studies would be the design and fabrication of a triaxial rhizometer [37] suitable for in situ XCT imaging, allowing kinematic studies of root-soil interaction under controlled confining pressures.  SRXCT studies in sampled cores rather than prepared soils may also be feasible, though collection of intact samples and repeatability of conditions is a substantial challenge at the ~ 5 mm scale required.

5.	Conclusions
This paper has shown the possibilities extended by a fusion of 4D SRXCT and kinematic quantification approaches, with both discrete and continuum implementations considered.  This has allowed dynamic, rhizosphere-scale soil-root interactions to be explicitly observed in vivo and in three dimensions.  The methods were suitable to quantify significant differences between soil kinematic parameters, considered with respect to distance from the root.  In a dry and compacted soil, roots without an intact root cap induced a significantly greater degree of soil deformation in the root growth direction than intact roots per unit of root elongation.  For decapped roots in both wet and dry compacted soils, deformation was significantly localised proximal to the root surface, with these differences being significant over a sub-millimetre distance away from the root.  The results produced by the continuum and discrete approaches did not differ significantly in this study.  The scientific findings corroborate the indications of the only extant study to compare phenotypes using digital correlation, showing that an intact root cap improves ingress to soils, inducing less soil deformation per unit of root growth, probably due to the effects of lubricating mucilage and border cells [9].  The results of the present work show this effect to be particularly important in drier soil.  The direct visualisation and quantification of root-scale interactions will be valuable to the development and validation of explicit models of root-soil mechanical interaction, which can be up-scaled to the root-system scale.  
6.	Authors Contributions
SDK co-designed the study, fabricated the growth environments, conducted the imaging, wrote the analysis scripts, analysed the data and co-wrote the manuscript.
TR co-designed the study, analysed and interpreted the data, and prepared the manuscript.
NK assisted with sample preparation and imaging, contributed to the development of the image analysis protocol, analysed and interpreted the data, and assisted with preparation of the manuscript.
DMMF assisted with imaging, contributed to the development of the kinematic analysis protocol, analysed and interpreted the data, and assisted with preparation of the manuscript
CS assisted with imaging, analysed and interpreted the data, and assisted with preparation of the manuscript.
LC assisted with imaging, analysed and interpreted the data, and assisted with preparation of the manuscript.
SD assisted with statistical testing, analysed and interpreted the data.
KD analysed and interpreted the data, and assisted with preparation of the manuscript.
IS assisted with preparation of the manuscript and provided technical assistance with the image analysis and kinematic protocols.
AvV analysed and interpreted the data.
7.	Acknowledgements
We gratefully acknowledge the provision of beam-time at the TOMCAT beamline of the Swiss Light Source (PSI, Switzerland), and the assistance of Dr Goran Lovric and Dr Federica Marone-Welford of PSI in data acquisition and reconstruction.  We acknowledge the assistance of the µVIS Centre for Computed Tomography in providing computing resources, and thank µVIS fellow Dr Richard Boardman for helpful discussions on kinematic analysis and data processing.
8.	Funding Statement
SDK, TR, KD, AvV, DMMF and CS are funded by ERC Consolidator grant 646809 (Data Intensive Modelling of the Rhizosphere Processes).  NK and LC are funded by BBSRC grant BB/L026058/1 (Rhizosphere by design: breeding to select root traits that physically manipulate soil). SK is also part supported by a University of Southampton New Frontiers Fellowship.
9.	Data statement
The data used for the analyses in this study are accessible under a creative commons licence (CC BY-NC 4.0) under DOI 10.5258/SOTON/D0176, and the analysis codes are accessible under a creative commons licence (CC BY-NC 4.0) under DOI 10.5258/SOTON/D0175.
10.	Competing Interests
We have no competing interests.
11.	References
1. 	Bengough AG, McKenzie BM, Hallett PD, Valentine TA. Root elongation, water stress, and mechanical impedance: a review of limiting stresses and beneficial root tip traits. J Exp Bot  [Internet]. 2011 Jan 1;62(1):59–68. Available from: http://jxb.oxfordjournals.org/content/62/1/59.abstract
2. 	Gewin V. Food: an underground revolution. Nat News. 2010;466(7306):552–3. 
3. 	Hettiaratchi DRP, Goss MJ, Harris JA, Nye PH, Smith KA. Soil compaction and plant root growth. Philos Trans R Soc London B Biol Sci. 1990;329(1255):343–55. 
4. 	Whalley WR, Leeds-Harrison PB, Clark LJ, Gowing DJG. Use of effective stress to predict the penetrometer resistance of unsaturated agricultural soils. Soil Tillage Res. 2005;84:18–27. 
5. 	Young IM, Crawford JW. Interactions and self-organization in the soil-microbe complex. Science. 2004;304:1634–7. 
6. 	Pierret A, Kirby M, Moran C. Simultaneous X-ray imaging of plant root growth and water uptake in thin-slab systems. Plant Soil. 2003;255:361–73. 
7. 	Peth S, Nellesen J, Fischer G, Horn R. Non-invasive 3D analysis of local soil deformation under mechanical and hydraulic stresses by μCT and digital image correlation. Soil Tillage Res. 2010;111(1):3–18. 
8. 	Keyes SD, Gillard F, Soper N, Mavrogordato MN, Sinclair I, Roose T. Mapping soil deformation around plant roots using in vivo 4D X-ray Computed Tomography and Digital Volume Correlation. J Biomech. 2016;49(9):1802–11. 
9. 	Vollsnes A V., Futsaether CM, Bengough AG. Quantifying rhizosphere particle movement around mutant maize roots using time-lapse imaging and particle image velocimetry. Eur J Soil Sci. 2010;61:926–39. 
10. 	Andò E, Hall SA, Viggiani G, Desrues J, Bésuelle P. Grain-scale experimental investigation of localised deformation in sand: a discrete particle tracking approach. Acta Geotech. 2012;7(1):1–13. 
11. 	Vardoulakis I, Sulem J. Bifurcation analysis in geomaterials. Blackie Acad Prof Glas. 1995; 
12. 	Keyes SD, Daly KR, Gostling NJ, Jones DL, Talboys P, Pinzer BR, et al. Rapid report High resolution synchrotron imaging of wheat root hairs growing in soil and image based modelling of phosphate uptake. New Phytol. 2013;198:1023–9. 
13. 	Daly KR, Keyes SD, Masum S, Roose T. Image-based modelling of nutrient movement in and around the rhizosphere. J Exp Bot. 2016;erv544. 
14. 	Nestler J, Keyes SD, Wissuwa M. Root hair formation in rice (Oryza sativa L.) differs between root types and is altered in artificial growth conditions. J Exp Bot. 2016;erw115. 
15. 	Iijima M, Higuchi T, Barlow PW, Bengough AG. Root cap removal increases root penetration resistance in maize (Zea mays L.). J Exp Bot. 2003;54(390):2105–9. 
16. 	Bengough AG, McKenzie BM. Sloughing of root cap cells decreases the frictional resistance to maize (Zea mays L.) root growth. J Exp Bot. 1997;48(4):885–93. 
17. 	Masum SA, Kirk GJD, Daly KR, Roose T. The effect of non‐uniform microscale distribution of sorption sites on solute diffusion in soil. Eur J Soil Sci. 2016;67(4):514–22. 
18. 	Gregory PJ, Hutchison DJ, Read DB, Jenneson PM, Gilboy WB, Morton EJ. Non-invasive imaging of roots with high-resolution X-ray micro-tomography. Plant Soil. 2003;255:351–9. 
19. 	Tracy SR, Daly KR, Sturrock CJ, Crout NMJ, Mooney SJ, Roose T. Three‐dimensional quantification of soil hydraulic properties using X‐ray Computed Tomography and image‐based modeling. Water Resour Res. 2015;51(2):1006–22. 
20. 	Moore R, Miller I. Cellular differentiation in root caps of Zea mays that do not secrete mucilage. Plant Cell Environ. 1993;16(8):1003–9. 
21. 	Paganin D, Mayo SC, Gureyev TE, Miller PR, Wilkins SW. Simultaneous phase and amplitude extraction from a single defocused image of a homogeneous object. J Microsc. 2002;206(1):33–40. 
22. 	Marone F, Stampanoni M. Regridding reconstruction algorithm for real-time tomographic imaging. J Synchrotron Radiat. 2012;19(6):1029–37. 
23. 	Dowd BA, Campbell GH, Marr RB, Nagarkar V V, Tipnis S V, Axe L, et al. Developments in synchrotron X-ray computed microtomography at the National Synchrotron Light Source. In: SPIE’s International Symposium on Optical Science, Engineering, and Instrumentation. International Society for Optics and Photonics; 1999. p. 224–36. 
24. 	Schindelin J, Arganda-Carreras I, Frise E, Kaynig V, Longair M, Pietzsch T, et al. Fiji: an open-source platform for biological-image analysis. Nat Methods. 2012;9(7):676–82. 
25. 	Ollion J, Cochennec J, Loll F, Escudé C, Boudier T. TANGO: A generic tool for high-throughput 3D image analysis for studying nuclear organization. Bioinformatics. 2013;29(14):1840–1. 
26. 	Legland D, Arganda-Carreras I, Andrey P. MorphoLibJ: integrated library and plugins for mathematical morphology with ImageJ. Bioinformatics. 2016 Nov;32(22):3532–4. 
27. 	Doube M, Kłosowski MM, Arganda-Carreras I, Cordelières FP, Dougherty RP, Jackson JS, et al. BoneJ: Free and extensible bone image analysis in ImageJ. Bone. 2010;47(6):1076–9. 
28. 	Gillard F, Boardman R, Mavrogordato M, Hollis D, Sinclair I, Pierron F, et al. The application of digital volume correlation (DVC) to study the microstructural behaviour of trabecular bone during compression. J Mech Behav Biomed Mater. 2014;29:480–99. 
29. 	Ahrens J, Geveci B, Law C. Paraview: An end user tool for large data visualization. Vis Handbook Ed by CD Hansen CR Johnson Elsevier. 2005; 
30. 	Bengough AG, Bransby MF, Hans J, McKenna SJ, Roberts TJ, Valentine TA. Root responses to soil physical conditions; growth dynamics from field to cell. J Exp Bot. 2005;57(2):437–47. 
31. 	Stewart JB, Moran CJ, Wood JT. Macropore sheath: quantification of plant root and soil macropore association. Plant Soil. 1999;211:59–67. 
32. 	Ehlers W, Köpke U, Hesse F, Böhm W. Penetration resistance and root growth of oats in tilled and untilled loess soil. Soil Tillage Res. 1983;3(3):261–75. 
33. 	Roose T, Keyes SD, Daly KR, Carminati A, Otten W, Vetterlein D, et al. Challenges in imaging and predictive modeling of rhizosphere processes. 2016; 
34. 	Zappala S, Helliwell JR, Tracy SR, Mairhofer S, Sturrock CJ, Pridmore T, et al. Effects of X-ray dose on rhizosphere studies using X-ray computed tomography. PLoS One. 2013;8(6):e67250. 
35. 	Tracy SR, Black CR, Roberts JA, Mooney SJ. Exploring the interacting effect of soil texture and bulk density on root system development in tomato (Solanum lycopersicum L.). Environ Exp Bot. 2013;91:38–47. 
36. 	Daly KR, Mooney SJ, Bennett MJ, Crout NMJ, Roose T, Tracy SR. Assessing the influence of the rhizosphere on soil hydraulic properties using X-ray computed tomography and numerical modelling. J Exp Bot. 2015;66(8):2305–14. 
37. 	Abdalla AM, Hettiaratchi DRP, Reece AR. The mechanics of root growth in granular media. J Agric Eng Res. 1969;14(3):236–48. 
12.	Figure Legends 
Figure 1 – Schematic of the root-soil system, as visualised by SRXCT imaging.  a) The system can be described according to the following definitions, with representative regions enclosed in hatched white lines: (i) the root; (ii) the primary mineral grains, which are rigid crystalline grains of quartz and similar phases (i.e. feldspars); (iii) the textural phase, comprising particles of clay, silt and organic matter whose discrete structure is largely irresolvable at the imaging resolution; (iv) fluid-filled pore spaces; (v) gas-filled pore spaces. b) The phases can be extracted as discrete representations through appropriate segmentation.  c) The DVC approach (left) considers displacement of grey-level texture within cubic sub-volumes (represented by the red hatched line), and the discrete approach (right) considers the displacement of discrete primary mineral grains.  Red crosses indicate centroid locations in each case. d) The kinematics of each grain and sub-volume in response to a root elongation vector  is described by vectors  and  respectively (solid white arrows), having lateral and vertical components ( , ) and ( , ), represented by the hatched white arrows.
Figure 2 – Schematic of experimental setup.  (a) Seedlings were affixed to the root guide using Plaster of Paris. (b)  Samples were affixed such that a 2-3 mm length of root protruded from the base of the root guide (indicated by arrow). (c) An assembled assay immediately prior to imaging.  The seed is surrounded by moist sand in region (i), and the root is free to elongate into the working volume of soil in region (ii).  Scale bars are 15 mm in all instances.
Figure 3 – Example segmentations are shown for wet loose soil (a), wet compacted soil (b) and dry compacted soil (c).  The left column shows the raw SRXCT data, and the right column shows the segmented representation, with discrete grey-levels representing primary mineral phase (white), the textural phase (grey), and the pore phase (black).  The red boxes indicate positions of sub-regions shown in Supplementary Figure S1.
Figure 4 – The process for extracting kinematics from raw SRXCT data is shown in schematic form. The original data (a), with the root outlined in white, is spatially down-sampled to a computationally tractable size.  b) For grain tracking, grains above a certain volume threshold are extracted (i), and for DVC, the domain is gridded into sub-volumes (ii).  In the DVC case, the array of centroid locations is illustrated with the grid of red dots, and a representative sub-volume (to scale) is shown in black. c) The prepared data are passed through the respective tracking algorithms. d) The results for each grain or sub-volume are combined with the Euclidean distance transform (i) and the distance transform of vertical column depth (ii), to quantify the spatial location of each sub-volume or grain centroid.  e) Example 3D vector fields are shown for the DVC (i) and grain tracking approaches (ii).  f) Plots show displacement intensity () of sub-volumes (i) and grains (ii), as a function of distance () from the root (the colourbar indicates depth down the soil column).  A decay with increasing distance is seen.  g) Plots show the displacement intensity () for sub-volumes (i) and grains (ii) as a function of the angle () between the displacement vector of each grain/ sub-volume () and the root growth vector () (the colourbar indicates distance () from the root surface).  A maxima at ~0.4 radians is seen.
Figure 5 – Digital sections along the root axis taken at the start and end of SRXCT time-series. (a) An intact root extending into dry consolidated soil.  The white arrow shows root elongation over the entire growth period.  Two grains are highlighted in red, whose displacement is essentially perpendicular to the root growth orientation. (b) An intact root extending into wet consolidated soil.  A series of grains initially situated lateral and basal to the original position of the root tip are translated to form a consolidated network via large displacements essentially aligned to the root growth orientation. (c) An intact root extending into dry consolidated soil.  A large macropore (highlighted in red) situated beneath the initial root tip position allows the root to elongate with minimal soil disturbance.
Figure 6 – The averaged elongation rate of the three biological replicates of each experimental condition.  Error bars show one standard deviation in the mean.  Elongation rates were higher for intact versus decapped roots.  In decapped samples, elongation was higher in wet consolidated versus dry consolidated soil, and higher in loose wet soil than consolidated wet soil.  In intact samples, the differences between soil treatments were not so pronounced 
Figure 7– Averaged results from the grain tracking quantification: (a) vertical displacement intensity (), (b) lateral displacement intensity (), (c) overall displacement intensity (), and (d) angle between root and grain displacement orientation (), all with respect to distance from the root surface ().  Values are averages of the three biological replicates over all time-steps, and the error bars show one standard deviation either side of the mean.
Figure 8– Averaged results from the DVC quantification: (a) vertical displacement intensity (), (b) lateral displacement intensity (), (c) overall displacement intensity (), and (d) angle between root and sub-volume displacement orientation (), all with respect to distance from the root surface ().  Values are averages of the three biological replicates over all time-steps, and the error bars show one standard deviation either side of the mean.
Figure 9 – 3D particle displacements shown for initial the initial growth step (6 minutes) of 6 samples, showing variation in shape of the displacement field. (a) A decapped sample in wet consolidated soil, with a uniformly distributed displacement field local to the root. (b) A decapped sample in wet loose soil, inducing extensive displacement across the majority of the domain. (c) An intact sample in wet consolidated soil, inducing localised and highly oriented displacement of material into a large macropore (seen on the left of the view in column 2). (d) An intact sample in wet loose soil, inducing localised and highly oriented displacement in the vicinity of a large grain (seen on the left of the view in column 1). (e) An intact sample in wet loose soil, inducing localised and highly oriented displacement with a strongly perpendicular component. (f) A decapped sample in wet loose soil inducing a localised and highly oriented displacement beneath the root tip, in the vicinity of a large macropore (seen towards the middle right of the view in column 1).