Grapevine Leaf Epiphytic Fungal and Bacterial Communities Are Influenced More by Spatial and Temporal Factors than by Powdery Mildew Fungicide Spray Programs

The plant phyllosphere comprises all aboveground parts of plants, such as leaves, stems, flowers, and fruits (Leveau 2019). In recent years, phyllosphere-associated microbiota have received increased attention for their crucial roles in plant disease tolerance, abiotic stress tolerance, host growth, and fruit quality (Gong and Xin 2021; Lindow and Brandl 2003; Rastogi et al. 2013; Vorholt 2012). In contrast to the soil rhizosphere microbiome that has been the primary focus of plant microbiome research, phyllosphere microorganisms exposed to the open environment encounter more challenges in surviving, resulting in their unique functions in metabolic potential, biocontrol, and climate change adaption (Bashir et al. 2022; Zhu et al. 2022). The plant leaf surface is a major component of the phyllosphere. Although the leaf is considered a relatively low nutrient habitat for microorganisms, the foliar microbial communities exhibit dynamic diversity and composition in response to the fluctuating environment (Pandiyan et al. 2021), leading to intricate microbe–microbe and plant–microbe interactions that can either benefit or harm plant health (Chaudhry et al. 2021).

Leaf epiphytic microorganisms are highly affected by environmental conditions such as temperature, moisture, and UV radiation (Leveau 2019). In addition to environmental factors, chemicals such as fungicides can alter the microenvironment on the leaf surface and further influence phyllosphere microbial communities (Leveau 2019). Although chemical applications are crucial to manage diseases and reduce crop yield loss, they may impact nontarget and potentially beneficial microorganisms that occupy the same niche in the plant habitat. Studies have shown that commonly used fungicides can negatively affect the phyllosphere fungal community diversity and/or richness of different crops, such as almond, soybean, and maize (Noel et al. 2022; Schaeffer et al. 2017), but have mixed effects on the fungal communities of wheat in that several saprotroph fungi rather than common wheat pathogenic fungi were significantly affected by fungicide treatment in their relative abundance (Karlsson et al. 2014; Knorr et al. 2019). It is, therefore, essential to understand the effects of chemical applications on the nontarget microorganisms in the perennial host grape (Vitis vinifera).

California is the largest producer of grapes in the United States, with more than 90% harvested from vineyards in the Central Valley (Alston et al. 2020; Geisseler and Horwath 2016). Among major grapevine diseases, powdery mildew is the most prevalent, causing heavy yield losses if not managed properly (Bois et al. 2017). The causal agent, Erysiphe necator (Schwein.), infects the aboveground plant tissues, including leaves, rachises, canes, and berries. The pathogen has season-long development as long as weather conditions favor its growth. Fungicide applications to control the pathogen include 10 to 14 annual sprays of systemic and contact chemistries with single and/or multisite modes of action (Bettiga et al. 2013; FRAC 2022; Kunova et al. 2021). How or if these commonly used fungicides for powdery mildew have nontarget effects on other grapevine leaf fungal and bacterial community compositions remains to be further investigated, particularly for the pathogen Botrytis cinerea, which causes gray mold and Botrytis bunch rot in grapes, and in which fungicide resistance development has been widely reported (Mlikota Gabler et al. 2003; Saito et al. 2019). Another nontarget fungal genus of interest is Alternaria, which was commonly found to be the dominant taxon on grapevine leaves (Gobbi et al. 2020; Pinto et al. 2014). As the impact of powdery mildew fungicides on Botrytis and Alternaria has not been fully assessed, we aimed to examine the compositions of these fungal genera as well as the entire fungal community under the fungicide treatments.

In recent years, the number of studies on grapevine-associated microbiomes has increased with the fast development of high-throughput next-generation sequencing and bioinformatic tools (Morgan et al. 2017). Studies in multiple vineyards have shown that microbial community diversity shifts, and the richness often decreases over the growing season, possibly due to the plant health conditions affected by the temporal fluctuations, chemical applications, and the production of phytoalexins (Fort et al. 2016; Pinto et al. 2014). In addition, the epiphytic and endophytic microbial compositions are influenced by plant phenological stages such as flowering, fruit set, veraison, and harvest (Carmichael et al. 2017; Liu and Howell 2021; Pinto et al. 2014), as well as plant tissue types such as root, leaf, fruit, cane, and sap (Deyett and Rolshausen 2020; Leveau and Tech 2011; Liu and Howell 2021; Zhang et al. 2017). Grape host genotype has also exhibited a role in microbiome differences, which could affect fruit quality and disease resistance (Singh et al. 2018a, b; Vionnet et al. 2018). Zhang et al. (2019) tested nine wine grape varieties with different sizes, growth patterns, and skin types in berries, and they found that the microbial community diversity was significantly higher in red wine than in white wine grape varieties.

On the other hand, Nerva et al. (2019) showed that grape cultivar Moscato was less susceptible to powdery mildew than cultivar Nebbiolo, and the two cultivars displayed different phyllosphere compositions under an elicitor treatment, suggesting that the host defense mechanisms could be genotype dependent. Moreover, the regional effect, often termed “terroir” in grapes, which could be attributed to soil and weather conditions, may play a more important role in shaping the phyllosphere microbiome (Carmichael et al. 2017; Fort et al. 2016; Liu et al. 2019, 2020; Singh et al. 2018b). A study showed that the geographical regions of vineyards distinct in temperature and terrain exerted a greater influence on berry fungal communities than grape cultivar or harvest time. In colder and higher altitude regions, fungal diversity and abundance were observed to be lower than in warmer and flat regions (Kioroglou et al. 2019).

Agricultural management systems also contribute to variances in the microbial community structure in grapevine leaf and berry sites (Castañeda et al. 2018; Kecskeméti et al. 2016; Morrison-Whittle et al. 2017; Setati et al. 2015; Varanda et al. 2016). However, studies on the effect of specific pesticide treatments on grape microbiomes are inconclusive in their results. For example, Čadež et al. (2010) tested iprodione (FRAC 2), pyrimethanil (FRAC 9), and cyprodinil + fludioxonil (FRAC 9 and 12) on grape berries and found that yeast communities differed in composition among treatments and had larger population sizes in the chemical-treated grape berries than in the nontreated control, but, overall, the treatment impact was not significant. In addition, the synthetic fungicide penconazole (FRAC 3) and inorganic fungicide copper (FRAC M1) have been compared with a biocontrol agent for their effects on grapevine leaf epiphytic microbiota, but both studies showed a minor impact of the treatments on nontarget bacterial and fungal communities (Gobbi et al. 2020; Perazzolli et al. 2014). Conversely, the application of a fungicide containing the active ingredient fenhexamid (FRAC 17) was shown to reduce the richness and diversity of wine juice microbial communities (Escribano-Viana et al. 2018).

