A low-cost greenhouse-based high-throughput phenotyping platform for genetic studies: a case study in maize under inoculation with plant growth-promoting bacteria

Greenhouse-based high-throughput phenotyping (HTP) presents a useful approach for studying novel plant growth-promoting bacteria (PGPB). Despite the potential of this approach to leverage genetic variability for breeding new maize cultivars exhibiting highly stable symbiosis with PGPB, greenhouse-based HTP platforms are not yet widely used because they are highly expensive; hence, it is challenging to perform HTP studies under a limited budget. In this study, we built a low-cost greenhouse-based HTP platform to collect growth-related image-derived phenotypes. We assessed 360 inbred maize lines with or without PGPB inoculation under nitrogen-limited conditions. Plant height, canopy coverage, and canopy volume obtained from photogrammetry were evaluated five times during early maize development. A plant biomass index was constructed as a function of plant height and canopy coverage. Inoculation with PGPB promoted plant growth. Phenotypic correlations between the image-derived phenotypes and manual measurements were at least 0.6. The genomic heritability estimates of the image-derived phenotypes ranged from 0.23 to 0.54. Moderate-to-strong genomic correlations between the plant biomass index and shoot dry mass (0.24–0.47) and between HTP-based plant height and manually measured plant height (0.55–0.68) across the developmental stages showed the utility of our HTP platform. Collectively, our results demonstrate the usefulness of the low-cost HTP platform for large-scale genetic and management studies to capture plant growth. Core ideas A low-cost greenhouse-based HTP platform was developed. Image-derived phenotypes presented moderate to high genomic heritabilities and correlations. Plant growth-promoting bacteria can improve plant resilience under nitrogen-limited conditions.


115
Low-cost high-throughput phenotyping platform 116 A low-cost greenhouse HTP platform was built, wherein the camera was positioned in a 117 way that it obtained images from directly above the plants. The system was built in a 118 conventional greenhouse with dimensions of 3.5 × 11 × 6 m height, width, and length, 119 respectively. A cooling wall and ventilation were used to maintain the desired temperature, 120 and additional luminosity was supplied using LED lamps.

121
The image capture system was inspired by the UAV flight plans. It consists of two 122 fixed parallel tracks (9 m) and one mobile perpendicular track (5 m). They were positioned 123 2.5 m above the ground. The two parallel tracks were fixed to the greenhouse roof, as 124 well as two support tracks to ensure stability and alignment. The parallel tracks move the 125 perpendicular track along the x-axis, whereas the perpendicular track moves the sensors Image processing and data extraction Multispectral images were processed by assembling orthomosaics and the dense point cloud 138 using Agisoft Metashape software (Agisoft LLC, St. Petersburg, Russia). The images were 139 imported, aligned, and optimized using GCP. This was followed by the calculation of the 140 dense point clouds and the stitching of orthomosaics.

141
The orthomosaics were analyzed using QGIS software (QGIS Development Team, 2021) 142 to obtain a shapefile for each plot. The plots were manually identified, and a geometry point

154
Dense cloud points were used to estimate plant height (PH HTP ) and canopy volume 155 (CV). Each point from the dense cloud point was composed of GPS coordinates (latitude, 156 longitude, and altitude in the universal transverse mercator). The dense cloud point data from the product of PH HTP (cm) and CC (cm 2 ) (Li et al., 2020). For CV, the dense cloud 163 points were filtered by colors using the "Select Points by Color" function in the Agisoft 164 Metashape software to remove the background. Plants were then reconstructed from the 165 point cloud data, and the CV was estimated using the α-shape algorithm (Lafarge and   166 Pateiro-Lopez, 2020). The algorithm requires an α value that controls the tightness of the 167 3-D reconstruction of the points. The optimal value of α that yielded the greatest correlation 168 with manual measurements was 0.01 (Moreno et al., 2020  A total of 13,826 single-nucleotide polymorphisms (SNPs) were available for the maize 189 inbred lines using a genotyping-by-sequencing method following the two-enzymes (PstI and 190 MseI) protocol (Sim et al., 2012;Poland et al., 2012). DNA was extracted using the 191 cetyltrimethylammonium bromide method (Doyle and Doyle, 1987). SNP calling was per-192 formed using the TASSEL 5.0 software (Bradbury et al., 2007) with B73 (B73-RefGen v4) 193 as the reference genome. The SNP markers were filtered if the call rate was less than 90%, 194 non-biallelic, and the minor allele frequency was less than 5%. Missing marker codes were  The image-derived phenotypes were collected over time to capture plant growth, as pre-206 viously described. For each replication, measurements were made at six time points defined 207 by the number of expanded leaves: 0 (before germination), 2, 3, 4, 5, and 6 (Hanway, 1966).

208
Since the genotypes presented expected inconsistencies in growth stages, it was determined 209 as the mode of the population at a given time. A time point before the germination step 210 was used to obtain the PH HTP . Heat accumulation was calculated from the growing degree and CV) and manually measured phenotypes (PH and SDM).

