Background: Image-based high-throughput phenotyping technologies have been rapidly developed in plant science recently, and they provide a great potential to gain more valuable information than traditionally destructive methods. Predicting plant biomass is regarded as a key purpose for plant breeders and ecologists. However, it is a great challenge to find a predictive biomass model across experiments. Results: In the present study, we constructed 4 predictive models to examine the quantitative relationship between image-based features and plant biomass accumulation. Our methodology has been applied to 3 consecutive barley (Hordeum vulgare) experiments with control and stress treatments. The results proved that plant biomass can be accurately predicted from image-based parameters using a random forest model. The high prediction accuracy based on this model will contribute to relieving the phenotyping bottleneck in biomass measurement in breeding applications. The prediction performance is still relatively high across experiments under similar conditions. The relative contribution of individual features for predicting biomass was further quantified, revealing new insights into the phenotypic determinants of the plant biomass outcome. Furthermore, methods could also be used to determine the most important image-based features related to plant biomass accumulation, which would be promising for subsequent genetic mapping to uncover the genetic basis of biomass. Conclusions: We have developed quantitative models to accurately predict plant biomass accumulation from image data. We anticipate that the analysis results will be useful to advance our views of the phenotypic determinants of plant biomass outcome, and the statistical methods can be broadly used for other plant species. Keywords: barley; high-throughput phenotyping; phenomics; biomass; modeling Received: 5 September 2017; Revised: 11 November 2017; Accepted: 9 January 2018 The Author(s) 2018. Published by Oxford University Press. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 2 Chen et al. predictive models has not been well evaluated. Further, it is still Introduction challenging to apply these models across experiments that are Biomass accumulation is an important indicator of crop final performed in different environmental conditions or with differ- product and plant performance. It is thus considered a key ent treatments due to a lack of datasets for assessment. In this trait in plant breeding, agriculture improvement, and ecologi- study, we present a general framework for investigating the rela- cal applications. The conventional approach of measuring plant tionships between plant biomass (referred to as shoot biomass biomass is very time-consuming and labor-intensive as plants hereafter) and image-derived parameters. We applied a multi- need to be harvested destructively to obtain the fresh or dry tude of supervised and unsupervised statistical methods to in- weight . Moreover, the destructive method makes multiple vestigate different aspects of biomass determinants by a list of measurements of the same plant over time impossible. With the representative phenotypic traits in 3 consecutive experiments in development of new technology, digital image analysis has been barley. The results showed that image-based features can accu- used more broadly in many fields, as well as in plant research rately predict plant biomass output and collectively reflect large [2–4]. It allows faster and more accurate plant phenotyping and proportions of the variation in biomass accumulation. We elu- has been proposed as an alternative way to infer plant biomass cidated the relative importance of different feature categories [2, 3, 5]. and of individual features in the prediction of biomass accumu- In recent years, plant biomass has been subject to intensive lation. The differences in the contribution of the image-based investigation by using high-throughput phenotyping (HTP) ap- features for the prediction of the 2 types of biomass measure- proaches in both controlled growth chambers [2, 3, 6–11]and ments, fresh weight and dry weight, were compared as well. Fur- field environments [ 5, 12–17], demonstrating that the ability of thermore, our models were tested for the possibility of predict- imaging-based methods to infer plant biomass accumulation. ing plant biomass in different experiments with different treat- For example, significant genotypic and environmental effects ments. on plant biomass in Setaria were revealed by the Bellwether Phenotyping Platform in a controlled environmental condi- tion . Yang et al.  showed that predicted rice biomass Results (including shoot fresh and dry weight) based on image-derived Development of statistical models for modeling plant morphological and texture features provided a relatively more biomass accumulation using image-based features complete representation than manual measurements in dis- secting its genetic architecture. In this regard, optimized mod- In the previous studies [19, 20], we have shown that a single els plus image-derived features from HTP systems will improve phenotypic trait—3-dimensional digital volume, which is a de- the power of dissecting genetic architecture of complex traits. rived feature from projected side and top areas—can be reason- Although there are some developed models for predicting plant ably predictive to estimate plant biomass accumulation. We ex- biomass, most of them have certain limitations. For example, pect that the predictive power could be improved when mul- Golzarian et al. (2011) modeled the plant biomass (dry weight) in tiple phenotypic traits are combined in a prediction model as wheat (Triticum aestivum L.) as a linear function of projected area, plant biomass is determined not only by structural features assuming plant density was constant. However, this method un- but also by density (physiological properties). To further inves- derestimated the dry weight of salt-stressed plants and over- tigate the relationship between image-derived parameters and estimated that of control plants. Even though the authors ar- plant biomass accumulation, deep phenotyping data that con- gued that the bias was largely related to plant age and the model tain both structural (e.g., geometric traits) and physiological might be improved by including the factor of plant age , the traits (e.g., plant moisture content as reflected by near-infrared differences in plant density between stressed and control plants [NIR]-related traits) were analyzed (Fig. 1A and B). Pot weights of may have been caused by different physiological properties of the plants were not included for the analysis, although they were plants rather than plant age. In another study, Busemeyer et weighed regularly. It might reflect the growth tendency of the al. (2013) developed a calibrated biomass determination model whole plant (shoots and roots), where herein we focused mainly for triticale (x Triticosecale Wittmack L.) under field conditions on shoots. based on multiple linear regression analysis of a diverse set Models were constructed to quantify the ability of imaging- of parameters, considering both the volume of the plants and based features to statistically predict the biomass accumulation. their density. Indeed, this model largely improved the predic- The models were developed by using 4 widely used machine- tion accuracy of the calibration models based on a single type learning methods (Fig. 1C): multivariate linear regression (MLR), of parameters and can precisely predict biomass accumulation multivariate adaptive regression splines (MARS), random for- across environments . However, Buesmeyer et al. (2013) used est (RF), and support vector regression (SVR), which have ex- very limited traits for the model and question whether it could tensively been used in accurate prediction of gene expression be applied broadly in other cases. As mentioned by Yang et al. [21–25] and DNA methylation levels [26–29]. We combined the (2014), noticeable improvement was achieved by adding mor- biomass measurements (fresh weight [FW] and/or dry weight phological features or texture features to the biomass-predicting [DW]) with image-based features and then divided them into a model . This suggests that adding more information/traits training dataset and a test dataset. A model was trained on the could improve the predictive performance of models. Therefore, training dataset and has then been applied to the test dataset a more effective and powerful model is needed to overcome to predict the plant biomass. The relationship between plant these limitations and to allow better utilization of the image- biomass accumulation and image-based features was assessed based plant features, which are obtained from noninvasive phe- based on the criterion of the Pearson correlation coefficient ( r) notyping approaches. between the predicted values and the actual values, or the co- Individual studies have recently shown that the predic- efficient of determination ( R ; the percentage of variance of tion accuracy of plant biomass based on image-derived fea- biomass explained by the model) (Fig. 1D). tures is relatively high even using the simplest linear regres- Our methodology was applied to 3 consecutive experiments sion models [3, 10, 18]. However, the performance of nonlinear (Fig. 2A; Supplemental Table S1 and Data S1), which were Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 Image-derived parameters predicting plant biomass 3 Figure 1: Modeling pipeline for predicting plant biomass accumulation based on image-derived parameters. A, Input data, including high-throughput image data and manually measured biomass data. Plants were phenotyped using various cameras such as visible (or color), fluorescence, and near-infrared sensors. I mage analysis was performed with IAP software  for feature extraction. The same plants were harvested and measured at the end of growth. Generally, 2 types of biomass were measured: fresh weight and dry weight. B, Trait processing. All the phenotypic traits were grouped into 4 categories: geometric, color-related, FLUO-related, and NIR- related traits. Phenotypic data were subjected to quality check to remove low-quality data. C, Each plant was described by a list of traits, resulting in a predictor matrix whose rows represent plants and columns represent image-based traits. This matrix was used to predict plant biomass accumulation by MLR, MARS, RF, and SVR models. The right panel represents the schema of model validation. In the first schema, a dataset (Dataset 1) was divided into training set and testing s et in a 10-fold cross-validation manner. In the second schema, the whole of 1 dataset (Dataset 1) was used for training and another dataset (Dataset 2) was used for testing. D, Model selection, evaluation, and result interpretation. The correlation of the predicted values and measured values was used to assess the overall performance of the model. designed to investigate vegetative biomass accumulation in re- image processing pipeline IAP . A representative list of traits sponse to 2 different watering regimes under semicontrolled for each plant in the last growth day were selected to test their greenhouse conditions in a core set of barley cultivars by nonin- ability to predict plant biomass. vasive phenotyping [20, 30]. There were 312 plants with 18 geno- types for each experiment. Plants were monitored using 3 types Coordinated patterns of plant image–based profiles and of sensors (visible, fluorescence [FLUO], and NIR) in a LemnaTec- their relation to plant biomass Scanalyzer 3D imaging system. An extensive list of phenotypic traits ranging from geometric (shape descriptors) to physiologi- We extracted a list of representative and nonredundant pheno- cal properties (i.e., color-, FLUO-, and NIR-related traits) could be typic traits for each plant from image datasets for each experi- extracted from the image data (Supplemental Data S1) using our ment (see Materials and Methods) (Fig. 1B). In common for these Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 4 Chen et al. Figure 2: Predictability of image-based traits to plant biomass. A, Schema depicting 3 consecutive high-throughput phenotyping experiments in barley. Plants in each experiment were harvested for biomass measurements: fresh weight (for all experiments) and dry weight (only for experiment 1). B, Scatter plots showing projections of the top 4 PCs based on PCA of image-based data. The component scores (shown in points) are colored and shaped according to the experiments (as legend listed in the box). The component loading vectors (represented in lines) of all traits (as colored according to their categories) were superimposed proportionally to their contribution. C, Boxplot showing the distribution of FW across different experiments. D, A dendrogram from cluster analysis based on the means of FW data over genotypes. E, Pearson’s correlation (mean values in the 3 datasets) between image-based traits and FW. Traits with the largest mean correlations values are labeled: 1: sum of leaf length (side view); 2: sum of FLUO intensity (side); 3: plant area border length (side); 4: sum of NIR intensity (top); 5: sum of FLUO intensity (top); 6: projected area (top); 7: projected area (side); and 8: digital volume. experiments, overall 36 high-quality traits that describe plant (e.g., FW) of plants from different experiments with different growth status in the last growth day were obtained. As a result, treatments was significantly different (2-way ANOVA, P < 2e- each dataset was assigned a matrix whose elements were the 16) (Fig. 2C). The relationship was reflected by a dendrogram signals of different features in different plants (Fig. 1C). Prin- from cluster analysis based on the means of FW over genotypes cipal component analysis (PCA) (Fig. 2B) was applied to these (Fig. 2D). Furthermore, the overall phenotypic patterns of these datasets. We found that plants from different experiments with plants were similar to their biomass output (Fig. 2B–D), revealing different treatments showed clearly distinct patterns of pheno- that these image-based features were potential factors reflect- typic profiles. For instance, stressed plants and control plants ing the accumulation of plant biomass. We thus explored the were separated using PCA by their first principal component relationship between the signals of these image-based features (PC1) and also by the top clusters obtained in hierarchical clus- and the level of plant biomass output. We calculated the cor- tering analysis (HCA), while plants from different experiments relation coefficients for each dataset. The correlation patterns were distinguished by PC2 and PC3 in PCA or subordinate clus- were consistent for different datasets, and more than half of the ters in HCA. Accordingly, it could be observed that the biomass features revealed high correlation coefficients ( r > 0.5) (Fig. 2E). Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 Image-derived parameters predicting plant biomass 5 Interestingly, both structural features (such as digital volume, features in a full model (Fig. 4A and B). Strikingly, the pre- projected area, and the length of the projected plant area border) dictability of other types of features (such as color-related and and density-related features (such as NIR and FLUO intensities) FLUO-based traits) was substantial, indicating that these traits were involved in the top ranked features. may act as unforeseen factors in biomass prediction. In addi- tion, the NIR-based features showed higher predictive capability for FW than for DW in control and stressed plants, revealing Relating image-based signals to plant biomass output that NIR signals were important factors in determining FW The above analyses suggest that plant biomass can be at least accumulation. partially inferred from image-based features. To examine which Next, we investigated the relative importance (RI) of each fea- model has the best performance and to select an appropriate ture for predicting biomass using a full model in the whole set model for biomass prediction, we then applied our regression of plants (i.e., “control + stressed plants”) (Fig. 4C and D, upper models (Fig. 1C) to predict plant biomass using image-based fea- panels). In an RF model, the RI of a feature is calculated as the tures. Our analyses were focused on the first experiment (i.e., increase of prediction error (%IncMSE) when phenotypic data for Exp 1), as the phenotypic traits of the corresponding dataset this feature are permuted , and thus indicates the contri- have been intensively investigated in our previous study . In bution of the feature after considering its intercorrelation in a this experiment, plant biomass was quantified in 2 forms: FW model. We found that the top 10 most important features in the and DW. We selected a collection of 45 image-derived parame- full model for predicting FW and DW included both structure- ters from this dataset that were nonredundant and highly rep- and density-related traits. As expected, projected area (from resentative. side or top view) and digital volume were the top ranked fea- We next triedtopredict FWandDWbased on this set tures, which have individually been considered proxies of shoot of image-derived features using 4 different regression models biomass in previous studies [3, 20, 18, 32–37]. However, several (MLR, RF, SVR, and MARS) (Fig. 3). The models were tested on geometric and color-related features that are top ranked in the control plants, stressed plants, and the whole set of plants, re- prediction have not been used in biomass predictions in previ- spectively (Fig. 3A and C). The prediction accuracy of our models ous analysis, although they are widely available among pheno- (the correlation coefficients between the predicted biomass and typing platforms. the actual biomass) was first compared with the ability of indi- In principle, we would expect that highly important features vidual features to predict biomass. It was found that our models in the full model would be related to a high predictive power generally showed better prediction power than the single digital in a degenerate model. Surprisingly, there was no clear corre- volume-based prediction (Fig. 3B and D), indicating that addi- lation observed between the feature importance and its predic- tional features improved the predictive power. Then the perfor- tive power (Fig. 4C and D). For example, several color-related and mance of these models was compared and evaluated. Overall, NIR-based features that were in the top 10 list of the most im- the performance of all the tested models was roughly similar portant features revealed insubstantial predictive power in indi- for the prediction of both FW (Fig. 3B) and DW (Fig. 3D) under vidual models. This observation implies that the relation of the stressed conditions. The prediction accuracy of our models is underlying biomass determinants is extremely complex and not still comparable to the results from previous studies [3, 6, 18] a linear combination of the investigated features. based on MLR models, even though many more features were Furthermore, we compared the relative importance of each considered in our study. The RF model slightly outperformed feature in predicting FW and DW (Fig. 4E). Although a positive other models in predicting biomass of control plants, account- correlation (r = 0.88) between the feature importance for FW and 2 2 ing for the most variance ( R = 0.85 for FW and R = 0.62 for DW could be observed, several features showed large differences DW) (Fig. 3B and D, left panels) and showed the best prediction in their ability to interpret FW or DW, including “nir.intensity” accuracy (Pearson’s correlation r = 0.93 for FW and r = 0.80 for (derived from side view images), “compactness.01” (top), DW) (Fig. 3B and D, middle panels). Of note, RF is the only model “hull.pc1” (top), “leaf.count” (side), “hsv.h.average” (top), and showing better performance than single digital volume-based “lab.a.mean” (top). For instance, NIR intensity and plant com- prediction (Fig. 3D). In this study, we focused on the results from pactness (top view) may be important for predicting FW but not the RF method in the rest of analysis, although results from dif- for DW. We also performed the above analyses by using only ferent methods were highly consistent and led to the same con- control (Supplemental Fig. S1) or stressed plants (Supplemental clusions. Fig. S2), respectively. We found that the patterns of feature im- portance were distinct between these 2 groups of plants. For ex- ample, NIR intensity was ranked as a top 5 feature for predict- Relative importance of different image-based features ing FW for stressed plants but was not substantially important for predicting plant biomass for control plants. These findings suggest that there are differ- As mentioned above, the image-based features could be classi- ences in underlying plant biomass determinants in these kinds fied broadly into 4 categories: plant structure properties, color- of treatment situations that are also reflected by their image- related features, NIR signals, and FLUO-based traits (Fig. 1B). based phenotypic traits. The last 3 types of features reflect plant physiological properties and can be considered plant density–related traits and are thus Image-based features are predictive of plant biomass related to their fresh or dry matter content. For each individual across experiments with similar conditions or feature or each type of features, we constructed a degenerate treatments model for biomass prediction using the corresponding feature(s) as the predictor(s). We compared the capability of each individ- In order to explore whether our models were generalizable ual or type of feature for predicting biomass accumulation in the across different experiments, we applied our models trained in first experiment (i.e., experiment 1). Geometric features showed 1 experiment to predict biomass (herein FW) in other experi- the most predictive power among the 4 categories for prediction ments using a common set of features. Examples of such cross- of both FW and DW, but were slightly less predictive than all experiment predictions are shown in Fig. 5A. We tested and il- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 6 Chen et al. Figure 3: Quantitative relationship between image-based features and plant biomass. A and C, Scatter plots of manually measured plant biomass (fresh weight and dry weight) vs predicted biomass values using 4 prediction models: MLR, MARS, RF, and SVR. The red line indicates the expected prediction (y = x ). The quantitative relationship between image-based features and biomass was evaluated by Pearson’s correlation coefficient (PCC r and its corresponding P-value), RMSRE , and the percentage of variance explained by the models (the coefficient of determination R ). B and D, Summary of the predictive power of each regression model. The results were based on 10-fold cross-validation with 10 trials. Models were evaluated based on control plants, stressed plants, and the whole set of plants. The solid lines represent the predictive performance based on the single “digital volume” feature. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 Image-derived parameters predicting plant biomass 7 Figure 4: The relative importance of image-based features in the prediction of plant biomass. The capabilities of different types of image-based features to predict plant biomass based on evaluation of either fresh weight (A) or dry weight (B). The overall predictive accuracies of each type of features are indicated. Grey bars denote the predictive accuracy using all features. The relative importance of each feature in the RF model (upper panel) and the predictive accuracy of each individual feature as the single predictor (lower panel) based on investigation of either FW (C) or DW (D). The calculation was based on the whole set of plants (control and stressed plants). Note that feature labels are shared in the upper and lower panels. Features are shown in numbers as ordered by their names. The 3 features highlighted in the red dash box are digital volume, projected side area, and projected top area. E, Comparison of the relative importance of features in prediction of FW and DW. The top 6 most different features are highlighted and labeled. lustrated all possibilities for cross-prediction using the whole set generally captured the relationships among the various image- of plants in the corresponding experiment. In general, the pre- based features. However, the third experiment had relatively diction accuracy within individual experiments remained high weaker correlations than the other 2 experiments for predict- 2 2 (r > 0.97 and R > 0.93 for all 3 experiments) (Fig. 5B), revealing ing biomass (with r > 0.81 and R > 0.65 ) (Fig. 5A). This might that our models were effectively predicting plant biomass based be mainly due to seasonal (temperature and illumination) dif- on image-derived feature signals among different experiments. ferences that caused different plants behaviors, namely lower Moreover, the prediction accuracy for cross-experiment predic- biomass for both control and stressed plants in experiment 3 tion, especially between the first 2 experiments ( r > 0.97 and . This suggests that different plant growth conditions might R > 0.94 ), was still relatively high, implying that our models cause some variation for cross-experiment prediction. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 8 Chen et al. Figure 5: Comparison of prediction accuracy across different experiments. A, Biomass prediction across experiments. Models were trained using data from 1 experiment and were applied to another experiment for prediction. The whole set of plants (i.e., “control + stressed” plant) was used in the analysis. Brown triangles denote stressed plants, and green circles control plants. The red box indicates that the prediction accuracy is relatively high between experiments 1 and 2. B, Boxplots of coefficient determination (R , left), Pearson’s correlation coefficients ( r, middle), and the root mean squared relative error (right) for different comparisons. “Within” denotes a model trained and tested on data from the same dataset with specific treatments (control, stress, or both), and “Cross” represents a model trained on 1 d ataset and tested on another dataset. “Control → stress” denotes a model trained on data with control treatment and tested on data with stress treatment, and vice versa for “stress → control.” The number of possible analyses for each category is shown above the boxes. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 Image-derived parameters predicting plant biomass 9 Figure 6: Comparison of prediction accuracy across different treatments. Refer to Fig. 5A for legend. The analysis was performed for control and stressed plants separately. At the same time, we tested cross-predictability of our mod- less area out of range. Therefore, image quality would be an im- els using treatment-specific data in the experiments (Fig. 6). portant variation source for our modeling and should be taken Similar results were obtained as above using the whole dataset into consideration for any application. (Fig. 5B). The weak predictive power for cross-prediction involv- ing control plants from the third experiment was most clearly Discussion observable in the low accuracy in the biomass prediction of this particular subset of plants. Generally, control and stressed Biomass is a complex but important trait in functional ecol- plants were found to have very weak predictive power when ogy and agronomy for studying plant growth, crop produc- related to each other (Fig. 6), as also supported by the distinct ing potential, and plant regeneration capabilities. Many differ- patterns of relative feature importance between these 2 plant ent techniques, either destructive or nondestructive, have been groups (Supplemental Figs S1 and S2). For each experiment, the used to estimate biomass [1–3, 5–17]. Compared with the tra- prediction accuracy was higher for stressed plants compared ditional destructive methods for measuring biomass, nonde- with control plants. This might have resulted from the imaging structive imaging methods provide a faster, more accurate ap- analysis process. Relatively small plants, stressed plants in this proach for plant phenotyping. In recent years, more and more case, would gain more clear images due to less overlapping or high-throughput plant phenotyping platforms have been set Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 10 Chen et al. up and applied worldwide. Accordingly, it becomes a current parameter (such as projected area) or several geometric pa- challenge to establish models utilizing the big datasets gained rameters, our analyses extended these studies by incorporat- from high-throughput imaging systems. Accurately predicting ing more representative features that cover both structural biomass from image data requires efficient mathematical mod- and physiological-related properties into a more sophisticated els as well as representative image-derived features. model. Although the predictive power of our model is roughly In this study, we have presented a systematic analysis of re- higher than that of single feature-based prediction, such as dig- lationships between plant biomass accumulation and image- ital volume (Fig. 3), our model also reveals the relative con- derived signals to confirm the assumption that biomass can be tribution of individual features in the prediction of biomass. The accurately predicted from image-based parameters. We built a information regarding the importance of each feature will offer random forest model of biomass accumulation using a com- new insights into the phenotypic determinants of plant biomass prehensive list of representative image-based features. The outcome. Interestingly, we found that several top ranked fea- comparison between a random forest model and alternative tures, such as digital volume and NIR intensity, showed genetic regression models indicated that the RF model outperforms correlations with biomass of fresh weight (Fig. 4C) , implying other models in terms of (1) better predictive power—especially that these top-ranked features may represent the main “phe- in comparison with the linear model, confirming the com- notypic components” of biomass outcome and that they can be plex phenotypic architecture of biomass, (2) better outperfor- further used to dissect genetic components underlying biomass mance than a single feature prediction model—arguing the accumulation. As image-based high-throughput phenotyping in complex phenotypic makeup of biomass, and (3) feasible bio- plants developed mainly in recent years, and therefore few cor- logical interpretability—the ability to readily extract informa- responding modeling studies have been performed, we believe tion about the importance of each feature in prediction. The that our model could be further improved when new types of high prediction accuracy based on this model, in particular the cameras and/or newly defined features become available. cross-experiment performance, is promising to relieve the phe- In summary, we have developed a quantitative model for notyping bottleneck in biomass measurement in breeding appli- dissecting the phenotypic components of biomass accumula- cations. For example, based on the established small reference tion based on image data. Apart from predicting biomass out- dataset that is used to train an RF model, it is possible to pre- come, the methods can be used to determine the most impor- dict biomass in several large plant populations within 1 experi- tant image-based features related to plant biomass accumula- ment or across several experiments using image data by taking tion, which are promising for subsequent genetic mapping to advantage of high-throughput phenotyping technologies. Alter- uncover the genetic basis of biomass. natively, the model can be trained from a much larger reference panel of plants that are grown in diverse environmental condi- Potential implications tions, which is then applied to a diverse set of experiments. The first evidence for this notion is the observation that our model As high-throughput plant phenotyping is a technique that is be- showed more predictive power in plants with 2 treatments than coming more and more widely used for automated phenotype in with a single treatment (Fig. 3, B and D). Indeed, when applying plant research, especially in plant breeding, we anticipate that our model to the combined dataset from all 3 experiments, we the methodologies proposed in this work will have various po- found that the prediction accuracy remains very high (R = 0.96 tential applications. We anticipate that the analysis results will and r = 0.98, average values from 10 times the 10-fold cross- be useful to advance our views of the phenotypic determinants validation). To keep the high prediction accuracy in other ap- of plant biomass outcome, and the statistical methods can be plications, there are some points that should be taken into con- broadly used for other plant species and therefore assist plant sideration. Considering the environmental effects on biomass breeding in the context of phenomics. accumulation, the application of our model will require the test- ing experiments to show similar conducted conditions as those of the reference experiments. This means that the plant cul- Materials and Methods tivation conditions should be standardized and any noise that Germplasm and experiments might lower image quality should be avoided. Another approach to improve applicability of models, which could not be tested in Barley plant image data were obtained as described previously this study, would be to improve the database for the training by [20, 30]. Briefly, a core set of 16 2-rowed spring barley culti- acquiring data from additional environmental sensors. Temper- vars (Hordeum vulgare L.) and 2 parental cultivars of a double ature, humidity, and illumination data would certainly help to haploid (DH) were monitored for vegetative biomass accumula- explain differences in the growth patterns among experiments, tion. Three independent experiments with identical setup were performed in different growth seasons. To this end, we expect performed in a (semi-)controlled greenhouse at IPK by using that our approach is extensible by incorporating such sensor the automated phenotyping and imaging platform LemnaTec- data in the data matrices. Furthermore, our results can provide Scanalyzer 3D. Experiments were performed consecutively from suggestive hints for biologists to set up phenotyping infrastruc- May to November 2011 over a period of 58 days each (Supple- tures for investigation of plant biomass. For instance, a visible mental Table S1). The greenhouse setup enabled sowing for the light imaging system would be sufficient to accurately predict next experiment 2 days before the old experiment ended. For fresh weight based on the observation that geometric features this, new pots were placed in the middle of the greenhouse, alone show high prediction accuracy (Fig. 4A). However, to in- while the old experiment was still on the conveyer belts. vestigate dry weight, it would be helpful to include an additional Each experiment consisted of 2 treatments: well-watered near-infrared camera system under normal growth conditions (control treatment) and water-limited (drought stress treat- and an additional fluorescence camera system under drought ment). In each treatment, 9 plants per core set cultivar as well as stress conditions (Fig. 4B). 6 plants per DH parent were tested. This resulted in a total of 312 In contrast to previous studies [2, 3, 6, 7, 18, 32–37], in which plants per experiment, corresponding to the maximal capacity biomass was investigated using only a single image-derived of the phenotyping platform. Watering and imaging were per- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 Image-derived parameters predicting plant biomass 11 formed daily. Drought stress was imposed by intercepting water Models for predicting plant biomass supply from 27 days after sowing (DAS 27) to DAS 44. Stressed To understand the underlying relationship between image- plants were rewatered at DAS 45. In total, for each experiment, derived parameters and the accumulated biomass (such as FW about 100 GB of raw (image) data was accumulated. At the end of and DW), we constructed predictive models based on 4 different experiments (DAS 58), plants were harvested to measure above- machine-learning methods: MLR, MARS, RF, and SVR. In these ground biomass in the form of plant fresh weight (for all exper- models, the normalized phenotypic profile matrices X for n×m iments) and/or dry weight (for experiment 1). a representative list of phenotypic traits were used as predic- tors (explanatory variables) and the measured DW/FW as the re- Image analysis sponse variable Y . All these models were implemented in R (release 2.15.2) . Image datasets were processed by the barley analysis pipelines To assess the relative contribution of each phenotypic trait to in the IAP software (version 1.1.2) . Analyzed results were ex- predicting the biomass, we also calculated the relative feature ported in the csv file format via IAP functionalities, which can importance for each model. Specifically, for the MLR model, we be used for further data inspection. The result table includes used the “lm” function in the base installation packages. The columns for different phenotypic traits and rows as plants are relative importance of predictor variables in the MLR model was imaged over time. The corresponding metadata are included in estimated by a heuristic method , which decomposes the the resulting table as well. proportionate contribution of each predictor variable to R .For Each plant was characterized by a set of phenotypic traits MARS, we used the “earth” function in the earth Rpackage.The also referred to as features, which were grouped into 4 cate- “number of subsets (nsubsets)” criterion (counting the number gories: geometric features, fluorescence-related features, color- of model subsets that include the variable) was used to calculate related features, and near-infrared-related features. These traits the variables’ feature importance, which is implemented in the were defined by considering image information from differ- “evimp” function. For the RF model, we used the randomForest R ent cameras (visible light, fluorescence, and near infrared) and package, which implements Breiman’s random forest algorithm imaging views (side and top views). See the IAP online documen- . We chose the “%IncMSE” (increase of mean squared error) to tation  for details about trait definition. represent the criteria of relative importance measure. For SVR, we utilized the e1071 R package, which provides functionalities Feature selection to use the libsvm library . The absolute values of the coef- ficients of the normal vector to the “optimal” hyperplane can Feature selection was performed with the same procedure as be considered the relative importance of each predictor variable described in Chen et al. . We applied the feature selection contributing to regression [42, 43]. technique to each dataset. Generally, we captured almost identi- cal subset features from different datasets. We manually added Evaluation of the prediction models several representative traits due to removal by variance infla- tion factors. For example, the digital volume and projected area To evaluate the performance of the predictive models, we are highly correlated with each other, but we kept both of them, adopted a 10-fold cross-validation strategy to check the predic- because we would investigate the predictive power of both fea- tion power of each regression model. Specifically, each dataset tures. Moreover, the regression models we used are insensitive was randomly divided into a training set (90% of plants) and a to collinear features. We thus kept as much of the representative testing set (10% of plants). We trained a model on the training features as possible. To apply the prediction models among dif- data and then applied it to predict biomass for the testing data. ferent datasets, a common set of features supported by all the Afterwards, the predicted biomass in the testing set was com- datasets was used. pared with the manually measured biomass. The predictive ac- curacy of the model can be measured by Data transformation (1) the Pearson correlation coefficient ( r) between the predicted Each plant can be presented by a representative list of pheno- values and the observed values; typic traits, resulting in a matrix X for each experiment, n×m (2) the coefficient of determination ( R ), which equals the fraction where n is the number of plants and m is the number of pheno- of variance of biomass explained by the model, defined as typic traits. Missing values were filled by mean values of other replicated plants. To make the image-derived parameters from n 2 SS (y − yˆ ) res i i 2 i = 1 diverse sources comparable, we normalized the columns of X R = 1 − = 1 − , n 2 SS tot (y − y¯) i = 1 by dividing the values with the maximum value of each column across all plants. Plants with empty values of manual measure- ments (FW and DW) were discarded for analysis. These trans- where SS and SS are the sum of squares for residuals res tot formed datasets were subjected to regression models. and the total sum of squares, respectively, yˆ is the predicted biomass of the i th plant, y is the observed biomass of the i Hierarchical clustering analysis and PCA th plant, y¯ is the mean value of the observed biomass, and Hierarchical clustering analysis and principle component analy- (1) the root mean squared relative error of cross-validation is de- fined as sis were performed on the transformed data matrix X in the n×m same way as described in Chen et al. . We also performed HCA using the genotype-level mean value of FW data to check s y −yˆ i i i =1 y the similarity of overall plant growth patterns in different exper- RMSRE = , iments. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 12 Chen et al. where s denotes the sample size of the testing dataset. infrared; PCA: principal component analysis; PCC: Pearson We repeated the cross-validation procedure 10 times. The correlation coefficient; RF: random forest; RMSRE: root mean mean and standard deviation of the resulting R and RMSRE squared relative error; SVR: support vector regression. values were calculated across runs. To evaluate the applicability of our methods across seasons Competing interests (thus different growth environments) and treatments (e.g., con- The authors declare that they have no competing interests. trol vs drought stress) in the same season, we applied the mod- els in different contexts with cohort validation. Specifically, we Funding trained the biomass prediction models under 1 specific context and predicted biomass in another different context and vice This work was supported by the Leibniz Institute of Plant versa. The predictive accuracy of the model was evaluated based Genetics and Crop Plant Research (IPK), the Robert Bosch on the measures R and RMSRE, as described above. Further- Stiftung (32.5.8003.0116.0) and the Federal Agency for Agricul- more, the predictive power was reflected by the bias μ between ture and Food (BEL; 15/12-13, 530-06.01-BiKo CHN), and the Fed- the predicted and observed values, defined as eral Ministry of Education and Research (BMBF; 0315958A and 031A053B). This research was furthermore enabled with support 1 yˆ − y i i of the European Plant Phenotyping Network (EPPN; grant agree- μ = · , i =1 n y ment No. 284 443), funded by the FP7 Research Infrastructures Programme of the European Union. where n denotes the sample size of the dataset. This bias indi- cates over- (μ> 0 ) or underestimation (μ< 0 ) of biomass. Author contributions D.C. designed the research. C.K. and M.C. supervised the project. Availability of source code and requirements K.N. and G.A. performed the LemnaTec experiments. D.A. cre- ated the ISA-Tab formatted description and uploaded data Project name: Modeling of plant biomass accumulation with records in the PGP repository. J.M.P. and C.K. analyzed image HTP data data. D.C. implemented the methods, analyzed data, interpreted Project home page: https://github.com/htpmod/HTPmod the results, and wrote the manuscript, with contributions from Operating system(s): Windows, Linux, and Mac OS R.S. All authors read and approved the final version of the article. Programming language: R License: open source under GNU GPL v3.0 Acknowledgements We would like to thank Ingo Muc ¨ ke for his management of the Availability of supporting data and materials LemnaTec system operations. We thank Michael Ulrich for per- The raw image datasets, as well as analyzed data sup- forming software tests and helping in data analysis. We would porting the results of this article, are available in the like to thank Dr Malia Gehan and Dr Christian Fournier for their PGP repository under https://dx.doi.org/10.5447/IPK/2017/24, helpful comments and suggestions. https://dx.doi.org/10.5447/IPK/2017/25,and https://dx.doi.org/ 10.5447/IPK/2017/26, according to the ISA-Tab format and the References recommendations of the Minimum Information About a Plant Phenotyping Experiment (MIAPPE) standard . The selected 1. Catchpole WR, Wheeler CJ. Estimating plant biomass: a re- data for modeling are available in Supplemental Data S1. Sup- view of techniques. Austral Ecol 1992;17(2):121–31. porting data, including metadata tables, raw image files, and an 2. Tackenberg O. A new method for non-destructive measure- archival copy of HTPmod are also available via the GigaScience ment of biomass, growth rates, vertical biomass distribution repository, GigaDB . and dry matter content based on digital image analysis. Ann Bot 2007;99(4):777–83. Additional files 3. Golzarian MR, Frick RA, Rajendran K et al. Accurate infer- ence of shoot biomass from high-throughput images of ce- Supplemental Figure S1. The relative importance of image-based real plants. Plant Methods 2011;7(1):2. features in the prediction of biomass in control plants. Refer to 4. Fahlgren N, Gehan MA, Baxter I. Lights, camera, action: high- Fig. 4 for legend. The calculation was based on control plants. throughput plant phenotyping is ready for a close-up. Curr Supplemental Figure S2. The relative importance of image- Opin Plant Biol 2015;24:93–99. based features in the prediction of biomass in stressed plants. 5. Montes JM, Technow F, Dhillon BS et al. High-throughput Refer to Fig. 4 for legend. The calculation was based on stressed non-destructive biomass determination during early plant plants. development in maize under field conditions. Field Crops Res Supplemental Table S1. Overview of 3 high-throughput phe- 121(2):268–73. notyping experiments in barley. 6. Feng H, Jiang N, Huang C et al. A hyperspectral imaging sys- Supplemental Data S1. Manual data and image-derived data tem for an accurate prediction of the above-ground biomass in the 3 experiments. of individual rice plants. Rev Sci Instrum 2013;84(9):095107. 7. Neumann K, Zhao Y, Chu J et al. Genetic architecture and Abbreviations temporal patterns of biomass accumulation in spring barley revealed by image analysis. BMC Plant Biol 2017;17(1):137. DAS: days after sowing; DW: dry weight; FLUO: fluorescence; 8. Zhang X, Huang C, Wu D et al. High-throughput phenotyping FW: fresh weight; HCA: hierarchical clustering analysis; HTP: and QTL mapping reveals the genetic architecture of maize high-throughput phenotyping; MARS: multivariate adaptive re- plant growth. Plant Physiol 2017;173(3):1554–64. gression splines; MLR: multivariate linear regression; NIR: near- Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018 Image-derived parameters predicting plant biomass 13 9. Muraya MM, Chu J, Zhao Y et al. Genetic variation of growth 27. Zhang W, Spector T, Deloukas P et al. Predicting genome- dynamics in maize (Zea mays L.) revealed through automated wide DNA methylation using methylation marks, ge- non-invasive phenotyping. Plant J 2017;89(2):366–80. nomic position, and DNA regulatory elements. Genome Biol 10. Fahlgren N, Feldman M, Gehan MA et al. A versatile pheno- 2015;16(1):14–14. typing system and analytics platform reveals diverse tem- 28. Das R, Dimitrova N, Xuan Z et al. Computational prediction poral responses to water availability in Setaria. Mol Plant of methylation status in human genomic sequences. Proc 2015;8(10):1520–35. Natl Acad Sci U S A 2006;103(28):10713–6. 11. Yang WN, Zilong Guo ZL, Huang CL et al. Combining 29. Zheng H, Wu H, Li J et al. CpGIMethPred: computational high-throughput phenotyping and genome-wide associa- model for predicting methylation status of CpG islands in tion studies to reveal natural genetic variation in rice. Nat human genome. BMC Med Genomics 2013;6(Suppl 1):S13. Commun 2014;5:5087. 30. Neumann K, Klukas C, Friedel S et al. Dissecting spatiotem- 12. Ehlert D, Horn H-J, Adamek R. Measuring crop biomass poral biomass accumulation in barley under different wa- density by laser triangulation. Comput Electron Agric ter regimes using high-throughput image analysis. Plant Cell 2008;61(2):117–25. Environ 2015;38(10):1980–96. 13. Ehlert D, Heisig M, Adamek R. Suitability of a laser 31. Breiman L. Random forests. Machine Learning 2001;45(1): rangefinder to characterize winter wheat. Precis Agric 5–32. 2010;11(6):650–63. 32. Dietz H, Steinlein T. Determination of plant species cover by 14. Erdle K, Mistele B, Schmidhalter U. Comparison of active and means of image analysis. J Veg Sci 1996;7(1):131–6. passive spectral sensors in discriminating biomass parame- 33. Leister D, Varotto C, Pesaresi P et al. Large-scale evaluation ters and nitrogen status in wheat cultivars. Field Crops Res of plant growth in Arabidopsis thaliana by non-invasive image 2011;124(1):74–84. analysis. Plant Physiol Biochem 1999;37(9):671–8. 15. Busemeyer L, Ruckelshausen A, Moller K et al. Precision phe- 34. Paruelo JM, Lauenroth WK, Roset PA. Estimating above- notyping of biomass accumulation in triticale reveals tem- ground plant biomass using a photographic technique. J poral genetic patterns of regulation. Sci Rep 2013;3(1):2442. Range Manag 2000;(2):190–3. 16. Cao Q, Miao Y, Wang H et al. Non-destructive estimation of 35. Walter A, Scharr H, Gilmer F et al. Dynamics of seedling rice plant nitrogen status with crop circle multispectral ac- growth acclimation towards altered light conditions can tive canopy sensor. Field Crops Res 2013;154:133–44. be quantified via GROWSCREEN: a setup and procedure 17. Fernandez MGS, Bao Y, Tang L et al. High-throughput designed for rapid optical phenotyping of different plant phenotyping for biomass crops. Plant Physiol 2017; species. New Phytol 2007;174(2):447–55. 10.1104/pp.17.00707. 36. Arvidsson S, Perez-Rodriguez P, Mueller-Roeber B. A growth 18. Neilson EH, Edwards AM, Blomstedt CK et al. Utilization of a phenotyping pipeline for Arabidopsis thaliana integrating im- high-throughput shoot imaging system to examine the dy- age analysis and rosette area modeling for robust quantifi- namic phenotypic responses of a C4 cereal crop plant to ni- cation of genotype effects. New Phytol 2011;191(3):895–907. trogen and water deficiency over time. J Exp Bot 2015;(7). 37. Hairmansis A, Berger B, Tester M et al. Image-based pheno- 19. Klukas C, Chen D, Pape JM. Integrated Analysis Platform: an typing for non-destructive screening of different salinity tol- open-source information system for high-throughput plant erance traits in rice. Rice 2014;7(1):16–16. phenotyping. Plant Physiol 2014;165(2):506–18. 38. Klukas C. The Integrated Analysis Platform, 2014. http:// 20. Chen D, Neumann K, Friedel S et al. Dissecting the phe- iapg2p.sourceforge.net/documentation.pdf. notypic components of crop plant growth and drought re- 39. R Core Team, The R Project for Statistical Computing. sponses based on high-throughput image analysis. Plant Cell http://www.r-project.org/. 2014;26(12):4636–55. 40. Johnson JW. A heuristic method for estimating the relative 21. Cheng C, Alexander R, Min R et al. Understanding transcrip- weight of predictor variables in multiple regression. Multivar tional regulation by integrative analysis of transcription fac- Behav Res 2000;35(1):1–19. tor binding data. Genome Res 2012;22(9):1658–67. 41. Chang C-C, Lin C-J. LIBSVM: a library for support vector ma- 22. Cheng C, Gerstein M. Modeling the relative relationship chines. ACM Transact Intell Syst Technol 2011;2(3):27. of transcription factor binding and histone modifications 42. Loo LH, Wu LF, Altschuler SJ. Image-based multivariate pro- to gene expression levels in mouse embryonic stem cells. filing of drug responses from single cells. Nat Methods Nucleic Acids Res 2012;40(2):553–68. 2007;4(5):445–53. 23. Cheng C, Yan K-K, Yip KY et al. A statistical framework 43. Iyer-Pascuzzi AS, Symonova O, Mileyko Y et al. Imaging and for modeling gene expression using chromatin features analysis platform for automatic phenotyping and trait rank- and application to modENCODE datasets. Genome Biol ing of plant root systems. Plant Physiol 2010;152(3):1148–57. 2011;12(2):R15. 44. Arend D, Junker A, Scholz U et al. PGP repository: a plant 24. Dong X, Greven MC, Kundaje A et al. Modeling gene expres- phenomics and genomics data publication infrastructure. sion using chromatin features in various cellular contexts. Database (Oxford) 2016; doi:10.1093/database/baw033. Genome Biol 2012;13(9):R53. 45. Cwiek-Kupczynska H, Altmann T, Arend D et al. Measures 25. Karlic´ R, Chung H-R, Lasserre J et al. Histone modification for interoperability of phenotypic data: minimum informa- levels are predictive for gene expression. Proc Natl Acad Sci tion requirements and formatting. Plant Methods 2016;12:44. U S A 2010;107(7):2926–31. 46. Chen D, Shi R, Pape JM et al. Supporting data for “Pre- 26. Ma B, Wilker EH, Willis-Owen SAG et al. Predicting DNA dicting plant biomass accumulation from image-derived methylation level across human tissues. Nucleic Acids Res parameters.” GigaScience Database 2018. http://dx.doi.org 2014;42(6):3515–28. /10.5524/100392. Downloaded from https://academic.oup.com/gigascience/article-abstract/7/2/1/4810759 by Ed 'DeepDyve' Gillespie user on 16 March 2018
GigaScience – Oxford University Press
Published: Feb 1, 2018
It’s your single place to instantly
discover and read the research
that matters to you.
Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.
All for just $49/month
Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly
Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.
All the latest content is available, no embargo periods.
“Whoa! It’s like Spotify but for academic articles.”@Phil_Robichaud