In the present study, we examined the impact of powdery mildew spray applications in vineyards on foliar fungal and bacterial communities, including pathogenic, beneficial, and nonpathogenic microorganisms. Our specific aims were to (i) investigate the richness and diversity of fungal and bacterial communities under two different fungicide rotational programs; (ii) compare the composition and abundance of leaf epiphytic fungal and bacterial communities in response to the site, season, and fungicide treatment; and (iii) explore shifts in the relative abundances of fungal genera of interest, including Alternaria, Botrytis, and Erysiphe.

Materials and Methods

Vineyards, fungicide spray programs, and sampling of leaves

Two established commercial vineyard sites in Tulare County in Central Valley California planted with the table grape cultivar Autumn King were treated with the same two powdery mildew spray programs in 2020. The two locations were designated as L1 (located near 36°N, 119°W) and L2 (located near 35°N, 119°W). The L1 vineyard had a T-shaped trellis system with a weedy ground cover, whereas the L2 vineyard had a Y-shaped trellis system with no ground cover (Supplementary Fig. S1). According to the Web Soil Survey of the U.S. Department of Agriculture, the L1 vineyard is on Calgro sandy loam soil, classified as a coarse-loamy, mixed, superactive, thermic Typic Durixerepts, whereas the L2 vineyard is on a Colpien loam, a fine-loamy, mixed, superactive, thermic Calcic Pachic Haploxerolls. The nearest California Irrigation Management Information System stations from L1 and L2 recorded average daily temperature ranges of 59.4 to 87.3°F (15.2 to 30.7°C) and 60.3 to 88.4°F (15.7 to 31.3°C) and average daily relative humidity of 26 to 73% and 29 to 65%, respectively, for the duration of the sampling.

Two rotational spray programs (treatments) were applied at both vineyards with similar timelines according to the plant phenological stage (Table 1). Spray programs comprised five applications of site-specific synthetic fungicides from the pre-bloom stage until the grape cluster closure stage to control powdery mildew. The common names of the active ingredient in these fungicides included myclobutanil (DMI; FRAC 3), tebuconazole (DMI; FRAC 3), fluoropyram (SDHI; FRAC 7), trifloxystrobin (QoI; FRAC 11), quinoxyfen (quinoline; FRAC 13), metrafenone (aryl-phenyl-ketone; FRAC 50), and cyflufenamid (phenyl-acetamide; FRAC U6) (Table 1). Each of the two treatments was performed in four blocks arranged in a randomized complete block design in the vineyard. Each block represented a replicate and included four rows.

TABLE 1 Powdery mildew seasonal spray programs conducted in two vineyards in the Central Valley of California in 2020^a

Variable	S0	First spray	S1	Second spray	S2	Third spray	S3	Fourth spray	S4	Fifth spray	S5
Phenological stage		Prebloom		Bloom		Fruit set		Pea size		Cluster closure
Treatment 1		Cyflufenamid (U6)		Fluopyram (7) + trifloxystrobin (11)		Cyflufenamid (U6)		Fluopyram (7) + trifloxystrobin (11)		Myclobutanil (3)
Treatment 2		Metrafenone (50)		Quinoxyfen (13)		Cyflufenamid (U6)		Fluopyram (7) + tebuconazole (3)		Quinoxyfen (13)

^a The corresponding grapevine phenological stages and the sampling time points (S0 to S5) are shown. Active chemicals in the fungicide products are listed instead of the product names. The FRAC code for each fungicide is shown in the parentheses.

View as image HTML

Throughout the growing season, from the pre-bloom stage (May) to the preharvest stage (September), leaf samples were collected at six time points according to the application of the synthetic fungicide spray: S0 (before the first spray), S1 (after the first spray), S2 (after the second spray), S3 (after the third spray), S4 (after the fourth spray), and S5 (after the fifth spray). Generally, there was a 2-week interval between fungicide treatments, and an additional sulfur spray was applied during this time. The samples were collected around the seventh day after the fungicide spray and before the sulfur spray. For each replicate block, per treatment, per sampling time point, 16 to 24 grapevine leaf samples were collected randomly within the inner edge of the second row (between row 2 and row 3) to minimize exposure to other treatments in neighboring rows. The sampled leaves were collected at the L3 to L4 growth stage (Chitwood et al. 2014) and were not selected based on disease symptoms but were randomly collected from each treatment. The collected leaves were bagged by replicate and put in a cooler with ice packs before being transferred to the lab. A total of 96 samples were collected, with 48 samples from each vineyard.

Processing of leaf samples

Leaf samples were stored at 4°C in the lab and processed within 2 days of the collection date. To collect epiphytes, we used a previously described leaf washing method (Singh et al. 2018a, b), with slight modifications. Briefly, one to five leaves, depending on the leaf size, were loosely set inside a 50-ml propylene tube with sterile tweezers to generate space among leaf surfaces. Tweezers were sanitized between replicates. Leaves in individual tubes were covered with a washing solution of 0.15 M sodium chloride solution containing 0.01% Tween 20. Next, the tubes were set in a horizontal shaker (New Brunswick Scientific, Edison, NJ, U.S.A.) at 150 rpm for 30 min, followed by another 30 min on a vertical tube rotator (Scilogex, Rocky Hill, CT, U.S.A.). After shaking, tubes were transferred to an ultrasonic bath (Fisher Scientific, Hampton, NH, U.S.A.) for 10 min. All of the above steps were performed at room temperature. Leaves were then removed with sterile tweezers, and the tubes were centrifuged (Beckman Coulter, Pasadena, CA, U.S.A.) at 4,500 × g for 30 min at 4°C. After centrifugation, ∼5 ml of supernatant was retained per tube, discarding the remainder, and re-suspended and pooled with samples from the same replicate. Merged samples were centrifuged again at 4,500 × g for 20 min at 4°C. The supernatants were discarded until only approximately 1.5 ml remained in each tube. The microbial pellets were re-suspended, transferred to 2-ml microcentrifuge tubes, and stored at −20°C until DNA extraction.

