Global analysis of the apple fruit microbiome: are all apples the same?

Microbial diversity associated with Royal Gala apple

After removal of low‐quality sequences and plant sequences, 6.117.315 16S and 48.528.735 ITS2 reads were obtained and assigned to 20.072 bacterial and 16.241 fungal ASVs, respectively. The ASVs corresponded to 25 bacterial and six fungal phyla, which in turn were assigned to 558 bacterial and 822 fungal genera. The apple fungal community across the investigated countries was dominated by Ascomycota (79.8%) and Basidiomycota (9.3%), although, other phyla such as Chytridiomycota, Entomophthoromycota, Mortierellomycota and Mucoromycota were also detected at a lower relative abundance (data not shown). Within the Ascomycota, genera such as Aureobasidium (29.00%), Cladosporium (16.60%) and unidentified groups of Capnodiales (8.80%) and Pleosporaceae (8.50%) represented more than 60% of the total fungal community (Table S2). Sporobolomyces (5.70%), Filobasidium (4.20%) and Vishniacozyma (1.60%) were the predominant Basidomycota. Regarding bacteria, Proteobacteria (65.1%), Firmicutes (15.8%), Actinobacteria (15.1%) and Bacteroidetes (2.3%) were the most prevalent bacterial phyla, representing 98.3% of the entire bacterial community. The abundance distribution of the bacterial phyla was consistent across countries, except in Turkey where Firmicutes were more abundant than Proteobacteria compared with the other countries. Sphingomonas (12.40%), Erwinia (11.30%), Pseudomonas (9.30%), Bacillus (7.10%), unidentified Oxalobacteraceae (6.80%), Methylobacterium (6.20%) and unidentified Microbacteriaceae (5.90%) were the most abundant bacterial genera. (Table S2). Results of the linear discriminant analysis (LEfSe) revealed 90 fungal and 57 bacterial taxa characterized each of the investigated countries (Fig. 1). Samples from Turkey had the highest number of fungal genera (25), which included Penicillium, Zasmidium and Pseudomicrostroma. In contrast, samples from Spain had the lowest number of fungal genera (5), which included Monilinia, Vishniacozyma and Bensingtonia (Fig. 1A). Samples from Israel and the western USA had the highest number of unique bacterial taxa, while only one bacterial taxon, identified as within the Paenibacillaceae was observed in samples collected in Uruguay (Fig. 1B).

The effect of growing region on the microbial diversity of apple fruit

The geographical location in which apples were sampled had a significant effect on the microbial diversity associated with the fruit (Table 1). For example, country of origin (including location within a country) had a significant effect on the diversity of fungi (P < 2 × 10⁻¹⁶) and bacteria (P > 2 × 10⁻¹⁶). Similarly, although to a lesser extent, the effect of orchard on fungi (P < 2 × 10⁻¹⁶) and bacteria (P = 1.09 × 10⁻⁸) was also statistically significantly. In addition, the two‐way interactions between country and tissue type as well as orchard and tissue type were significant for both fungal and bacterial diversity. Pairwise comparison between Shannon diversity of the investigated countries indicated that both fungal and bacterial diversity differed significantly between locations and orchards within a location (Table S3). Samples from Italy had the highest fungal diversity, followed by samples from Turkey and Israel (Fig. 2A). The highest bacterial diversity was observed in apples collected from Italy, the USA and Switzerland (Fig. 2B).

Spatial variation in the apple microbiome

