Under Attack: Fungal Landscape of Spring Barley From Seedling Emergence to Harvest During a Severe Puccinia Hordei Epidemic

All plant tissues from leaves, stems and roots are hosting a wide diversity of fungal species. Our understanding of the assembly of this diversity of fungi during the plant growth cycle is limited. Here, we characterized the mycobiome of three spring barley cultivars at weekly intervals from seedling emergence to senescence and seed maturity. A notable proportion of members of the fungal communities were shared among different plant organs, but community dynamics were tissue specic. A severe attack of Puccinia hordei occurring during the vegetative stage had profound effects on the mycobiome, and P. hordei biomass displaced that of other taxa. Plant tissue type was the most important factor determining the mycobiome, but also plant age was contributing signicantly. Using a random forest model, we found that specic members of the mycobiome were responding differently to plant age. A co-occurrence network analysis revealed complex interactions among fungal OTUs, and network connectivity was changing as per plant growth stage and plant tissue type. This study contributes to the understanding of assembly of fungal communities in cereals by providing a detailed atlas of fungal communities associated with barley. This knowledge will be vital for microbiome assisted plant health management and our study will serve as an important baseline for future efforts to harness microbiota in cereal health.


Introduction
The ability of plants to recruit optimal microbiomes is intimately linked to plant growth and health.
Assembly of the plant-associated microbiota is a dynamic process that involves tight associations between the environment and the plant during the plant growth cycle from seed germination to plant senescence. The soil is a major driver of root-associated microbial communities as the primary inoculum source for the rhizosphere microbiota [1]. The phyllosphere microbiota is likewise driven by soil microbial communities [2], but also by deposition of a signi cant number of microorganisms directly from the air [3]. Origins of air-borne spores from especially fungal pathogens have been subject to several studies, and long-distance movement of spores of plant pathogens across continents has been demonstrated [4].
The majority of diseases of signi cant economic importance affecting barley (Hordeum vulgare L.) are caused by fungal pathogens such as Fusarium spp., Puccinia hordei, Pyrenophora graminea, P. teres, Blumeria graminis f.sp. hordei, Rhynchosporium commune, Ramularia collo-cygni and Stagonospora nodorum that all infect aboveground parts of barley [5]. Barley roots are likewise hosts for a number of yield-reducing fungal pathogens such as Fusarium spp., Gauemannomyces graminis var. tritici and Bipolaris sorokiniana [5]. All of these pathogens are tissue-speci c and prefer different developmental stages of the plant and certain environmental conditions such as temperature and humidity, to succeed [5]. Despite the detailed knowledge available for most cereal fungal pathogens of economic importance, only limited information is available on the spatiotemporal dynamics of fungal communities associated with the crops and their effect on the pathogens [6]. The mycobiome is most likely varying across the barley life cycle from early development to late senescence, as different community members colonize the plant depending on their lifestyles and requirements. Therefore, baseline information of fungal plant colonizers during the entire host growth cycle is needed to fully understand interactions between pathogens and fungal communities for sustainable disease control.
The early stages of plant development are highly determinative for the structure of the plant-associated microbiota [1], while it stabilizes at the later growth stages [7]. In a recent study, bacterial communities were highly dynamic during the vegetative phase of rice root development but ultimately stabilized as the plants matured [8]. However, the entire process of microbial enrichment and structural changes of communities throughout the plant growth cycle has only been poorly studied, particularly for fungal communities despite their important roles in plant health and productivity. Formation of plant-associated microbiota is driven by host age and development, among other factors, in such diverse crops as Medicago, pea, maize, sugar beet and wheat [9][10][11]. However, root microbiomes of rice were mainly in uenced by plant developmental stage and only to a lesser extent by the age itself [8,12].
In our study we asked the following questions: a) what is the microbial community composition, diversity and structure in leaf, stem, head and roots of barley? b) what are the temporal dynamics of fungal communities during the barley growth cycle? c) which are the fungal taxa that contributes mostly to the variation of the fungal communities over time? To answer these questions, we used a high-resolution metabarcoding approach in a highly detailed study to elucidate spatiotemporal dynamics of the barley mycobiota in different host tissues in three different commercial cultivars during a full growth season.
We sampled root, leaf, stem and head material weekly from seedling emergence to senescence to establish an atlas of fungal colonization of barley during all its developmental stages. At the same time, we visually assessed the presence of major barley diseases in the crop. A brown rust (Puccinia hordei) epidemic was peaking in July, and less severe attacks of other barley pathogens were observed during the growth phase of the crop. Our results highlighted major fungal species present in barley, and further described shifts in fungal composition both during barley development and the Puccinia hordei pathogen attack. This is the rst study to comprehensively highlight the fungal landscape of barley that will be a valuable contribution to the understanding of pathogens in a microbiome context.

