Skip to main content
Download PDF

Four-dimensional measurement of root system development using time-series three-dimensional volumetric data analysis by backward prediction

Plant Methods volume 18, Article number: 133 (2022) Cite this article

Abstract

Background

Root system architecture (RSA) is an essential characteristic for efficient water and nutrient absorption in terrestrial plants; its plasticity enables plants to respond to different soil environments. Better understanding of root plasticity is important in developing stress-tolerant crops. Non-invasive techniques that can measure roots in soils nondestructively, such as X-ray computed tomography (CT), are useful to evaluate RSA plasticity. However, although RSA plasticity can be measured by tracking individual root growth, only a few methods are available for tracking individual roots from time-series three-dimensional (3D) images.

Results

We developed a semi-automatic workflow that tracks individual root growth by vectorizing RSA from time-series 3D images via two major steps. The first step involves 3D alignment of the time-series RSA images by iterative closest point registration with point clouds generated by high-intensity particles in potted soils. This alignment ensures that the time-series RSA images overlap. The second step consists of backward prediction of vectorization, which is based on the phenomenon that the root length of the RSA vector at the earlier time point is shorter than that at the last time point. In other words, when CT scanning is performed at time point A and again at time point B for the same pot, the CT data and RSA vectors at time points A and B will almost overlap, but not where the roots have grown. We assumed that given a manually created RSA vector at the last time point of the time series, all RSA vectors except those at the last time point could be automatically predicted by referring to the corresponding RSA images. Using 21 time-series CT volumes of a potted plant of upland rice (Oryza sativa), this workflow revealed that the root elongation speed increased with age. Compared with a workflow that does not use backward prediction, the workflow with backward prediction reduced the manual labor time by 95%.

Conclusions

We developed a workflow to efficiently generate time-series RSA vectors from time-series X-ray CT volumes. We named this workflow 'RSAtrace4D' and are confident that it can be applied to the time-series analysis of RSA development and plasticity.

Background

Plant roots are essential for water and nutrient uptake from the soil. Root growth pattern through the soil volume is a major determinant of the ability of plants to absorb water and nutrients [1]. For example, by growing roots in the soil area with fertilizer or water, roots can absorb water and/or nutrients efficiently [2,3,4,5]. Conversely, by avoiding root growth in polluted soils, roots do not absorb substances detrimental to plant growth [6]. The resulting specific root growth pattern is known as root system architecture (RSA) [1] and this RSA plasticity is a promising breeding target for making crops resilient to stressful soils [7]. However, there are several issues to be addressed when measuring RSA plasticity in soils.

Measurements of RSA plasticity in soils must be nondestructive, because destructive measurements disrupt the three-dimensional (3D) growth pattern of the roots in the soil and make it impossible to observe changes in RSA over time. For example, shovelomics [8], which uses shovels to excavate the roots from the soil to measure RSA traits, or the monolith method [9], which uses boxes or cylinders with an open bottom to collect soil blocks by driving the monolith into the ground, are simple and commonly used destructive methods for RSA analysis; however, neither of them meets the conditions for measuring RSA plasticity. On the other hand, widely used nondestructive measurement techniques, such as the rhizotron, the minirhizotron, or the root box, involve installing transparent plexiglass plates or tubes in the soil to capture images of the roots growing adjacent to the plexiglass surface [10,11,12,13,14]. These nondestructive methods can be used to evaluate RSA development over time by acquiring images of the glass surface and isolating root segments in the images using an image processing approach [15,16,17]. However, these nondestructive methods are limited to the observation of RSA development occurring on the plexiglass surface, and the installation of plexiglass affects root growth.

Other nondestructive techniques for in vivo root measurement, use 3D imaging with X-ray computed tomography (CT) and magnetic resonance imaging (MRI) to scan the soil profile and create 3D images that include root segments [18]. In most cases, such root segments are extracted using image processing approaches to measure root distribution in the soil as an RSA trait [19,20,21]. Thus, X-ray CT and MRI effectively overcome the shortcomings mentioned above. However, most studies have only analyzed overall parameters, such as root distribution in the soil but not specific parameters, such as individual root growth rates, as it is difficult to track all individual roots in 3D images. Nonetheless, specific parameters are most important for estimating how each root responds to the particular soil environment in which it is trying to grow, and which varies greatly even over very short distances.Numerous studies have measured the specific parameters of individual roots in two-dimensional images. Thus, for example, in an experiment using a root box with soil substrates, the root elongation rate was measured by tracking rice (Oryza sativa) and sorghum (Sorghum bicolor) root tips [22]. Similarly, in Arabidopsis (Arabidopsis thaliana), the root elongation rate was measured using solid medium plates [23,24,25]. In both cases, cameras were placed at the root tips to measure root elongation. However, using these methods, only a limited number of roots could be observed, whereby, elongation of the entire RSA could not be measured. Furthermore, even without limiting the dimension, most studies that analyze the entire RSA focus on root distribution [26,27,28], as the root segmentation and skeletonization required for root distribution measurements are relatively easy. However, the problem is more complicated when a large number of roots is to be tracked. Currently, there are very few methods for measuring specific parameters of all the roots that make up the RSA.