218
Likelihood-ratio and Wald tests 219 The following model was used to test the effects of genotype, management (B+ and B-), and their interaction.
where y is the vector of phenotypes; 1 is the vector of ones; X 1 , X 2 , and X 3 are the incidence 220 matrices for the fixed effects; Z 1 and Z 2 are the incidence matrices for the random effects; µ 221 is the overall mean; r, b, and m are the fixed effects for replication, block within replication, 222 and management (B+ and B-), respectively; g ∼ N (0, Gσ 2 g ) is the vector of random effect 223 of genotype; gm ∼ N (0, G ⊗ Iσ 2 gm ) is the vector of random effects of the interaction between 224 genotype and management; and ∼ N (0, Iσ 2 ) is the random residual effect. Here G is the 225 additive genomic relationship matrix (VanRaden, 2008); I is the identity matrix; σ 2 g is the 226 additive genomic variance; σ 2 gm is the genotype-by-management interaction variance; and 227 σ 2 is the residual variance. The significance of random and fixed effects was assessed using 228 the Wald and likelihood-ratio tests, respectively. The analysis was performed using the R 229 package ASReml-R (Butler et al., 2017).

230
Bayesian genomic best linear unbiased prediction 231 Univariate and bivariate Bayesian genomic best linear unbiased prediction (GBLUP) models were used to estimate genomic heritability and genomic correlation separately for B+ and B-. These Bayesian models were the same as those used for the Wald and likelihood-ratio tests, but the management (m) and genotype-by-management interaction terms (gm) were dropped. For the univariate model, a flat prior was assigned to r and b. The variance components, σ 2 g and σ 2 , were drawn from a scaled inverse χ 2 distribution. For the bivariate model, y is the vector of phenotypes of two responses; g ∼ N (0, Σ g ⊗ G) is the vector of genotypes; ∼ N (0, Σ ⊗ I) is the residual; ⊗ is the Kronecker product; and Σ g and Σ are the variance-covariance matrices for additive genomic and residual effects taking the forms where subscripts 1 and 2 refer to the first and second phenotypes. An inverse Wishart 232 distribution was assigned to Σ g and Σ with degrees of freedom ν = 4 and scale matrix 233 S such that the prior means of Σ g and Σ equal half of the phenotypic variance. All the 234 Bayesian GBLUP models were fitted using 60,000 Markov chain Monte Carlo samples, 10,000 235 burn-in, and a thinning rate of 60 implemented in JWAS software (Cheng et al., 2018a,b).

236
Model convergence was assessed using trace plots of the posterior distributions of the variance 237 components.

238
Heritability and genomic correlation 239 The variance components obtained from the univariate Bayesian GBLUP were used to estimate genomic heritability using the following formula.
where n r is the number of replications (2) Table S5). Figure 3 shows the growth patterns of the image- correlations were observed between the HTP and manually measured phenotypes (Table 2).  ligule in the HTP platform is a challenging task because PH HTP is based on plant height 353 projection ( Figure S1).

354
There are several greenhouse-based HTP platforms available that differ in terms of pre-