Host-mediated selection experimental design. (A) Representative of a single selection line of 15 rice plants. 1) Enrichment Generation (EG): rice are transplanted into 1:1 mixtures of field soil and sterilized calcined clay grow for 40 days in well-watered conditions. 2) Roots and adhering soil are harvested, minced, and ground, then mixed into sterilized calcined clay to inoculate the first selection generation (SG1). 3) SG plants grow for 30 days in well-watered conditions. 4) Drought treatment is initiated; plants are imaged daily to track phenotypes for 10 days. 5) The three best-performing plants are identified using Normalized Difference Vegetative Index (NDVI) and biomass data. 6) Roots and soils from best-performing plants are harvested and diluted with sterilized calcined clay to inoculate the subsequent selection generation (SG2). Steps 3-6 are then repeated iteratively for additional selection generations. (B) An overview of the type and number of selection and control lines for an individual soil treatment. Selection Lines (SL 1-6) consist of 15 plants per generation each and were maintained as independent lines across generations (i.e. transfers of selected inocula occurred within SLs). Sterile- and Live-Inoculated control lines (SI and LI), on the other hand, each included 30 plants per generation and were inoculated by an amalgam of material from all selection lines within their soil treatment; SI inocula were autoclave-sterilized prior to each generation. Like their SL counterparts, SI plants were well-watered throughout the enrichment generation and subjected to drought in selection generations; LI plants were well-watered in all generations. (C) Rice Field (green), Desert (red), and Serpentine Seep (blue) soil treatments—inoculated by field soils—each included six selection lines (SL 1-6) of 15 plants per generation. Calcined Clay (yellow) did not receive any inoculum in EG, instead assembling microbiomes composed of environmental taxa, and included just three SLs of 15 plants per generation. The number and length of columns are proportional to the number of lines and selection generations, respectively, associated with each soil treatment.

Field sites, soils, and starting diversity of source inocula. (A) Collection sites for Rice Field, Desert, and Serpentine Seep. Calcined Clay “source inocula” was simply calcined clay (Profile Greens Grade, Cat.No. FZ-TP) that was rinsed and autoclave-sterilized before use. (B) Field soils differed in physico-chemical characteristics (Supplemental Table 1). Rice Field was visibly clay-dominant, Desert predominantly sandy, and Serpentine Seep high in organic matter, but less clay-like compared to Rice Field. (C) Bacterial diversity and community compositions were distinct for each soil, even at broad taxonomic levels (Supplemental Figure 1).

Rice phenotypes across host-mediated selection generations. (A) Representative images of Rice Field selection line plants in the first (SG1) and last (SG6) selection generations at 1, 5, and 10 days after drought. NDVI values are shown in the bottom left corner of each image. (B) Drought tolerance, quantified as the area under the curve of NDVI values (AUC NDVI), significantly improved for all soil treatments except Serpentine Seep. Here, we show data for all selection lines combined (gray bars) and assessed the response of AUC NDVI to generation via linear mixed effects models (formula = AUC NDVI ∼ generation + (1|selection line); Rice Field: β=0.17, t = 28.8, p < 2e-16; Desert: β=0.19, t = 14.1, p < 2e-16; Serpentine Seep: β=-0.003, t = - 0.3, p = NS; Calcined Clay: β=0.17, t = 8.7, p < 0.001). Filled bars represent the mean AUC NDVI value of all selection line (SL) plants +/- one standard error; white, unfilled bars are data for sterile-inoculated (SI) plants. Because selection lines were independent of one another, we also tested them individually via OLS regression: all selection lines of Rice Field (6/6), Desert (6/6), and Calcined Clay (3/3) significantly improved over time; Serpentine Seep included selection lines that both significantly improved (2/6) and worsened (2/6) (Supplemental Figure 3). (C) Shoot dry weight biomass also significantly increased across selection generations for all soil treatments, except Serpentine Seep (Rice Field: β=0.02, t =7.0, p < 2e-16; Desert: β=0.06, t =11.6, p < 2e-16; Serpentine Seep: β=0.006, t =1.2, p = NS; Calcined Clay: β=0.05, t =10.3, p < 2e-16). Tested individually, all selection lines of Rice Field (6/6) and Calcined Clay (3/3) significantly improved over time; most Desert selection lines (5/6) significantly improved; Serpentine Seep included one selection line that significantly increased (1/6) and one that significantly decreased (1/6) (Supplemental Figure 3).