DNA extraction and amplicon sequencing preparation

Genomic DNA extractions were carried out using a ZymoBiomics DNA MicroPrep Kit (Zymo Research, Irvine, CA, U.S.A.). The samples were first centrifuged at 15,000 × g for 10 min at 4°C to concentrate microbial pellets in ∼250 ul, and the extraction proceeded following the manufacturer's manual. The resulting 96 DNA samples were submitted to Microbiome Insights (Vancouver, BC, Canada) to complete the following processes. DNA extracts were purified with a Genomic DNA Clean & Concentrator kit (Zymo Research) and quantified with a Qubit dsDNA HS Assay kit (Invitrogen, Thermo Fisher Scientific, Waltham, MA, U.S.A.). PCR and library preparation were conducted targeting either the prokaryotic 16S rRNA V4 region with 515F/806R primers (515F: 5′-GTGYCAGCMGCCGCGGTAA-3′; 806R: 5′-GGACTACNVGGGTWTCTAAT-3′) (Apprill et al. 2015; Parada et al. 2016) or fungal ITS2 rRNA region with fITS7/ITS4 primers (fITS7: 5′-CCGTGARTCATCGAATCTTTG-3′; ITS4: 5′-CCTCCGCTTATTGATATGC-3′) (Ihrmark et al. 2012; White et al. 1990). Amplicon sequencing was performed using Illumina 2 × 250 bp MiSeq (Illumina, San Diego, CA, U.S.A.).

Bioinformatic and statistical analyses

The amplicon sequence data were first examined and processed using the QIIME2 pipeline (v 2021.11) (Bolyen et al. 2019). Paired-end reads were demultiplexed with the q2-demux function. The reads were denoised, trimmed/truncated, and merged with DADA2 (Callahan et al. 2016) via the q2-dada2 function to filter out low-quality (quality score < 18) and chimeric reads. The parameters were as follows: The first five nucleotides were trimmed in both ITS and 16S sequences. In ITS, the forward and reverse reads were truncated to 243 and 152 bp, respectively. In 16S, the forward and reverse reads were truncated to 249 and 184 bp, respectively. The retained amplicon sequence variants (ASVs) were aligned, and the phylogenetic trees were established with MAFFT (Katoh and Standley 2013) and FastTree2 (Price et al. 2010) via the q2-phylogeny align-to-tree-mafft-fasttree function. ASVs were clustered into operational taxonomic units (OTUs) via the q2-vsearch function (Rognes et al. 2016). For ITS, an open-reference clustering method was used with 99% identity against the UNITE reference database (v 8.3) (Abarenkov et al. 2021; Nilsson et al. 2019), and for 16S, a de novo clustering method was used with a 97% identity threshold. The taxonomic classifiers were trained with the q2-feature-classifier fit-classifier-naive-bayes function using the UNITE (v 8.3) and RDP (v Taxonomy 18) reference databases for ITS and 16S, respectively (Abarenkov et al. 2021; Cole et al. 2014; Nilsson et al. 2019; Wang et al. 2007). The resulting OTU taxonomic table, abundance table, representative sequences, and rooted tree data from QIIME2 were then imported to R (v 4.1.3) (R Core Team 2022) using R package qiime2R (v 0.99.6) (Bisanz 2018).

Data processing, rarefaction curve, and alpha/beta diversity analyses were conducted using R packages phyloseq (v 1.40.0) (McMurdie and Holmes 2013), vegan (v 2.6-2) (Oksanen et al. 2022), and microeco (v 0.13.0) (Liu et al. 2021). Table data output and data visualization were performed using R packages ranacapa (v 0.1.0) (Kandlikar et al. 2018) and ggplot2 (Wickham 2016), respectively. Before rarefaction, the 16S OTUs taxonomically identified as “mitochondria” or “chloroplast” were removed from the dataset. Rarefaction was performed using the minimum sampling depth (25,265 for ITS and 5,270 for 16S) to ensure even coverage among all the remaining samples. Alpha diversity for the richness and evenness of local microbial communities was measured by the Shannon diversity index (Shannon 1948). The Shapiro-Wilk test was conducted to test the normality at α = 0.05 for the alpha index data, and the Kruskal-Wallis test (α = 0.05) was performed for non-normal distributed data to calculate the statistical differences. When significance was observed (P < 0.05) among sampling time points, a multiple pairwise comparison was performed with the Wilcoxon rank sum test. Beta diversity was assessed using a Bray-Curtis distance matrix (Bray and Curtis 1957) to test the dissimilarity of microbial communities between groups in locations, sampling time points, and treatments. Principal coordinates analysis (PCoA) ordination (Gower 1966) was conducted to visualize the patterns and identify microbial community groupings along with 95% confidence intervals. Permutational multivariate analysis of variance (adonis test) (Anderson 2001) with 999 permutations based on the Bray-Curtis distance matrix was calculated to test whether the compared groups were compositionally distinct (α = 0.05).

Comparisons across taxonomic composition and relative abundances were analyzed using the heatmap function in R package ampvis2 (Andersen et al. 2018) at the taxonomic levels of class/genus for fungal communities and phylum/class for bacterial/archaeal communities. To visualize the flow of taxa changes over sampling time points, alluvial plots were conducted using R package ggalluvial (Brunson and Read 2020), showing the 10 most abundant fungal genera and bacterial/archaeal classes. We also tested the abundance differences for fungal genera of interest (i.e., Alternaria, Botrytis, and Erysiphe) across sampling time points and between treatments in the two locations. One-way ANOVA (α = 0.05) was conducted using JMP software (v 13) (SAS Institute, Cary, NC, U.S.A.), followed by the post-hoc Tukey HSD test for multiple pairwise comparison among the sampling time points when significance was observed.

To gain insights into the ecological roles of the fungal communities, we assigned functional guilds to fungal genera using the FungalTraits database (Põlme et al. 2020). To test the abundance differences of each functional guild among the sampling time points for each individual vineyard, we used one-way ANOVAs (α = 0.05) and post-hoc Tukey HSD tests using JMP software (v 13) (SAS Institute).

