An automated, high-throughput method for standardizing image color profiles to improve image-based plant phenotyping

Software

Image analysis was performed with the C++ source file included with this manuscript and must be compiled against OpenCV (only tested on 3.1). Statistical analyses and graphics were done using R version 3.4.4 ( R Core Team, 2018) with the following packages: ggplot 2.2.1, ggthemes 3.4.0, reshape2 1.4.3, plyr 1.8.4, grid 3.4.4, gridextra 2.2.1, vegan 2.4-6. T-test, and Wilcoxon signed-rank test were done with base R functions. Constrained analysis of principal coordinates was done with vegan.

Plant growth and imaging conditions

The plant growth phenotyping data used in this experiment is a part of a large dataset. Plant growth conditions were based on a previous experiment (Veley et al., 2017) with adjustments. Details are as follows: Two sorghum genotypes (Sorghumbicolor (L.) Moench, genotypes BTx623 and China 17) were germinated between pieces of damp filter paper in petri dishes at 30 °C for 3 days prior to planting. The genotype information was collapsed throughout this analysis. Square pots (7.5 cm wide, 20 cm high) fitted with drainage trays were pre-filled with Profile Field & Fairway™ calcined clay mixture (Hummert International, Earth City, Missouri). On day 4 (4 DAP), germinated seeds were transplanted into pre-moistened pots, barcoded (including genotype identification, treatment group, and a unique pot identification number), randomized and scanned onto the Bellwether Phenotyping System (Conviron, day/night temperature: 32 °C/22 °C, day/night humidity: 40%/50%, day length: 14 hr (ZT 0–14 = hr 0,700–2,100), night length: 10 hr (ZT 14–24 = hr 2,100–0,700), light source: metal halide and high pressure sodium, light intensity: 400 µmol/m²/s). Plants received 40 mL fertilizer treatments (see below) daily, and were watered to weight (640 g starting dry weight, 1,070 g target weight DAP 4–11, 960 g DAP 12–25, not including weight of phenotyping system carrier) daily by the system using distilled water. Fertilizer treatment information was collapsed for all analysis except where indicated, and were as follows:

Fertilizer treatments

Control fertilizer treatment (100% N): 6.5 mM KNO₃, 4.0 mM Ca(NO₃)₂4H₂O, 1.0 mM NH₄H₂PO₄, 2.0 mM MgSO₄7H₂O, micronutrients, pH  4.6.

Low nitrogen fertilizer treatment (10% N): 0.65 mM KNO₃, 4.95 mM KCl, 0.4 mM Ca(NO₃)₂4H₂O, 3.6 mM CaCl₂2H₂O, 0.1 mM NH₄H₂PO₄, 0.9 mM KH₂PO₄, 2.0 mM MgSO₄7H₂O, micronutrients, pH  5.0.

The same micronutrients were used for both treatments: 4.6 µM H₃BO₃, 0.5 µM MnCl₂4H₂O, 0.2 µM ZnSO₄7H₂O, 0.1 µM (NH₄)₆Mo₇O₂₄4H₂O, 0.2 µM CuSO₄5H₂O, 71.4 µM Fe-EDTA.

Standardization method

The image standardization method used here is based on a color transfer method that can adjust the colors in an image to match a target image color profile (Gong, Finlayson & Fisher, 2016). The goal is to create a transform such that when applied to the values of every pixel in a source image, it returns values mapped to a target image profile. It was shown that the color transfer is a single homography from source to target (Gong, Finlayson & Fisher, 2016). Hue-preserving methods have been characterized in detail and may be sufficient for some datasets (Mackiewicz, Andersen & Finlayson, 2016). However, we aimed to change the hue profile to eliminate hue bias introduced by light source quality batch effects. Color data from the source and the target image are used to define the homography and the color data must be sampled from homologous points, so a reference is included in the images from which the values of R, G, and B can be measured. We included a ColorChecker Passport Photo (X-Rite, Inc.), which has a panel of 24 industry standard color reference chips, within each image. A target image (or reference) is declared and the process of computing the homography and applying it to a given source image is as follows:

1. Let T and S be matrices containing measurements for the R, G, and B components of each of the ColorChecker reference chips pictured in the target image and source image respectively. (1)${\displaystyle T=\left[\begin{array}{ccc}\hfill T_R\hfill & \hfill T_G\hfill & \hfill T_B\hfill \end{array}\right]=\left[\begin{array}{ccc}\hfill t_{r_1}\hfill & \hfill t_{g_1}\hfill & \hfill t_{b_1}\hfill \\ \hfill t_{r_2}\hfill & \hfill t_{g_2}\hfill & \hfill t_{b_2}\hfill \\ \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill \\ \hfill t_{r_{24}}\hfill & \hfill t_{g_{24}}\hfill & \hfill t_{b_{24}}\hfill \end{array}\right]\text{and}S=\left[\begin{array}{ccc}\hfill S_R\hfill & \hfill S_G\hfill & \hfill S_B\hfill \end{array}\right]=\left[\begin{array}{ccc}\hfill s_{r_1}\hfill & \hfill s_{g_1}\hfill & \hfill s_{b_1}\hfill \\ \hfill s_{r_2}\hfill & \hfill s_{g_2}\hfill & \hfill s_{b_2}\hfill \\ \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill \\ \hfill s_{r_{24}}\hfill & \hfill s_{g_{24}}\hfill & \hfill s_{b_{24}}\hfill \end{array}\right]}$

