Defining Composition and Function of the Rhizosphere Microbiota of Barley Genotypes Exposed to Growth-Limiting Nitrogen Supplies

ABSTRACT The microbiota populating the rhizosphere, the interface between roots and soil, can modulate plant growth, development, and health. These microbial communities are not stochastically assembled from the surrounding soil, but their composition and putative function are controlled, at least partially, by the host plant. Here, we use the staple cereal barley as a model to gain novel insights into the impact of differential applications of nitrogen, a rate-limiting step for global crop production, on the host genetic control of the rhizosphere microbiota. Using a high-throughput amplicon sequencing survey, we determined that nitrogen availability for plant uptake is a factor promoting the selective enrichment of individual taxa in the rhizosphere of wild and domesticated barley genotypes. Shotgun sequencing and metagenome-assembled genomes revealed that this taxonomic diversification is mirrored by a functional specialization, manifested by the differential enrichment of multiple Gene Ontology terms, of the microbiota of plants exposed to nitrogen conditions limiting barley growth. Finally, a plant soil feedback experiment revealed that host control of the barley microbiota underpins the assembly of a phylogenetically diverse group of bacteria putatively required to sustain plant performance under nitrogen-limiting supplies. Taken together, our observations indicate that under nitrogen conditions limiting plant growth, host-microbe and microbe-microbe interactions fine-tune the host genetic selection of the barley microbiota at both taxonomic and functional levels. The disruption of these recruitment cues negatively impacts plant growth. IMPORTANCE The microbiota inhabiting the rhizosphere, the thin layer of soil surrounding plant roots, can promote the growth, development, and health of their host plants. Previous research indicated that differences in the genetic composition of the host plant coincide with variations in the composition of the rhizosphere microbiota. This is particularly evident when looking at the microbiota associated with input-demanding modern cultivated varieties and their wild relatives, which have evolved under marginal conditions. However, the functional significance of these differences remains to be fully elucidated. We investigated the rhizosphere microbiota of wild and cultivated genotypes of the global crop barley and determined that nutrient conditions limiting plant growth amplify the host control on microbes at the root-soil interface. This is reflected in a plant- and genotype-dependent functional specialization of the rhizosphere microbiota, which appears to be required for optimal plant growth. These findings provide novel insights into the significance of the rhizosphere microbiota for plant growth and sustainable agriculture.

limiting plant growth amplify the host control on microbes at the root-soil interface. This is reflected in a plant-and genotype-dependent functional specialization of the rhizosphere microbiota, which appears to be required for optimal plant growth. These findings provide novel insights into the significance of the rhizosphere microbiota for plant growth and sustainable agriculture. KEYWORDS barley, metagenomics, nitrogen, rhizosphere-inhabiting microbes T o sustainably enhance global food security, innovative strategies to increase crop production while preserving natural resources are required (1)(2)(3). Capitalizing on the microbial communities thriving in association with plants, collectively referred to as the plant microbiota (4,5), has been identified as one of these innovative strategies (6). For instance, members of the microbiota populating the rhizosphere, the interface between roots and soil, can provide their plant host with access to mineral nutrients and protection against abiotic and biotic stresses (7). Thus, applications of the plant microbiota have the potential to integrate and progressively replace nonrenewable inputs in crop production (8).
This potential is of particular interest for alternatives to nitrogen (N) applications to staple crops, such as cereals, as approximately 50% of applied fertilizers are lost either to the atmosphere or in groundwater (9,10), largely as a consequence of microbial denitrification and nitrification processes. Soil microbes can contribute to the release of nitrogen from soil organic matter (SOM) for plant uptake (11). These mineralization processes are estimated to contribute more than 50% of crop nitrogen (12), even in intensively fertilized systems, and are the fundamental basis of sustained plant productivity in uncultivated soils, as typically more than 90% of soil N is present in organic forms (13). The importance of the plant in influencing these microbial mineralization processes has been highlighted by the phenomenon of rhizosphere priming effects (14), where root release of organic compounds impacts rates of SOM decomposition and nitrogen mobilization (15). Therefore, elucidating the relationships between rhizosphere microbiota composition and nitrogen availability for plant uptake can be a key toward sustainable crop production (16).
The composition of the rhizosphere microbiota is driven, at least in part, by the genetics of its host plants (4,17). In turn, the processes of domestication and breeding selection, which progressively differentiated wild ancestors from modern, "elite" cultivated varieties (18), modulated the plant's capacity of shaping the microbiota thriving at the root-soil interface (19,20). As crop wild relatives have evolved in marginal soils (i.e., not exposed to synthetic fertilizers), their microbiota may be equipped with beneficial functions for sustainable agriculture (7,21). Despite that the impact of plant domestication on the rhizosphere microbiota has been studied in multiple plant species (22)(23)(24)(25)(26)(27), the significance of microbial diversification between wild and cultivated plant genotypes of the same species remains to be fully elucidated (21,28).
Barley (Hordeum vulgare), the fourth most cultivated cereal worldwide (29), represents an attractive model to investigate host genetic control of the rhizosphere microbiota within a framework of plant domestication. For instance, we previously demonstrated that domesticated (H. vulgare subsp. vulgare) and wild (H. vulgare subsp. spontaneum) barley genotypes host microbiotas of contrasting composition (30,31). More recently, we gathered novel insights into the genetic basis of this host-mediated microbiota diversification (32)(33)(34). In parallel, investigations targeting specific microbial genes indicated that barley plants may exert a control on microbes underpinning the nitrogen biogeochemical cycle (35) and that this effect is dependent on community composition (36). However, it is unclear how genetic differences between wild and domesticated genotypes may impact the composition and function of the rhizosphere microbiota of plants exposed to contrasting nitrogen supplies, in particular the ones limiting plant growth.
To address this knowledge gap, in this investigation, we used barley as an experimental model and state-of-the art sequencing approaches to test three interconnected hypotheses. First, we hypothesized that the host control on rhizosphere bacteria is modulated by, and responds to, nitrogen availability for plant uptake. Specifically, we anticipated that differences in microbiota composition among barley genotypes are maximal under limiting nitrogen supplies, when plants rely on their microbiota for Ncycling processes to support optimal growth. We further hypothesized that, under conditions limiting barley growth, the plant's reliance on the rhizosphere microbes will be manifested by a functional diversification mediated, at least partially, by the host genotype. Finally, we hypothesized that that these distinct structural and functional configurations of the microbiota contributed to differential plant growth responses.