Alpha and beta diversity of bacterial communities over time. (A) Shannon-Index diversity and (B) number of observed ASVs for soil treatments over time. Each point represents a bacterial community from one of several sample types, including original field soils (FS), enrichment generation (EG) inocula (FS::EG), root-associated selection line (SL), or root associated sterile-inoculated (SI); sample types are indicated by shape. EG inocula were simply field soil diluted by sterilized calcined clay. Sample sizes and numbers of selection generations differ between soil treatments; having observed increasing similarities in microbiome composition between soil treatments after four selection generations, we decided to focus effort and resources on Rice Field, which we carried forward an additional two selection generations. (C) Bray-Curtis (BC) dissimilarity between soil treatment selection lines plotted in two-dimensional space. Each point represents an individual sample, where color indicates soil treatment and size indicates generation. BC values are the average of 1000 permutations where each sample was individually subset to 5000 reads. (D-G) Similar to (C), but panels are individual soil treatments and include both SL (filled) and SI (unfilled) samples. Despite forming distinct clusters, SL and SI microbiomes generally follow parallel compositional trajectories through time. For ordinations in panels C and D, Rice Field samples were subset to match the sequencing effort of other soil treatments prior to calculating BC values.

Dispersal homogenizes microbiomes over time. Using both original field soils (FS) and enrichment generation inocula (FS::EG, field soil mixed with sterile calcined clay), we assigned ASVs as originally belonging to one or multiple soil origins and is indicated by color. For example, if an ASV was observed in field soil or inocula samples for Rice Field, but not Desert or Serpentine Seep, it was designated a Rice Field ASV; if observed in Rice Field and Serpentine Seep, but not Desert, it was designated a “Rice Field + Serpentine Seep” ASV. ASVs that were not observed in either the field soil samples or source inocula were assigned “unk/env” to indicate their origin as unknown, but likely from the growth chamber environment. Notably, sterile-inoculated (SI) lines of each soil treatment—which acquired microbiomes from the environment—were each biased towards their respective soil treatments. Calcined Clay, which received no initial inoculum in the enrichment generation was more evenly represented by ASVs from other soil treatments.

Microbiomes become overrepresented by fast-growing taxa over time. We determined ASV doubling times by matching ASVs to the Estimated Growth rates from gRodon Online (EGGO) database. We then performed 1000 bootstraps randomly subsetting 100 ASVs that were present in each generation for selection line (SL) and sterile-inoculated (SI) microbiomes. Each point represents the mean predicted doubling time of the ASV subset in each bootstrap. By setting generation as a numeric factor and performing OLS regression, we found that predicted doubling time significantly decreased across selection generations in both SL (t = -90.9, p<2e-16) and SI (t = -19.6, p<2e-16) samples.

Identifying phenotype-associated ASVs. (A) A phylogeny of Rice Field ASVs remaining after four selection generations with branches colored by phylum; black bars and corresponding labels highlight notable clades at lower taxonomic levels. The inner (AUC NDVI) and outer (Shoot biomass) tiled rings indicate whether an ASV showed a strong positive (red) or negative (blue) response to phenotype, or lack thereof (white). To determine whether a response was strong or not, we ran two mixed effects linear models for each ASV-phenotype combination. One model included selection line as a random factor, while the other included selection generation as a random factor nested within selection line; doing so allowed us to better hone in on ASVs that were both enriched over time and consistently correlated with phenotypes within selection generations. Responses were considered strong if the t-value resulting from one or both models was two or more standard deviations away from the mean. Gray bars extending outward from the tiled rings indicate the median relative abundance of ASVs in the final selection generation (SG6). Symbols beyond the gray bars denote whether ASVs were included in the numerator (red) or denominator (blue) of the balances shown in panel C. (B) The mean, combined relative abundance of drought enriched ASVs over time; error bars (present, but small) are +/- one standard error of the mean. Drought enriched ASVs were identified using the effect and t-test functions in the ALDEx2 R package by comparing CLR-transformed abundances between SL (droughted) and LI (not droughted) plants at each selection generation. (C) Balances successfully predict NDVI and biomass phenotypes. In brief, balances are compositionally-aware log contrasts in which the numerator and denominator represent two sets of ASVs whose ratio can be used to predict a continuous phenotype (Rivera-Pinto et al. 2018). We composed distinct balances for each phenotype and performed OLS regressions on individual lines; in each case, SI lines included, balances were significant predictors of plant phenotype. Statistics and black trendlines correspond to the overall regression; colored trendlines show the relationship between balance values and phenotype within generations.

