Automated high-resolution structure analysis of plant root with a morphological image filtering algorithm

Background : Research on rice (Oryza sativa) roots requires the automatic analysis of root structure using image processing. It is challenging for a digital filter to identify the roots from the obscured and cluttered background, and to separate the primary roots from lateral roots. The original Frangi filter (FF), presented by Alejandro F. Frangi in 1998, is a low-pass filter dedicated to blood vessel image enhancement. Considering the similarity between vessels and roots, the FF is applied to identify the roots. However, the original FF only enhances the tube-like primary roots but erases the lateral roots. Hence, a new method is developed to meet the demands by simultaneously maintaining the primary and lateral morphological structure of roots. Result : In this work, a crucial part of the FF, Gaussian filtering, is redesigned to discriminate against the primary and lateral roots in a 2-D image. Inspired by the structure-awareness of the FF, an Improved Frangi Filtering Algorithm (IFFA) designed for plant roots is proposed. First, Multilevel image thresholding, connected-components labeling and width correction are used to optimize the resultant binary image. Then, to enhance the local structure, the truncated Gaussian kernel is modified resulting in more discernible lateral roots. Compared to the original FF and the Automatic Root Image Analysis (ARIA), a commercial software, IFFA is a faster and more accurate algorithm achieving an identification accuracy of 97.48%. Conclusion : IFFA is an effective 2-D filtering approach to enhance the roots of rice (Oryza sativa) for segmentation and further biological research. IFFA is faster than ARIA and the original FF, and IFFA’s accuracy outperforms its counterparts as per the Intersection Over Union (IOU) and Dice Similarity Coefficient (DSC) criteria.


Background
In the field of biology and genetic breeding, plant phenotypic refers to the external traits of organisms determined by the genotype and environment of the plant, such as shape, structure, size, and color. Researches on the aboveground parts of the plants have made great progress. While the study of underground parts, especially the roots of plants, encounters difficulties in observing and analyzing a large number of plants.
Root system structure shows a great importance in studying interspecific interactions and genetic improvement of crops [1,2]. To acquire a high-quality root-system image and analyze the root structure, much research had been done. In Gernot Bodner's research, a data processing pipeline was developed for automatic root segmentation [3]. Pfeifer developed a comprehensive segmentation method suitable for comparatively large columns sampled in situ [4]. Christian and Anjali both established a platform to image the roots automatically and non-destructively [5,6]. Morphological filters possess effective noise suppression with reduced geometrical feature blurring [7]. The main features were extracted using a morphological filter algorithm to reconstruct fine roots with high fidelity [8]. WinRHIZO, EZ-RHIZO, and ROOTEDGE are software designed for the measurement of root system structure [9][10][11]. In addition to 2D image processing and X-ray micro-tomography [12], Suqin focused on 3D reconstruction and dynamic modeling [13]. Compared to the former methods, a costeffective and efficient way to automatically extract precise root-system structure from high-resolution root images is presented in this paper. Since rice is the main food source, this work would take rice roots as a research object.
Usually, plant roots need to be observed in a nondestructive manner. Dragging roots out leads to an unrecoverable injure to the structure of the root system. To solve these problems, Randy and Anjali used transparent media to observe the root structure. [14,6] In this work, transparent plastic bags with the nutrient solution are used to cultivate the plants. CCD cameras are used to gather pictures of roots. A typical image is shown in Figure 1.

Figure 1 A sample picture of roots
Difficulties arise from the fact that the root color is similar to the background color, besides, the main roots and lateral roots overlap each other. Furthermore, lateral roots are thin and make it harder to differentiate them. Frangi proposed an effective filter to enhance vessel structures with the eventual goal of vessel segmentation [15]. The FF smooths the image using Gaussian filtering, find the Hessian matrix, and the tubular object enhanced. Gaussian filtering is essentially low-pass filtering. If the window width is the scale factor σ, the filter output is the largest and the effect is best when it matches the actual width of the root. As shown in Figure 2, the original image has a lot of noise, which is close to the lateral roots on the scale. In the background, some disturbing indentations are similar to the color of the roots. When using the original FF, if the noise filtering is to be ensured, a large number of side roots are filtered out ( Figure 3). When the side roots are to be retained, there will be more noise and adhesion and blurring between different roots ( Figure 4). In general, the output is not satisfactory or suitable for research purposes. Some improvements should be made before applying the FF on roots enhancement.