The effect of tissue types on fungal diversity (Shannon) was not statistically significant when tissue samples from all countries were grouped together (F = 2.45, P = 0.0875). The interaction between country and tissue type, as well as between orchard and tissue type, however, were significant (Table 1). In the majority of the orchards, calyx‐end tissue exhibited a higher fungal diversity, followed by peel and stem‐end tissues, with a few exceptions observed in samples collected from Uruguay, Turkey and Spain (Fig. 4A). In contrast, tissue type had a significant effect on bacterial diversity, regardless of the sampling location (F = 68.794, P = 2 × 10⁻¹⁶), as well as in the interaction between country and tissue, as well as orchard and tissue (Table 1). Stem‐end tissues harboured the highest bacterial diversity relative to fruit peel and calyx‐end tissues, except in the New Brunswick, Canada samples (Fig. 4B). PERMANOVA analysis indicated that tissue type, as well as the interaction between tissue type and country, and tissue type and orchard, had a significant effect on fungal community composition (Table 2). This effect was observed in all of the investigated orchards in all countries, based on the results of the PCoA analysis where samples collected from apple calyx‐end, stem‐end and peel, tissues clustered separately from each other (Fig. 5A). Similar results were also found for the bacterial community which differed significantly in all orchards, (Fig. 5B).

The core microbiome of Royal Gala apple

The global core of the apple microbiome, defined at taxa present in at least 75% of the samples, consisted of six fungal genera, namely: Aureobasidium, Cladosporium, Alternaria, Filobasidium, Vishniacozyma and Sporobolomyces and two bacterial genera namely: Sphingomonas and Methylobacterium. While none of the bacterial genera was found to be prevalent in 90% of the samples, the fungal genera Aureobasidium and Cladosporium were found in up to 96% of the samples. Interestingly, the community composition of Sphingomonas was sufficient to distinguish between most of the investigated countries and showed niche specialization within the fruit i.e. stem‐end, calyx‐end and peel tissues harboured different Sphingomonas communities (Fig. S2). Similar results were also observed for Aureobasidium, a core fungal genus, however, species variability was limited, and differences were attributed to niche specialization in the different tissue‐types (data not shown).

In order to detect potential interactions between core and non‐core groups we depicted co‐occurrences by constructing a correlation matrix based on normalized distribution patterns of bacterial and fungal genera (Fig. 6A). Clustering pattern indicates that genera can be divided into five key groups of co‐occurring species mixing bacterial and fungal genera. Core species are distributed in two clusters, each hosting one of the two most dominant Ascomycota genera ‐ Aureobasidium (green) and Cladosporium (purple). Microbiomes with a high relative abundance of Aureobasidium and a low abundance of Cladosporium were characterized in Switzerland, USA and Canada; alternatively, high numbers of Cladosporium and low numbers of Aureobasidium were described in Israel and Turkey (Fig. S1). Considering the significant negative and positive interactions between genera, core species were found to have a significantly higher number of interactions in comparison to non‐core species with an average node degree of 19.125 neighbours in comparison to 12.23 in none core species (Fig. S2). A network formed by the interactions of core genera with core and non‐core groups is consistent of 142 edges and connects 8 and 60 core and non‐core genera, respectively (Fig. 6B). The highest number of interactions −30– was recorded for one core genus – Sphingomonas. Using the network, we could identify potentially useful relationships among and between core and non‐core genera within the microbial community (Fig. 6B). For example, the core genera Methylobacterium is positively associated with Burkholderiales—a group that includes reported biocontrol agents (Angeli et al., 2019), and a negative association with a known apple pathogen Podosphaera. These co‐occurrence associations can be indicative of cooperative and competitive interactions, respectively, and can serve the design of experiments to assess these interactions in vitro and on the fruit.

Geographical location