Once roots have extended into different portions of the soil profile, their spatial position remains unchanged unless pressure is applied to the soil or soil volume changes owing to changes in soil moisture content [29, 30]. In other words, provided a root does not change its spatial coordinates once it has elongated, then roots developing at different times will necessarily overlap among themselves, such that matching underground sections at different plant growth stages is relatively easy [27]. In this case, RSA development can be easily calculated from the overlapping sections of the soil using the RSA vector that represents the RSA skeleton as a set of sequences of coordinate points [31, 32]. Therefore, if the RSA vector for all sampling time points is available, we can measure RSA development from the overlapping portions [33]. Nonetheless, the issue with this methodology is how to efficiently generate an RSA vector for all sampling time points. A possible solution is to predict a vector at one time point, given a vector at another time point.

There are two types of predictions: forward and backward. Forward prediction predicts the future after a certain time point, whereas backward prediction predicts the past [34]. Both types of prediction have been used with sequential digital data, such as satellite images at different time points [35] and video frames [34, 36]. We hypothesized that these predictions may be possible for sequential vector data as well. Focusing on the root tip, the forward prediction needs to estimate the direction of root growth, whereas the backward prediction does not need to estimate the direction, because the one-time RSA vector perfectly overlaps with the RSA vector at a previous time point. Therefore, backward prediction is suitable for predicting time-series RSA vectors.

In this study, we developed a workflow that semi-automatically created time-series RSA vectors from X-ray CT images captured every day and calculated the time-series local RSA parameters. This was achieved by using a soil substrate whose volume barely changes with water content [19] and backward prediction to compute the RSA vector at all time points in the time series based on the RSA vector at the last time point of the time series.

Results

Workflow without backward prediction

When backward prediction is not used, an RSA vector must be created for all the time points. Figure 1 shows a schematic for measuring RSA traits from time-series X-ray CT volumes without backward prediction. In this example, six X-ray CT scans of a single potted rice plant resulted in six X-ray CT volumes. Each X-ray CT volume was converted to an RSA-segmented volume using the RSA visualization software RSAvis3D [19], and then to the RSA vector using the RSA vectorization software RSAtrace3D [31]. Then, RSA traits were calculated based on the six RSA vectors. Segmentation is fully automated; however, vectorization requires human intervention [19, 31]. Given \(n\) number of time points, the effort required is \(n\) times greater than when the number of time points is one.

Fig. 1
figure 1

Workflow without backward prediction. Segmentation: CT volumes were converted to RSA-segmented volumes by RSAvis3D. Projection views of 3D RSA-segmented volumes are shown. Vectorization: RSA-segmented volumes were vectorized to RSA vector. Projection views of 3D RSA-segmented volumes merged with 3D vectors. RSA vector is drawn in green

Full size image

Workflow with backward prediction

Figure 2 shows a schematic diagram for measuring RSA traits from time-series X-ray CT volumes using backward prediction. To perform backward prediction, the time-series CT volumes should be aligned in 3D because the position and angle of the pot on the turntable in the X-ray CT machine are different every day. After alignment, each X-ray CT volume was converted to an RSA-segmented volume. The vectorization step differs from that without backward prediction. We must create an RSA vector only at the last time point of the time series. The remaining RSA vectors at other time points were automatically created by backward prediction. Compared with the workflow without backward prediction, this workflow reduces the amount of work requiring human intervention to \(1/n\).

Fig. 2
figure 2

Workflow with backward prediction. Segmentation: CT volumes were converted to RSA-segmented volumes by RSAvis3D. Projection views of 3D RSA-segmented volumes are shown. Vectorization: RSA-segmented volumes were vectorized to an RSA vector. Projection views of 3D RSA-segmented volumes merged with 3D vectors. RSA vector is drawn in green. Backward prediction: The ochre arrows indicate the direction of the backward prediction

Full size image

Volume alignment by registration

Generally, 3D volume alignment is performed by extracting feature points from a 3D volume and performing registration with these points [37]. However, the feature points of the roots change as roots grow. Therefore, extraction of feature points from the soil substrate is desirable. To this purpose, we focused on particles such as minerals, which show high-intensity signals in the CT images of soil substrates. Minerals absorb X-rays well; therefore, they are represented as white dots in the CT volume (Fig. 3a). A 3D point cloud from the dots in the CT volume was created and converted to a 2D point cloud because the vertical positions did not need to be aligned (Fig. 3b). These point clouds were designated as SBI (Soil Block Identifier). If a soil substrate with minimal volume changes is used, the time-series SBI should be identical. We performed ICP (iterative closest point) registration [38] with the time-series SBI to align the time-series CT volumes (Fig. 3c). ICP registration using SBI was designated SBI-ICP registration. Subsequently, we used 27 CT volumes of rice RSA from 7 to 27 days after sowing (DAS) as sample materials (Teramoto et al., 2020). Figure 3d, e, show the overlaid images of 21 top-view projections of rice RSA from 7 to 27 DAS with and without SBI-ICP registration, respectively. Without SBI-ICP registration, the overlaid image was blurred (Fig. 3d), whereas, it became clear with SBI-ICP registration (Fig. 3e). The sequential animations are shown in Additional file 1: Movie S1. These results indicate that the SBI-ICP registration fully aligned CT volumes.

Fig. 3
figure 3

SBI-ICP registration. a A x–y slice of an X-ray CT image. Image size is 307.2 mm × 307.2 mm. The yellow arrowhead indicates high intensity signals derived from minerals. b A point cloud based on high intensity signals derived from minerals. The point cloud was designated as SBI (soil block identifier). c ICP (iterative closest point) registration with time-series X-ray CT volumes. d Overlaid top-view projections of 21 RSA-segmented volumes without and e with ICP registration