Co-occurrence network analysis was conducted to investigate the potential connection and correlation within the bacterial/archaeal or fungal communities in individual vineyards. The iNAP (Feng et al. 2022) web-based pipeline was used to construct and analyze networks with the Molecular Ecological Network Analyses framework (Deng et al. 2012). Briefly, OTUs detected in fewer than 80% of the samples were first filtered out. The selected majority was calculated using Spearman's rank correlation among OTUs to create the adjacent correlation matrix. The random matrix theory approach (Luo et al. 2007) was applied to obtain a cutoff threshold value for the network. For ITS, the cutoff values were 0.65 and 0.68 for vineyards L1 and L2, respectively. For 16S, the cutoff values were 0.65 and 0.61 for vineyards L1 and L2, respectively. The network topological attributes were calculated by the iNAP pipeline and then visualized using Cytoscape software (v 3.9.1) (Shannon et al. 2003). Network nodes were labeled at the genus level for ITS OTUs and order level for 16S OTUs to look for potential keystone taxa in each local network.

The amplicon sequencing raw data were submitted to the National Center for Biotechnology Information Sequence Read Archive database. The BioProject ID is PRJNA899094.

Results

Sequence analysis

Illumina paired-end sequences resulted in a total of 5,923,803 reads for ITS and 2,714,623 reads for 16S rRNA genes across all samples. After quality filtering, 4,267,669 reads remained for ITS with 4,204 unique ASVs, and for 16S rRNA, 1,820,759 reads remained with 2,390 unique ASVs. The median number of sequences per sample was 42,478.5 and 18,699.5 in ITS and 16S rRNA, respectively. A total of 95 samples passed quality checks and were kept; one sample from location L2, sampling time point S4, under treatment 2 contained a low number of quality reads and was excluded from analyses. After OTU clustering and removing mitochondrial and chloroplast sequences, another sample, “Location L1/Sampling time point S4/Treatment 1,” containing mostly chloroplast reads in the 16S dataset, was excluded specifically for 16S data analysis. In total, 3,033 and 1,303 OTUs for ITS and 16S rRNA were obtained, respectively. The rarefaction curves for the observed species showed that all samples had reached a plateau phase at the sampling depth of 25,265 reads for ITS and 5,270 reads for 16S, indicating the sufficient sampling size (Supplementary Fig. S2). After rarefaction, a total of 2,877 OTUs in ITS and 1,218 OTUs in 16S rRNA were applied for the downstream diversity analyses.

Diversities and dissimilarities of microbial communities

The Shannon diversity indices, representing the richness and the evenness of communities, were not significantly different based on location or by spray treatments in either fungal (Supplementary Fig. S3A and B) or bacterial/archaeal (Supplementary Fig. S3D and E) communities. However, the sampling time point had a significant effect on the alpha diversity of both communities (P < 0.001) (Fig. 1). In fungal communities, the earliest sampling time point, S0, showed the lowest diversity, whereas the latest sampling time point, S5, displayed higher diversity compared with other sampling points (Fig. 1A). On the other hand, bacterial/archaeal communities had a relatively reverse trend, where S0 displayed the highest alpha diversity and significantly differed from other time points (except S3) (Fig. 1C). The results suggest that the grapevine epiphytic microbial community diversity shifted temporally with bacterial and fungal communities demonstrating contrasting trends.

Download as PowerPoint

A permutational multivariate analysis of variance test for beta diversity based on Bray-Curtis distances revealed that location and sampling time point factors significantly influenced the compositional variances for both fungal and bacterial/archaeal communities (P < 0.001), whereas treatments (i.e., fungicide spray programs) did not significantly impact the structure of either fungal or bacterial/archaeal community composition (Tables 2 and 3, respectively). For fungal communities, PCoA explained 30 and 14.6% of the fungal community variances across axis 1 and 2, respectively (Fig. 1B). Fungal community dissimilarities based on location and sampling time points were distinct, with communities in the locations diverging more with time (Fig. 1B). Conversely, samples did not distinctly cluster based on treatment (Supplementary Fig. S3C). For bacterial communities, the PCoA plots described 23.9 and 14% of the variability in axis 1 and axis 2, respectively (Fig. 1D). Similar to the fungal community, we found that location and time point both impacted bacterial/archaeal community dissimilarity (Fig. 1D), but spray treatment did not (Supplementary Fig. S3F).

TABLE 2 Permutational multivariate analysis of variance (PERMANOVA) (adonis test) with 999 permutations using Bray-Curtis dissimilarity metric in comparison of the effect of location, sampling time, and treatment on fungal community

Variable	df	SS	R2	F	P valuea
Location	1	1.1546	0.1341	23.7977	0.001***
Sampling time	5	3.1613	0.3673	13.0316	0.001***
Treatment	1	0.0701	0.0082	1.4444	0.161
Residual	87	4.2210	0.4904
Total	94	8.6070	1.0000

^a *** indicates significance at the P < 0.001 level.

View as image HTML

TABLE 3 Permutational multivariate analysis of variance (PERMANOVA) (adonis test) with 999 permutations using Bray-Curtis dissimilarity metric in comparison of the effect of location, sampling time, and treatment on bacterial community composition

Variable	df	SS	R2	F	P valuea
Location	1	0.6797	0.1004	14.1549	0.001***
Sampling time	5	1.9357	0.2858	8.0624	0.001***
Treatment	1	0.0278	0.0041	0.5796	0.925
Residual	86	4.1295	0.6097
Total	93	6.7727	1.0000

^a *** indicates significance at the P < 0.001 level.

View as image HTML

Microbial community taxonomic compositions and dynamics