Apples represent a major item of export and are shipped globally. Therefore, it is of importance to determine if the structure of the apple fruit microbiome is relatively uniform regardless of where the fruit is produced. Rather than the presence of a uniform microbiome, the present study revealed that geographical location is a principle factor determining the structure of the apple fruit microbiome. Fungal communities, however, were more affected by geographical location (country and site within a country) than bacterial communities. The stability of the fruit‐associated bacterial community, relative to their fungal counterparts, has been previously reported in stored apples (Wassermann et al., 2019b; Abdelfattah et al., 2020). The higher level of variation in the fungal community may be potentially attributed to the fact that fungal assemblages appear to be derived from regional fungal pools with limited dispersal capability (Lumibao et al., 2019). In addition, we observed that as the variation in the microbial communities among sites was positively correlated with the distance between those locations, especially for fungi. For example, variations in fungal and bacterial communities associated with apple fruit were lower at a local scale, i.e., among orchards within the same geographical location, sites within a country, e.g., eastern and western USA and Canada and increased at the country level. Furthermore, a continental pattern can be drawn especially for the bacterial community which in one hand indicates adaptation of the apple microbiome to local environments, and on the other hand, it may be explained by the metacommunity theory. A metacommunity is defined as a set of local communities that are linked by dispersal of multiple potentially interacting species (Leibold et al., 2004). However, the present study had an insufficient distribution of samples to evaluate this premise. Nevertheless, the geographical location has been previously reported to be one of the most important determinants of the structure of the plant microbiome (Mezzasalma et al., 2018; Lin et al., 2020). A study of the maize rhizosphere found that location had a higher impact on the plant microbiome than genotype (Peiffer et al., 2013). Similarly, a study of the global citrus rhizosphere microbiome reported large variations in community structure that were attributed to geographical location (samples collected in different countries; Xu et al., 2018). The large‐scale variations between countries, together with the similarity observed among apple microbial communities within a country or region within a country, suggests that the structure of the microbial community on apple fruit is locally adapted to local environmental conditions that influence microbial diversity and composition (Hoostal et al., 2008). In this regard, it is also commonly recognized that the humid, wet conditions present in the eastern portions of the USA and Canada, present a much greater disease and pest challenge than the dry conditions present in the western USA and Canada. This is especially supported by the differences in diversity levels between these two contrasting environments, although more evident for the fungal community (e.g., Figure 3).

Tissue type

Plants tissues provide a variety of niches that can harbour distinct microbial communities. Plant roots, leaves, flowers, fruit as well as other organs, represent different microhabitats, each with specific features that favour the growth of specific microorganisms in these organs. Different tissue types within the same organ have been previously reported to exhibit spatial variations in microbial community structure. For example, the upper and lower leaf sides, as well as the peel and pulp of various fruits, including apple, have been reported to exhibit differences in microbial community structure (Abdelfattah et al., 2016; Vionnet et al., 2018; Abdelfattah et al., 2020; Piombo et al., 2020). The experimental design used in the present study was selected to determine if spatial variations within a fruit are global, i.e., will be evident regardless of geographical location and the variety of environmental conditions present in the different sites. Results indicated that the effect of fruit tissue‐type on the composition of the microbial community was rather limited, R ² = 0.0266 for fungi and R ² = 0.05191 for bacteria, yet significant, i.e., P = 0.001. A larger effect was observed, however, when individual orchards were analysed separately (Figs 5 and 6). Spatial variations in fungal and bacterial community composition and Shannon diversity due to tissue‐type was consistently observed in all of the investigated orchards. These results, along with previous studies, confirms that spatial variation in the structure of the microbial community exist between tissue‐types (calyx‐end, stem‐end and peel) at a global level. Since geographical location, is the main factor shaping the structure of the apple microbiome, however, the effect of tissue‐type is greatly reduced when samples of tissue‐types are pooled across countries. Notably, the association of a distinct microbiome with such a small environmental niche (tissue‐type) suggests specialized adaptation and function to those microhabitats. We suggest that the conditions (morphological, nutrient and environmental) present in each of these microhabitats (tissue‐types) could play an important role in determining community structure. For instance, the calyx‐end is an open site that may create special niche for specialized fungi such as Alternaria and other fungal pathogens which can cause internal rots. Interestingly, Erwinia species were found at higher abundance in the Calyx‐end tissue compared with the other tissue types, especially in Canadian apples. This can be explained by the fact that the calyx contains floral residues which are most affected by Erwinia amylvora, the cause of fire blight disease of pome fruit.

Core microbiome