Study sites, treatments and plant material
The study site was located at an experimental research eld at Flakkebjerg, Denmark (55.321547 N, 11.385 E). Commercial cultivars of spring barley (Cumulus, Quench and Milford) were sown in spring 2015 in 22 m 2 plots in four replicates in a split-plot design. The plots were not treated with fungicides.
Disease severities of some of the most important barley leaf diseases caused by Puccinia hordei, Rhychosporium commune, Pyrenophora teres, Blumeria graminis, and Ramularia collo-cygni were visually assessed by estimating the percentage leaf coverage of speci c symptoms on ag leaves according to EPPO guidelines [13].
Ten plants from each of 4 replicated plots were randomly collected regardless of symptoms on the leaves at weekly intervals starting from emergence (April 24th ) and ending at seed maturity (August 21st ) by uprooting the plants. Until week 6, only root and leaf samples were collected but later roots, leaves, stem and heads were sampled. At week 17 and 18, only heads were sampled. Roots were gently washed in sterile water before processing. The leaf just below the ag leaf was taken for analysis. Tissue from the ten plant samples from each replicate were pooled to constitute leaf, root, stem, and head samples from each cultivar, time-point and replicate plot, resulting in 624 samples. All samples were lyophilized for 48 hours, ground using sterile steel beads in a Geno/Grinder 2000 (OPS Diagnostics, Bridgewater, NJ, USA) at 1500 rpm for 3 x 30s, lyophilized again for 48 hours and then ground using a mortar and pestle to obtain a ne powder. An overview of the samples can be seen in Table ST1.
Library preparation and amplicon sequencing 100 mg of ground material was used for DNA extraction using KingFisher TM mL (Thermo Fisher Scienti c Inc., Waltham, MA, USA) with a sbeadex TM kit (LGC Genomics, Middlesex, UK). Next, we prepared sequencing libraries from the 624 DNA samples. For ampli cation of the fungal ITS2 region, fITS7 and ITS4 primers were used [14]. PCR was performed in a reaction mixture of 25 µl consisting of 1 × PCR buffer, 1.5 mM MgCl2, 0.2 mM dNTPs, 1 µM primers, 1 U GoTaq Flexi polymerase (Promega Corporation, Madison, USA) and 1µl DNA template. PCR was conducted in a GeneAmp PCR System 9700 thermal cycler (Thermo Fisher Scienti c) using 94°C for 5 min, followed by 25 cycles at 94°C for 30 s, 57°C for 30 s, 72°C for 1 min, and a step at 72 •C for 10 min. Dual indexing was used to allow pooling of 96 samples. For indexing, primers including tags were used in a PCR for 10 cycles, with the thermal cycler program described above. In addition to dual indexing, tags were added to the forward primer as internal barcode [15]. Primer sequences are listed in Table ST2. PCR products were pooled, precipitated and re-eluted as described earlier [16]. To extract amplicons (300-450 bp) the pooled DNA was separated on an agarose gel and the band of expected size was excised and extracted using QIAquick Gel Extraction Kit (Qiagen), and the concentration evaluated using a spectrophotometer. All samples were shipped to Euro ns MWG (Ebersberg, Germany) for sequencing on their Illumina MiSeq platform using dual indexing strategy. All sequence les from this study were deposited in the NCBI sequence read archive (PRJNA756450).