Of the 10 most dominant fungal genera detected across the growing season in the two vineyards, seven belonged to Ascomycota, including genera Alternaria, Ascochyta, Stemphylium, Aureobasidium, Botrytis, Taphrina, and Preussia; three belonged to Basidiomycota, including genera Vishniacozyma, Filobasidium, and Sporobolomyces, all basidiomycetous yeast or yeast-like fungi (Fig. 2A). The most abundant fungal genus was Alternaria, accounting for 9.6 to 24% of the fungal composition in vineyard L1 and 14 to 22.3% in vineyard L2, and its relative abundances tended to increase at the later stage in both vineyards (Fig. 2A; Supplementary Fig. S4A). The genus Botrytis, containing a major grapevine pathogen, was detected at higher levels in L2, revealing a marked increase in relative abundance across the season, from 1.2% at S0 to 12.3% at S5. In contrast, L1 had a much lower relative abundance of Botrytis at every sampling time point, ranging from 0 to 0.2% (Fig. 2A; Supplementary Fig. S4A). For other plant-pathogenic genera, Ascochyta and Taphrina, the relative abundances decreased in both vineyards with time, whereas the relative abundance of Stemphylium increased in L1 and slightly decreased in L2 (Fig. 2A). Preussia, a saprophytic and endophytic fungal genus, slightly increased in the later season (Fig. 2A). The yeast-like fungi Aureobasidium showed a steadier presence over time, accounting for 1 to 3% of community abundance in L1 and 1.8 to 4% in L2 (Fig. 2A; Supplementary Fig. S4A). The three yeast fungi Vishniacozyma, Filobasidium, and Sporobolomyces displayed a similar patten with high relative abundance in the early season, which dropped later in the season (Fig. 2A). Overall, vineyard L1 had a greater relative abundance of yeast fungi than L2. The decline in the relative abundance of members in the genera Vishniacozyma and Sporobolomyces in L1 over time was especially noticeable as the abundance of Vishniacozyma decreased from 12.2% at S0 to 4.2% at S5, and for Sporobolomyces, from 12% at S0 to 1.6% at S5 (Fig. 2A; Supplementary Fig. S4A).

Download as PowerPoint

For bacteria/archaea, the most abundant class Actinobacteria accounted for 20.8 to 33.7% in L1 and 32 to 42% in L2 (Supplementary Fig. S4B). Bacilli, the second most abundant class, was more abundant in L1 (15 to 37.8%) than L2 (12.6 to 25.7%). In addition, Bacilli and the other major Firmicutes bacteria, Clostridia and Erysipelotrichia, revealed a trend of increasing abundance over the growing season (Fig. 2B; Supplementary Fig. S4B). On the contrary, the major classes Alphaproteobacteria, Gammaproteobacteria, and Betaproteobacteria, all belonging to the phylum Proteobacteria, decreased over time, except for Gammaproteobacteria, which fluctuated with time in L1 (Fig. 2B; Supplementary Fig. S4B). The declining shift was also found in the three Bacteroidetes classes, i.e., Cytophagia, Bacteroidia, and Flavobacteriia. Although these Bacteroidetes were the major taxa observed, they accounted for less than 2% abundance in both locations (Fig. 2B; Supplementary Fig. S4B).

When evaluating disparities between abundant taxa, the clear shifts observed for location and time point are not apparent when communities are contrasted based on treatment, consistent with whether bacterial/archaeal or fungal communities were assessed (Supplementary Fig. S5). This suggests that the two powdery mildew spray programs played a minor role in driving the composition of the dominant microbial taxa at the class and genus levels.

Comparison of the abundance of Alternaria, Botrytis, and Erysiphe over time and treatment

Three fungal genera Alternaria, Botrytis, and Erysiphe, known to include common plant pathogens, were investigated for changes in their abundance over sampling time points and treatments. When compared across sampling time points in the two vineyards, only Botrytis showed a significant difference in abundance between the two locations (P < 0.001). Vineyard L2 had a greater abundance of Botrytis compared with vineyard L1, where Botrytis counts were low throughout the growing season (Fig. 3A). In L2, Botrytis showed a significant increase in abundance, especially in the later season (sampling time point S5) (P < 0.001) (Fig. 3A). Alternaria was the most prevalent genus in the fungal microbiota and showed no significant differences in abundance between the two vineyards. However, its growth exhibited a similar pattern in both L1 and L2, where relative abundance was greatest during the middle season (sampling time point S2 to S4) compared with the earlier (S0 and S1) and later season time points (S5) (P < 0.001) (Fig. 3A). Erysiphe (grape powdery mildew genus) was the target fungal pathogen of the fungicide spray programs implemented in the vineyards, and it continuously displayed low abundance across sampling time points in both locations (Fig. 3A). When comparing the abundance of each fungal genus between the two spray program treatments within each location, we found no significant difference across all three genera (Fig. 3B). The abundance of each genus was compared between the two treatments (Trt1 and Trt2) at individual time points within each location. We found that all three genera displayed similar trends in both Trt1 and Trt2 over time, and no significant differences were found between the treatments in each comparison. The only exception was that a greater abundance of Alternaria was detected in Trt2 compared with Trt1 at S4 in L1 (P < 0.05) (Supplementary Fig. S6). This was consistent with the prior observation that individual spray programs did not exert a spray-specific effect on fungal community diversity.

Download as PowerPoint

Functional roles of fungal communities

In addition to the fungal compositional structure, we explored the role that these fungi might play in the leaf microenvironment based on their functional guilds (Fig. 4; Supplementary Table S1). Among functional guilds, the largest portion of fungal microbiota from the grapevine leaves belonged to plant pathogens, and this proportion increased with sampling time. In vineyard L2, plant-pathogenic fungi had a significantly higher relative abundance at S5 compared with other sampling time points (P < 0.001) (Fig. 4). Several saprotroph guilds were identified within the fungal community, with soil saprotrophs accounting for the highest relative abundance in this category, followed by nectar/tap saprotroph. Nectar/tap saprotroph was the most dominant guild in vineyard L1, especially in the early season (S0 and S1), compared with later time points (P < 0.001) and guilds. Fungi categorized as litter saprotroph, wood saprotroph, and pollen saprotroph were also found in both vineyards across the season, although they were a smaller cohort of the fungal community, with only mild shifts in their relative abundances (Fig. 4). Foliar fungal endophyte was another relatively large portion in our collection. Although they proportionally declined across the season in L2, no specific pattern was observed in L1 (Fig. 4). Fungal epiphytes revealed low amounts throughout the growing season in both vineyards (Fig. 4).

Download as PowerPoint

The abundances by Trt1 and Trt2 were also compared at each sampling time point for each location to evaluate the fungicide spray effect on the functional guild distribution over time. The result showed that both treatments in each location displayed a similar pattern in the abundance shift, as shown in Figure 5, and the abundances were not different between Trt1 and Trt2 in most fungal guilds at each sampling time point. Only the category of foliar fungal endophyte had a greater abundance in Trt1 than Trt2 at S3 in L1 (P < 0.05), litter saprotroph was more abundant in Trt2 than Trt1 at S5 in L2 (P < 0.05), and pollen saprotroph was more abundant in Trt2 than Trt1 at S4 in both locations (P < 0.05 for both) (Supplementary Fig. S7).