2. Extend S to include the square and cube of each element. (2)${\displaystyle S=\left[\begin{array}{ccccccccc}\hfill s_{r_1}\hfill & \hfill s_{g_1}\hfill & \hfill s_{b_1}\hfill & \hfill s_{r_1}^2\hfill & \hfill s_{g_1}^2\hfill & \hfill s_{b_1}^2\hfill & \hfill s_{r_1}^3\hfill & \hfill s_{g_1}^3\hfill & \hfill s_{b_1}^3\hfill \\ \hfill s_{r_2}\hfill & \hfill s_{g_2}\hfill & \hfill s_{b_2}\hfill & \hfill s_{r_2}^2\hfill & \hfill s_{g_2}^2\hfill & \hfill s_{b_2}^2\hfill & \hfill s_{r_2}^3\hfill & \hfill s_{g_2}^3\hfill & \hfill s_{b_2}^3\hfill \\ \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill \\ \hfill s_{r_{24}}\hfill & \hfill s_{g_{24}}\hfill & \hfill s_{b_{24}}\hfill & \hfill s_{r_{24}}^2\hfill & \hfill s_{g_{24}}^2\hfill & \hfill s_{b_{24}}^2\hfill & \hfill s_{r_{24}}^3\hfill & \hfill s_{g_{24}}^3\hfill & \hfill s_{b_{24}}^3\hfill \end{array}\right]}$

3. Calculate M, the Moore–Penrose inverse matrix of S (* denotes transpose). (3)${\displaystyle M={\left(S^{\ast }S\right)}^{-1}{S}^{\ast }}$

4. Estimate standardization vectors of each R, G, and B channel by multiplying M with each column of T. (4)${\displaystyle \left[\begin{array}{c}\hfill R_h\hfill \\ \hfill G_h\hfill \\ \hfill B_h\hfill \end{array}\right]=\left[\begin{array}{c}\hfill M{T}_R\hfill \\ \hfill M{T}_G\hfill \\ \hfill M{T}_B\hfill \end{array}\right]=\left[\begin{array}{ccccccccc}\hfill {\hat{r}}_r\hfill & \hfill {\hat{r}}_g\hfill & \hfill {\hat{r}}_b\hfill & \hfill {\hat{r}}_{r^2}\hfill & \hfill {\hat{r}}_{g^2}\hfill & \hfill {\hat{r}}_{b^2}\hfill & \hfill {\hat{r}}_{r^3}\hfill & \hfill {\hat{r}}_{g^3}\hfill & \hfill {\hat{r}}_{b^3}\hfill \\ \hfill {\hat{g}}_r\hfill & \hfill {\hat{g}}_g\hfill & \hfill {\hat{g}}_b\hfill & \hfill {\hat{g}}_{r^2}\hfill & \hfill {\hat{g}}_{g^2}\hfill & \hfill {\hat{g}}_{b^2}\hfill & \hfill {\hat{g}}_{r^3}\hfill & \hfill {\hat{g}}_{g^3}\hfill & \hfill {\hat{g}}_{b^3}\hfill \\ \hfill {\hat{b}}_r\hfill & \hfill {\hat{b}}_g\hfill & \hfill {\hat{b}}_b\hfill & \hfill {\hat{b}}_{r^2}\hfill & \hfill {\hat{b}}_{g^2}\hfill & \hfill {\hat{b}}_{b^2}\hfill & \hfill {\hat{b}}_{r^3}\hfill & \hfill {\hat{b}}_{g^3}\hfill & \hfill {\hat{b}}_{b^3}\hfill \end{array}\right]}$

5. To standardize the R, G, and B components of each pixel, i, in the source image, apply each standardization vector to S_i, the linear, quadratic, and cubic RGB components of i. ${\displaystyle S_{i,standardized}=\left(\begin{array}{ccc}\hfill S_i{R}_h\hfill & \hfill S_i{G}_h\hfill & \hfill S_i{B}_h\hfill \end{array}\right)}$ (5)${\displaystyle \text{where}{S}_i=\left(\begin{array}{ccccccccc}\hfill r_i\hfill & \hfill g_i\hfill & \hfill b_i\hfill & \hfill r_i^2\hfill & \hfill g_i^2\hfill & \hfill b_i^2\hfill & \hfill r_i^3\hfill & \hfill g_i^3\hfill & \hfill b_i^3\hfill \end{array}\right).}$