RESULTS
Nitrogen conditions limiting plant growth amplify the host effect on the barley rhizosphere microbiota. To gain insights into the role played by nitrogen availability for plant uptake on the composition of the barley bacterial microbiota, we selected one reference barley cultivar, Morex (here "Elite") and two wild genotypes from the B1K collection (37), B1K-12 and B1K-31 (here "Desert" and "North", respectively). The rationale for this choice was two pronged. First, we previously characterized these genotypes for their capacity to recruit distinct microbiotas and genetic relatedness (30,31). Second, the wild genotypes are representative of the two main barley ecotypes identified in the Southern Levant as drivers of plant adaptation to the environment (38,39). Consequently, and despite the limited number, these genotypes may capture the "extremes" of the evolutionary pressure on the host recruitment cues of the barley microbiota. Plants were grown under glasshouse conditions in an agricultural soil previously used for microbiota investigations and designated "Quarryfield" (31,33,34,40). Pots containing the individual genotypes and unplanted soil controls (here "Bulk") were supplemented with three modified Hoagland's solution preparations (41) containing all essential macronutrients and micronutrients and three levels of mineral nitrogen (Table S1 in the supplemental material), the optimum required for barley growth (N100%), a quarter dose (N25%), or no nitrogen (N0%). At early stem elongation (Fig. S1), which represents the onset of maximum nitrogen uptake for small grain cereals (42), plants were harvested, and total DNA preparations were obtained from rhizosphere and unplanted soil specimens. In parallel, we determined aboveground plant biomass, plant nitrogen content in leaves, and concentrations of ammonium (NH 4 1 ) and nitrate (NO 3 -) in rhizosphere and unplanted soil samples.
We observed that plant performance was affected by the N application; aboveground biomass and plant nitrogen content were significantly lower at N0% than at N100%, with N25% yielding intermediate values (Kruskal-Wallis test followed by Dunn post hoc test, individual P values of ,0.05, false-discovery rate [FDR] corrected; Fig. 1), compatible with a nitrogen deficiency status for barley growth. Likewise, the residual nitrogen in the rhizosphere at the completion of the experiments, measured as a concentration of ammonium and nitrate, respectively, displayed a significant decrease in the values recorded for N100% to N25% and from the latter to N0%. (Kruskal-Wallis test followed by Dunn post hoc test, individual P values of ,0.05, FDR corrected; Fig. 1).
In parallel, we generated a 16S rRNA gene amplicon sequencing library from the obtained rhizosphere and unplanted soil controls and identified 26,411 individual amplicon sequencing variants (ASVs) accruing from 6,097,808 sequencing reads. After pruning in silico ASVs representing either host (i.e., plastid-or mitochondrial-derived sequences) or environmental contaminations and low-count ASVs, 5,035,790 reads were retained, representing over 82% of the initial data set. Canonical analysis of principal coordinates (CAP) differentiated bulk soil from rhizosphere profiles, as evidenced by a segregation of either class of samples along the axis accounting for the largest source of variation ( Fig. 2A). Furthermore, we observed a "gradient" along the axis accounting for the second source of variation aligned with the treatment effect, in particular for rhizosphere samples ( Fig. 2A). The sample effect (i.e., either bulk soils or the rhizosphere of the individual genotypes) exerted the primary impact on the bacterial communities thriving at the root-soil interface (permutational multivariate analysis of  Fig. 2A) followed by the nitrogen treatment effect (PERMANOVA, R 2 = 0.105, P = 0.0004, 5,000 permutations; Fig. 2A) and their interaction term (PERMANOVA, R 2 = 0.098, P = 0.0380, 5,000 permutations; Fig. 2A).
To further examine the impact of treatment on the abundance of individual prokaryotic ASVs underpinning host-mediated diversification, we performed a set of pairwise comparisons between barley genotypes at the different N levels. We observed that N0% was associated with the largest number of differentially recruited ASVs, while higher N levels progressively obliterated recruitment differences among genotypes (Wald test, individual P values of ,0.05, FDR corrected; Fig. 2B). Of note, and congruent with previous experiments conducted in the same soil (31), the pair Elite-Desert yielded the highest number of differentially recruited ASVs (Fig. 2B).
Taken together, these observations indicate that nitrogen availability for plant uptake is a factor in (i) modulating the microhabitat-and genotype-dependent recruitment cues of the barley bacterial microbiota by (ii) promoting the selective enrichment of individual taxa in the rhizosphere and (iii) whose magnitude is maximized when no nitrogen is applied to the system.
The metabolic potential of the rhizosphere microbiota exposed to nitrogen conditions limiting barley growth. We generated over 412 million paired-end metagenomic reads from 12 additional samples to gain insights into the functional significance of microbiota diversification in plants exposed to nitrogen conditions limiting barley growth. These represented three biological replicates each of Bulk soil and the rhizospheres of Elite, North, and Desert exposed to the N0% treatment. Upon in silico removal of low-quality sequences and sequences matching the barley genome, likely representing "host contaminations" (Fig. S2), taxonomic classification of the sequencing reads at kingdom level revealed that Bacteria outnumbered Fungi by 2 orders of magnitude, regardless of the sample investigated (Fig. 3A). Closer inspection of the data classified within the kingdom Fungi revealed no significant differences among samples for sequences assigned to the class Glomeromycetes, which we used as a proxy   Fig. 3B). Although the separation between replicates of the same genotype, in particular the Elite-Desert pair, was manifested exclusively when looking at bacteria, we identified a comparable effect of the sample type on composition of both bacterial and fungal communities. For instance, the R 2 value computed for normalized relative abundances returned values between 0.66 and 0.68 for the bacterial and fungal component, respectively (PERMANOVA, 5,000 permutations, individual P values of ,0.01; Fig. 3C and D). Next, we mined the metagenomic data set for sequencing reads associated with known genes underpinning the nitrogen biogeochemical cycle. We were able to identify genes implicated in processes as diverse as nitrification, denitrification, and nitrate reduction and synthesis and degradation of nitrogen-containing organic compounds, Individual shapes in the plot denote individual biological replicates whose color and shape depict sample type and nitrogen treatment, respectively, as indicated in the bottom part of the figure. Numbers in the plots depict the proportion of variance (R 2 ) explained by the factors "Sample," "Treatment," and their interactions, respectively. Asterisks associated with the R 2 value denote its significance; Sample, P = 0.0002; Treatment, P = 0.0004; Sample * Treatment, P = 0.0380; Adonis test, 5,000 permutations. (B) Horizontal blue bars denote the number of ASVs differentially enriched (Wald test, individual P values of ,0.05, FDR corrected) between the Elite and two wild barley genotypes at different nitrogen treatments as recapitulated by the shape and color scheme. Vertical bars depict the number of differentially enriched ASVs unique for or shared among two or more pairwise comparisons highlighted by the interconnected dots underneath the vertical bars. although the abundances of genes associated with the individual process did not discriminate between barley genotypes (Wald test, individual P values of .0.05, FDR corrected; Fig. S3). This suggests that, under the tested conditions, host control of the nitrogen biogeochemical cycle does not represent the main driver of the functional diversification of the barley rhizosphere microbiota. This motivated us to further discern the metabolic capacity of barley-associated communities by assembling metagenomic reads and predicting their encoded proteins (Table 1). Predicted proteins were clustered, resulting in 10,554,104 representative sequences. The representative protein sequences were subjected to functional enrichment analysis to identify Gene Ontology (GO) categories differentially enriched in the barley rhizosphere. We observed a consistent "rhizosphere effect" in the functional potential of the barley microbiota manifested by a spatial separation of plant-associated communities from Bulk soil in an ordination (Fig. 4A) sustained by a differential enrichment of multiple GO categories (Fig. 4B). Closer inspection of these categories revealed a significant enrichment of multiple GO terms in each of the tested genotypes  . In C and D, cluster dendrograms constructed using Bray-Curtis dissimilarity matrices of the metagenomic sequencing reads (counts per million) assigned to the family level in bacteria and fungi, respectively. Individual shapes denote individual biological replicates whose color depicts sample type as indicated in the bottom part of the figure. Numbers associated with each dendrogram depict the proportion of variance (R 2 ) explained by the factor "Sample" in bacteria or fungi, respectively. Asterisks associated with the R 2 value denote its significance; Sample bacteria, P = 0.0012; Sample fungi, P = 0.0004; Adonis test, 5,000 permutations.  Table 2). In particular, the microbiota associated with Desert, North, and Elite genotypes was enriched for GO terms implicated in carbohydrate metabolic processing, cell adhesion, pathogenesis, response to abiotic stimulus, responses to chemical, protein-containing complex assembly, and bacterial-type flagellum-dependent cell motility. These enrichments appear congruent with the adaptation of polymicrobial communities to a host capable of providing substrates for microbial growth. Conversely, Bulk soil specimens were enriched predominantly for functions implicated in photosynthesis and sporulation, which are congruent with microbial adaptation to a lack of organic resources, such as the case in unplanted soils.
To gain a finer view of the functional diversification of the Bulk and rhizosphere microbiotas, we performed a cluster analysis of individual GO terms on the top 10 clusters differentiating between samples (Fig. S4). For each cluster, we determined the significance of individual terms in pairwise comparisons between Bulk soil and rhizosphere samples and, within the latter, between genotypes (Wald test, individual P values of ,0.05, FDR corrected, Data Set S1 at https://doi.org/10.5281/zenodo.7119900). This allowed us to implicate nitrate transporters with functions putatively underpinning multitrophic interactions, such as response to reactive oxygen species and the type VI secretion system. These two functions were also significantly enriched in and differentiated between Elite and Desert communities (Wald test, individual P values of ,0.05, FDR corrected, Data Set S1 at https://doi.org/10.5281/zenodo.7119900, cluster 6). Conversely, ammonium transporters were identified as a depleted function in rhizosphere communities (Wald test, individual P values of ,0.05, FDR corrected, Data Set S1 at https://doi.org/10.5281/zenodo.7119900, clusters 5 and 8) as were functions implicated in phosphate metabolism, including "cellular  phosphate homeostasis", "negative regulation of phosphate metabolic process," and "phosphate ion transmembrane transport" (Wald test, individual P values of ,0.05, FDR corrected, Data Set S1 at https://doi.org/10.5281/zenodo.7119900, clusters 5 and 8). The overarching picture emerging in this investigation was that, at the metagenomic resolution we obtained, the major effect on the functional potential of the microbiota is exerted by the microhabitat (i.e., Bulk versus rhizosphere). Conversely, the effect of the host genotype appears confined to a limited number of individual GO terms and, congruent with the 16S rRNA gene survey, manifested predominantly in the comparison between Elite and Desert genotypes. Genome reconstruction of the bacteria populating the barley root-soil interface. As a first step toward linking structural and functional diversification of the barley microbiota, we attempted to reconstruct genomes of individual bacteria proliferating at the root-soil interface. We assembled the generated metagenomic reads and combined contigs with similar nucleotide composition and differential abundance across samples. This resulted in the reconstruction of 67 metagenome-assembled genomes (MAGs) with a completion of .50% according to the presence of a minimal set of essential genes and a proportion of contamination less than 10% (Methods). These MAGs were taxonomically affiliated with 14 different bacterial classes, and their genomes were systematically mined for the top 10 GO terms significantly enriched in the rhizosphere samples compared to Bulk soil controls (Wald test, individual P values of ,0.05, FDR corrected; Fig. 5). Next, we determined co-occurrence patterns between these terms and identified two clusters. One of those linking genomes coding for "photosynthesis," "carbohydrate metabolic process," and "iron-sulfur cluster assembly," with another one linking "cellular component organization or biogenesis," "response to chemical," "bacterial-type flagellum-dependent cell motility," and "protein-containing complex assembly" (Pearson correlation, individual P value of ,0.05; Fig. 6). When we interpolated the results of these two analyses, we observed that this second cluster is predominantly represented by MAGs classified as Proteobacteria, while the "carbohydrate metabolic process" defining the first cluster was preferentially associated with MAGs classified as Bacteroidia. For 11 of the 15 MAGs assigned to this class, the presence of "carbohydrate metabolic process" predicted a significant enrichment of given MAGs in the microbiota of the Elite variety (Fig. S5, Wald test, individual P value of ,0.05, FDR corrected). Conversely, among Proteobacteria MAGs, the selected GO terms failed to predict enrichment patterns in a given plant genotype (Fig. S5, Wald test, individual P value of ,0.05, FDR corrected).
A distinct bacterial consortium is associated with optimum barley growth under nitrogen-limiting supply. To establish a causal relationship between structural and functional configurations of the rhizosphere microbiota and plant growth, we performed a plant-soil feedback experiment by growing the Elite variety in soils previously used for the growth of either domesticated or wild genotypes amended with an N0% solution (here "conditioned soil"). For this analysis, we focused on the pair Elite-Desert, as these genotypes displayed the most contrasting microbiota ( Fig. 2 and 3). The conditioned soils were used either in their "native" form or subjected to a heat treatment, which we hypothesized would lead to a disruption of the taxonomic and functional configurations of the microbiota (Fig. 7A).
Plants grown in the "heat-treated" soil displayed a growth deficit compared to their native counterpart, although these differences were significant only for the Desert-conditioned soil (one-way analysis of variance  a Desert-mediated impact on the bacterial but not the fungal communities populating the conditioned soils (Kruskal-Wallis test followed by a Dunn post hoc test, P = 0.039, FDR corrected; Fig. S7). This observation motivated us to gain insights into the taxonomic composition of the bacterial communities inhabiting the conditioned soil. A 16S rRNA gene amplicon library constructed from the samples subjected to the feedback experiments   Fig. S8). Likewise, we observed a compositional shift between heattreated and native samples (PERMANOVA, R 2 Treatment = 0.289, P = 0.0002, 5,000 permutations; Fig. 7C). At the end of cultivation, Bulk soil communities could be separated along the axis accounting for the second source of variation, while the taxonomic composition of rhizosphere samples appeared "to converge" to common profiles. Congruently, when we inspected the individual bacteria underpinning this diversification, we identified a phylogenetically diverse group of 10 ASVs, whose enrichment in native samples, accounting for ;7% of the sequencing reads, discriminated between their microbiotas and both unplanted soil and heat-treated samples (Wald test, P value of ,0.05, FDR corrected; Fig. 7D).
Taken together, our data suggest that the heat treatment of the soil substrate led to a scenario comparable to a dysbiosis (43) of the rhizosphere microbiota, defined by a low-diversity community composition associated with reduced growth of the plant host. As the median aboveground dry weight recorded in the feedback experiment is comparable to that recorded for barley plants grown for the first time in Quarryfield soil (compare data of Fig. 7B with the data of Fig. 1), the heat treatment appears to disrupt the capacity of Elite varieties to assemble a taxonomically diverse bacterial consortium associated with optimum barley growth.

DISCUSSION
Our investigation revealed that nitrogen availability for plant uptake impacts the magnitude of the host "genotype effect" on the barley rhizosphere, measured as number of ASVs differentially recruited between genotypes. This effect was maximal at N0%, when measurements of residual N in the rhizosphere were 0 mg kg 21 , while it was minimal at N100%, with ;200 mg kg 21 residual NO 3 in the rhizosphere. This is reminiscent of the observation that in Medicago truncatula, a model legume, the host controls the microbiota in a nitrogen-and genotype-dependent manner (44). Although the modulation of the rhizosphere microbiota in legumes has been associated with plant genes implicated in the establishment of symbiosis with nitrogen-fixing bacteria rather than nitrogen nutritional status per se (45,46), it is conceivable that this latter aspect impacts, at least in part, host-microbe interactions in the barley rhizosphere. This would be congruent with observations gathered from rice, a cereal phylogenetically related to barley, in which the nitrate transporter NRT1.1B emerged as a critical regulator of both nitrogen use efficiency and microbiota recruitment (47). Likewise, a regulator of lateral root development, designated LRT1, emerged as a determinant of microbiota recruitment in maize plants exposed to limiting nitrogen supplies (48). Rhizodeposition (i.e., the plant-mediated release of organic compounds in the rhizosphere) may represent the nexus between plant adaptation to limiting nitrogen supply and microbiota recruitment (16). Consistently, the availability of organic carbon in barley rhizodeposits is inversely correlated with the amount of nitrate concentration in shoots (49). As wild and domesticated barley plants display differential responses to nitrogen fertilization (50) and a genotype-dependent control of rhizodeposition (51), the characterization of primary and secondary metabolites released in the barley rhizosphere may provide mechanistic insights into microbiota diversification in barley. However, as recent investigations revealed that the host genetic control of the rhizosphere microbiota in wild and domesticated barley display a quantitative inheritance (32,34), additional experiments with dedicated genetic material are required to untangle the molecular mechanisms linking nitrogen availability with microbiota diversification.
The observation that the magnitude of host control on the microbiota was greater when plants were exposed to a nitrogen supply limiting barley growth motivated us to embark on a metagenomic survey of this condition. This approach revealed that the microbial communities proliferating at the barley root-soil interface are largely dominated by bacteria; more than 98% of the annotated sequences at the phylum level were classified as bacteria. This is strikingly similar to a previous investigation conducted in a soil with different physical and chemical characteristics (30). The dominance of bacterial sequences over other members of the microbiota is not unusual in soil metagenomes (52), although both the protocol used for microbial DNA preparation and the databases used for sequencing annotation (53) can artificially inflate the proportion of bacteria among the analyzed metagenomes. Despite this potential caveat, we demonstrated with an independent quantification that bacterial gene copy number exceeded that from fungal sources by several orders of magnitude, as previously reported (54,55). It is however important to consider that PCR-based approaches may fail to provide an accurate estimation of fungal biomass due to nuclear exchanges on these filamentous microorganisms (56). A closer inspection of the bacterial and fungal abundances revealed no significant differences among samples for arbuscular mycorrhizal fungi (AMF) symbionts of barley. This apparent discrepancy, in view that AMF can mobilize nitrogen for plant uptake (57), is, however, congruent with a previous investigation conducted with field-grown barley plants, where no genotype effect on AMF root colonization was observed regardless of the nitrogen regimen (58). Furthermore, this observation is also consistent with the fact that our N0% treatment was associated with a replete amount of phosphorus for barley growth, a condition known to suppress AMF colonization (59). Finally, it is important to mention that the microbiota inhabiting soils with a pH below 7, such as Quarryfield, is less conducive to AMF activity than the microbiota inhabiting neutral to alkaline substrates (60).
For these reasons, we decided to focus on the functional characterization of the bacterial component of the microbiota of plants exposed to limiting nitrogen supplies. This allowed us to identify three main GO categories enriched in and differentiating rhizosphere samples from Bulk soil, namely, "carbohydrate metabolic process," "response to chemical," and "pathogenesis." Of interest is the GO category "carbohydrate metabolic process," whose enrichment emerged as both microhabitat and genotype dependent, which is congruent with previous observations that root-derived dissolved organic carbon and carbohydrate utilization by soil microbes display a host genetic component in wild and domesticated barley genotypes (32,51). Likewise, "response to chemical" may mirror the adaptation of rhizosphere communities to plant secondary metabolites, released through rhizodeposition, capable of selectively impacting microbial proliferation, as observed in barley (33) and other cereals (61,62).
Conversely, the GO category "pathogenesis" appears difficult to reconcile with the fact that no obvious symptoms of disease were observed in our samples. However, studies conducted with the model plant Arabidopsis thaliana revealed that components of the host immune system are required for the establishment of a diverse and functional microbiota at the root-soil interface (63); this suggests that the endogenous barley microbiota has evolved the capacity of modulating host immune responses to colonize the rhizosphere. This scenario appears further corroborated by the enrichment of the GO category "bacterial-type flagellum-dependent cell motility"; despite the fact that molecular components of this machinery have been considered a paradigmatic epitope of the plant immune system (64), it is now emerging that their recognition by host plants contributes to signal modulation and microbiota establishment (65). Similarly, enrichment of the GO category "response to reactive oxygen species" in the microbiota of the Elite plants may further explain the role of this class of compounds in modulating plant-associated bacterial communities (66). A prediction of these observations is that components of the barley immune system may act as a "checkpoint" for the taxonomic and functional composition of the rhizosphere microbiota.
Substrate availability and interorganismal relationships appear to be determinants also for the Bulk soil communities, as mirrored by the significant enrichment of the GO terms "photosynthesis," "antibiotic biosynthetic process," and "sporulation." The absence of a source of organic compounds such as rhizodeposition creates a niche for the proliferation of CO 2 -fixing microorganisms, which are ubiquitous in the soil ecosystem (67). Likewise, the enrichment of "antibiotic biosynthetic process" is congruent with what was observed for agricultural soils in a cross-microbiome survey (68), while sporulation underpins microbial adaptation to stressful soil conditions (69). Furthermore, unplanted soil communities display a differential enrichment for function implicated in phosphorous homeostasis. As the relative abundances of carbon, nitrogen, and phosphorus can be considered constrained in microbial biomass (70), this observation suggests that, although phosphorous was applied with the nutrient solution to all specimens, this element may act as a limiting factor predominantly for unplanted soil communities, where the lack of exudates reduces phosphorous solubility. Taken together, these observations provide mechanistic insights into the multistep selection process differentiating rhizosphere communities from Bulk soil communities (4,27), implicating the modulation of host immune responses as one of the requirements for bacterial establishment in the rhizosphere of plants exposed to limiting nitrogen supply. However, as these experiments were performed in a single soil type, caution is required in extrapolating the results as being indicative of general phenomena applicable across all soils. Further metagenomics investigations with plants exposed to replete nitrogen conditions, benefiting also from the latest development in sequencing technologies (71), will be required to accurately gauge the impact of this resource (or lack thereof) on the functional potential of the barley microbiota.
Despite the fact that the 67 MAGs generated in this work accounted for less than 10% of the metagenomic reads, these figures are aligned with what has been recently observed for the rhizosphere of sorghum (72), a cereal phylogenetically related to barley. This effort allowed us to identify genomes belonging to not only the dominant phyla of the plant microbiota (i.e., Actinobacteria, Bacteroidetes, Firmicutes, and Proteobacteria) but also members of additional classes, such as an individual member of the metabolically diverse and yet poorly characterized Zixibacteria phylum (73)(74)(75). Furthermore, mapping reads associated with the GO terms differentially enriched between microhabitat and genotype allowed us to gain novel insights into the relationships between taxonomic and functional composition of the barley microbiota. For instance, we observed an association between the GO category "carbohydrate metabolic process" and the enrichment of members of the phylum Bacteroidetes in the Elite rhizosphere. As cell wall features represent a recruitment cue for the plant microbiota (76), this enrichment may mirror the capacity to degrade complex polysaccharides encoded by members of this phylum (77,78). The observed genotype-specific enrichment may be further explained by polymorphisms of barley genes regulating carbohydrate composition in the cell wall (79,80).
The plant-soil feedback experiment we implemented suggested that a functional rhizosphere microbiota is required for optimal barley growth under nutrient-limiting conditions. Although not significantly different, mean values of aboveground biomass of Elite plants recorded in the Desert-conditioned soil were higher than those recorded from the soil conditioned with the same genotype. Despite that phylogenetic relatedness between condition and focal species in plant-soil feedback experiments appears unrelated to the strength of the feedback itself (81,82), compositional shifts between the conditioned and focal microbiota tend to be associated with enhanced plant growth (83). However, significant differences in growth were observed when Elite plants were exposed to heat-inactivated soils, which are associated with a reduction of alpha diversity indices in the rhizosphere, a condition that has been previously linked to stressful soil conditions (83). In turn, this effect could be due to the treatment on the microbes per se, the release of mineral nutrients, and/or the disruption-liable carbon compounds released through exudates (84) by conditioning plants capable of modulating individual members of the barley microbiota (33). Unlike recent observations gathered from plant-soil feedback experiments of maize plants exposed to limiting nitrogen conditions (48), what emerged from our study is the control exerted by the recipient genotype on the resulting bacterial communities. This was manifested by the microbiota of plants exposed to either Desert-conditioned or Elite-conditioned soil "converging" toward a phylogenetically conserved bacterial consortium. This is in accordance with data gathered from rice, using both soil feedback experiments (85) and synthetic communities (47), indicating the host genotype as a driver of a plant growthpromoting microbiota. Likewise, a recently developed indexed bacterial collection of the barley rhizosphere microbiota indicated a growth promotion potential for members of the phyla differentially recruited in the feedback experiment (86,87). Taken together, this suggests that the enriched bacteria represent a consortium of beneficial bacteria required for optimum barley growth whose recruitment is driven, at least in part, by the host genotype.
Conclusions. Our results point to nitrogen availability for plant uptake as inversely correlated with the magnitude of host genetic control on the taxonomic composition of the barley rhizosphere microbiota. When nitrogen supply limits barley growth, wild and domesticated genotypes retain specific functional signatures, which appear to be encoded by distinct bacterial members of the microbiota. Although we found evidence for nitrogen metabolism executed by these communities, adaptation to the plant immune system emerged as an additional recruitment cue for the barley microbiota. Plant-soil feedback experiments suggest that these distinct compositional and functional configurations of the microbiota can be "rewired" by the host genotype, leading to a recruitment of a consortium of bacteria putatively required for optimum plant growth. Thanks to recent insights into barley genes shaping the rhizosphere microbiota (32,34), these concepts can now be tested under laboratory and field conditions to expedite the development of plant varieties profiting from improved yields with reduced impacts of N fertilization on the environment.

MATERIALS AND METHODS
Experimental design. This investigation consists of three distinct but interconnected experiments. For each experiment, plants were maintained under controlled conditions in the same soil type designated Quarryfield (see Soil below), and individual samples were arranged in a completely randomized design. For the first experiment, we grew individual biological replicates (i.e., pots) Elite, Desert, and North and Bulk soil controls exposed to three different nitrogen treatments, designated N0%, N25%, and N100% (see Nitrogen treatments below), according to the following scheme and subjected all samples to 16S rRNA gene amplicon sequencing: N0% Desert = 5, N0% North = 3, N0% Elite = 4, N0% Bulk = 4, N25% Desert = 5, N25% North = 3, N25% Elite = 4, N25% Bulk = 4, N100% Desert = 5, N100% North = 4, N100% Elite = 3, and N100% Bulk = 4. Alongside these samples, we prepared two additional Bulk soil controls amended with a plug of the agar substrate used for seed germination. The total number of sequenced samples was 50. In the second experiment, we grew and subjected to shotgun metagenomic sequencing three individual biological replicates (i.e., individual plants in individual pots) of the genotypes Elite, Desert, and North and three Bulk soil controls exposed to N0% treatment. The total number of sequenced samples was 12. In the third and final experiment, we grew and subjected to 16S rRNA gene amplicon sequencing individual biological replicates (i.e., individual plants in individual pots) of the Elite genotype soil controls in Quarryfield soils that were previously conditioned (see Plant-soil feedback experiment below) with either the Elite or Desert genotype in native form or after heat treatment. For the former, we also contemplated Bulk soil control pots. After discarding pots with no detectable plant growth, the numbers of rhizosphere samples exposed to Elite-conditioned soil retained for sequencing were 11 for Elite-native rhizosphere and 7 for Elite-native Bulk , the numbers of rhizosphere samples exposed to Desert-conditioned soil were 14 for Desert-native rhizosphere and 9 for Desert-native Bulk . The numbers of sequenced rhizosphere samples exposed to heat-treated soil were 15 for Elite-treated rhizosphere and 15 for Desert-treated rhizosphere . The total number of sequenced samples was 71.
Soil. Soil was sampled from the agricultural research fields of the James Hutton Institute, Invergowrie, Scotland, UK, in the Quarryfield site (56°27'50N 3°4'290W). This is a sandy silt loam soil with a pH of 6.2 and 5% organic matter content. The nitrogen content of this soil was 1.8 mg kg 21 ammonium and 13.5 mg kg 21 nitrate. The site was left unplanted and unfertilized in the 3 years preceding the investigations.
Plant material and growth conditions. Barley seeds of the domesticated (Hordeum vulgare subsp. vulgare) and wild (Hordeum vulgare subsp. spontaneum) genotypes, the variety Morex (i.e., Elite), and the accessions B1K-12 (i.e., Desert) and B1K-31 (i.e., North) were surface sterilized as previously reported (88) and germinated on 0.5% agar plates at room temperature. Seedlings displaying comparable rootlet development were sown individually in 12-cm diameter pots containing approximately 500 g of the Quarryfield soil, from which stones and large debris were manually removed. Unplanted pots filled with the same soil (i.e., Bulk soil controls) were maintained in the same glasshouse and subjected to the same treatments as planted pots. One-week-old plantlets were transferred for 2 weeks to a growth room at 4°C for vernalization. Following the vernalization period, plants were maintained in a randomized design in a climatic-controlled glasshouse at a 18/14°C (day/night) temperature regimen with 16 h of daylight that was supplemented with artificial lighting to maintain a minimum light intensity of 200 mmol quanta m 22 s 21 until early stem elongation (Fig. S1 in the supplemental material). Watering was performed weekly as indicated (see Nitrogen treatments below). Pots were rotated on a weekly basis to minimize potential biases associated with given positions in the glasshouse. Nitrogen treatments. The nutrient solutions described in this study (i.e., N100%, N25%, and N0%) are reported in Table S1. Nutrient solutions were applied at a rate of 25 mL per kg of soil each week. Applications started 2 days after planting, were interrupted during vernalization, and were reinstated once the plants were transferred to the growing glasshouse and they reached early stem elongation. Fourteen treatments were applied with a total of 312.5 mg of NO 3 and 81.4 mg of NH 4 1 for the N100% solution, 78.1 mg of NO 3 and 20.0 mg of NH 4 1 for the N25% solution, and 0 mg of NO 3 and NH 4 1 for the N0% solution applied per pot.
Plant and soil nitrogen determination. To assess the N content of the plant, at the time of sampling, a newly expanded leaf was sectioned from every plant, freeze-dried, ball milled, and measured for N content in an Elemental Analyzer CE-440 (Exeter Analytical, Inc., UK). The soil from the pots was sieved through a 2 mm mesh sieve and mixed. Five grams of soil was added to 25 mL of 1 M KCl, and the resulting solution was mixed in a tube roller for 1 h at ;150 rpm. Supernatant was transferred to 50 mL Falcon tubes and centrifuged for 15 min at 5,000 rpm, then the supernatant was subjected to another round of centrifugation. The supernatant was transferred to a Falcon tube and analyzed with a Discrete Analyser Konelab Aqua 20 (Thermo Fisher, Waltham, USA) in the analytical services of The James Hutton Institute (Aberdeen, UK). In parallel, ;10 g from the sieved soil was oven dried at 70°C for 48 h, and dry weight was recorded to express the analytical results in NO 3 and NH 4 1 in mg N kg 21 of soil. Bulk soil and rhizosphere DNA preparation. At early stem elongation, plants were excavated from the soil, and the stems were separated from the roots. The uppermost 6 cm of the root system was detached from the rest of the root corpus and processed for further analysis. The sampled aboveground material was oven dried at 70°C for 48 h, and the dry weight was recorded. The roots were shaken manually to remove loosely attached soil. For each barley plant, the seminal root system and the attached soil layer was collected and placed in a sterile 50 mL Falcon tube containing 15 mL of phosphate-buffered saline (PBS). The rhizosphere was operationally defined, for these experiments, as the soil attached to this part of the roots and extracted through this procedure. The samples were then vortexed for 30 s and transferred to a second 50 mL Falcon containing 15 mL of PBS and vortexed again to ensure the dislodging and suspension of the rhizosphere. Then, the two Falcon tubes with the rhizosphere suspension were combined and centrifuged at 1,500 Â g for 20 min to precipitate the rhizosphere soil into a pellet, flash-frozen with liquid nitrogen, and stored at 280°C until further analysis. In addition, we incubated water agar plugs (;1 cm 3 ) into two unplanted soil pots, and we maintained them as control samples among the experimental pots to monitor the effect of this medium on the soil microbial communities. DNA was extracted from unplanted soil and rhizosphere samples using a FastDNA spin kit for soil (MP Biomedicals, Solon, USA) according to the manufacturer's recommendations and stored at 220°C.
Preparation of 16S rRNA gene amplicon pools. The hypervariable V4 region of the small subunit rRNA gene was the target of amplification using the PCR primer pair 515F (59-GTGCCAGCMGCCGCGGTAA-39) and 806R (59-GGACTACHVGGGTWTCTAAT-39). The PCR primers had incorporated an Illumina flow cell adapter at their 59 termini, and the reverse primers contained 12 bp of unique "barcode" for simultaneous sequencing of several samples (89). PCRs were performed using 50 ng of metagenomic DNA per sample using the Kapa HiFi HotStart PCR kit (Kapa Biosystems, Wilmington, USA). The individual PCRs were performed in 20 mL final volume and contained 4 mL of 5Â Kapa HiFi buffer, 10 mg of bovine serum albumin (BSA; Roche, Mannheim, Germany), 0.6 mL of a 10 mM Kapa dNTP solution, 0.6 mL of 10 mM solutions of the 515F and 806R PCR primers, and 0.25 mL of Kapa HiFi polymerase. The reactions were performed using the following program: 94°C (3 min), followed by 35 cycles of 98°C (30 s), 50°C (30 s), and 72°C (1 min), and a final step of 72°C (10 min). For each primer combination, a no-template control (NTC) was included in the reactions. To minimize amplification biases, PCRs were performed in triplicate, and at least two independent master mixes per barcode were generated (i.e., 6 reactions/sample). PCRs were pooled in a barcode-dependent manner, and an aliquot of each amplification product was inspected on a 1.5% agarose gel. Only samples whose NTCs yielded an undetectable PCR amplification were retained for further analysis. PCR purification was performed using an Agencourt AMPure XP kit (Beckman Coulter, Brea, USA) with 0.7 mL of AmPure XP beads per 1 mL of sample. Following purification, each sample was quantified using PicoGreen (Thermo Fisher Scientific, Watham, USA), and individual barcode samples were pooled in an equimolar ratio to generate amplicon libraries.
Illumina 16S rRNA gene amplicon sequencing. The pooled amplicon library was submitted to the Genomics facility at The James Hutton Institute (Invergowrie, UK) for quality control (Bioanalyzer, Agilent Technologies), processing, and sequencing. Amplicon libraries were supplemented with 15% PhiX control library (Illumina). The resulting high-quality libraries were run at 10 pM final concentration on an Illumina MiSeq system with paired-end 2 Â 150 bp reads (89) to generate sequencing FASTQ files.
Amplicon sequencing read processing. Sequence reads were subjected to quality assessment using FastQC (90). ASVs were then generated using DADA2 version 1.10 (91) and R 3.5.1 (92) following the basic methodology outlined in the DADA2 pipeline tutorial (93). Read filtering was performed using the DADA2 paired FastqFilter method, trimming 10 bp of sequence from the 59 end of each read using a truncQ parameter of 2 and maxEE of 2. The remainder of the reads was of high quality; consequently, no 39 trimming was deemed necessary. The dada2::learn_errors() method was run to determine the error model with a MAX_CONSIST parameter of 20, following which the error model converged after 9 and 12 rounds for the forward and reverse reads, respectively. The dada2::dada() method was then run with the resulting error model to denoise the reads using sample pooling, followed by read merging and chimera removal using the consensus method. Taxonomy assignment was performed using the Ribosomal Database Project (RDP) naive Bayesian classifier through the dada2::assignTaxonomy() method with the SILVA database (94) (version 138) using a minimum bootstrap confidence of 50. The DADA2 outputs were finally converted to a Phyloseq object (version 1.26.1) (95). The Phyloseq objects for both the nitrogen gradient and the plant-soil feedback experiments were initially merged. Next, sequences classified as either "chloroplast" or "mitochondria" were pruned in silico from the merged object. Likewise, ASVs matching a list of potential contaminants of the lab (96) were removed as well as ASVs lacking a taxonomic classification at the phylum level (i.e., "NA"). We further applied abundance filtering and retained ASVs occurring with at least 20 reads in 2% of the samples. Finally, the Phyloseq objects were rarefied at 25,000 reads per sample, as recommended for groups with large differences in library size (97), before downstream analyses.
Metagenome sequencing, annotation, and analysis. We generated a new set of 3 bulk soil and 3 rhizosphere DNA preparations from each of the three genotypes tested (i.e., Desert, North, and Elite) from specimens maintained in Quarryfield soil under N0% conditions as described above. These 12 new preparations were quantified and submitted to the LGC Genomics sequencing service (Berlin, Germany) where they were used to generate DNA shotgun libraries using the Ovation Rapid DR multiplex system 1-96 NuGEN (Leek, The Netherlands) kit following the manufacturer's recommendations. These libraries were run simultaneously into an individual Illumina NextSeq500 run following the manufacturer's recommendations with 2 Â 150 bp chemistry and generated a total of 412,385,413 read pairs. After sequencing, read pairs were demultiplexed according to the samples' barcodes using the Illumina bcl2fastq2.17.1.14 software.
Metagenome analysis was conducted according to the general approach of Hoyles et al. (98), using updated tools where appropriate. Sequence reads were quality assessed using FastQC and quality/ adapter trimmed using TrimGalore (99), using a quality cutoff of 20, a minimum sequence length of 75 bp, and removing terminal N bases. Taxonomic classification of the sequence reads was performed using Kraken 2.0.9 (100) with the Kraken PlusPFP database (101), which incorporates protozoa, fungi, and plants in addition to the archaea, bacteria, and viruses present in the standard database. Host contamination was removed by alignment against the Morex V2 barley genome sequence (102) using BWA MEM (103), and nonaligning reads were extracted from the resulting BAM files using SAMtools (104). Metagenome assembly was conducted using MegaHit version 1.2.9 (105) with the "meta-large" preset. Predicted proteins were produced from all assemblies using Prodigal version 2.6.3 (106), which were then clustered using MMseqs2 version 11.e1a1c (107) and the "easy_cluster" method. Abundance of predicted proteins in each sample was determined by alignment of sequence reads against the representative cDNA sequences of the clusters using the Burrows-Wheeler Aligner Maximal Exact Matches (BWA MEM) and determining the read counts associated with each sequence using a custom PySAM (108) script. Functional annotations of the protein sequences were performed using InterProScan 5-50-84.0 (109) and Interpro version 84.0. GO terms were enumerated using a custom Python script, which assessed the number of occurrences of each term in each sample based on the previously determined abundance of each annotated sequence. GO terms were mapped to the metagenomics GO slim subset dated 2020 March 23 (110) using the Map2Slim function of OWLtools (111). Functional enrichment analysis was performed using DESeq2 (112) version 1.26.0.
Nitrogen cycle gene analysis. Abundance of nitrogen cycle genes was determined using the NCycProfiler tool of NCycDB (117) with the diamond method. Pairwise t tests were performed between the samples of each group within each gene to identify combinations with statistical differences between samples (Benjamini-Hochberg corrected FDR of ,0.05).
Plant-soil feedback experiment. We grew the Desert and Elite genotypes in Quarryfield soil supplemented with a N0% nutrient solution under controlled conditions (see Plant growth conditions). We selected these two genotypes because, in the tested soils, they host a taxonomic and functionally distinct microbiota. At early stem elongation, we removed the plants from the soil, and we harvested the residual soil and kept it separated in a genotype-dependent manner. We reasoned that at the end of cultivation, the soils would have been enriched, at least partially, for specific microbial taxa and functions associated with either genotype. This residual soil, either in a "native form" (i.e., not further treated after sampling) or after being exposed to a heat treatment (126°C for 1 h, repeated twice at an interval of ;12 h), was used as a substrate for subsequent cultivation of a recipient Elite barley genotype. These plants were maintained under controlled conditions (see Plant growth conditions) and supplemented with an N25% solution to compensate for the near-complete depletion of this mineral in the previous cycle of cultivation (compare the NH 4 1 and NO 3 concentrations of rhizosphere specimens at N0% and N25% in Fig. 1). At early stem elongation, plants were harvested, and their aboveground biomass was determined after drying stems and leaves at 70°C for 48 h. At the end of each replicated experiment, the residual soil was collected and subjected to chemical and physical characterization (Yara United Kingdom, Ltd., Grimsby, United Kingdom). A quantitative real-time PCR assay was used to quantify the bacterial and fungal DNA fractions in samples from the conditioned soil experiment as follows. DNA samples were diluted to 10 ng/mL and successively diluted in a serial manner to a final concentration of 0.01 ng/mL. This final dilution was used for both the femto fungal DNA quantification kit and femto bacterial DNA quantification kit (Zymo Research), and quantification was conducted according to the manufacturer's protocol. Briefly, 2 mL of the 0.01 ng/mL dilution of each sample was used together with 18 mL of the corresponding fungal or bacterial master mix. Two microliters of the fungal or bacterial standards was also used to create the respective quantification curves. DNA samples from the conditioned soil experiment were randomized in the 96-well plates, using a minimum of 11 biological replicates per treatment. Quantification was performed in a StepOne thermocycler (Applied Biosystems by Life Technologies) following the cycling protocols of each of the above-mentioned bacterial and fungal kits.
Statistical analyses on the univariate data set and amplicon sequencing. Data analysis was performed in R software using a custom script with the following packages: Phyloseq (95) version 1.36.0 for preprocessing and alpha-and beta-diversity analyses, ggplot2 version 3.3.4 (118) for data visualization, vegan version 2.5-7 (119) for statistical analysis of beta-diversity, PMCMR version 4.3 (120) for nonparametric analysis of variance. For any univariate data set used (e.g., aboveground biomass), the normality of the data distribution was checked using a Shapiro-Wilk test. For data sets that were normally distributed, the significance of the imposed comparisons was assessed by an ANOVA test followed by a Tukey post hoc test. Nonparametric analysis of variance tests were performed by using a Kruskal-Wallis rank sum test followed by a Dunn's post hoc test with the functions kruskal.test and the posthoc.kruskal.dunn.test, respectively, from the package PMCMR. We used Spearman's rank correlation to determine the similarity between unplanted soil profiles and Bulk soil samples amended with water agar plugs (Table S2). Analysis of the differentially enriched ASVs was performed (i) between individual genotypes and Bulk soil samples to assess the sample effect and (ii) between the rhizosphere samples to assess the genotype effect. The genotype effect was further corrected for a microhabitat effect (i.e., for each genotype, only ASVs enriched against both unplanted soil and at least another barley genotype were retained for further analysis). The analysis was performed using the DESeq2 package (112) version 1.32.0 consisting of a moderated shrinkage estimation for dispersions and fold changes as an input for a pairwise Wald test. This method identifies the number of ASVs significantly enriched in pairwise comparisons with an adjusted P value (false-discovery rate [FDR] of ,0.05). This method was selected because it outperforms other hypothesis-testing approaches when data are not normally distributed and a limited number of individual replicates per condition are available (97). The scripts used to analyze the data and generate the figures of this study are available at https:// github.com/BulgarelliD-Lab/Barley-NT-2020.

SUPPLEMENTAL MATERIAL
Supplemental material is available online only.

ACKNOWLEDGMENTS
We are particularly grateful to Eyal Fridman (Agricultural Research Organisation, Bet Dagan, Israel) for providing us with the B1K seeds used in this study. We thank Lawrie Brown (The James Hutton Institute, Invergowrie) for advising us in the preparation and application of the nutrient solutions and Niamh Johnston (Nuffield Research Placements Student) and Jim Wilde (The James Hutton Institute, Invergowrie) for technical assistance during the feedback experiment. We thank LGC Genomics Gmbh (Berlin, Germany) for generating the shot gun metagenomics reads. For the purpose of open access, the authors have applied a CC BY public copyright license to any author accepted manuscript version arising from this submission. The experimental work presented in the manuscript was supported by a Royal Society of Edinburgh/Scottish Government Personal Research Fellowship cofunded by Marie Curie Actions, a Carnegie Trust for the Universities of Scotland Research Incentive grant (RIG007411), and a UK Research and Innovation grant (BB/S002871/1) awarded to D.B. R.A.T. was supported by a Scottish Food Security Alliance-Crops studentship, provided by the University of Dundee, the University of Aberdeen, and the James Hutton Institute. R.K. was partially supported by a British Society for Plant Pathology MSc/MRes Bursary Scheme. The metagenomic data analysis was supported by the H2020 Innovation Action "CIRCLES" (European Commission, grant agreement 818290) awarded to the University of Dundee. The Genomics facility at the James Hutton Institute is funded by the Rural & Environment Science & Analytical Services Division of the Scottish Government.