Bioinformatics and statistical analysis
Analysis of sequence reads was done using VSEARCH version 2.6 [17]. The samples were split based on indices and internal barcodes, and paired end reads were joined with an overlapping minimum read length of 30 base pairs followed by chimera detection, dereplication, ITS extraction and clustering [18]. Taxonomy assignments for the clustered operational taxonomic units (OTUs) was done using the UNITE v.8 database for fungal OTUs [19], and OTUs assigned as unknown at kingdom level were removed. OTU table, taxa assignments and metadata were arranged in R v 4.1.0 [20] using 'phyloseq' package [21] for data analysis and visualization.
Samples representing less than 500 reads were excluded and the OTU table was either transformed to relative abundance or rari ed before performing alpha and beta diversity based calculations. Shannon index and species richness were calculated using 'vegan' package [22] on rari ed OTU table 100 times, and calculating mean values. Diversity analysis, species richness and dissimilarity analysis were carried out using the 'vegan' package. For beta diversity calculations, the OTU table was transformed to relative abundances to calculate Bray-Curtis for fungal communities, and partitioning of variance was carried out for the PERMANOVA test using 'adonis' function from "vegan" package. In addition, the relationship of time and community structure was tested on Bray-Curtis dissimilarity using Mantel test. Indicator OTUs for each plant tissue were computed using "labsdv" package in R.
A Random forest regression model on relative abundances of OTUs was used to predict plant age, and OTUs responding to host age [23]. OTU table and taxa abundances were set as predictors, while host age was set as response in the model implemented in randomForest v 4.6 package in R [24]. The model was applied to the different plant tissues by splitting data. Rare OTUs (<50 reads; <3 samples) were removed.
For predicted vs actual plant age, we used 75% of a randomly selected dataset for training, and the rest for testing the model. We evaluated the model using mean squared error (MSE) percentage values to identify major age-associated OTUs.
Fungal OTU correlations were visualized using networks, and network properties were computed using the "igraph" package [25] as described earlier [26]. The OTU table was trimmed as described above and normalized as read abundance counts per million using the "edgeR" package [27]. Spearman correlations on all OTU pairs were carried out and highly correlated OTUs with Spearman's rank correlations > 0.6 for positive correlations and < − 0.6 for negative correlations, and p values < 0.001 were used for creating the network graphics. All correlated OTUs were visualized in a network, where OTUs were set as nodes, and correlation as edges. All networks were subjected to a Fruchterman-Reingold layout with 999 permutations.

Data structure
Sequencing of the fungal ITS2 from 624 wheat samples resulted in 12.2 million reads after quality control. These reads were clustered into 422 OTUs excluding singletons and OTUs that were occurring in less than three samples (Table ST1). Samples containing less than 500 reads were removed from the dataset leaving 573 samples for diversity-based calculations.

Disease assessments
The severity of barley diseases was visually assessed in the plots (Table ST3). Symptoms characteristic of P. hordei were observed in the plots of all three cultivars, with the highest severity in the susceptible cultivar Quench (ST3). The timing of appearance of rust symptoms was different in the three cultivars. Symptoms were visible already in week 8 in the susceptible Quench increasing until weeks 12 and 13. In Cumulus (moderately resistant) symptoms became visible in week 11 and peaking in week 14, while in the resistant Milford, symptoms were not visible until week 15. In contrast, relative sequence read abundances were strikingly similar in all three cultivars (OTU1 was assigned to P. hordei) (Figure 1). We observed less severe symptoms of Rhyncosporium commune, Pyrenophora teres, Blumeria graminis and Ramularia collo-cygni in the plots, except that R. commune symptoms were notable in weeks 11 to 14, and P. teres symptoms were notable in weeks 7 to 14 in Cumulus. B. graminis and R. collo-cygni symptoms were observed in Milford at the later growth stages.
Diversity of fungal communities varied across the season Root communities had signi cantly higher alpha diversities compared to communities from leaves, stems and heads (Wilcoxon test p <0.001 in 'observed' and 'Shannon index') ( Figure SF 1). Fungal alpha diversity varied across plant age. In the stems and heads, alpha diversity was declining over the season ( Figure 2), while fungal diversity in leaves was rapidly increasing during an initial recruitment period peaking in week 6 followed by a gradual decrease until week 11 when diversity again increased until maturity. The diversity in roots gradually increased until a peak around week 11 and after that alpha diversity was declining until harvest.