Download as PowerPoint

Network analysis

We investigated the structure and relationships among the ubiquitous taxa using co-occurrence network interaction patterns. For ITS, after filtering fungal OTUs based on distribution, 63 and 62 OTUs in L1 and L2 remained as network nodes. The L1 network contained 333 links with 134 negative and 199 positive correlations (R² of power-law = 0.7831; average degree = 10.5714; average clustering coefficient = 0.5177), and the L2 network contained 391 links with 230 negative and 161 positive correlations (R² of power-law = 0.6606; average degree = 12.6129; average clustering coefficient = 0.514) (Supplementary Tables S2 and S3). Among the nodes, 17 and 28 connector hubs were identified in L1 and L2, respectively. The most connected node in L1 was an OTU from the genus Vishniacozyma with the maximal degree of 36, whereas an OTU of the genus Papiliotrema with a maximal degree of 31 was the most connected in vineyard L2 (Fig. 5A; Supplementary Tables S2 and S3). In the L1 network, there were two crucial nodes. One node, from the genus Vishniacozyma, displayed the maximal betweenness centrality and closeness centrality, suggesting its importance in connecting other nodes, and the second, from the family Sporormiaceae (genus unidentified), had the highest eigenvector centrality, indicating its influence on other major nodes (Supplementary Table S2). On the other hand, four OTUs were found to play influential roles in the L2 network. An OTU from the genus Podosphaera was the node with maximal betweenness centrality and eigenvector centrality, whereas OTUs from the genera Vishniacozyma, Papiliotrema, Filobasidium, and Podosphaera had the greatest closeness centrality (Supplementary Tables S2 and S3).

Both locations were composed of seven modules (Fig. 5A). In the L1 network, five modules were more closely connected with all the connector hubs included, whereas two modules with only two nodes in each were separated from the larger networking cluster (Fig. 5A; Supplementary Table S2). In one of the small modules, Alternaria, the most abundant fungal genus, was connected only with another dominant genus, Ascochyta, with a positive correlation. Moreover, within the larger module, the powdery mildew pathogen Erysiphe showed a single negative connection with one of the module connector hubs, Papiliotrema (Supplementary Table S3). Erysiphe was not among the major nodes in L2. In the L2 network, four modules forming a larger networking cluster and three small individual modules were observed (Fig. 5A). Alternaria in L2 was still in the same group as Ascochyta but displayed more complex interactions with other nodes within the larger cluster. Botrytis, which only appeared in the L2 network, acted as a connector hub with a node degree of 25, including 14 positive and 11 negative correlations to other nodes. In addition, 20 of the links from Botrytis were connected to other network hubs. Podosphaera, another genus containing powdery mildew pathogen species, was shown in L2 only and acted as a hub node in the network. However, no correlation was found between Botrytis and Podosphaera (Fig. 5A; Supplementary Table S3). Furthermore, a total of six fungal genera were presented in both networks as connector hubs, including five yeast or yeast-like fungi, Cystofilobasidium, Filobasidium, Papiliotrema, and two Vishniacozyma spp., and one grapevine pathogen, Ganoderma (Adaskaveg and Gilbertson 1987) (Fig. 5A).

The overall connections for the 16S network were fewer than those found in fungal networks (Fig. 5B). Fifty and 49 nodes were obtained in L1 and L2, respectively. The L1 network contained 97 links with 33 negative and 64 positive correlations (R² of power-law = 0.8003; average degree = 3.88; average clustering coefficient = 0.312), and the L2 network contained 129 links with 76 negative and 53 positive correlations (R² of power-law = 0.8695; average degree = 5.2653; average clustering coefficient = 0.3683) (Supplementary Tables S2 and S3). Among the nodes, two connector hubs (out 1094 and OTU 697) and one module hub (OTU 551) were identified in L1. The most connected node in L1 was OTU 696 (order Micrococcales) with the maximal degree of 13, and it also displayed the highest eigenvector centrality. The OTU 1094 (order Bacillales) was another influential node with both the maximal betweenness centrality and closeness centrality. In the L2 16S network, six connect hubs were identified. OTU 1104 (order Erysipelotrichales) had the maximal degree of 20 and was the most important node in the network with the highest betweenness centrality, closeness centrality, and eigenvector centrality (Supplementary Tables S2 and S3).

Seven modules were found in both locations in the 16S network (Fig. 5B). In the L1 network, hub nodes OTU 1094 (order Bacillales), OTU 697 (order Micrococcales), and OTU 551 (order Sphingomonadales) were all in the same module. Overall, the two orders Bacillales and Micrococcales accounted for the most nodes in both L1 and L2 networks (Fig. 5B). Noticeably, OTU 1104 (order Erysipelotrichales) was the most influential node in the L2 network, and its class, Erysipelotrichia, was the fourth most abundant taxa in the bacterial/archaeal communities (Fig. 2B; Supplementary Fig. S4B).

Discussion

We investigated how grapevine leaf microbiota were affected under two seasonal powdery mildew fungicide spray programs in two vineyards located in the Central Valley of California using 16S rRNA and ITS2 to decipher bacterial/archaeal and fungal community assemblages, respectively. Our results demonstrated that sampling time points (temporal factor) and vineyard location (spatial factor) rather than fungicide treatments were the driving factors impacting microbial diversity, especially for fungal communities.

The temporal changes in fungal and bacterial/archaeal community compositions observed in both vineyards were likely driven by grapevine phenological stages and environmental conditions. Previous studies have demonstrated the impact of the host growth stage on the dynamic shift of microbial communities on grapevine leaves (Pinto et al. 2014) and berries (Carmichael et al. 2017). For example, Liu and Howell (2021) revealed significant fluctuations in fungal communities under different grapevine developmental stages, plant tissues, and climate conditions. Consistent with previous results, we observed that the sampling time point had the greatest effect on the grapevine leaf microbiome, with increasing microbial diversity and distinct microbial compositions as the season progressed.