The resulting image will have the same color profile, contrast, and brightness of the target image. To quantify the deviance, denoted as D, between the source and target images, extend T to include quadratic and cubic terms for each R, G, and B exactly as (2) and estimate standardizations for the six new columns in T.

1. Construct matrix H using all nine standardizations vectors. (6)${\displaystyle H_{9x9}=\left[\begin{array}{c}\hfill R_h\hfill \\ \hfill G_h\hfill \\ \hfill B_h\hfill \\ \hfill R_h^2\hfill \\ \hfill G_h^2\hfill \\ \hfill B_h^2\hfill \\ \hfill R_h^3\hfill \\ \hfill G_h^3\hfill \\ \hfill B_h^3\hfill \end{array}\right]=\left[\begin{array}{c}\hfill M{T}_R\hfill \\ \hfill M{T}_G\hfill \\ \hfill M{T}_B\hfill \\ \hfill M{T}_{R^2}\hfill \\ \hfill M{T}_{G^2}\hfill \\ \hfill M{T}_{B^2}\hfill \\ \hfill M{T}_{R^3}\hfill \\ \hfill M{T}_{G^3}\hfill \\ \hfill M{T}_{B^3}\hfill \end{array}\right].}$

2. Compute deviance from reference. (7)${\displaystyle D=1-\text{det}\left(H\right).}$

If the source image has an identical color profile as the target image, H is the identity matrix and the matrix determinant is 1. So that no change in profile is reported as a 0, one minus the determinant is the returned value.

Identifying previously unknown source of variance

The Bellwether Phenotyping Facility at the Donald Danforth Plant Science Center is a controlled-environment growth facility that utilizes a Scanalyzer 3D-HT plant-to-sensor system (LemnaTec GmbH) to move plants through stationary imaging cabinets (Fahlgren, Gehan & Baxter, 2015). The design of the imaging cabinets is intended to produce consistent images with uniform backgrounds. To test the consistency of images collected by the Bellwether system, we imaged sorghum plants for 11 days and included a ColorChecker Passport panel in the field of view so image color profiles could be compared. We designated an image taken at noon on the first day of the experiment to be the reference, and we calculated the deviance (D) for every other image. We observed that there was large variance in D in the image set (median = −0.370, mean = −8.756, variance = 530.991). By manually inspecting images that showed a large and small D, we observed that a large, negative D corresponded to a decrease image brightness. Furthermore, we observed that D had a periodic distribution with the time of day the image was taken. Collapsing all the images taken over the 11-day span to a 24-hour window shows that there is a clear time of day effect (Fig. 1A). Images taken during experimental nighttime hours had decreased brightness compared to images taken during daytime hours. (Wilcoxon rank-sum test, p-value < 0.001). We tested two hypotheses that could account for the observed diurnal pattern: (1) light leakage into the imaging cabinet through its vents during the day, or (2) the temperature of the room was influencing the brightness of the fluorescent bulbs that illuminate the cabinet. To test the former hypothesis, we took two pictures containing a ColorChecker Passport panel: one with the vents open and another with the vents covered. Setting one of images as the reference and standardizing the other to it resulted in D =  − 0.454, which is insufficient to explain the large absolute D values observed in the plant image set. To test the latter hypothesis, a two-day test was done with a 14 hr/10 hr (day/night) photoperiod. For the first 24 hrs the temperature cycled with the photoperiod (32 °C/22 °C, day/night). For the second 24 hrs temperature was held at a constant 32 °C. Five thermometers were placed throughout the imaging cabinet and recorded temperature every minute. Images containing a ColorChecker Passport panel were captured every hour (ZT0-ZT18; to mimic the period of time images were collected in the full experiment) and used to calculate D for each image. Temperatures within the chamber tracked closely with those in the room (Figs. 1B and 1C). We observed that under diurnal conditions, D varies with temperature, and in particular, large absolute values of D were observed for temperatures below 28 °C (Fig. 1D). Additionally, the variation of D was consistent with the variation observed for the plant image set both in frequency and magnitude (Fig. 1). Furthermore, when the ambient temperature remained constant, D did not change (Fig. 1E). The imaging cabinet is illuminated with fluorescent bulbs and it is known that light output from this type of bulb is affected by ambient temperature (Bleeker & Veenstra, 1990). Using the measure of deviance defined in Eq. (7), we were able to detect and identify a source of variance within our plant imaging facility.

Standardization of images improves shape and color measurements