A core microbiome is a set of microbes consistently present over time on a specific host and is likely to be critical to host development, health and functioning (Berg et al., 2020). Defining the core microbiome enables researchers to filter out transient associations and focus on stable taxa with a greater likelihood of influencing host phenotype and is therefore essential in exploring the potential for pre/probiotic treatments that support host health (Berg et al., 2020). In this study, the core microbiome of apple fruit was defined as fungal and bacterial taxa present in at least 75% of all samples. The core microbiome of apple fruit accounted for 13.6% of the bacterial and 63.4% of the fungal communities across all the investigated locations. We found two bacterial genera, namely Sphingomonas and Methylobacterium and six fungal genera i.e. Aureobasidium, Cladosporium, Alternaria, Filobasidium, Vishniacozyma and Sporobolomyces. This is a considerably low number of taxa, relative to other reported core microbiomes in plants (Pfeiffer et al., 2016). However, this can be attributed to the high number of samples in the present study; which lowers the probability that same taxon will be present in all samples and the evaluation of samples from different countries and tissue‐types.

Sphingomonas, a gram‐negative, non‐motile, aerobic bacterial genus, is known for its bioremediation of heavy metals and biodegradation of polycyclic aromatic hydrocarbons, and is associated with plant growth promotion through its ability to produce gibberellins and indole acetic acid in response to different abiotic stress conditions, such as drought, salinity and heavy metal stresses (Asaf et al., 2020). Interestingly, those phytohormones are also involved in fruit maturation, development and quality. For example, fruit‐set in tomato (Solanum lycopersicum) depends on gibberellins and auxins (Serrani et al., 2010; Liu et al., 2018b). Similarly, Methylobacterium is a gram‐negative, aerobic, motile bacterial genus with plant growth‐promoting properties (Krug et al., 2020). Sphingomonas and Methylobacterium have been previously reported as a component of the apple microbiome and as two of their predominate genera (Liu et al., 2018a; Wassermann et al., 2019a; Wassermann et al., 2019b; Abdelfattah et al., 2020), as well as a component of the core microbiome in several other plant species (Pirttilä et al., 2000; Delmotte et al., 2009; Mezzasalma et al., 2018; Trivedi et al., 2020). Aureobasidium and Cladosporium have also been reported as a common member of the microbiome of apple (Abdelfattah et al., 2016; Wassermann et al., 2019b) and other plants (Abdelfattah et al., 2015; Abdelfattah et al., 2018b; Abdelfattah et al., 2019). These taxa can be found as endophytes or epiphytes in association with various plant organs, e.g., leaves, flowers, fruit, seed and so forth. Although the core microbiome is typically considered to have a high level of specificity between species, the common reporting of these taxa suggests the possibility of a core microbiome that is shared between different plant species. This commonality is expected to exist at the level of genus and that some degree of species specificity may exist. In this regard, we found that the core bacterial genera, Sphingomonas and Methylobacterium accounted for a considerable fraction of the observed variation between the investigated locations, as well as tissue types. For example, the community composition of either Sphingomonas or Methylobacterium was sufficient to distinguish between most of the investigated countries. Similar results were also observed for Aureobasidium, a core fungal genus, however, species variability was limited, and differences were attributed to niche specialization in the different tissue‐types. Notably, both bacterial genera appeared to be distinct in tissue types. The geographical location demonstrated to be an important determinants of the Methylobacterium community composition in the plant phyllosphere (Knief et al.,  2010 ). The presence of distinct Sphingomonas community in different fruit tissue‐types suggests site‐specialization to these microhabitats. Interestingly, the majority of the fungal core microbiome was represented by yeasts with known antagonistic activity against pre‐ and postharvest pathogens. Despite being one of the most common fungi associated with apples, Penicillium, the causal agent of the most important apple postharvest disease, blue mould (Vero et al., 2002; Hocking, 2014; Ballester et al., 2017), was not found to be a component of the core microbiome. Penicillium species are able to grow and proliferate at low temperatures during cold storage, giving them an advantage over other fungal species (Abdelfattah et al., 2020). In this regard and considering samples in the present study were collected immediately after harvest, it can explain the low prevalence and the absence of Penicillium from the apple core microbiome. Information about the core microbiome can be further used to develop biological control strategies against apple diseases and disorders. Though core species, by definition, are detected across all samples, their relative abundance pattern vary and, in some cases, forms characteristic groups of microorganisms. Dissecting the microbiome into co‐occurrence modules can serve the construction of synthetic communities with distinct function (Vannier et al., 2019). For example, such associations can serve the design of multiple‐species synthetic communities for achieving an efficient biocontrol activity. Alternatively, dissecting the microbiome into microbial modules can allow limiting the search for a single efficient antagonist agent. In the context of the apple fruit microbiome, co‐occurrence patterns have stratified the fruit microbiome into five key groups with core genera located in two of them: one cluster with Aureobasidium, and the second with Cladosporium, the two most abundant Ascomycota genera. Though most of the significant interactions detected in the network are positive, some negative associations allow formulating predictions for potential biocontrol agents against pathogens. Based on the network view, experimental design of potential biocontrol agent could compare the activity of a single microorganism versus consortium representing a native co‐occurring module. Potential biocontrol strategies can hence benefit from the network view of microbiome interactions allow to go beyond the single biocontrol agent to the educated design of a biocontrol consortium.