Diversity due to different site locations has been reported to be another major influential factor in the shifts of grape microbial composition (Bokulich et al. 2014; Carmichael et al. 2017; Kioroglou et al. 2019). The spatial effect could originate from environmental factors such as soil, agricultural practices, and the surrounding microclimate. For example, the soil microbiome in the vineyard has been associated with changes in grapevine endophytic and epiphytic compositions, with some dominant soil taxa found in varying amounts in the phyllosphere (Liu et al. 2019; Martins et al. 2013; Zarraonaindia et al. 2015; Zhang et al. 2017). In our study, the two vineyards were located within the same county and were located approximately 87 kilometers apart. As the climate conditions (moisture, temperature, and UV index) were similar, the significant differences in microbial diversity and composition found between the two sites were likely associated with locale-specific aspects such as soil composition and production systems. Soil moisture and texture have been shown to shape soil microbial communities and functions (Wang et al. 2020; Xia et al. 2020). Overall, the L1 soil may have had higher soil salinity and/or a hard layer at about 1 m, with moderate drainage. The L2 site has comparatively more clay than the L1 site, is moderately to well drained, and does not have historical soil salinity concerns. These variances in soil texture and likely in soil moisture may impact soil microbial communities that could be the underlying factor driving site variance in the phyllosphere community. Furthermore, different ground-covering features and trellis systems were employed in the two vineyards. There is evidence that grapevine plant health could benefit from under-vine vegetation cover, which moderates the impact of precipitation and improves the soil ecosystem by increasing the diversity of soil microbiota (Chou et al. 2018; Vanden Heuvel and Centinari 2021). Moreover, different kinds of vineyard trellis systems resulted in different canopy architectures and sizes, which created variable grape microclimate environments influencing UV exposure and humidity in the canopy (Yu et al. 2022). The site environmental conditions and agricultural management practices, such as ground cover and trellis system, may have stronger effects than the spray program on the grapevine microbiome. However, further studies are needed to assess the impact of these site-specific factors.

Several grape microbiome studies have shown nontarget effects of fungicide applications, although many focused on wine grapes and the yeast fungi contributing to fermentation. Results were contradictory in terms of the microbial diversity, the composition of the microbial community, and whether the impact of the fungicide application was significant (Čadež et al. 2010; Grangeteau et al. 2017; Sumby et al. 2021). In our study, in contrast to time and location, the powdery mildew fungicide spray programs did not have a substantial influence on the grapevine leaf microbial diversity for fungal or bacterial/archaeal communities. Nevertheless, as we examined the total effect of the whole spray program, the result reflects the overall impact of the season-long program and not the individual effect of a fungicide spray. Moreover, the FRAC groups for the synthetic fungicides used in the spray program Trt1 included FRAC 3, 7, 11, and U6, whereas spray program Trt2 included 3, 7, 13, 50, and U6 (Table 1). Among them, FRAC 7 and 11 belong to the same mode of action that interferes with the cellular respiration. FRAC 3 has a mode of action for inhibiting sterol biosynthesis, FRAC 13 for signal transduction inhibition, FRAC 50 for microtubule assembly inhibition, and U6 for an unknown mechanism (FRAC 2022). The overlap in shared fungicides could have affected community composition, but this seems unlikely because major shifts after the application of FRAC 7, 11, or 3 were not observed. It is also possible that sulfur, which was applied between each fungicide application for both programs, may have a strong effect on the microbial community, effectively masking any program-specific effects. Future work is needed to decipher the effect of the individual fungicides from each FRAC group, as their impact on the grapevine microbiota has not been fully investigated.

Among the dominant fungal genera in our result, the filamentous fungus Alternaria and yeast-like fungus Aureobasidium were commonly reported as the major fungi in the grapevine phyllosphere (Castañeda et al. 2018; Morgan et al. 2017; Pinto et al. 2014; Singh et al. 2018a, b). Alternaria is a ubiquitous genus in the field that includes many species functioning as phytopathogens, saprophytes, and endophytes and was also found to be the main taxa in the grape leaf and berry endosphere (Varanda et al. 2016; Wijekoon and Quill 2021). It has been reported that the endophyte Alternaria could be a potential biocontrol agent for grape downy mildew (Musetti et al. 2006). In our study, the genus was prevalent in the leaf fungal microbiota across the growing season in both locations, with a significant dominance at the middle stage S2 to S4 corresponding to the plant bloom to pea-size berry phenological stage. Aureobasidium, another widely present fungus, displayed a relatively steady presence throughout the growing season. It has been explored for its antagonistic properties against plant pathogens and its production of compounds that benefit the winemaking process (Bozoudi and Tsaltas 2018; Pinto et al. 2018). Interestingly, the overall fungal composition profile obtained in our study was mostly in accordance with Singh et al. (2018a, b), where Ascochyta, Botrytis, Filobasidium, Sporobolomyces, Stemphylium, Taphrina, and Vishniacozyma were also reported as the core fungal taxa in grapevine leaf microbiota across several grape varieties in France. Despite the distinct continental regions and grapevine genotypes, the similar taxonomical compositions could be due to strong driving forces from the “agro-climatic zone” that combines general natural environment factors and human practices (Singh et al. 2018b).

Botrytis, the genus containing Botrytis cinerea, a pathogen causing grape bunch rot and gray mold, was found to substantially increase in vineyard L2 over time, whereas its abundance was consistently low in L1. These abundance differences were not affected by treatment, confirming that spatial and temporal effects were the major drivers for microbial composition changes instead of fungicide spray programs. However, Erysiphe, the target of the fungicide sprays, revealed low presence throughout the season, suggesting the effectiveness of the spray programs at reducing Erysiphe necator. Other major plant pathogens from the top 10 taxa included Ascochyta, which causes Ascochyta blight in chickpea, lentil, and pea crops (Kaiser 1997; Pande et al. 2005); Stemphylium, a causal agent of Stemphylium leaf blight with a broad host range (Hay et al. 2021); and Taphrina, which is known for host plant deformity disease, such as peach leaf curl (Tsai et al. 2014). All of these pathogens have been rarely reported in grapes, and their roles in grapevine health are still unknown. Preussia, a saprophyte and endophyte found in soil and plant debris, has been reported to promote plant growth by producing nitric oxide, indole acetic acid, and gibberellins (Al-Hosni et al. 2018; Khan et al. 2016; Mapperson et al. 2014). Whether it has benefits for grapevine growth remains to be investigated.