Full size image

Overview of backward prediction

A schematic of backward prediction is shown in Fig. 4. Given three root segments at three different time points (Fig. 4a), the root segments at earlier time points should be shorter than the root segment at the last time point (Fig. 4b). If we can determine how much shorter it is, backward prediction is possible. The RSA vector at the last time point will be overlaid completely on the root segment at the last time point, but incompletely on the root segment at the earlier time point (Fig. 4c). Therefore, with a discriminator that can judge whether the vector node overlaps with the root segment or not, we can estimate the length of the root segment and perform backward prediction.

Fig. 4
figure 4

Overview of backward prediction in this study. a Time course diagram of cultivation. There are three sampling time points, A, B, and C. The last data of the time series is obtained at sampling time point C. b Schematic diagram showing that the roots are shorter before a certain time point. The black segment indicates a root at time points A, B, and C. c Schematic diagram of overlaying a root segment drawn in (b) and a RSA vector. Green circles indicate vector nodes and green lines indicate a connection between the nodes. Each node was scored based on whether the node overlapped with the root segment. RSA vector at sampling time point C should be manually prepared. RSA vectors at sampling time points A and B are automatically predicted by backward prediction

Full size image

Discriminator for overlapping

We designed a convolutional neural network (CNN)-based simple discriminator to determine the overlap between the vector node and the root segment. This discriminator was trained using the RSA segments and vector data at the last time point. Volume blocks (17 × 17 × 17) were sampled from the RSA-segmented data as positive and negative training data, and positive and negative data were collected at coordinates where vector nodes were located and not located, respectively (Fig. 5a). The blocks were converted to scalars using the CNN, and the scalars were normalized to 0–1 using the sigmoid function as scores (Fig. 5b). A score of 1 indicated that the input data was an RSA segment and a score of 0 indicated that it was not. Because scoring of faint segments, such as a seminal root (indicated by arrowheads in Fig. 5a), by the discriminator trained with all vector data was difficult, sub-discriminators fine-tuned by roots were used for backward prediction (Fig. 5c).

Fig. 5
figure 5

Discriminator for overlapping. a Based on the vector node coordinates, blocks with size 17 × 17 × 17 voxels were cropped as positive and negative training data. The white arrowheads indicate faint root segments. b Convolutional neural network of the discriminator. Conv3D: 3D convolution, BN, batch normalization. c Fine-tuning of the discriminator model for each root

Full size image

Fully automated backward prediction