Libraries and sequencing, data processing, downstream and statistical analysis

Lyophilized samples were homogenized, and their DNA was extracted using DNeasy PowerLyzer PowerSoil Kit (Qiagen, Germantown, MD, USA). Initial tissue disruption of 250 mg was performed with a Qiagen PowerLyzer 24 Homogenizer (Qiagen, Germantown, MD, USA). DNA extractions were automated using a Qiagen QiaCube (Qiagen, Germantown, MS, USA), using the processing routine recommended by the manufacturer for the PowerSoil kit. Extracted DNA was used as the template for amplicon PCR reactions that amplified the bacterial 16S ribosomal region and the fungal internal transcribed spacer (ITS) region. The V4 region of 16S rRNA was amplified using the universal primers 515F (Parada et al., 2016) and 806R (Apprill et al., 2015) in conjunction with peptide nucleic acids (PNAs) (PNA Bio) added to inhibit amplification of ribosomal and mitochondrial sequences (Lundberg et al., 2013). ITS amplicons were amplified using ITS3/KYO2 (Toju et al., 2012) and ITS4 (White et al., 1990) primers along with a custom‐designed blocking oligo designed to inhibit amplification of the host DNA (5′ ATTGATATGCTTAAATTCAGCGGGTAACCCCGCCTGACCTGGGGTCGCGTT‐C3 spacer 3′). All primers were modified to include the necessary Illumina adapters (www.illumina.com) for subsequent PCR addition of Illumina indexes for multiplexing.

For bacteria, PCR reactions were conducted in a total volume of 25 μl containing 12.5 μl of KAPA HiFi HotStart ReadyMix (Kapa Biosystems), 1.0 μl of each primer (10 μM), 2.5 ul of mitochondrial PNA (5 uM), 2.5 μl of plastid PNA (5 uM), 2.5 μl of DNA template and 3 μl nuclease‐free water. Reactions were incubated in a T100 thermal cycler (BioRad) at 95°C for 5 min followed by 30 cycles of 95°C for 30 s, 78°C for 5 s, 55°C for 30 s, 72°C for 30 s and a final extension at 72°C for 5 min. For fungal (ITS) amplicon generation, 25 μl PCR reactions contained 12.5 μl of KAPA HiFi HotStart ReadyMix (Kapa Biosystems), 1.0 μl of each primer (10 μM), 1.0 μl of blocking oligo (10 uM), 2.5 μl of DNA template and 7 μl nuclease‐free water. Reactions were incubated in a T100 thermal cycler (BioRad) at 95°C for 5 min followed by 30 cycles of 95°C for 30 s, 55°C for 30 s, 72°C for 30 s a final extension at 72°C for 5 min. Library preparation following amplicon PCR was performed as specified in the Illumina 16S Metagenomic Sequencing Library Preparation guide precisely as outlined. Subsequent library size, quality and confirmation of the absence of adapter dimers were performed on an Agilent 2100 Bioanalyzer (Agilent). Paired‐end sequencing of amplicons was done on an Illumina MiSeq (Illumina) sequencer with a V3 600‐cycle Reagent Kit (Illumina).