Given that there was variable brightness in the plant image set, the commonly used fixed-threshold segmentation routines produced noisy results. Dark, low-quality images that were segmented using a fixed-threshold routine optimized for bright, high-quality images resulted in either background pixels being classified as object, object pixels being classified as background or both (Figs. 2A–2D). For example, segmentation of an uncorrected sample image take during the experimental night fails to remove a large amount of background near darker parts of the image and salt and pepper noise around the plant. In contrast, segmentation of the transformed image results in clean segmentation of the plant without any additional filtering needed.

To further explore how our standardization method affects phenotype analysis, a set of morphological shapes of each plant was measured before and after standardization. We calculated the difference between the measurements obtained from the original images and the standardized images and grouped the results into experimental day and night (Fig. 2E). The values of the shapes are on vastly different scales so normalization of those values by subtracting the grand mean and dividing by the variance was done for each image and for each shape. The effects in this scale are equivalent to the effects in the native scale for each shape. Many commonly reported measurements (area, convex hull area, width, perimeter, circularity, solidity, center of mass x-coordinate, height, ellipse minor axis, fractal dimension, roundness, eccentricity, center of mass y-coordinate, ellipse angle, convex hull vertices) are significantly different using an unequal variance t-test (p-value < 0.05) between images taken during daytime and nighttime temperatures. Other common measurements such as ellipse major axis, ellipse center x-coordinate, aspect ratio, and ellipse center y-coordinate were not significantly different (Fig. 2E). Plant area (cm²) was previously shown to correlate to fresh weight biomass on the Bellwether platform (Fahlgren, Gehan & Baxter, 2015), and we observed a similar relationship in the current experiment with a small subset of plants (n = 162) that were manually weighed on the last day of the experiment (Fig. 2F). Color profile standardization had no effect on fresh-weight biomass prediction as the slopes and intercepts of the linear models using the unstandardized and standardized datasets were not significantly different (p-value > 0.05, two-way ANOVA). Plant area responds to salt and pepper noise in a one-to-one relationship, so a relatively small amount of noise compared to total plant area will not produce large effects. In contrast, shapes like perimeter are more sensitive and produce larger differences because inclusion of one background pixel can add one or more pixels to the resulting measurement (Figs. 2G and 2H).

Next, we sought to understand how the color profile of the object was affected by our standardization method. Hue is a measure of the color in degrees around a circle of an object in an image, and we previously used hue to measure plant responses to abiotic stress (Moroney et al., 2002; Veley et al., 2017). We compared the hue profiles of plants from original and standardized images and observed an increase in signal above 120° (green) in the nighttime temperature unstandardized images relative to the daytime temperature unstandardized images (Figs. 3A and 3C). Additionally, we observed a large amount of noise below 90° (red-yellow) in the same comparison. Visual inspection of a subset of images suggested that the former is due to misinterpretation of the true object pixels’ shade of green and the latter is due to background being included in the hue histogram. The standardized image sets result in consistent color profiles between the nighttime and daytime images (Figs. 3B and 3D). To statistically test for differences between hue histograms, a two-sample Kolmogorov–Smirnov test was applied to each pairwise comparison. Hypothesis testing was done using the following: (8)${\displaystyle K>c\left(\alpha \right)\sqrt{\frac{n+m}{nm}}}$ (9)${\displaystyle c\left(\alpha \right)=\sqrt{-\frac{1}{2}\text{ln}\left(\frac{\alpha }{2}\right)}}$ (10)${\displaystyle K=\text{sup}\left|\begin{array}{c}\hfill F_x\left(\text{hue}\right)-F_y\left(\text{hue}\right)\hfill \end{array}\right|}$ F is a cumulative distribution of samples x and y, each with sample sizes n and m, respectively. In Eq. (8), if the left-hand side is greater than the right-hand side at α = 0.05, then we reject the null hypothesis that the two histograms are the same distribution. Using this test, we determined that the standardization procedure significantly alters the color profile of the images taken under nighttime temperatures but not the images taken under daytime temperatures, and the standardized image set are no longer significantly different between day and night temperatures (Fig. 3E). This suggests that the standardization method accurately maps the color space from the dark, low-quality images to the same color space as the bright, high-quality images. Changes in hue spectra have been shown to be associated to varying levels of nitrogen in sorghum and rice (Veley et al., 2017; Wang et al., 2014). We hypothesized that variation in color profiles due to differences in light quality would affect our ability to accurately measure nitrogen treatment effects in the present dataset. We used a constrained analysis of principal coordinates conditioned on time, in days, to determine whether treatment groups could be separated by hue values (Figs. 3F–3I). Consistent with previous results, clear separation of the nitrogen treatment groups was observed for both unstandardized and standardized images taken during the daytime (Figs. 3F and 3G; Veley et al., 2017; Wang et al., 2014). In contrast, unstandardized nighttime images appeared to have relatively similar hue spectra across nitrogen treatments, whereas standardized nighttime images show clear separation (Figs. 3H and 3I).