- © 2010 American Society of Plant Biologists
Abstract
The ability to nondestructively image and automatically phenotype complex root systems, like those of rice (Oryza sativa), is fundamental to identifying genes underlying root system architecture (RSA). Although root systems are central to plant fitness, identifying genes responsible for RSA remains an underexplored opportunity for crop improvement. Here we describe a nondestructive imaging and analysis system for automated phenotyping and trait ranking of RSA. Using this system, we image rice roots from 12 genotypes. We automatically estimate RSA traits previously identified as important to plant function. In addition, we expand the suite of features examined for RSA to include traits that more comprehensively describe monocot RSA but that are difficult to measure with traditional methods. Using 16 automatically acquired phenotypic traits for 2,297 images from 118 individuals, we observe (1) wide variation in phenotypes among the genotypes surveyed; and (2) greater intergenotype variance of RSA features than variance within a genotype. RSA trait values are integrated into a computational pipeline that utilizes supervised learning methods to determine which traits best separate two genotypes, and then ranks the traits according to their contribution to each pairwise comparison. This trait-ranking step identifies candidate traits for subsequent quantitative trait loci analysis and demonstrates that depth and average radius are key contributors to differences in rice RSA within our set of genotypes. Our results suggest a strong genetic component underlying rice RSA. This work enables the automatic phenotyping of RSA of individuals within mapping populations, providing an integrative framework for quantitative trait loci analysis of RSA.
Root systems are complex three-dimensional (3D) structures that provide functions central to plant fitness, such as water and nutrient acquisition. Plant health and survival are dependent on root system architecture (RSA), the spatial configuration of different types and ages of roots emerging from a single plant (Lynch, 1995). RSA varies dramatically within and between species, allowing for soil exploration in diverse conditions (Fitter, 2002). Age is also a major factor in RSA; young plants have less complex root systems, but as plants mature their root systems become correspondingly more complicated. Modification of RSA could contribute to improvements of desirable agronomic traits such as yield, drought tolerance, and resistance to nutrient deficiencies (Tuberosa et al., 2002; Beebe et al., 2006; Steele et al., 2006, 2007). However, largely due to the difficulty of noninvasively observing root architecture, knowledge regarding the genes underlying components of RSA is limited. Most information concerning RSA derives from the relatively simple dicot root system of laboratory-grown young Arabidopsis (Arabidopsis thaliana) seedlings (Armengaud et al., 2009; French et al., 2009), and comparatively little is known about the substantially more complicated RSA of monocots like rice (Oryza sativa) or maize (Zea mays; Hochholdinger and Tuberosa, 2009).
Traditional methods of observing RSA, such as excavation or washed soil cores, destroy the topology of the root system, but in recent years several newer techniques have been used to nondestructively image root systems. X-ray computed tomography (Gregory et al., 2003; Perret et al., 2007) has the advantage of observing roots grown in soil, but can be limited by cost, imaging time, resolution, or container size, depending on the type of x-rays and computed tomography scanner used. NMR (van der Weerd et al., 2001; Jahnke et al., 2009) and magnetic resonance imaging (Asseng et al., 2000; van As, 2007) are both costly and the equipment is not easily accessible. NMR is also very sensitive to the type of media used for plant growth. Laser scanning of root systems (Fang et al., 2009) provides precise measurements but requires longer imager times and can be expensive. Stereoscopic methods using plants grown in transparent media (Wulfsohn et al., 1999) are limited by the need to manually move the microscope stage to see different parts of the root system and to image roots horizontally. Further, the time to image increases as root system complexity increases. Additional nondestructive methods, such as flatbed scanning of petri plates (Armengaud et al., 2009; French et al., 2009), thin gel chambers (Hargreaves et al., 2009), or slim growth pouches (Hund et al., 2009) work well for simpler root systems, but obscure the more complex 3D RSA of monocot root systems.
Quantitative trait loci (QTL) studies suggest that many traits fundamental to RSA are controlled by numerous small-effect loci (de Dorlodot et al., 2007; Hochholdinger and Tuberosa, 2009). In rice, most QTL studies have relied on root features such as maximum root length, number, or thickness (Kamoshita et al., 2002; Li et al., 2005; MacMillan et al., 2006). However, these features are poor estimators of a root system's spatial distribution. The spatial distribution of roots is a key factor in the plant's ability to uptake water and nutrients, but a comprehensive description of RSA is currently lacking. RSA-associated QTL have yet to be cloned, due, in part, to the difficulty in accurately measuring the components of RSA in mapping populations.
Many different programs are available for the measurement of RSA. Some programs are specialized for use with certain types of images, such as those from a minirhizotron (Zeng et al., 2008; Rootfly [www.ces.clemson.edu/ ∼ stb/rootfly/], RootView [www.mv.helsinki.fi/aphalo/RootView.html], and RooTracker [www.biology.duke.edu/rootracker]). Others are not automated and require significant user input (Le Bot et al., 2010). Additional programs accurately characterize the simpler dicot root system of Arabidopsis (Armengaud et al., 2009; French et al., 2009), but are not as suitable for the intricate monocot root systems of rice and maize.
The identification of genes underlying the RSA of rice is thus a challenging task, as it requires three technological innovations: first, an imaging platform to nondestructively image rice root systems, second, an analysis system that automatically and accurately phenotypes RSA, and third, a statistical method to determine which comprehensive RSA traits differ between genotypes and hence, are likely QTL mapping candidates. We report here methods to address each of the above challenges, and apply these methods to a set of more than 2,000 images from 12 rice genotypes.
RESULTS
Imaging Platform to Nondestructively Observe RSA in Real Time
To address the first challenge, we developed a nondestructive imaging platform that allows us to view the entire root system from any angle for at least 2 weeks. Rice seeds were germinated in 0.3% agar in 2 L ungraduated cylinders (see “Materials and Methods” for cylinder dimensions). We tested several different types of agar, including Phytagel (Sigma), gellan gum (Sigma), granulated agar (DIFCO), and Gelzan CM (Sigma). We chose Gelzan CM as it provided the highest level of optical clarity.
For imaging, cylinders were placed on an automated turntable connected to a laptop PC and digital camera. The turntable rotates through a predetermined angle, stops, and the camera acquires an image. Using the same parameters (lighting, camera focus, distance to the cylinder) for all images, 20 angles around the cylinder were imaged (Fig. 1). Using our gel-based platform, we generated approximately 2,300 images of multiple individuals of 12 different rice genotypes at the 14th d after planting (dap; Fig. 1). For most genotypes, only plants that had germinated within 1 d of each other were used for analysis. One limitation of previous nondestructive methods is the length of time to image crop plants. Using this platform we can image one plant, regardless of size, in approximately 10 min, making this method amenable to imaging mapping populations in the first step of QTL analysis.
Gel-based growth platform. A, 93-11 growing in gel in a 2 L ungraduated cylinder. The cylinder is placed on a lightbox on the imaging turntable and lit from the side. B, Images were acquired by zooming in on the root system. C, Cropped images from multiple angles were used for analysis. Approximately 20 angles (every 18°) from each plant were used for analysis; four are shown here.
Trait Descriptor Development
Dicot root systems such as those found in young Arabidopsis plants can often be efficiently characterized by features such as lateral root number or length of specific zones within a single root (Armengaud et al., 2008; French et al., 2009). However, these relatively simple features are not sufficient to fully characterize the complexity of rice RSA. Therefore, we chose to phenotype traits that captured aspects of the spatial distribution of the root system. We included standard traits that are typically used to describe root systems, such as width (Price et al., 2002; de Dorlodot et al., 2007), total root length (Armengaud et al., 2009), specific root length (SRL; Eissenstat, 1991), and volume (Qu et al., 2008). However, our imaging and automatic phenotyping system allows us to include novel traits that more comprehensively describe the monocot root system.
The standard RSA features we estimated include median and maximum number of roots, average root radius, SRL, total surface area, total root length, depth, volume, length distribution, maximum horizontal width, and width-to-depth ratio. The novel traits are network perimeter, solidity, convex area, network area, and bushiness index. The definitions of each of these traits are as follows, where the abbreviation in parentheses denotes labels used in figures: (1) Median number of roots (MedR) is the result of a vertical line sweep in which the number of roots that crosses a horizontal line as the line moves vertically down the image was estimated, then the median of all values for the entire network was calculated. (2) Maximum number of roots (MaxR) is calculated by sorting the number of roots crossing a horizontal line in the vertical line sweep from smallest to largest, the maximum number is considered to be the 84th percentile value (1 sd). MedR and MaxR computed using the line sweep approach give a good estimation when roots are well separated, but it can underestimate the number of roots when roots occlude each other. (3) Average root radius (Av. Radius) is the mean value of the root radius estimation computed for all pixels of the medial axis of the entire root system. (4) The total root length (Total Length) is the number of pixels in the medial axis of the root system. Note that estimating the length of the medial axis using Euclidean distance between pixels is a more rigorous calculation of the total root length, however our estimation computed as the number of medial axis pixels is faster and is highly correlated with the length of the medial axis. (5) Volume is estimated as the sum of cross-sectional areas for all pixels of the medial axis of the root system. More precisely, if the medial axis of a root system consists of N pixels, and ri is a radius of a cross section corresponding to the pixel i in the medial axis, then:
(6) SRL is the total root length divided by root system volume. (7) Total surface area (SA) is estimated as a sum of circumferences of cross sections corresponding to all pixels in the medial axis of the root system:
(8) Perimeter (Perim) is the total number of network pixels connected to a background pixel (using a 8-nearest neighbor neighborhood). (9) Convex area (ConvA) is the area of the convex hull that encompasses the image. Here, the convex hull is the smallest convex set of pixels that contains all other pixels in the root system. (10) Network area is the number of network pixels in the root system. (11) Solidity is the fraction equal to the network area divided by the convex area. (12) Bushiness is the ratio of the maximum to the median number of roots. (13) Maximum width (Max Width) is a maximum horizontal width of the whole RSA. (14) Depth is a maximum vertical distance reached by the root system. (15) Width-to-depth ratio (W/D) is a ratio of the maximum width to depth. (16) Length distribution (Length Distr) is a ratio of the total root length in the upper one-third of the root depth to the total root length located in the lower two-thirds of the root depth. Note that many of these traits depend on calculation of the medial axis and the radius of a root. The medial axis is defined as the points within the root system that have more than one equidistant point on the root system's boundary. In practice, the medial axis is a series of curves comprised of the union of all center lines of the branching root system. The radius of a root at a medial pixel is estimated as a distance from this pixel to the closest boundary pixel. We used bwmorph and bwdist Matlab functions to compute the medial axis and to estimate the distance from the medial axis to the closest boundary points. Note that occlusions or crossings between roots can distort the medial axis when analyzing individual two-dimensional (2D) images, an issue that will be resolved, in part, by transitioning to 3D analysis in future work.
High-Throughput Analysis of Root Systems from 12 Rice Varieties Using an Automatic Image Analysis Pipeline
Next, we created an automatic phenotyping system to accompany the gel-based imaging platform. This analysis system is comprised of a pipeline that preprocesses each image, calculates root features for each image, and then combines all of the phenotyping information into a comprehensive trait-ranking step (Supplemental Fig. S1). In the first step, the original images were cropped and converted to binary images using an adaptive thresholding method that accounts for variations in background image intensity. Next, an automatic phenotyping system calculated the 16 different RSA traits discussed above for each 2D image.
We tested our pipeline and trait analysis using 12 genotypes of rice (93-11, Caipo, Basmati 217, Teqing, Moroberekan, Jefferson, Nipponbare, Lemont, IR64, Bala, Azucena, and Carolina Gold). Genotypes were chosen for their parentage in available QTL mapping populations and their diversity within rice. RSA traits for most of these varieties are not well known. Using the above features, we automatically calculated RSA features of these genotypes at the 14th dap. We found a wide diversity of root system structures within the 12 genotypes (Fig. 2; Supplemental Table S1). Some root systems have characteristic shapes, for example Caipo and Azucena individuals displayed a long root system with a small width-to-depth ratio, whereas Basmati 217 individuals had a dense system of shorter roots. We observed a trend such that indica varieties tended to have more shallow root systems with wider RSA (higher width-to-depth ratio) compared to the longer root systems with lower width-to-depth ratios of many japonica varieties. Additionally, RSA was consistent with conditions in which a genotype is typically grown. For example, the longer roots of Caipo and Moroberekan are well suited for the upland soils of Brazil and West Africa, respectively, while the shallower RSA of IR64 and Teqing are sufficient for irrigated conditions. Most genotypes exhibited some degree of branching, but not in a way that could be easily categorized in a qualitative fashion.
Root architecture variation is higher among genotypes than within a genotype. Each set of three images represents a different rice genotype used in the analysis. Within each set are three individual plants of that genotype.
Our automatic phenotyping also illustrated the remarkable similarity of RSA within a genotype (Fig. 2), despite the highly plastic nature of RSA traits (Malamy, 2005). We examined the mean RSA traits of all 12 genotypes 14 dap (Supplemental Table S1). Approximately 20 images, taken every 18 degrees, were used from each individual plant. For each image, we calculated the 16 RSA features described above. Statistical analyses of the difference between mean traits of distinct genotypes demonstrated that RSA varied significantly less within a genotype than between genotypes (Figs. 2 and 3; see “Materials and Methods”). We chose to examine pair-wise comparisons of all 12 genotypes, as representations of the types of root systems within rice. Of the 66 possible pair-wise comparisons between 12 genotypes, 64 had at least one trait that varied significantly between a pair of genotypes even when accounting for the effect of performing multiple comparisons (see Supplemental Table S2 and “Materials and Methods”). Notably, this was true even when analyzing plants that had germinated at slightly different times. These results suggest a strong genetic component underlying rice varietal mean RSA traits.
SVM analysis identifies distinguishing traits in root systems. A to C, SVM comparison of Basmati 217 and Jefferson. D to F, SVM comparison of Lemont and Nipponbare. A and D, Images from eight individuals each of Basmati 217 and Jefferson (A) and Lemont and Nipponbare (D). B and E, Feature maps of root architecture. Each row represents a unique image of a root, approximately 20 images per individual. Columns represent features whose values are normalized to have zero mean and unit sd. Images are ordered according to their distance from the separating hyperplane. The black dotted line indicates the classification boundary between the two genotypes. The bar on the left indicates the genotype. Misclassified images are seen on the left bar as blue lines above the dotted line or red lines below. C and F, Normalized mean feature values show the intra-and intervarietal variation for the 12 varieties examined. The points represent average normalized values, and the small lines represent ses. The two vertical lines represent a visual aid to follow the difference between two specific varieties.
Trait Parameter Validation
To validate our trait measurements, we created a wire model of a root system and used our pipeline for imaging and analysis. We then compared these data to manual measurements of the wire model not in the gel. Because of strong reflection from the model on the glass when lit from the sides, we used a UV lightbox positioned behind the cylinder to light the model in the gel (Supplemental Fig. S2). We chose four traits that could reliably be measured by hand to compare to the imaging and computational analysis: Depth, Max_Width, Av. Radius, and Total Length. Hand measurements differed from 0.67% to 6.5% compared to the average measurements computed for 20 images, indicating the expected level of accuracy of our imaging and computational pipeline in quantifying RSA (Supplemental Table S3). Previous reports examining the accuracy of a nondestructive technique have shown differences within 10% for root length (Gregory et al., 2003), 13% to 22% for average root segment length and number of segments (Perret et al., 2007), and 21% to 41% for root length per volume (Heeraman et al., 1997). Since many of the other traits we describe, such as surface area, volume, and SRL, involve the computation of the four traits validated here, the other traits in our pipeline are also likely to be reliable.
RSA Trait Ranking
Our ultimate goal is to identify genes that underlie RSA. To this end, we wanted a simple method to identify which set of traits best distinguishes RSA between two genotypes, for example, between the parents of a mapping population. We can then use these ranked traits as a prioritized list of QTL mapping candidates. We chose to rank traits according to the degree to which they can be used to separate a pair of genotypes. We ranked traits using a support vector machine (SVM) analysis based on the set of 16 calculated features for each genotype. We tried several different statistical methods, including principal component analysis, but none of these provided the accuracy of SVM. SVM is a supervised learning method for pattern recognition given labeled data (Vapnik, 1998). Several different types of SVMs are available. We chose the standard approach of utilizing linear SVM with binary classification as this method consistently performed well. SVMs have been used in biological systems to classify microarray expression data and predict phenotypes (Furey et al., 2000; Ramaswamy et al., 2001). Most recently, SVMs have been used to analyze maize roots (Zhong et al., 2009). Here, the class of each individual is the rice genotype, e.g. Basmati 217 or Jefferson. Our goal in using SVM analysis was to determine which traits contributed the most to distinguishing between examples from two distinct genotypic classes.
The first step in the SVM analysis was to determine whether root systems of different genotypes could be separated. We used an SVM classifier to separate images of each possible pair of genotypes (see “Materials and Methods” for details). We visualized the results of SVM using a heatmap (Fig. 3; Supplemental Fig. S3). As an example, we present the SVM-derived heatmaps for comparing Basmati 217 and Jefferson along with Lemont and Nipponbare in Figure 3. For all 66 of the total pairwise classifications (Supplemental Fig. S3), the SVM classifier successfully determined the genotype based on single RSA images, the vast majority with greater than 95% accuracy (Supplemental Tables S3 and S4). High classification accuracy ensures that feature rankings obtained from the SVM analysis are meaningful.
Feature ranking via SVM allowed us to determine the relative importance of each feature in each pairwise classification. Figure 4 depicts the normalized feature scores of RSA images from all 66 pairwise classifications (see “Materials and Methods” for explanation of how rankings are calculated). The numerical values of feature ranks for each pairwise comparison between individuals are presented in Supplemental Table S4. Importantly no single feature was the top-ranked feature for all comparisons; depth was the most frequently top-ranked feature (11/66 comparisons) followed by average root radius (top ranked in 9/66 comparisons). The diversity of top-ranked features demonstrates the breadth of variation within the rice germplasm chosen and suggests that features useful for QTL analysis are dependent on the mapping population assayed.
To ensure that the simultaneous analysis of 16 features did not lead to successful classification of genotypes by chance (and therefore spurious trait ranking), we assessed the statistical significance of the accuracy of each pairwise classification via a reshuffling test (see “Materials and Methods” and two illustrative examples in Supplemental Figs. S4 and S5). The statistical test demonstrated that the classifications are unlikely to be the result of chance correlations. For image-based classification, all 66/66 pairwise comparisons were significant (P < 0.001). However, because of correlations in the RSA of images from the same individual, we also developed an individual-based classification and permutation scheme (see “Materials and Methods”). Using the individual-based analysis, we find that all 66 pairwise comparisons are significant (P < 0.05), which confirms that the classification accuracies are unlikely to be the result of chance correlations.
DISCUSSION
We have described a noninvasive imaging and accompanying analysis platform to automatically phenotype and rank traits in complex root systems. Our imaging platform overcomes several limitations of previous noninvasive platforms. It allows observation of unconstrained root growth in three dimensions for at least 2 weeks, while maintaining short imaging times. This is invaluable for QTL experiments, in which hundreds of plants must be imaged. This system is also highly flexible. If a longer growth time is required, the size of the container used for growth can be increased. Further, strength or nutrient composition of the gellan gum can be altered without compromising image quality.
Traditional measurements of root architecture, while useful for describing aspects of the rice root system, do not account for the spatial distribution and the complexity of rice RSA. We therefore chose to phenotype both traditional and novel traits for our analysis. Several of the traditional features we estimate, such as SRL and average root radius, are known to be important for plant fitness (Eissenstat, 1991; Casper and Jackson, 1997; Hodge et al., 1999; Robinson et al., 1999; Hodge, 2009). It is probable that our novel traits, which capture the spatial distribution of the root system, are also contributors to plant fitness. For example, solidity and bushiness are likely to be functionally important indicators of root foraging, but have been inaccessible due to prior imaging and analysis limitations (Hodge, 2004; Kembel et al., 2008). Our method thus allows a more comprehensive picture of RSA and expands the number of RSA QTL that can be identified. However, as both traditional and novel traits are proxies for root system functions, further work is needed to determine their functional importance to RSA.
Our phenotyping demonstrated that RSA varies less within varieties than between varieties, suggesting a strong genetic component to RSA. In rice, key contributors to this variation are the depth and average root radius, as these traits were the top-ranked traits in 20 of 66 pairwise comparisons. Trait ranking showed that the features that best differentiate two genotypes depend on the genotypes compared. For example, between Nipponbare and Jefferson, perimeter and width-to-depth ratio are the top-ranked features (Fig. 4; Supplemental Table S4), while between Nipponbare and Azucena, the volume and SRL are top ranked (Fig. 4; Supplemental Table S4). Some features, such as depth and average root radius, were consistently among the top ranked, while other features, such as the bushiness index, convex area, surface area, maximum width, and width-to-depth ratio, were never among the top-ranked features. However, calculation of the rank of each feature in the classification showed that all these features contributed considerably to specific comparisons (Supplemental Table S4). Notably, even with two highly similar root systems (Fig. 3, D–F) we were able to identify traits that discriminate between the genotypes. This was important, as it suggests that our analysis can identify traits useful in mapping populations, which typically have subtle gradations of phenotypes and that cannot be distinguished by eye.
Graphical depiction of the relative importance of each feature to each of the 66 pairwise comparisons for the 12 rice varieties. Pairs of varieties in each comparison are listed at left. Each row represents the relative importance of every feature (contained in columns), where the contribution of each feature is coded according to the legend color bar at right. Feature rankings are the average of multiple SVM comparisons, each of which has been normalized such that the sum of squared rankings is equal to 1. Actual numerical values for each feature for each comparison are listed in Supplemental Table S4.
Currently, we use multiple images in 2D as a proxy for 3D representation, as 3D reconstruction is not possible with the present version of our setup. Some of our features are 2D traits, such as perimeter, convex area, network area, and bushiness index. However, for other traits like average root radius, SRL, depth, and solidity, we estimate a 3D feature using multiple 2D images. Although our wire model validation demonstrates that multiple 2D images and their analyses provide good estimates of some RSA parameters, ideally, 3D images are needed to provide trait measurements of highest accuracy. We are currently developing the next generation of technology for acquisition, reconstruction, and processing of 3D images of plant roots.
One limitation of our platform is the artificial growth media, since artificial systems like gellan gum differ from soil. While in most cases the link between field traits and RSA has still to be demonstrated, a recent report showed high correlation between root architecture parameters for phytagel-grown soybean (Glycine max) and biomass and nutrient content for field-grown soybean (Fang et al., 2009). Further, QTL for monocot root traits in hydroponics and grain yield in the field have been shown to overlap (Tuberosa et al., 2002). Additional work has identified QTL for root traits in hydroponics in the same chromosomal region with those previously identified in soil-grown rice (Cui et al., 2008), while others have demonstrated strong correlations between hydroponic root growth and field drought resistance (Price et al., 1997), though many of these traits were QTLs of large effect.
Our platform overcomes limitations of previous imaging methods by rapid acquisition of high-quality images. These images, combined with our automated feature calculation, make it possible to phenotype both conventional traits and those that were formerly more difficult to estimate but thought to be crucial to resource uptake by roots (Hodge, 2004; Kembel et al., 2008). Our method thus expands the number of functionally important RSA QTL that can be identified, and in turn, will lead to broader knowledge regarding the genetic mechanisms underlying RSA. In addition, our approach is flexible in terms of the species, traits, and conditions assayed. It is also relatively low cost, and the technology is readily accessible. The feature ranking is automated and rapid, permitting comprehensive phenotyping, and reducing researcher bias. The methodology presented here will contribute to the goal of identifying genes underlying RSA, enabling advances in crop productivity.
MATERIALS AND METHODS
Plant Growth
Seeds were dehulled and sterilized with 10% peroxide for 10 min, followed by 70% ethanol for 1 min, and rinsing three times with sterile water. Sterilized seeds were planted at approximately 1 cm below the surface of the gel. Plants were grown in 2 L ungraduated borosilicate cylinders (Fisher) filled with 750 mL Yoshida's nutrient solution (Yoshida et al., 1976) and 0.3% Gelzan (the highest grade of gellan gum available; Sigma). Cylinders are approximately 520 mm high, with an o.d. of 82.5 mm. To prevent the agar from moving during imaging, six to eight glass pipette tips were bonded to the bottom of each cylinder with Sylgard 184 (Robert McKeown Co.). Plants were grown for 14 d at 12-h day/night, 28°C day, 25°C night. Images were taken 14 dap. Supplemental Fig. S6 details the growth protocol. Twelve rice (Oryza sativa) genotypes were used: 93-11 (9311; indica), Caipo (Cai; tropical japonica), Basmati 217 (Ba217; Basmati), Jefferson (Jef; tropical japonica), Teqing (Teqing; indica), Lemont (Lemont; tropical japonica), Moroberekan (Moro; tropical japonica), Nipponbare (Nipp; temperate japonica), IR64 (IR64; indica), Carolina Gold (Gold; japonica), Bala (indica), and Azucena (tropical japonica). Genotypes were chosen based on parentage of available QTL mapping populations and to represent both of the two major subpopulations of rice (indica and japonica). Analyses were performed on six to 18 individuals of each genotype.
Imaging Platform
Plants were imaged using a PhotoCapture360 turntable and software from Ortery Technologies connected to a Canon PowerShot G7 digital camera and Dell Latitude 620 laptop computer. Cylinders were lit from the sides with fluorescent light bulbs, and from below with a UV light box (VWR); lighting conditions were similar for all experiments. Images from 20 angles per plant per day were acquired. Each plant required approximately 10 min to image. Blurred images were discarded, and the rest were cropped to remove the seed, aerial structures of the plant, and the sides of the cylinder using a batch cropping function in Adobe Photoshop. Images were converted to binary format using an adaptive thresholding method coded in MATLAB. In less than 5% of the images, small blemishes on the cylinders or in the gels that could impact the thresholding process were removed in Adobe Photoshop (no alterations were made to the root system).
Image Preprocessing
A gray-scale image X is given as input in jpeg format to a Matlab routine that performs an adaptive thresholding filter. The filter step is as follows. First, the entire image is broken up into smaller square arrays of length w, defaulted to be 100 × 100 pixels. Next, within each subimage, the mean pixel intensity, I0, is calculated. The values of all other pixels that are a certain fraction C, defaulted to 5%, above I0, are considered to be part of the root network, all others are considered to be part of the background (and are set to 1 and 0, respectively). Finally, all connected clusters of 1-s (i.e. part of the network) are found, and those clusters that are smaller than some critical size, b defaulted to be 5% of the window size, are removed. Because thresholding is local, this method is better suited for nonstationary background intensity as is commonly found in images of root systems. The transformed black and white image will be noted as Y.
Quantification of RSA Traits
The black and white image Y from the image preprocessing step is then provided, as input, to the routine that quantifies network features. The following traits are calculated: perimeter, solidity, convex area, network area, median number of roots, maximum number of roots, bushiness, average root radius, SRL, total length, surface area, maximum width, depth, volume, length distribution, and width-to-depth ratio. The result of this step is a transformation:
where vi denotes a vector of quantitative RSA features for each image Yi.
SVM Analysis
The previous steps described the transformation of grayscale images Xi into black and white images Yi and then feature vectors vi. Each image has a set of labels, Li, which contains information on the individual, including its rice genotype. For pairwise classification, each feature vector can be labeled as 1 (for individuals of genotype 1) or 2 (for individuals of genotype 2).
The SVM method involves a training and a testing step. In training, some fraction (generally half) of the vectors vi are selected and presented to the classifier. We trained a classifier on a randomly chosen subset of half of the images from a pair of genotypes from the 14th dap (the training set) and then used the classifier to classify the other half of the images (the test set). The objective of the classifier is to find a hyperplane (in this case a 15-dimensional hyperplane in a 16-dimensional space) that separates the vectors into two groups—those of type 1 and those of type 2—such that they are maximally separated on opposite sides of the hyperplane (known as the maximum-margin hyperplane). We use a linear SVM classifier that returns a normal vector to the hyperplane. The absolute values of the normalized coordinates of this vector indicate the importance of directions for separation and can be considered as weights of contribution of each feature to the classification. The maximum-margin hyperplane is then used to classify images based on their position with respect to the hyperplane. Figure 4 and Supplemental Table S4 show feature ranks (averaged over multiple SVM comparisons) that are normalized components of the normal vector to the separating hyperplane used for pairwise comparisons.
The accuracy of SVM classification depends on the training set. Since the training set is picked randomly, each SVM analysis results in a slightly different accuracy. Percentage accuracies in the text and figures are the median of all the SVM analyses for a pairwise comparison. A 100% accuracy corresponds to achieving 100% accuracy in half or more of the many SVM trials run for a pairwise comparison.
Statistical Analysis
Intra- versus Intervarietal Feature Variance
The dataset of RSA images has both intra- and intervarietal variance for the 16 features described above. Each genotype is represented by six to 18 individuals with approximately 20 images for each individual. Hence the intravariation of the 16 features can be explained by the variation of RSA among individuals of the same variety as well as among images of the same individual taken from different angles. Intervariance is caused by genetic differences between 12 rice genotypes. We observe that intergenotype variance of RSA features is greater than the variance within a genotype (Figs. 2 and 3). To verify this, we performed a permutation test. We decided to use the permutation test since it does not make any assumptions as to the distribution of the observed values. This test estimates the significance of the difference between the mean values (DI) of a root feature of two genotypes without any assumptions on data distribution. For each pair of genotypes and for each root feature we shuffled the labels of individuals (and hence shuffled the labels of all images of those individuals) 1,000 times and for each shuffle we computed the difference, DIi, i = 1 to 1,000, of the means of feature values given reshuffled labels. The P value of the measured difference DI was estimated as a fraction of the differences that were greater than the measured difference, #(DIi ≥ DI)/1,000. The difference between mean values was considered significant if the corresponding P value was less than 0.0031, which corresponds to 0.05 level of significance with a Bonferroni correction for performing multiple comparisons of 16 features.
The test demonstrated that for each trait there is a pairwise classification for which the mean value of the trait varies significantly (see Supplemental Table S2). In addition, for 64 pairs of genotypes there is at least one, and on average eight, features that vary significantly. Notice that high P values correspond to close mean values plotted on Figure 3.
The P values can be used to obtain alternative feature rankings. However, these rankings take into account only one feature, disregarding its relation to the other 16 features. Nevertheless, features that are near the top of SVM-based rankings have significant P values for most pairwise comparisons (see Supplemental Tables S2 and S4). For example, for a given pairwise comparison, the top-ranked feature using SVM image-based analysis exhibits a statistically significant difference in mean RSA value using a permutation test 42/66 times. By contrast, the bottom-ranked feature using SVM image-based analysis exhibits a statistically significant difference in mean RSA value using a permutation test 18/66 times.
Classification Accuracy
The pairwise classification of images based on traits was tested to see if classification accuracies could have been the result of chance correlations. We randomly shuffled the image labels of all images within a pairwise comparison and redid the SVM analysis 1,000 times (once for each reshuffling). The mean accuracy of image classification for pairwise comparisons given random shuffling of image genotype labels was generally between 40% to 70% (Supplemental Figs. S4 and S5). In all cases the actual accuracy of SVM pairwise classification (Supplemental Table S5) far exceeded the randomization result, demonstrating that the classifications are not the result of chance correlations (P < 0.001).
The above results may be affected by correlations among images of the same individual. Therefore, we have performed the same procedure for blocks of images, where each block represented a different individual. This means that SVM training, classification, and permutations were done using individual-based schemes (classification of a block split by the hyperplane was resolved using the majority rule). A permutation test using 1,000 reshufflings at the individual level showed that the accuracy of all 66 pairwise individual-based comparisons (Supplemental Table S6) were statistically significant with P value less than 0.05. This further strengthens our finding that classifications are unlikely to be the result of chance correlations.
The difference in significance levels between the image-based and the individual-based schemes is caused by the lack of data for the latter scheme (on average 20 individuals versus an average of 390 images in each pairwise comparison), which also resulted in a drop in accuracy from an average of 99% to an average of 92% (Supplemental Table S6). We expect that the accuracy of individual-based classification can be improved by modifications to image preprocessing steps that will remove large-scale noise present in some images.
Supplemental Data
The following materials are available in the online version of this article.
Supplemental Figure S1. Image analysis pipeline.
Supplemental Figure S2. Wire model.
Supplemental Figure S3. SVM heatmaps of the 66 possible pairwise comparisons.
Supplemental Figure S4. Randomization of the root system images of the pairwise comparison for Basmati217 and Jefferson.
Supplemental Figure S5. Randomization of the root system images of the pairwise comparison for Lemont and Nipponbare.
Supplemental Figure S6. Growth and imaging protocol.
Supplemental Table S1. Trait mean values and sds for all 16 traits and 12 genotypes.
Supplemental Table S2. P values of the differences between the mean values of root features for all 66 pairwise classifications.
Supplemental Table S3. Comparison of four traits computed in image analysis pipeline with hand measurements for the wire model.
Supplemental Table S4. Normalized feature weights in the normal of the separating hyperplane for each pairwise comparison for the 12 rice varieties.
Supplemental Table S5. Median accuracy of each image-based pairwise comparison for the 12 rice varieties.
Supplemental Table S6. Median accuracy of each individual-based pairwise comparison for the 12 rice varieties.
Acknowledgments
We thank S. McCouch for seeds, R. Clark for technical assistance, and J. Lynch, L. Kochian, T. Mitchell-Olds, T. Ehrlich, and members of the Benfey lab for comments on the manuscript.
Footnotes
↵1 This work was supported by the National Science Foundation (Plant Genome Research Program grant nos. 0606873 and 0820624 to P.N.B., J.H., and J.S.W.) and Defense Advanced Research Projects Agency (grant no. HR0011–05–1–0057 to J.H., J.S.W., and P.N.B.). J.S.W. holds a Career Award at the Scientific Interface from the Burroughs Wellcome Fund. A.I.P. is supported by a National Institutes of Health Ruth Kirchstein National Research Service Award postdoctoral fellowship.
↵2 These authors contributed equally to the article.
↵3 These authors contributed equally to the article.
The author responsible for distribution of materials integral to the findings presented in this article in accordance with the policy described in the Instructions for Authors (www.plantphysiol.org) is: Philip N. Benfey (philip.benfey{at}duke.edu).
↵[W] The online version of this article contains Web-only data.
↵[OA] Open Access articles can be viewed online without a subscription.
- Received November 9, 2009.
- Accepted January 20, 2010.
- Published January 27, 2010.
REFERENCES