Plant tissue and plant age are important determinants of the barley mycobiome
Using a multivariate PERMANOVA test performed on a Bray-Curtis dissimilarity matrix, we found that a notable part of the variation in the dataset was explained by plant tissue type (17%), plant age (2%) and plant tissue and age interactions (4%), whereas host genotype (cultivar) had a negligible effect on the fungal communities (Table 1). After splitting the dataset into the different plant tissue types, plant age was the strongest factor shaping fungal communities explaining 44% of the variation in the leaves followed by heads (36%), stems (23%) and roots (7%). Ordinate analysis showed that the communities were loosely clustered in the early growth stages, while stabilizing and forming tight clusters in all tissue types at the later growth stages, implying a stronger selection imposed on the fungal communities at the later growth stages (Figure 3). Stem and leaf originating communities were similar, albeit forming their own separate clusters in the PCoA plots. Root fungal communities stabilized at later growth stages and separated clearly from communities in above-ground plant tissues. A mantel test using Bray-Curtis dissimilarities revealed time-decay relationships of fungal communities in the different plant tissues. The strongest effect of plant age on the community structure was observed in leaves (Mantel statistic R: 0.57, p < 0.001) followed by heads (Mantel statistic R: 0.43, p < 0.001), stems (Mantel statistic R: 0.26, p < 0.001) and roots (Mantel statistic R: 0.22, p < 0.001). In addition, temporal slopes based on linear regression showed the highest value for heads (0.04) followed by shoots (0.03), stems (0.02) and roots (0.01) (Figure SF 2).
Members of the fungal community uctuate with plant age Relative abundances of the most common genera showed that root fungal communities were highly different from the communities in above-ground plant organs. While comparing communities from leaves, heads and stems, the communities were overlapping, although relative abundances was different in the different tissues (Figure 4), and we found several fungal OTUs that were shared among roots and above-ground tissues (Figure SF 3). The temporal uctuations of fungal read abundances at genus level in the different tissue types was highly dynamic (Figure 4). Con rming visual assessments (ST3), relative read abundance of Puccinia (OTU_1) was peaking between weeks 8 and 13 ( Figure 1). Cladosporium, Alternaria and Vishniacomyza were found in all tissues, but Cladosporium dominated in the leaves and heads, while Alternaria was dominating in the later growth stages in aerial parts. At senescence, roots were colonized by members of the genus Falciphora, a genus of dark septate fungi that may act as growth promoters [28]. Rhyncosporium secalis which is the causal organism of leaf blotch on barley [5] was dominant on stems in weeks 11 and 12. The basidiomycete yeast Sporobolomyces, which is often found in cereals [6] was highly abundant in the heads in weeks 12-14.
Next, we used a Random Forest (RF) approach to model plant age as a function of fungal composition (ST 4). The association between plant age and fungal composition was strongest in the shoots explaining 72% of the variation, followed by 61% in roots, 63 % in heads and 51% in stems. We then used 75% of the data as the training dataset, and 25% for the test dataset, for making actual vs prediction age scatter plots. Our RF model successfully predicted plant age based on the mycobiomes in the different plant parts, with the strongest predictions in shoots (R 2 =0.97) ( Figure SF 4). In roots, the fungal genera Olpidium and Articulospora were the most important taxa contributing to the age-related variation, while Dioszegia and Phaeosphaeriaceae scored highest in the shoots ( Figure 5). Similarly, Pyrenophora and Alternaria were the most important in stems, and nally, Filobasidium and Epicoccum were the most important features in heads. Feature OTUs for plant age were characteristically varying in relative abundance over time. For instance, the relative abundance of Olpidium (OTU_54) in the roots were gradually increasing until week 4 and then decreasing. On the contrary, a steep increase in the relative abundance of Falciphora (OTU_10), and Periconia (OTU_21) was observed in roots ( Figure 5).

Fungal networks in shoots and roots
To elucidate fungal co-occurrences across the plant growth cycle, we visualized co-occurrences in network graphs using weekly data for the leaf ( Figure 6A) and root samples (SF 5). In general, a higher number of co-occurrences were found in roots throughout the plant growth stages. In leaves, cooccurrence networks were weak early in the season but increased dramatically at week 12 ( Figure 6B). In the late stages, there was a notable relatively constant percentage of negative co-occurrences (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30) in the leaves. In the roots, on the contrary, numbers of fungal co-occurrences, both positive and negative were found to uctuate through the season. Generally, highly connected OTUs we different from week to week. However, Alternaria (OTU_12, OTU_435), Periconia (OTU 21), Myrmecridium (OTU_25) and Capnodiales (OTU_318) were found to be highly connected in at least three time points (ST 5).