Data analysis

Demultiplexing, merging forward and reverse reads, quality filtering and trimming and ASVs generation were done using DADA2 as integrated in Qiime2 ADDIN EN.CITE (Callahan et al., 2016; Bolyen et al., 2019). Taxonomic assignment of ASVs was done using a similarity threshold of 97% against the Greengenes database for 16S reads and against the UNITE database for ITS reads (Abarenkov et al., 2010; DeSantis et al., 2006). Plant related sequences e.g. chloroplast and mitochondria, were filtered out before the following analyses. Rarefaction to an even sequencing depth of 1000 reads per sample was used to normalize ITS and 300 reads for the 16S features tables which were used to calculate Shannon diversity. To evaluate the effect of the geographical location, orchard and tissue type (i.e., stem‐end, calyx‐end and peel) on the fungal and bacterial diversity, Analysis of Variance Model (a wrapper for fitting linear models) was used. Here, we modelled Shannon index as a response variable and country, orchard, and tissue as fixed effects. In addition, pairwise comparisons between Shannon index of the fungal and bacterial diversity among the sampling locations (countries) was done using Wilcox test and the p values were corrected using FDR method. For the subsequent analysis, the unrarefied ASV table was normalized using MetagenomeSeq's Cumulative Sum Scaling (CSS; Paulson et al., 2013). Differences in community composition between the investigated countries, orchards and tissue types were tested using Permutational Multivariate Analysis of Variance Using adonis (~PERMANOVA) in vegan R with 999 permutations. Here, we modelled the distance among countries in their community composition, based on Bray–Curtis dissimilarity index as a response variable and country, orchard and tissue as fixed effects. Hierarchical clustering of the community composition based on Bray–Curtis dissimilarity distances with ‘average’ as the clustering method was performed using hclust in R Package ‘stats” version 4.0.1 (Oksanen et al., 2007). These results were visualized in dendrograms using fviz_dend function in the R package factoextra version 1.0.7 (Kassambara and Mundt, 2020).

The core microbiome was calculated based on genera present in at least 75% of the investigated samples using core function in Microbiome package (Salonen et al., 2012). Interactions between core and non‐core species were limited to genera whose normalized relative abundance >0.1% (average across replicas) in at least a single sample. Co‐occurrences were described based on Spearman's rho correlation coefficients between the distribution patterns of the genera joining the normalized bacterial and fungal tables. Scores were calculated using ‘Pandas.corr’ python package v1.1.0. Correlation matrix and visualized using ‘seaborn.clustermap’ python package v0.10.1. Co‐occurrence networks were generated based on correlation scores. Network visualization and the positioning of the nodes were calculated according to the force‐directed Fruchterman & Reingold algorithm used for calculating layouts of simple undirected graphs (Buchfink et al., 2015). The algorithm was implemented using nx.draw function via the ‘pos’ parameter in the ‘NetworkX’ python package v1.11. Node degree was calculated using the nx.degree function. Visualization was generated using ‘Plotly’ python package v4.9.0. Linear discriminant analysis effect size (LEfSe; Segata et al., 2011) was used for biomarker discovery to determine a list of taxa that best characterize each geographical location (Country). Higher LEfSe score indicate higher consistency of differences in relative abundance between taxa of each country.