Method
Automated platform for image acquiring. A multi-angle camera layout is adopted to improve the efficiency of capturing images from the bag seedlings. To avoid the interference of ambient light, the imaging equipment is mounted in a dark-room. As shown in Figure 5, four cameras are installed in different angles to capture four root images each time. The dark-room consists of an aluminum-profile frame wrapped with metal sheets to block the ambient light. This procedure is controlled by a PLC controller, which receives the signals generated by the sensors mounted on the equipment. We used four CCD industrial cameras with a resolution of 4608×3288 pixels to ensure the detection of the lateral root from different perspectives. (1) where ( * ) stands for the convolution operation, H stands for the Hessian matrix. Before calculating the Hessian matrix, Gaussian smoothing is applied to remove noise by reducing the deviation of the second-order partial derivative. Gaussian smoothing is described as follows, where σ is standard deviation.
Since H is a symmetric matrix, its eigenvalues could be used to construct the enhancement filter. The two eigenvalues, 1 , 2 , can be calculated by When the scale factor σ matches the actual width of the vessel, the output of the filter is maximized. So as a spatial scale factor, iterations can obtain different scales of output.
The half-width of the rectangular window of the local characteristic analysis is generally 3σ. When the diameter of the root is smaller than the width and height of the corresponding rectangular window of the current scale, the eigenvalues of the Hessian matrix of the tubular root satisfy , S with 1 and 2 as The responsivity of p in vessel areas is where the background of an image would be cut off, is a parameter controlling the difference between thread objects and block objects. ∈ [10 −5 , 10 −6 ] and ∈ [3,6] are parameters controlling the smoothness of thread objects.
The half-width of the rectangular window of each pixel is generally taken to be 3 . When the window corresponding to different scales matches the diameter of the root, the maximum response occurs. Changing the value of , and record the maximum response of each pixel and finally reconstruct the image, the trace of roots would be enhanced.
Target improving aspects and processing procedure. Drawn from engineering experience, a successful image filtering method should meet the following requirements, removing noise, reserving contents, fast process, and strong robustness. To achieve these goals, the processing procedure is designed as shown in Figure 6. In the following sections, each of the tasks will be detailed.
Truncated Gaussian kernel. In the original FF procedure, Gaussian filtering and calculation of the Hessian matrix is completed in the same step, the width of the Gaussian kernel is given as, h = 2 ⊔ (3σ) + 1 (15) where σ is the standard deviation of the Gaussian kernel. A large kernel enhances the ability to remove noise, however, the blurring effect may cause the loss of meaningful high-frequency components. In this case, when σ ranges from 1 to 10, h ranges from 7 to 70. The Gaussian kernel is significantly too large to reserve the lateral roots. Thus, the truncated Gaussian kernel is proposed. Its width is changed to h = ⊔ (2σ) (16) From Equation 16, we can understand the physical meanings of the filter kernel, i.e. σ responses to the high-frequency components and limits the thread-shape lowfrequency components. Designed in this way, the filter is sensitive to the object morphology. In its results, the outlier noises can be removed while the lateral roots are reserved. Figure 7 shows the comparison of the original FF and the FF with truncated Gaussian kernel. Adjustment of variance. In the original FF algorithm, variance σ varies from 1 to 10, be a step of 0.5. It results in 19 scans for processing one image. If the FF is used to filter the image of high quality or large quantity, the time consumed is unacceptable.
Finding a way to minimize the range of variance would benefit the whole procedure. It is also a self-adaptive way to handle images of different scales. The intersection-over-union (IOU) method is introduced. It is a measurement of the similarity between two images, widely used in machine vision. During one iteration of the FF, if the IOU of the original grayscale image and current output is more than 0.9, the iteration will terminate. Variance σ ranges from 1 to 5 is obtained and saving over half of the processing time.
Multilevel image thresholding (Otsu). The Otsu method is an adaptive threshold segmentation method. A threshold t is defined to divide a grayscale image into an object and background [16]. The ratio of all pixels with a grayscale of i to the total number of pixels defined as . The possibility of object O and background B could be calculated as follows The average grayscale of and is as follows The average grayscale of the whole image is as follows The variance is (t) = ω ( − ) 2 + ( − ) 2 (20) By changing the value of t to maximize d, the best threshold for segmentation can be obtained.
The multi-threshold Otsu method is an alternative formulation for Otsu's method. Regardless of the number of classes being considered during the thresholding process, the sum of the cumulative probability functions of M classes equals one, and the mean of the image brightness is equal to the sum of the means of M classes weighted by their cumulative probabilities [17]. In practice, the MATLAB function "multithresh" is used to calculate the threshold needed. The number of classes could be controlled in the function, usually 3 or 4. Figure 7 shows some different output images of the Multithreshold Otsu method. "3-C-1" means divide the point set by "three" threshold and remove "one" part nearest to black (0/255). Taking image details, sharpness and processing time into consideration, "3-C-1" is chosen to be applied to the whole dataset.
Label connected components. This step can be regarded as the final step to remove noise and unnecessary parts in images. A watershed algorithm is used to color the connecting component, making the image clearer for observation. Figure 8 shows the original and final output to illustrate the effect of IFFA.
Width correction. Width correction is proposed to meet different requirements. By using IFFA, the width of roots would be magnified. If the width of roots is most important for users, the real form should be reserved. Set the IOU is limited in step c (adjustment of variance) to 0.7, the true width of roots would be reserved. Figure 9 shows an example of width correction.
Label connected components. This step can be regarded as the final step to remove noise and unnecessary parts in images. A watershed algorithm is used to color the connecting component, making the image clearer for observation. Figure 8 shows the original and final output to illustrate the effect of IFFA.
(a) 3-C-1 (b) 3-C-2 (c) 4-C-1 (d) 4-C-2 Figure 8 Output images of the Multi-threshold Otsu method "3-C1" means divide the point set by "three" threshold and remove "one" part nearest to black (0/255) Width correction. Width correction is proposed to meet different requirements. By using IFFA, the width of roots would be magnified. If the width of roots is most important for users, the real form should be reserved. Set the IOU is limited in step c (adjustment of variance) to 0.7, the true width of roots would be reserved. Figure 9 shows an example of width correction.