Bacterial communities associated with increased plant biomass and plant vigor during drought are enriched for both common and distinct functions. By comparing the metagenome-assembled genomes (MAGs) of taxa from balance group numerators with all other MAGs in the dataset, we identified high level functional categories that were significantly enriched in the genomes of these groups. Selection/Background refers to the ratio of gene counts among phenotype-associated taxa relative to all taxa for each functional category. The size of points indicates the background gene count, or the total number of genes within a functional category represented by all MAGs. We determined p-values with hypergeometric tests and applied multiple testing correction before determining significance. Functional categories that were significantly enriched in balance numerator groups were then further investigated to identify significantly enriched individual COGs (Supplemental Figures 10-11)

Bacterial diversity of field soils collected as candidate source inocula for host-mediated selection. (A) Relative abundance barplots colored by Phylum show large differences between soil bacterial communities even at coarse taxonomic levels. Each bar represents an individual technical replicate. (B) Shannon-Index diversity for each soil. The mean number of observed ASVs are shown above each boxplot. (C) NMDS-ordinated Bray-Curtis dissimilarity between soils. (D-F) Phenotypic data from a pilot study where rice plants were either inoculated with the field soils shown in A-C or mock inoculated.

(A) Simplified schematic of lightbox design depicting full spectrum (white) and 730nm (red) LED light strips as well as the position of the standard (white) and modified (red) cameras. The interior of the box was painted with an ultra light-absorbing black paint. A shared shutter trigger ensured images from both cameras were captured simultaneously; up to five plants were imaged at the same time. (B) A representative drawing of the portions of visible and near-infrared light each camera’s red, green, and blue channels were sensitive to. Beneath, black bars indicate the range of emission spectra of light sources as well as wavelengths expected to be absorbed (ABS) and fluoresced (FLO) by chlorophyll. (C) Comparisons of the standard and modified cameras’ median values for red, green, and blue channels individually; each dot represents an individual plant. (D) NDVI responds to rice drought stress status over time. Generally, the interval between two time points with the largest difference corresponds to the visible onset of leaf curling symptoms (data not shown). (E) At harvest, NDVI is highly predictive of rice percent water content.

AUC NDVI and biomass phenotypes for individual selection lines. (A) AUC NDVI scores over time. Statistics and trendlines result from OLS regression. Each point represents an individual plant. Filled circles correspond to selection line (SL) plants; open circles correspond to sterile-inoculated (SI) plants. (B) Shoot dry weight biomass over time.

Rice growth conditions. (A) Rice were grown in spatially discrete conetainers and water reservoirs to minimize admixture of microbiomes. A small ball of polyester fiber was placed at the bottom of conetainers to prevent substrate from spilling out. Reservoirs were filled with milliQ water which calcined clay substrate was able to wick to the top surface of the container. (B) Host-mediated selection was performed in a walk-in growth chamber set to 16-hour days, during which daytime temperature and humidity were held at 28°C and 80%, respectively, and relaxed to 24°C and 75% humidity at night. Each soil treatment was grown on its own shelving unit, but plants were rotated within that unit every three days to ensure each plant experienced the same overall environment and ballast exposure. The lightbox with a rack of plants ready to be imaged can be seen in the bottom right corner.