Discussion
Recent studies of plant microbiomes have deepened our understanding of the development of microbial diversity in plants [29]. While plant genotype, inoculum sources, and environmental factors are important factors shaping the plant microbiome, the assembly of microbial communities in different tissues during an entire plant growth cycle has not been given much attention. We monitored weekly changes in fungal communities during an entire life cycle from seedling emergence to grain maturity in roots as well as leaves, stems and heads of barley. Pathogens are prominent members of the plant microbiome, and our study highlighted the dynamics of important barley pathogens and their relation to the associated microbial communities. Despite the growing awareness of the role of microbial communities in plant health and productivity, the effects of these communities on disease progression are not well understood. It is well known, however, that some microorganisms act as antagonists while others promote pathogens or act as secondary pathogens [30]. Therefore, not only the pathogens themselves should be the subject of study of plant disease biology, but also the plant-associated microbiota. As a rst step toward a holistic view of pathogens existing and acting in a highly diverse and competitive microbial community, we characterized the barley microbiome and created a baseline atlas of fungi associated with different tissues of barley during an entire growth season and during a severe brown rust epidemic.
Fungal disease progression during the barley growth-cycle. In the sampling season 2015, a severe epidemic of brown rust, P. hordei, was recorded in Denmark [31] and was also observed in our experimental plots starting from week 8 in the susceptible cultivar Quench, week in Cumuilus and from week 13 in the resistant Milford. These observations were supported by a dramatic increase in sequence read numbers identi ed as belonging to P. hordei already from week 8 in the leaves of all three cultivars, suggesting that metabarcoding is more sensitive for detecting early pathogen infections compared to visual assessments of symptoms. The cultivar differences of observed brown rust symptoms were not re ected in similar differences in read abundances, indicating that cultivar differences were not dependent on fungal biomass, but only in symptom expression. The massive infections observed by P. hordei caused a decrease in the relative read abundances and diversity of other fungi, which could be caused by a displacement of other species in host tissues by P. hordei, although technical artefacts caused by the relative read abundance measurements likely also contributed to the decrease in diversity [32]. We observed an increase in alpha-diversity indices in the leaves after shoot emergence until week 7 coinciding with the start of the P. hordei epidemic. This could be caused by initial colonization of the native leaf surface by deposition of spores from the air [33]. After the epidemic, diversity increased until senescence, likely caused by the displacement of P. hordei by re-colonization by other taxa.
Plant tissue type forms the mycobiome. Alpha diversity was highest in roots compared to the aboveground tissues, as also observed in other crops [26,34]. A large number of fungal OTUs were shared among root and shoot tissues, which is consistent with earlier studies [35]. We further demonstrated that this overlap between below-and above-ground tissues exists from the early crop development and sustains until late plant maturity.
We found that barley genotype had negligible effects on fungal communities, explaining only 1 -3 % of the variation in the different tissues. The role of plant genotype at ecotype/cultivar level is disputed and is likely associated with the level of genetic heterogeneity among genotypes. We used high-yielding commercial cultivars adapted to Northern European conditions, and thus genotypes were highly similar. Detailed studies of host genetic control of microbiome assembly has only recently been initiated [36,37]. For instance, in a recent study in wheat using cultivars having different pathogen resistance genes, it was found that genotype was important for determining the leaf-associated microbiome [38]. In switchgrass, on the contrary, plant ecotype was only a minor factor determining foliar fungal endophyte communities [39], as was also the case in the wild perennial plant Boechera stricta [35] and in maize [40]. However, these different effects could be highly dependent on the presence of e.g. pathogen resistance genes.
Plant tissue type was explaining 17% of the variation in the total dataset, and in PCoA plots communities clearly clustered according to tissue type. Communities from roots were highly different from the above ground tissues and this effect was particularly strong at later growth stages. Although microbial communities assemble at the early stages of growth already before emergence of seedlings [1,7], they have not stabilized yet, as could also be observed by the loose clustering in the rst weeks of the microbiomes in the PCoA plots. More tightly clustered communities formed at later stages, suggesting that initial assembly of communities is stochastic, while deterministic formation of the communities is occurring at later stages of plant development determined by tissue-speci c factors [41]. Temporal development of fungal communities in root and shoot tissue has not yet been reported in detail yet although several studies show clear differences in roots vs shoots communities at certain growth stages of crop development [11,26]. The type of habitat in leaves, stems and heads determines community structures, and fungal members that have specialized habitat preferences depending on growth habits [42], for instance, Epicoccum and Alternaria evolved to colonize certain tissues [43].
Temporal dynamics of mycobiota in barley. As the plant develops and differentiates into various plant tissues, its suitability as a host for microorganisms changes and, thus, niches will be created while others will disappear [29]. In line with this, our results demonstrated highly plastic communities during plant growth in the different plant tissues. For instance, root alpha diversity increased from week 1 to 11, as early growth stages are essential for microbial community assembly [7] [44,45]. However, in our experiments we only observed minor changes at the earliest sampling weeks most likely caused by assembly already before seedling emergence. Community structure analysis showed temporal decay of dissimilarity of the fungal communities in barley. The time-decay in fungal communities may overlap with the plant tissue specialization, and formation of microbe-speci c niches in the plant [46]. We found stronger time-decay relationships in the above-ground tissues meaning that the rate of community turnover was higher in those tissues. The constant deposition of fungal spores from air and the considerable variation in temperature, wind and humidity may explain this [47,48]. Further, the rate of time decay (slope of the regression line) was highest in the heads. Cereal head development is known for the rapid tissue development comprising cascades of events from emergence, owering to grain lling which could explain the rapid turnover rate of fungal communities. Weaker but signi cant time decay rates were seen in roots indicating that root mycobiota were quickly stabilizing and preventing new microbial members to become established [49]. To pinpoint OTUs that were the most responsive to plant growth and age, we used a machine learning approach. Machine learning methods such as Random Forest modelling is widely being used to gain insight into environmental or host characteristics in both human [50], and plant microbiome studies [8]. Here we used an RF model for age prediction of OTUs and we identi ed age-predicting OTUs speci c for each tissue type with the strongest prediction in leaves. A reason for this could be that leaves are exposed to air-borne inoculum in which concentrations are subject to e.g. plant cover differences in the landscape during the season. In contrast, soil represents a more constant microbial pool during the growth season.
Microbial co-occurrence networks. We assessed and visualized fungal inter-relations in co-occurrence networks. Microbial networks are widely used to infer putative microbial interactions, and to predict the most highly interacting members of the community. Our results demonstrated co-occurrences of fungi in all the tissues and networks that were highly dynamic according depending on plant age. Not surprisingly, the networks in the root-associated microbiomes had a higher number of co-occurrences, caused by the higher number of taxa associated with the roots. In the leaves, co-occurrence networks were relatively weak at the early growth stages becoming stronger at the later growth stages. Remarkably, there was a higher number of negative co-occurrences at the later stages of growth. We speculate that this could be caused by increasing competition among fungi for declining food resources in senescing tissues [51] [6,52]. In the roots, on the contrary, fungal co-occurrence networks were relatively constant across barley growth stages, probably again re ecting the stable soil environment.

Conclusion
Mycobiome development in cereal crops is a critical component to enhance our understanding of crop health and productivity. This study demonstrated that fungal communities during barley growth and development are highly dynamic with a strong time-decay relationship. The rate of temporal decay was highest in above ground parts, which was further supported by the fact that community turnover is higher in shoots compared to roots. We further showed that pathogen infection had profound effects on fungal communities. Finally, using random forest modelling, we pinpointed tissue-speci c fungal taxa that were playing important roles for community turnover across the temporal scale. To our knowledge, this is the most comprehensive study of the barley mycobiome to map temporal shifts from seedling emergence to seed maturity and will thus serve as an important baseline for future studies on the cereal mycobiome.   Alpha diversity estimated using observed richness and Shannon diversity of fungal communities in roots, leaves, stems and heads from seedling emergence (week 1) to harvest (week 18). The alpha diversity of each sample is indicated with a dots and the regression curve is indicated with the line and shade showing 95% con dence interval.   The mean squared error percentage for random forest model trained to predict plant age and top 10 the most in uencing OTUs for plant age in a) root, b) shoot, c) stem and d) heads are shown.