Yeast and yeast-like fungi have been another focus when it comes to the grape microbiome due to their role as biocontrol agents in plant health and metabolic activity in wine fermentation (Di Canito et al. 2021; Setati et al. 2012). For example, some Aureobasidium, Candida, and Metschnikowia species were used to control Botrytis and Penicillium (Spadaro and Droby 2016). In our study, yeast Vishniacozyma and Sporobolomyces and yeast-like Aureobasidium and Filobasidium were identified as the main taxa. These genera have been reported as antagonistic or potential biocontrols against Botrytis cinerea (Carmichael et al. 2019; Gramisci et al. 2018). In both vineyards L1 and L2, we observed a general decreasing trend of Vishniacozyma, Filobasidium, and Sporobolomyces and an increasing trend in plant pathogens Botrytis in L2 and Stemphylium in L1 as the season progressed. Whether the abundance of these yeast/yeast-like fungi has a direct association with the abundance of plant pathogens in the vineyard needs to be further examined and confirmed with culture-dependent assays.

Ecological guilds represent groups performing similar functional roles in an ecological community and, therefore, can be useful for inferring and predicting the functional characteristics of the microbial community (Djemiel et al. 2022; Simberloff and Dayan 1991). With the FungalTraits analysis, we found that, overall, the majority of the fungal taxa were assigned to the functional categories of plant pathogens and saprotrophs. However, the guild assigned for each fungal taxa was based on its primary ecological role, and many taxa serve more than one role in the ecosystem (Põlme et al. 2020). Regional differences were noticeable, as vineyard L2 had a significant increase in plant pathogen guild over time, which corresponded to a drastic increase in Botrytis in the later season, whereas this trend was not observed in L1. Among the several saprotrophs we listed, soil saprotroph, nectar/tap saprotroph, and litter saprotroph were the major functional groups, whereas wood saprotroph and pollen saprotroph were the minor groups. The presence of soil saprotrophic fungi (including Vishniacozyma) on grapevine leaves could be due to the fungi trafficking from soil to leaf. Following the guilds of the plant pathogen and soil saprotroph, foliar endophyte was the third major guild in the fungal community. This group mainly comprised two genera, Trichomerium and Septorioides, which have been reported to be dominant endophytes in Chinese chestnut (LaBonte et al. 2018) and loblolly pine tree (Oono et al. 2015), respectively. However, little information about these two genera related to grapevine is available, and their roles remain unknown. The guild of epiphyte was found to contain only a few taxa in our collection, suggesting that most fungi dwelling on the grapevine leaf function in one or more other major ecological roles.

Network analysis explores the connectedness of microbial interactions and offers insight into the dynamics and associations between taxa within an ecological community. Finding highly connected hubs and keystone taxa with network analysis can identify the driving elements of the dynamic interactions shaping their abundance (Banerjee et al. 2018; van der Heijden and Hartmann 2016). In our study, although the mycobiota co-occurrence network in the two vineyards showed equal modules, L2 had more connector hubs than L1, suggesting higher connectivity and more complex interactions in L2. Some of the connector hubs were also the dominant fungal taxa, for example Filobasidium, Vishniacozyma, and Stemphylium in L1 and Filobasidium, Vishniacozyma, Ascochyta, Botrytis, and Taphrina in L2. Noticeably, our results also showed that the abundance of individual taxa does not necessarily reflect the level of connectivity in the network. For example, although it was the most abundant taxa, Alternaria was not a hub in either location. In vineyard L1, Alternaria was only connected with Ascochyta, with a negative interaction, suggesting that these two dominant fungi might compete for the same resources. On the other hand, many taxa served as connector hubs despite their infrequent relative abundance. For example, Podosphaera, a hub genus only found in L2, is a genus that contains pathogen species that cause powdery mildew in apples and peaches. As peach is also a common crop in the Central Valley in California, the existence of Podosphaera could be expected, but how it interacts with other fungi or affects grapevine health is unknown. Ascochyta, Botrytis, and Taphrina were possibly the most important genera containing pathogens in L2, considering their abundance and network role as connector hubs. Interestingly, Botrytis had negative interactions with other yeast or yeast-like fungi, such as Vishniacozyma (top taxa), Filobasidium (top taxa), Cystofilobasidium, and Papiliotrema, suggesting possible antagonistic ability against Botrytis. Botrytis was also negatively associated with other dominant pathogenic genera such as Alternaria, Ascochyta, and Taphrina, indicating a possible competition relationship among these potential plant pathogens. In the L1 network, the powdery mildew pathogen genus Erysiphe had a negative association with the yeast-like fungi Papiliotrema. Although useful to generate new hypotheses, co-occurrence network analysis is based on OTU correlations, which can have biases and limitations arising from the compositionality of microbiome data, the presence of rare taxa, and the various environmental factors that influence the microbial interactions (Matchado et al. 2021). Thus, further experiments are needed to solidify the associations among the taxa.

Among the most abundant bacterial taxa, three classes from the phylum Firmicutes composed the largest portion of the community, including Bacilli and Erysipelotrichia, and OTUs of these classes were found to play crucial roles in the co-occurrence network. Interestingly, Bacilli was reported as an effective bio-control agent for some grape pathogens, including Botrytis cinerea (Aziz et al. 2016; Bruisson et al. 2019; Calvo-Garrido et al. 2019; Pacifico et al. 2019). Other dominant taxa belonging to Actinobacteriota, Bacteroidetes, and Proteobacteria were also typically reported in the grapevine phyllosphere microbiome (Morgan et al. 2017).

Microbiome studies have opened up new insights into agricultural management, aiming to promote plant health and resilience under biotic and abiotic stressors (Berg et al. 2017; French et al. 2021). Our study focused on changes in the grapevine leaf microbiome in California vineyards under seasonal fungicide spray programs to provide a general picture of the microbial communities and evaluate the impact of fungicide applications. Our results confirmed the major influence of spatial and temporal factors on microbial diversity and their dynamic shifts compared with the minimal effects of specific fungicide spray programs. Moreover, the core grapevine leaf microbiota was consistent with several previous studies and included potential biocontrol agents from fungi, yeast, and bacteria. With functional guilds and network analysis, we provide further insight into the ecological roles of the predominant grapevine leaf mycobiota and their potential interactions on the leaf surface. Future studies exploring the associations of these functional guilds and the role that individual fungicide sprays play in vineyard ecology are needed to further understand these complex relationships.

Acknowledgments

We thank Gabriel Torres for conducting the fungicide application experiments and Karen Vasquez for collecting leaf samples.