We vectorized the RSA-segmented volume of the last day (27 DAS) using RSAtrace3D [31] and performed a backward prediction. An example of backward prediction for one root at 10 DAS is shown in Fig. 6a. Approximately 120 vector nodes were used to represent the roots. The score for each node was calculated using a discriminator, and then plotted. We found a border between scores 0 and 1, indicating that the root length was shortened to the border. For all roots, all borders from 7 to 27 DAS were estimated to calculate the relative length (RL) from root tip (Fig. 6b); 0.0 indicated the full-length and 1.0 indicated roots with length 0. There were 26 roots, but one root (Root #02) did not show a change in RL (Additional file 2: Movie S2) because this root was full-length at 7 DAS. Thus, 25 trajectories were obtained. The RL of all roots was 0.0 at 27 DAS because the roots at 27 DAS were the longest. Through backward tracing, the RL of each root approached 1.0. Most of the trajectories were represented as lines bent at two points. An example of this is shown in Fig. 6c. The slope of the line indicates that RL changed. In other words, the roots elongated by spending the days it took RL to change from 0.0 to 1.0. The two points of bending can be regarded as the beginning and termination of root elongation. Given that the days of the beginning and ending of root elongation are \({x}_{1}\) and \({x}_{2}\), respectively, and that the starting point is \({P}_{1}\left({x}_{1}, 1\right)\) while the end point is \({P}_{2}\left({x}_{2}, 0\right)\), the two points can be calculated by approximating the trajectory using Eq. (1):

Fig. 6
figure 6

A test case of backward prediction. Time-series X-ray CT volumes from 7 to 27 days after sowing was used. a An example of backward prediction of a root. Node indexes (the number of nodes counted from the root tip) and scores calculated by the discriminator were plotted. The red dashed line is represented by the equation in the figure. b An example of backward prediction of all 25 roots. Node index was converted into relative length (RL) from the root tip; 0.0 indicates the full-length and 1.0 indicates the root with length 0. Each trajectory corresponded to each root. c An example of fitting for basic trajectories. A trajectory was approximated with a line having two bending points. The red dashed line is represented by the equation in the figure. d An example of fitting for overlapped trajectories. e An example of fitting for incomplete trajectories. f All 25 trajectories after approximation

Full size image
$$\begin{array}{c}y=min\left\{max\left[\frac{1}{\left({x}_{1}-{x}_{2}\right)}\left(x-{x}_{2}\right), 0\right],1\right\}\end{array}$$
(1)

However, some trajectories, such as those shown in Fig. 6d and e, could not be represented as lines that bent at two points. As shown in Fig. 6d, RL stopped changing approximately at 21 DAS, and approached 0.0 again at approximately 12 DAS. This type of trajectory was created by overlapping with another trajectory from 7 to 21 DAS. Therefore, the slope from 21 to 24 DAS reflected the original root growth. When the period from 7 to 12 DAS was ignored, the trajectory in Fig. 6d could be regarded as a line that bent into two points. Given that \(ythr\) is RL when the change in RL stopped, the two points can be calculated by approximating the trajectory using Eq. (2):

$$\begin{array}{c}y=min\left\{max\left[\frac{1}{\left({x}_{1}-{x}_{2}\right)}\left(x-{x}_{2}\right), 0\right],ythr\right\}\end{array}$$
(2)

As shown in Fig. 6e, RL had already started at 7 DAS, and the slope was interrupted before it reached 1.0. In this case, assuming the slope extends to 0 DAS, the two points can be calculated by approximating the trajectory using Eq. (1). All trajectories after approximation are plotted in Fig. 6f, which represents the elongation pattern of each root.

Reconstruction of RSA vector and quantification of root growth

Using backward-prediction results, RSA vectors at all time points were reconstructed (Additional file 3: Movie S3). This was achieved by shortening the vector at 27 DAS using RL. Vectors almost identical to those created without backward prediction were created. This time-series of vector data was used to evaluate the transition of root number and elongation rate. The root number increased linearly (Fig. 7a), and the elongation rate was significantly higher in late emerging roots (Fig. 7b). These results suggest that the elongation rate of the roots depends on the timing of root emergence.

Fig. 7
figure 7

Calculated root growth parameters in rice. a Changes in the root number over time. b Relationship between days after sowing (DAS) when roots emerge, and root elongation rate. * indicates significant Pearson’s correlation (significance level is 5%)

Full size image

Improving efficiency through backward prediction

We evaluated the amount of effort that could be saved by using backward prediction (Table 1). The program was run in a 64-bit Ubuntu 20.04 LTS (CPU: Intel® Core™ i7-8700 CPU@3.20 GHz, memory:32 GB, GPU: NVIDIA GeForce RTX 2080 Ti). Without backward prediction, RSA segmentation and vectorization required 231 min. Instead, with backward prediction, the segmentation time was the same but the vectorization time, which is the only process that requires human intervention, was reduced to 1/21. It takes some additional time for registration and backward prediction, but all the analyses were concluded in 46 min.

Table 1 Time required to quantify time-series X-ray CT images for 21 days.
Full size table

Discussion

An efficient method for creating RSA vectors from time-series 3D volumetric data is required for the measurement of RSA plasticity. In this study, we developed a semi-automatic workflow for this purpose, which consists of SBI-ICP registration based on the distribution of strong signals in soils and backward prediction that automatically creates an RSA vector before a given reference point. Using the time-series RSA vector created, we successfully measured RSA development. Thus, for example, we measured the elongation rate of individual roots. We propose that this workflow is applicable to the study of RSA plasticity. To the best of our knowledge, this is the first study to quantify RSA development by backward prediction of RSA vectors with time-series 3D volumetric data. The implementation of this workflow, which is specified for rice, was named RSAtrace4D and is available at the GitHub repository (https://github.com/st707311g/RSAtrace4D).

RSAtrace4D enables a previously difficult analysis of the growth of the individual roots that make up RSA. Using 21 CT images of an upland rice variety at 21 time points, we demonstrated that the elongation rate of roots was higher at later developmental stages (Fig. 7). If such data can be easily obtained, it can be used not only to evaluate varietal differences but also to evaluate RSA plasticity in response to soil environment. For example, RSAtrace4D could be used for developmental analysis of RSA in response to soil fertility level [26, 39,40,41,42], local soil compaction conditions [43,44,45], or global stress factors, such as drought [46], or high temperature [47]. Therefore, we believe that RSAtrace4D is an effective tool for elucidating the mechanisms of root development in heterogeneous environments and during environmental stress.

Using 21 CT images, we manually vectorized RSA at the last time point of the time series, and the RSA vectors at all time points were automatically created using backward prediction. Given \(n\) number of time points, the required labor was calculated as \(1/n\), indicating that the greater the number of time-series data, the more labor-saving it is. Additional processing time is required, as noted in Table 1; however, these processes are fully automated and have no impact on labor. The factors contributing to the success of this method are the fixation of roots in soils once they have grown, and alignment of CT volumes with strong signal in the soil (SBI-ICP registration). We used calcined clay Profile® Greens Grade™ (PROFILE Products, Buffalo, Illinois, USA) as a soil-like substrate [19]. Regardless of soil moisture conditions, its features showed little change in volume. Furthermore, the strong signal particles were moderately contaminated (Fig. 3a), making them easy to use for SBI-ICP registration. However, particularly near the base, roots are dense, and there is a risk that the position of high-signal particles near the base may shift as the roots grow. Because ICP registration uses a point cloud and not a single point, some shifting with growth occurred which was not a problem for the month-long cultivation in this study. If SBI-ICP registration is performed using natural soil, the point cloud may be altered due to changes in soil volume during cultivation, and SBI-ICP registration may not work optimally. In such cases, the field soil should be adjusted to mix with a soil substrate whose volume does not change, such as calcined clay, so that the point cloud is not altered. Alternatively, plants should be grown in an environment where the soil volume is not altered.

The method described herein that uses backward prediction has the potential for application to developmental analysis for aboveground parts of RSA, although there are some obstacles. In the current time-series analysis of aboveground traits, organ-level development can be tracked over time by plant part matching at different growth stages [48, 49]. This was achieved by matching the plant parts in 3D mesh data by aligning the 3D mesh data and regarding the overlapping parts as the same parts [49], or by matching the plant parts in 3D point-cloud data by isolating key points encoding both semantic and topology of different growth stages and performing registration of each point [48]. In both cases, there are certain requirements for time-series analysis. The first is to acquire time-series data at a high frequency to facilitate matching. A long interval between successive data acquisition events will result in the inability to perform matching due to changes in plant structure as the plant grows. Second, the matching feature points must be designated. Changing the feature points during growth makes matching more difficult. Therefore, we assumed that backward prediction for the aboveground parts would be possible if scanning intervals were as short as possible, and the shift of feature points should be compensated.

For multiple plants in one volume, we assumed that this workflow would work properly although it is necessary to vectorize each individual plant. Therefore, it is expected to be used to study plant-plant interactions. Furthermore, this workflow is considered applicable even if the type of data changes because the discriminator is retrained for each pot. The registration algorithm may need to be modified when the soil type or other parameters are changed, but the source code is available on GitHub; therefore, researchers are free to modify it. The only step that requires human intervention is the creation of vector data at the final sampling time point. At present, full automation of RSA vectorization is one of the issues that should be solved to accelerate RSA research [50]. Therefore, automatically creating RSA vectors at the last time point is a challenge that needs to be solved to enable fully automatic RSA measurements of time-series image data. Because analysis over time is required for both above- and below-ground structures, such analysis will become more widespread based on this method.

Conclusions

We developed a semi-automatic workflow for the measurement of 3D root system development from time-series X-ray computed tomography volumes using backward prediction. We observed changes in root elongation rate with growth, which indicated that this workflow is a useful tool for root growth measurement, and that it can be applied to the study of RSA plasticity responses to varying soil environments.

Methods

We used X-ray CT data as previously described, together which details of the plant materials, growth conditions, and X-ray CT scanning conditions [19]. A brief description of the process is provided below.

Plant materials and growth conditions

The upland rice variety ‘Kinandang Patong’ (IRGC #23,364) was used in the experiments described herein. Plants were cultivated in a growth chamber for 28 days in 25 cm deep and 20 cm diameter pots (TSP2530P, Tecs, Itako, Ibaraki, Japan) filled with Profile® (Greens Grade™, PROFILE Products, Buffalo, Illinois, USA), an inorganic soil amendment commonly used in cultivation systems for X-ray CT scanning. A hydroponic solution (pH 5.5) consisting of 1.23 mM NO3−, 0.41 mM NH4+, 0.18 mM H2PO4−, 1.00 mM SO42−, 1.78 mM K+, 0.55 mM Mg2+, 0.37 mM Ca2+, and 8.9 μM Fe3+ was used to saturate the growth media and ion-exchanged water was supplied from the pot bottom during cultivation. The growth conditions were 14 h photoperiod regime, temperature from 25 to 30 °C, and humidity from 50 to 60%.

X-ray CT scanning and root segmentation

Rice roots were imaged daily from 7 to 21 DAS using the X-ray CT system inspeXio SMX-225CT FPD HR (Shimadzu Corporation, Nakagyo-ku, Kyoto, Japan). The scanning conditions were as follows: tube voltage was set at 225 kV, tube current at 500 μA, 1.0-mm Cu (copper) filter, 1200 projections, using a signal averaging two frames over 360° without binning (pixel detector resolution: 3000 × 3000), at 4.0 fps. The final spatial resolution was 0.3 mm, corresponding to a total volume of 30.72 × 30.72 × 25.8 cm3. Root segments were isolated using the segmentation software RSAvis3D [19]. To eliminate the roots that changed their growth direction by touching the pot wall, an area 18 cm in diameter was isolated.

SBI-ICP registration

Python version 3.8.12 was used for SBI-ICP registration [51]. The CT volumes were loaded as NumPy arrays [52]. The regions with high intensity were isolated by thresholding, and their centroids were calculated by the ‘regionprops’ function from scikit-image package version 0.18.3 [53]. The centroids were converted into a 3D point cloud by Open3D package [54] and ICP registration was performed by the ‘registration_icp’ function in Open3D package (0.13.0).

Backward prediction

Again, Python version 3.8.12 was used for this purpose. The registered CT volumes were converted into RSA-segmented volumes using RSAvis3D. Twenty six roots in the RSA-segmented volume at 27 DAS were vectorized using RSAtrace3D. The 2835 node points that make up the vector were extracted as positive nodes. Negative nodes of the same number were randomly selected such as not to overlap with positive nodes (Fig. 5a). A 17 × 17 × 17 block was extracted from each node. Using 5670 nodes, the U-Net model was trained for 100 epochs. The structure of the model is shown in Fig. 5b. The loss was defined as the mean squared error (MSE), and the Adam optimizer was used at a learning rate of 1e-5. Training was stopped when the loss fell below 0.05 five times consecutively. By fine-tuning the trained model, the submodel for each root was trained again (Fig. 5c).

Using the trained models, the RSA vectors from 7 to 26 DAS were predicted. A score was calculated for each node, and the length of the vector was determined. Then, a graph was created with the x-axis as the node index (the number of nodes counted from the root tip) and the y-axis as the score (Fig. 6a). By approximating this plot with Eq. (3), the borderline between scores of 0 and 1, \(xthr\), was calculated.

$$ y = f\left( x \right) = \left\{ {\begin{array}{*{20}l} {0,} \hfill & {x;xthr} \hfill \\ {1,} \hfill & {x \ge xthr} \hfill \\ \end{array} } \right. $$
(3)

Approximation was performed using brute force, and \(xthr\), whose MSE was the smallest, was selected. Thus, \(xthr\) was converted into RL, and a graph was created with the x-axis as DAS and the y-axis as RL (Fig. 6b). DAS at the time root elongation began (\({x}_{1}\)), and when it was concluded (\({x}_{2}\)), were calculated by approximating each line on the graph with Eq. (4).

$$\begin{array}{c}y=f\left(x\right)=min\left\{max\left[\frac{1}{\left({x}_{1}-{x}_{2}\right)}\left(x-{x}_{2}\right), 0\right], ythr\right\}\end{array}$$
(4)

where \(ythr\) is the y value at which changes in RL are stalled (Fig. 6d). Approximation was performed using brute force, and \({x}_{1}\) and \({x}_{2}\), for which MSE was smallest.

Calculation of root number and root elongation rate

Root length was calculated from the root length at 27 DAS and RL. Given that the root length on the last day is \(Length\) and the root length at n DAS is \({Length}_{n}\), this value was calculated using Eq. (5).

$$\begin{array}{c}{Length}_{n}=Length\times \left(1-\mathrm{RL}\right)\end{array}$$
(5)

If \({Length}_{n}\) is zero or RL is one, the root has not yet initiated at \(n\) DAS. The number of roots for which \({Length}_{n}\) was greater than zero was defined as the root number at \(n\) DAS. Pearson's correlation test was performed using the ‘cor.test’ function of R, version 3.6.3.

Availability of data and materials

The datasets used in this study are available at the GitHub repository (https://github.com/st707311g/) and the project homepage (https://rootomics.dna.affrc.go.jp/en/). If not, they are available from the corresponding author on reasonable request.

Abbreviations

3D:

Three-dimensional

CNN:

Convolutional neural network

CT:

Computed tomography

DAS:

Days after sowing

ICP:

Iterative closest point

MRI:

Magnetic resonance imaging

RL:

Relative length

RSA:

Root system architecture

SBI:

Soil block identifier

References

  1. Lynch J. Root architecture and plant productivity. Plant Physiol. 1995;109:7–13.

    Article  CAS  Google Scholar 

  2. Gowariker V, Krishnamurthy VN, Gowariker S, Dhanorkar M, Paranjape K. The fertilizer encyclopedia. Hoboken: John Wiley & Sons; 2009.

    Google Scholar 

  3. Uga Y, Sugimoto K, Ogawa S, Rane J, Ishitani M, Hara N, et al. Control of root system architecture by deeper rooting 1 increases rice yield under drought conditions. Nat Genet. 2013;45:1097–102.

    Article  CAS  Google Scholar 

  4. Kitomi Y, Hanzawa E, Kuya N, Inoue H, Hara N, Kawai S, et al. Root angle modifications by the DRO1 homolog improve rice yields in saline paddy fields. Proc Natl Acad Sci U S A. 2020;117:21242–50.

    Article  CAS  Google Scholar 

  5. Oo AZ, Tsujimoto Y, Mukai M, Nishigaki T, Takai T, Uga Y. Synergy between a shallow root system with a DRO1 homologue and localized P application improves P uptake of lowland rice. Sci Rep. 2021;11:9484.

    Article  CAS  Google Scholar 

  6. Khare D, Mitsuda N, Lee S, Song WY, Hwang D, Ohme-Takagi M, et al. Root avoidance of toxic metals requires the GeBP-LIKE 4 transcription factor in Arabidopsis thaliana. New Phytol. 2017;213:1257–73.

    Article  CAS  Google Scholar 

  7. Uga Y. Challenges to design-oriented breeding of root system architecture adapted to climate change. Breed Sci. 2021;71:3–12.

    Article  CAS  Google Scholar 

  8. Trachsel S, Kaeppler SM, Brown KM, Lynch JP. Shovelomics: high throughput phenotyping of maize (Zea mays L.) root architecture in the field. Plant Soil. 2011;341:75–87.

    Article  CAS  Google Scholar 

  9. Teramoto S, Kitomi Y, Nishijima R, Takayasu S, Maruyama N, Uga Y. Backhoe-assisted monolith method for plant root phenotyping under upland conditions. Breed Sci. 2019;69:508–13.

    Article  Google Scholar 

  10. Cheng W, Coleman DC, Box JE. Measuring root turnover using the minirhizotron technique. Agric Ecosyst Environ. 1991;34:261–7.

    Article  Google Scholar 

  11. Satomura T, Fukuzawa K, Horikoshi T. Considerations in the study of tree fine-root turnover with minirhizotrons. Plant Root. 2007;1:34–45.

    Article  Google Scholar 

  12. Eshel A, Beeckman T. Plant roots: the hidden half. Florida: CRC Press; 2013.

    Book  Google Scholar 

  13. Huck MG, Taylor HM. The rhizotron as a tool for root research. Adv Agron. 1982;35:1–35.

    Article  Google Scholar 

  14. Neufeld HS, Durall DM, Rich PM, Tingey DT. A rootbox for quantitative observations on intact entire root systems. Plant Soil. 1989;117:295–8.

    Article  Google Scholar 

  15. Wang T, Rostamza M, Song Z, Wang L, McNickle G, Iyer-Pascuzzi AS, et al. SegRoot: a high throughput segmentation method for root image analysis. Comput Electron Agric. 2019;162:845–54.

    Article  Google Scholar 

  16. Xu W, Yu G, Zare A, Zurweller B, Rowland DL, Reyes-Cabrera J, et al. Overcoming small minirhizotron datasets using transfer learning. Comput Electron Agric. 2020;175:105466.

    Article  Google Scholar 

  17. Smith AG, Petersen J, Selvan R, Rasmussen CR. Segmentation of roots in soil with U-Net. Plant Methods. 2020;16:13.

    Article  Google Scholar 

  18. de Dorlodot S, Forster B, Pagès L, Price A, Tuberosa R, Draye X. Root system architecture: opportunities and constraints for genetic improvement of crops. Trends Plant Sci. 2007;12:474–81.

    Article  Google Scholar 

  19. Teramoto S, Takayasu S, Kitomi Y, Arai-Sanoh Y, Tanabata T, Uga Y. High-throughput three-dimensional visualization of root system architecture of rice using X-ray computed tomography. Plant Methods. 2020;16:66.

    Article  Google Scholar 

  20. Gao W, Schlüter S, Blaser SRGA, Shen J, Vetterlein D. A shape-based method for automatic and rapid segmentation of roots in soil from X-ray computed tomography images: Rootine. Plant Soil Plant and Soil. 2019;441:643–55.

    Article  CAS  Google Scholar 

  21. Mairhofer S, Zappala S, Tracy SR, Sturrock C, Bennett M, Mooney SJ, et al. RooTrak: automated recovery of three-dimensional plant root architecture in soil from X-Ray microcomputed tomography images using visual tracking. Plant Physiol. 2012;158:561–9.

    Article  CAS  Google Scholar 

  22. Iijima M, Oribe Y, Horibe Y, Kono Y. Time lapse analysis of root elongation rates of rice and sorghum during the day and night. Ann Bot. 1998;81:603–7.

    Article  Google Scholar 

  23. Yazdanbakhsh N, Fisahn J. Analysis of Arabidopsis thaliana root growth kinetics with high temporal and spatial resolution. Ann Bot. 2010;105:783–91.

    Article  Google Scholar 

  24. Yazdanbakhsh N, Fisahn J. Stable diurnal growth rhythms modulate root elongation of Arabidopsis thaliana. Plant Root. 2011;5:17–23.

    Article  Google Scholar 

  25. Fisahn J, Yazdanbakhsh N, Klingele E, Barlow P. Arabidopsis thaliana root growth kinetics and lunisolar tidal acceleration. New Phytol. 2012;195:346–55.

    Article  Google Scholar 

  26. Gao W, Blaser SRGA, Schlüter S, Shen J, Vetterlein D. Effect of localised phosphorus application on root growth and soil nutrient dynamics in situ—comparison of maize (Zea mays) and faba bean (Vicia faba) at the seedling stage. Plant Soil Plant and Soil. 2019;441:469–83.

    Article  CAS  Google Scholar 

  27. Rellán-Álvarez R, Lobet G, Lindner H, Pradier PL, Sebastian J, Yee MC, et al. GLO-Roots: an imaging platform enabling multidimensional characterization of soil-grown root systems. Elife. 2015;4:1–26.

    Article  Google Scholar 

  28. Bontpart T, Concha C, Giuffrida MV, Robertson I, Admkie K, Degefu T, et al. Affordable and robust phenotyping framework to analyse root system architecture of soil-grown plants. Plant J. 2020;103:2330–43.

    Article  CAS  Google Scholar 

  29. Parker JC, Amos DF, Kaster DL. An evaluation of several methods of estimating soil volume change. Soil Sci Soc Am J. 1977;41:1059–64.

    Article  Google Scholar 

  30. Giraldez JV, Sposito G, Delgado C. A general soil volume change equation: I. The two-parameter model. Soil Sci Soc Am J. 1983;47:419–22.

    Article  Google Scholar 

  31. Teramoto S, Tanabata T, Uga Y. RSAtrace3D: robust vectorization software for measuring monocot root system architecture. BMC Plant Biol. 2021;21:398.

    Article  Google Scholar 

  32. Lobet G, Pagès L, Draye X. A novel image-analysis toolbox enabling quantitative analysis of root system architecture. Plant Physiol. 2011;157:29–39.

    Article  CAS  Google Scholar 

  33. Möller B, Chen H, Schmidt T, Zieschank A, Patzak R, Türke M, et al. rhizoTrak: a flexible open source Fiji plugin for user-friendly manual annotation of time-series images from minirhizotrons. Plant Soil. 2019;444:519–34.

    Article  Google Scholar 

  34. Bull DR. Communicating pictures: a course in image and video coding. Massachusetts: Academic Press; 2014.

    Book  Google Scholar 

  35. Jia D, Song C, Cheng C, Shen S, Ning L, Hui C. A novel deep learning-based spatiotemporal fusion method for combining satellite images with different resolutions using a two-stream convolutional neural network. Remote Sens. 2020;12:698.

    Article  Google Scholar 

  36. Bovik AC. The essential guide to video processing. Massachusetts: Academic Press; 2009.

    Google Scholar 

  37. Zitova B, Flusser J. Image registration methods: a survey. Image Vis Comput Elsevier. 2003;21:977–1000.

    Article  Google Scholar 

  38. Sharp GC, Lee SW, Wehe DK. ICP registration using invariant features. IEEE Trans Pattern Anal Mach Intell. 2002;24:90–102.

    Article  Google Scholar 

  39. Yamazaki K, Ohmori Y, Fujiwara T. A positive tropism of rice roots toward a nutrient source. Plant Cell Physiol. 2020;61:546–53.

    Article  CAS  Google Scholar 

  40. Hodge A. The plastic plant: root responses to heterogeneous supplies of nutrients. New Phytol. 2004;162:9–24.

    Article  Google Scholar 

  41. Flavel RJ, Guppy CN, Tighe MK, Watt M, Young IM. Quantifying the response of wheat (Triticum aestivum L) root system architecture to phosphorus in an Oxisol. Plant Soil. 2014;385:303–10.

    Article  CAS  Google Scholar 

  42. Bardhan K, York LM, Hasanuzzaman M, Parekh V, Jena S, Pandya MN. Can smart nutrient applications optimize the plant’s hidden half to improve drought resistance? Physiol Plant. 2021;172:1007–15.

    Article  CAS  Google Scholar 

  43. Vanhees DJ, Schneider HM, Sidhu JS, Loades KW, Bengough AG, Bennett MJ, et al. Soil penetration by maize roots is negatively related to ethylene-induced thickening. Plant Cell Environ. 2022;45:789–804.

    Article  CAS  Google Scholar 

  44. Mcneill A, Kolesik P. X-ray CT investigations of intact soil cores with and without living crop roots. SuperSoil 2004 3rd Aust New Zel Soils Conf. 2004.

  45. Correa J, Postma JA, Watt M, Wojciechowski T. Soil compaction and the architectural plasticity of root systems. J Exp Bot. 2019;70:6019–34.

    Article  CAS  Google Scholar 

  46. Numajiri Y, Yoshino K, Teramoto S, Hayashi A, Nishijima R, Tanaka T, et al. iPOTs: Internet of Things-based pot system controlling optional treatment of soil water condition for plant phenotyping under drought stress. Plant J. 2021;107:1569–80.

    Article  CAS  Google Scholar 

  47. Luo H, Xu H, Chu C, He F, Fang S. High temperature can change root system architecture and intensify root interactions of plant seedlings. Front Plant Sci. 2020;11:160.

    Article  Google Scholar 

  48. Magistri F, Chebrolu N, Stachniss C. Segmentation-based 4D registration of plants point clouds for phenotyping. IEEE Int Conf Intell Robot Syst. 2020. 2433–9.

  49. Paproki A, Sirault X, Berry S, Furbank R, Fripp J. A novel mesh processing based technique for 3D plant analysis. BMC Plant Biol. 2012;12:63.

    Article  Google Scholar 

  50. Teramoto S, Uga Y. Improving the efficiency of plant root system phenotyping through digitization and automation. Breed Sci. 2022;72:48–55.

    Article  Google Scholar 

  51. Van Rossum G, Drake FL. Python 3 reference manual. Scotts Valley: CreateSpace; 2009.

    Google Scholar 

  52. Harris CR, Millman KJ, van der Walt SJ, Gommers R, Virtanen P, Cournapeau D, et al. Array programming with NumPy. Nature. 2020;585:357–62.

    Article  CAS  Google Scholar 

  53. der Walt S, Schönberger JL, Nunez-Iglesias J, Boulogne F, Warner JD, Yager N, et al. scikit-image: image processing in Python. PeerJ. 2014;2:e453.

    Article  Google Scholar 

  54. Zhou Q-Y, Park J, Koltun V. Open3D: A modern library for 3D data processing. arXiv. 2018. https://doi.org/10.48550/arXiv.1801.09847.

    Article  Google Scholar 

Download references

Acknowledgements

Not applicable.

Funding

This work was supported by the JST CREST Grant Number JPMJCR17O1, Japan.

Author information

Authors and Affiliations

  1. Institute of Crop Sciences, National Agriculture & Food Research Organization, Tsukuba, Ibaraki, 305-8602, Japan

    Shota Teramoto & Yusaku Uga

Authors
  1. Shota Teramoto
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Yusaku Uga
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

ST designed the study, developed, and evaluated the software, and wrote the manuscript. YU coordinated the project, designed the study, and revised the manuscript. Both authors have read and approved the final manuscript.

Corresponding author

Correspondence to Yusaku Uga.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Movie S1.

Animations of 21 top-view projections of rice RSA from 7 to 27 DAS with and without SBI-ICP registration.

Additional file 2: Movie S2.

Animations of backward prediction results from 7 to 27 DAS.

Additional file 3: Movie S3.

Animations of 21 top-view projections of rice RSA computed using RSA vectors generated by backward prediction. The color scale indicates the depth.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Teramoto, S., Uga, Y. Four-dimensional measurement of root system development using time-series three-dimensional volumetric data analysis by backward prediction. Plant Methods 18, 133 (2022). https://doi.org/10.1186/s13007-022-00968-x

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13007-022-00968-x

Keywords

Plant Methods

ISSN: 1746-4811