Bray-Curtis dissimilarity (A-B) and Aitchison distance (C-D) between and within soil treatments across selection generations. Panels A and C show dissimilarity/distances between soil treatment selection lines by BC dissimilarity and Aitchison distance, respectively. Panels B and E are similar, except dissimilarity/distance is between selection line and steril-inoculated microbiomes within soil treatments. BC values are the average of 1000 permutations where each sample was individually subset to 5000 reads. Aitchison distance was calculated as the Euclidean distance between samples after read counts were transformed via robust centered log ratio transformation. Letters above boxplots are the result of a post-hoc Tukey HSD with 95% confidence intervals (p<0.05) to detect significant differences between groups. Rice Field samples were subset to match the sequencing effort of other soil treatments when generating data for panels A and C.

The number and abundance of ASVs shared between soil treatments increased over time. (A) Venn diagrams for the numbers of ASVs shared between soil treatments at each generation. Venn components are colored by soil treatment: Rice Field (green), Desert (red), Serpentine Seep (blue), and Calcined Clay (yellow). (B) The combined relative abundance of the top 20 ASVs commonly enriched in all soil treatments over time. Each point represents an individual sample. Patterns of enrichment were similar for both selection line (filled circles) and sterile-inoculated (open circles) microbiomes.

Neutral community models for Rice Field selection line microbiomes at each selection generation. (A) We fit ASV relative abundance data to neutral models using the previous generation’s community to parameterize migration and replacement rates. Gray points fall within the 95% confidence interval for abundance-occupancy values (x and y axes, respectively) expected if the community were assembled via neutral processes alone. Blue and red points indicate ASVs outside this prediction and therefore more likely experiencing positive and purifying selection, respectively. R-squared values for each plot indicate the overall model fit; blue and red bars in the bottom right corner of each plot show the total number of ASVs above and below neutral model predictions. (B) Given the improved fit of neutral models in later selection generations, neutral processes appear to become more important to microbiome assembly over time. Additionally, neutral models consistently outperformed binomial ones, which produced similar abundance-occupancy curves but without parameterization.

Most abundant ASVs within Rice Field microbiomes. One ASV in particular—identified as belonging to the genus Ideonella—rose to high abundance in Rice Field microbiomes, sometimes accounting for more than 50% of the relative abundance of the community. Often, the combined relative abundance of next top ten most abundant taxa was still less than that of Ideonella sp. Regardless, both showed patterns of enrichment over time in all selection lines as well as sterile- and live-inoculated microbiomes. Bars represent mean relative abundance +/- one standard error.

ASVs identifies as significantly and consistently enriched in either selection line (drought) or live-inoculated (water) microbiomes. After filtering for FDR-corrected significant enrichment and an effect size greater than 0.5 (indicating significant and strong enrichment signal) we further subset to ASVs meeting these criteria in at least three selection generations. In total, 30 drought and 30 water indicators remained.

Enrichments of individual COGs within significantly enriched categories associated with AUC NDVI balance numerator taxa. As in Figure 7, Selection/Background refers to the ratio of gene counts among phenotype-associated taxa relative to all taxa, but for individual COGs rather than COG functional categories. All COGs shown were significantly enriched in balance taxa, even after applying multiple testing corrections scaled by the total number of COGs within each functional category. Each column in the heatmap represents an individual MAG; columns are grouped by phylum. Numbers within heatmap cells indicate the COG gene count for each MAG.

Enrichments of individual COGs within significantly enriched categories associated with shoot biomass balance numerator taxa. As in Figure 7, Selection/Background refers to the ratio of gene counts among phenotype-associated taxa relative to all taxa, but for individual COGs rather than COG functional categories. All COGs shown were significantly enriched in balance taxa, even after applying multiple testing corrections scaled by the total number of COGs within each functional category. Each column in the heatmap represents an individual MAG; columns are grouped by phylum. Numbers within heatmap cells indicate the COG gene count for each MAG.

NDVI values for plants at 10 days after drought (DAD 10) decrease as rice shoot dry weight biomass increases. Data shown here are for Rice Field selection line plants only, but this pattern held true for other soil treatments and SI plants as well.