Analysis of binary image.
As long as the binary image is extracted, several indexes of the root system could be calculated. The width-depth ratio (WDR) is defined in the algorithm as follows, WDR = ℎ ℎ ℎ (21) Figure 9 Final output of the IFFA Furthermore, WDR is an initial indicator to study the drought resistance of plants.
Assuming that every pixel stands for identical mass, the center of mass (CoM) is defined as follows, CoM = min(∑ ⃗⃗⃗⃗ − 0 ) (22) By changing the value of σ, the number of main roots and lateral roots can be calculated. Thus, the ratio of lateral roots to the main roots and the density of lateral roots are obtained. The other indexes, such as primary root length, total root length, can also be determined.

Results
Comparison with the original Frangi filter. A dataset label by humans is introduced. Let IMG stands for the output image and MASK stands for the labeled unfiltered image. Intersection over Union (IOU) is defined as follows, The higher the IOU value is, the more accurate the object is reserved. Pixel accuracy (PA) is a standard measuring the ratio of correct pixels to all pixels. The Dice Similarity Coefficient (DSC) is defined as follows, (a) original gray image (b) σ up to 5, after Otsu (c) width correction, σ up to 2, after Otsu Figure 10 Example of width correction Similar to the IOU standard, DSC indicates the overlap measurement between the ground truth root regions and the segmentation results [18]. Table 1 shows the average processing time and IOU of FF and IFFA to the labeled image and Figure 11 illustrates the FF and IFFA results for the same labeled raw image.  Comparison with ARIA. ARIA as commercial software is extensively used in analyzing data of roots [19]. However, after image filtering, some main roots and most of the lateral roots would be removed. The starting point of roots needs to be manually labeled before the calculation and it could handle an RGB image. Figure 12 shows the interface of ARIA with optimal filtering parameter settings. Table 2 shows the time comparison between ARIA and IFFA.

Figure 12
Parameter settings and filtered result of ARIA

Discussion
We have shown that the IFFA method can achieve outstanding performance in the enhancement of images of rice roots which grow in transparent plastic bags. Replacing the Gaussian filter with a truncated Gaussian filter, the lateral roots will be protected from being erased. Applying Otsu multi-threshold segmentation and width correction to the images, we enhance the root structure which is pure and sharp. The results show an effective way to extract the rice root structure.
Compared to the original FF, IFFA is faster in operation. The IOU standard helps to provide a satisfactory enhancement. The Ostu method removes all the gray noise pixels, and we could choose to retain more lateral roots (and more noise) or less. ARIA seems to remove nearly all the lateral roots and cause distortion, so it could not be used to make valuable comparisons.
We observe that the pixel-level accuracy of IFFA is over 97%. Nevertheless, one in ten lateral roots might be removed during the filtering process, and overlapping at the primary root scale may cause distortions. However, these imperfections can be improved by accumulating more sample cases to determine an optimal filter parameter, such as the truncated width.

Conclusions
A high-performance and well-defined automatic root image processing tool is important for the agronomic study and plant physiological investigation. Based on the original Frangi filter, a 2-D image processing algorithm for plant root structure analysis is improved, which effectively filters out the in-bag root background, enhances the primary roots while keeping the lateral root almost intact. The proposed Improved Frangi Filtering Algorithm outperforms other root image enhancement algorithms in both pixel resolution and image signal integrity, and it can be used for automatic rice root phenotyping.

Declarations
We declare no conflict of interest.