Is the Distribution of Two Rare Orchis Sister Species Limited by Their Main Mycobiont?

Abstract

As orchids rely on their mycorrhizal fungi for nutrient supply, their spatial range is dependent on the distribution of orchid mycorrhizal (OM) fungi. We addressed possible correlations between mycorrhizal specificity and the geographic distribution of orchids and OM fungi in three populations of the rare sister species Orchis patens and O. canariensis. Metabarcoding of the fungal ITS2 region indicated that, although adult plants of either species were colonized by several ceratobasidioid, tulasnelloid, sebacinoid and serendipitoid fungi, the mycobiont spectra were dominated by Tulasnella helicospora (which occurred in 100% of examined plants with high read numbers), which is a globally distributed fungus. In vitro assays with a T. helicospora isolate obtained from O. patens indicated the effectiveness of this OM fungus at germinating seeds of its native host. At a local scale, higher read numbers for T. helicospora were found in soil samples collected underneath O. patens roots than at locations unoccupied by the orchid. Although these findings suggest that the geographical pattern of the main fungal symbiont does not limit the distribution of O. patens and O. canariensis at this scale, the actual causal link between orchid and OM fungal occurrence/abundance still needs to be better understood.

1. Introduction

Biotic interactions are usually only measured on local (population) scales, so they are often seen as not playing an important role in determining the larger, broad-scale geographical distributions of taxa, which have traditionally been assumed as being determined mainly by abiotic variables, such as climatic and edaphic factors [1]. However, it is increasingly recognized that the current observed range limits of a plant species can also be affected by its interactions with other organisms, including predators, competitors, and mutualists such as pollinators and mycorrhizal fungi [2,3,4,5]. Orchids have associations with root-associating (orchid mycorrhizal, OM) fungi which they depend on to germinate from seed to a degree unparalleled in the plant kingdom [6]. Dependence on OM fungi is indeed extreme in orchids, which rely on their mycobionts for seed germination and seedling establishment, and in most cases, they maintain associations with these fungi into adulthood as well [7]. OM associations exhibit a broad range of specificities and levels of dependency on fungi [8,9,10]. Some orchid species maintain an association with a single mycobiont, whereas others go through a facultative or obligatory mycobiont switch at the transition to the adult stage or during environmental perturbations [11], or expand their association in adulthood. Specificity appears to be generally narrower in seedlings than in adult plants [12,13], and seedling mycobionts are usually a subset of the fungi colonizing roots of the adults [14]. Mycorrhizal fungi may therefore be key drivers of the size and distribution of orchid populations [15,16]. Since many orchids are currently threatened or endangered [6,17], a deeper appreciation of the factors influencing their species range might be critical to their long-term protection, by allowing us to predict whether orchids can colonize new habitats or become threatened in others under conditions of rapid environmental change [1].

Studies in both Orchidoideae and Cypripedioideae have indicated that the evolutionary history and therefore phylogenetic relatedness have a strong imprint on OM fungal spectra [13,18]. Broad interactions, in which species interact with multiple partners, can feature varying evolutionary specialization, and for a given focal species, the extent of evolutionary specialization on a partner species likely depends on the degree to which interacting with that partner fulfills a requirement of the focal species’ niche, as well as on whether other potential partners overlap in their ability to fill that niche [18]. Furthermore, the actual host range of a species depends not only on the specialization degree between the partners, but also on the chance they have to interact, i.e., in their biogeographic range limits [18,19]. However, currently little is known about the geographical distribution of OM fungi [1,20].

Various investigations have addressed potential connections between orchid specificity for symbiotic fungi and rarity [10,15,16]. It has been suggested that low specificity evolved to overcome limitations posed by a rare fungal associate. Orchid species that are generalist in their OM fungal preferences, or selective with one fungus that is widespread, can have a much larger distribution than species associating with few and/or rare fungi [21]. A comparison of very rare, uncommon, and common orchids (rarity being also positively related to the degree of habitat specificity) revealed no significant difference in the number of fungi they associated with. However, in the case of multiple fungal symbionts, the degree of (co)dominance of the different mycobionts should also be taken into account. Furthermore, an orchid requiring a specific but common fungus could even have greater opportunities for germination than a species associating with many fungi, each of which are very uncommon [15]. Indeed, the surveys that have attempted to sample the large-scale distribution of mycorrhizal fungi associating with a particular orchid species have shown that the wide distribution of some orchid species may, to some extent, be explained by the widespread occurrence of their mycorrhizal associates [20,21,22,23,24,25,26,27,28], or that, even if OM fungal communities differ largely across wide geographical areas, the distribution of orchids is not necessarily limited by the distribution of specific OM fungi [9,29]. Indeed, widespread orchids are often mycorrhizal generalists featuring flexibility in the OM fungi with which they associate at a geographical scale [28]. Highly restricted orchid species, such as Platanthera (Piperia) yadonii and Liparis loeselli, may also adopt a generalist or opportunistic approach in selecting a wide array of OM fungi [9,30]. By contrast, the rarity of the fungal symbiont(s) may constrain the distribution of the orchids featuring extreme specificity in their mycorrhizal partnerships, as recently shown for Platanthera praeclara [31].

Orchis patens sensu lato is a rare orchid, closely related to the Orchis spitzelii group [32], which exhibits a highly disjoint distribution (occurring in northern Africa—Algeria and Tunisia—northern Italy—a few spots around Genoa in Liguria—and the Canary Islands—all islands except Lanzarote and Fuertevenura) and is included in the European IUCN Red List as Endangered (EN) [33]. Orchis patens s.l. was in fact recently split into O. patens sensu stricto (occurring in northern Africa and northern Italy) and O. canariensis [34]. These two orchid sister species, therefore, represent an opportunity to study possible correlations between phylogenetic relatedness, mycorrhizal specificity and the geographic distribution of both orchids and OM fungi. We hypothesized that (i) due to the different geographical range of the two species, the identity of the mycorrhizal associate(s) varies between O. patens and O. canariensis and (ii) the narrow orchid distribution corresponds to the rarity of the OM symbiont(s).

2. Materials and Methods

2.1. Plant Species and Material

Sampling was carried out during spring (March–May) 2016, for a preliminary analysis by using the cloning technique from Orchis patens (from Italy) and O. canariensis (from Tenerife). For metabarcoding use, during the spring of 2017 and 2018, O. canariensis roots were collected from eight and ten plants at Gran Canaria (Figure 1a) and La Gomera from the National Park of Garajonay (Figure 1b), respectively, while for O. patens (Figure 1c), eight root samples were collected, five from Santa Margherita Ligure (Park of Portofino) and three from Breccanecca, two localities (at <20 km distance) in Liguria (Northwestern Italy), the only Italian region hosting this species. Roots, after collection, were surface cleaned and dried with silica gel until DNA extraction. In addition, three roots sampled from as many Orchis spitzelii plants collected in the Park of Majella (Abruzzo, Central Italy) were analyzed to identify their main OM fungal symbiont(s). Three fresh roots of Italian O. patens of about 10 cm were also collected from the Park of Portofino (locality Madonna della Neve) for fungal isolation.

2.2. Molecular Analysis of Mycorrhizal Fungi

2.2.1. DNA Extraction from Roots

DNA was extracted from about 20 mg of dried roots with the DNeasy Plant Mini Kit (Qiagen, Hilden, Germany) according to the manufacturer’s instructions after disruption with TissueLyser. The quality and quantity of DNA samples were assessed by spectrophotometry (ND-1000 Spectrophotometer NanoDropH; Thermo Scientific, Wilmington, Germany) and on a 1% agarose gel.

2.2.2. Cloning and Fungal Amplification from Roots

The entire ITS (ITS1-5.8S-ITS2) region was PCR-amplified with the ITS1-OFa, ITS1-OFb and ITS4-OF primers (OF), and primers ITS1 and ITS4tul (Tul) specifically designed for OM fungi [35]. For OF primers, the amplifications were performed with the following PCR conditions: 2 min at 96 °C, 35 cycles of 30 s at 94 °C, 40 s at 58 °C and 45 s at 72 °C, with an extension of 10 min at 72 °C; while for Tul primers, conditions were the following: 5 min at 94°C, 35 cycles of 30 s at 94 °C, 45 s at 54 °C and 1 min at 72 °C, with an extension of 10 min at 72 °C. Five samples (two belonging to the Italian O. patens and three to O. canariensis from Tenerife) were cloned for a preliminary analysis by using pGem-T Easy Vector kit (Promega, Madison, WI, USA) according to the manufacturer’s instructions. Sequencing of the purified 96 ITS clones was performed by Eurofins Genomics by Sanger sequencing. The nuclear ribosomal ITS2 region was amplified, for fungal metabarcoding, from the previously amplified ITS products by means of semi-nested PCR. For the nested PCR, ITS3mod and ITS4 [36,37] tagged primers were used with the following PCR conditions: 30 s at 98 °C, 30 cycles of 10 s at 98 °C, 30 s at 64 °C and 20 s at 72 °C, with an extension of 10 min at 72 °C. All amplifications were carried out by mixing 0.4 U of Phusion High Fidelity DNA polymerase (Thermo Fisher Scientific, Waltham, MA, USA), 1× Phusion High Fidelity buffer, 0.5 µM of primers, 0.2 mM of each dNTPs and 1 µL of genomic DNA (20 ng), in a final volume of 25 µL. Each DNA extract was amplified in triplicate.

2.2.3. Sampling, DNA Extraction and Amplification from Soil

DNA was also extracted from 31 soil samples (approximately 500 mg each) collected in seven areas (Madonna della Neve, Gave, Frate, Mulino, Kulm, Felciara, Sella Porcile) within the Natural Park of Portofino (Liguria, Northern Italy), at sites where O. patens was either present (Madonna della Neve and Kulm, 8 samples) or not (Gave, Frate, Mulino, Felciara and Sella Porcile, 23 samples) at sampling time (orchid presence data at these sites are available from 2014 until present). Two soil samples, in particular (sample MN3 and sample K4 from Madonna della Neve and Kulm, respectively), were collected beneath O. patens roots. After DNA extraction with the Dneasy PowerSoil Kit (Qiagen, Hilden, Germany), following the manufacturer’s instruction, samples were amplified by using the above-mentioned primer pairs, with the same PCR conditions. Amplifications were performed in three replicates.

2.2.4. Amplicon Purification and Illumina Sequencing

All PCR products were checked on 1% agarose gel, replicates pooled together and then purified with the Wizard SV Gel and PCR Clean-Up System (Promega, Madison, WI, USA). Quantification was carried out with Qubit 2.0 (Thermo Fisher Scientific, Waltham, MA, USA), by means of the Qubit™ dsDNA HS Assay Kit. Sequencing libraries were then prepared by mixing equimolar amounts of purified PCR products. The libraries were paired-end sequenced using the Illumina MiSeq (2 × 250 bp) by IGA Technology Services Srl (Udine, Italy).

2.3. Fungal Isolation and Germination Trials

Cultures of potential OM fungi were obtained from intracellular hyphal pelotons isolated from surface sterilized pieces of O. patens roots collected from Madonna della Neve (Natural Park of Portofino, Italy). Roots were firstly washed under running tap water for 30 min to remove soil and adhering contaminants and then surface sterilized in a 1% sodium hypochlorite solution for five minutes. Roots were then rinsed three times with sterile deionized water, as described by Swarts and Dixon [6], and cut into tiny transversal sections, each one divided at least in two. Tissue blocks were then put onto Malt Extract Agar (MEA) and Potato Dextrose Agar (PDA) media enriched with 5 mL of streptomycin sulphate (1 g/70 mL). Hyphae growing from these pelotons were sub-cultured onto MEA medium to obtain pure cultures for long-term storage at room temperature (20 °C) in the dark. After a first screening based on morphological observations, Rhizoctonia-like fungi were used for Sanger sequencing following amplification of their ITS region by means of the primer pairs ITS1F-ITS4 [36,38] with the following PCR conditions: 5 min at 95 °C, 35 cycles of 30 s at 94 °C, 45 s at 54 °C and 1 min at 72 °C, followed by an extension of 10 min at 72 °C.

After molecular identification, fungal strains were used for germination trials and cryopreserved at the University of Turin. Seeds of the sympatric taxa O. provincialis and O. patens, as well as their hybrid Orchis X fallax (which also occurs at sites where O. patens is no longer present), were used for the germination tests. Four mature capsules (45 days after pollination) were collected from four plants at Breccanecca; dried seeds were pooled together (for each taxon) and stored at 4 °C until use. Seeds were sterilized in a 1% sodium hypochlorite solution for 20 min, and then washed three times in sterile water. Seeds were then sowed in four Petri dishes (for each taxon) containing sterile oat agar medium (OA) with the fungal inoculum as described by Ercole et al. [39]. Seeds were also sowed in four Petri dishes (for each taxon) with Malmgren Modified medium (M551, Phytotechnology) as the asymbiotic control. Petri dishes were kept in complete darkness at room temperature for two months and checked every two weeks under a stereomicroscope for germination (i.e., the development of protocorms bearing a shoot protomeristem). The germination percentage of each Petri dish was calculated as the average of the four quadrants (considering 100 embryonized seeds per quadrant), while the total germination of each species was calculated as the average of germination percentages of individual Petri dishes ± SE (standard error).

2.4. Data analysis and Bioinformatics

2.4.1. Bioinformatic Analysis

Paired-end reads from each library were initially merged using PEAR v.0.9.2 [40], setting the quality score threshold at 28 for trimming the low-quality part of reads and the minimum length of reads after trimming at 200 bp. The assembled reads were then processed by means of the Quantitative Insights into Microbial Ecology (QIIME) v.1.8 software package [41]. Following Voyron et al. [37], sequence processing and sample assignment were performed with a minimum sequence length cutoff of 200 bp, minimum Phred quality score of 28, calculated over a sliding window of 50 bp and allowing a maximum mismatch of 3 bp over the forward and reverse primers. When necessary, sequences were re-orientated to 5’ to 3’. De-multiplexing was carried out based on primers and tags. Chimeric sequences were discarded following a de novo (abundance-based) detection with USEARCH61 [42], as implemented in the QIIME pipeline. Operational taxonomic units (OTUs) were defined with an open reference-based clustering strategy, by means of the USEARCH61 method, at 98% similarity; OTUs encompassing less than 10 sequences were removed [37]. Taxonomy assignment and OTU picking were performed by using the full “UNITE+INSD”; dataset database version 7.2 for QIIME as a reference [43,44]; the BLAST algorithm [45] was used for taxonomy assignment, with an e-value of 1e^-5 as the threshold. The OTU representative sequences (i.e., the most abundant sequence within each OTU) were submitted to GenBank and recorded under the following string of accession numbers: MT487864-MT488022 and MT489316 (Tulasnella helicospora isolated strain).

2.4.2. Phylogenetic Analysis

The representative sequences of the tulasnelloid, ceratobasidioid, sebacinoid and serendipitoid OTUs were used to perform maximum likelihood (ML) analyses. Sequences used in the ML analyses comprised sequences from cloning, best BLAST hits and fungal sequences from terrestrial orchids, including the focal species, from different continents and environments, as well as from non-orchid plants, fungal strains and fruiting bodies. Sequences were aligned using the program MUSCLE [46] with default conditions for gap opening and gap extension penalty in MEGA v.7.0.26 [47] and edited for manual adjustment. ML estimation was performed with RAxML v.8 [48] through 1000 bootstrap replicates [49] using the GTR + GAMMA algorithm to perform a tree inference and search for a good topology on CIPRES Science Gataway [50]. Bootstrap values were mapped on the best tree using the –f option of RAxML and -x 12345 as a random seed. Nodes receiving a bootstrap support ≥70% were considered as well supported.

2.4.3. Statistical Analyses

To allow for comparisons among the different root samples, a rarefaction was performed (5178 sequencing per samples) on the total OTU table by means of the rarefy_even_depth function in the R (v 3.6.0) package phyloseq (v 1.22.3) [51]. Then, a subset of 37 OMF OTUs was selected as the final OTU table of orchid roots. To compare the different soil samples, a subsampling, with 19772 sequences per soil sample was performed. After rarefaction, a subset of only OM fungi was selected for the final OTU table of soils. Indicator species analysis, a classification-based method to measure associations between species and groups of sites [52], was carried out using the multipatt function in the indicspecies (v. 1.7.6) R package, with 999 permutations [53], in order to assess whether fungi (and, if so, which ones) were significantly associated with a particular soil group. The betadisper and permutest (with 999 permutations) functions in the R package vegan were first used to assess the multivariate homogeneity of group dispersions [54] and then the effect of genotype and/or site was evaluated using a two-way permutational multivariate analysis of variance (two-way PERMANOVA with contrasts, 999 permutations), as implemented in the adonis function of the vegan package (v 2.5-4) of R [54,55].

3. Results

3.1. Identification of OMF in O. patens and O. canariensis

Seventy-nine transformed clones out of 96 were correctly sequenced (82.3% success rate), 66 of which were from O. patens and 13 from O. canariensis. These sequences, which matched GenBank sequences of tulasnelloid fungi (a group of fungi known as OM symbionts, see below), exhibited a pairwise percentage identity of 97.6% and clustered together (RAxML analysis, data not shown), indicating that both the Orchis patens (from Italy) and Orchis canariensis (from Tenerife) analyzed plants were associated to the same tulasnelloid mycobiont.

An average of 25,082 reads per root sample were obtained from Illumina sequencing after the amplification of the ITS2 region from O. patens and O. canariensis roots. The bioinformatic analysis of the Illumina MiSeq reads yielded 37 OTUs assigned to ceratobasidioid, tulasnelloid, sebacinoid and serendipitoid fungi (15, 3, 8 and 11 OTUs, respectively).

The ceratobasidioid, sebacinoid and serendipitoid fungi found in O. patens and O. canariensis were closely related either to fungi amplified from the roots of other orchid species, the roots of non-orchid plants or soil (Supplementary Information Figures S1 and S2).

The three tulasnelloid OTUs also clustered with sequences from orchid roots (Figure 2). OTU6, which could be assigned to Tulasnella helicospora, exhibited 96.3–99.7% pairwise identity with the 79 preliminary ITS clones. Only a single fungal isolate from O. patens roots was confirmed by Sanger sequencing to be an OM fungus, belonging to Tulasnella. In the phylogenetic tree, the corresponding sequence clustered with OTU6 (and the cloned ITS from O. patens and O. canariensis), with which it had 100% pairwise identity.

3.2. OM Fungal Diversity in the Different Orchid Populations

OM fungal assemblages were compared in the three orchid populations sampled (O. patens from Liguria, O. canariensis from Gran Canaria and La Gomera). OM fungal communities differed significantly both between the two orchid species and between the two O. canariensis populations (PERMANOVA based on Bray–Curtis index, p < 0.001 and p = 0.0045, respectively).

The three populations shared six OTUs (the tulasnelloid OTU6, the ceratobasidioid OTUs 83, 160, 280 and 336, and the serendipitoid OTU436), corresponding to 16% of total OM OTUs. Another set of eight OTUs (the ceratobasidioid OTUs 230, 273, 274, 302 and 303, and the serendipitoid OTUs 460, 471 and 473; 22% tot. OMF OTUs) was exclusively shared by the two Spanish populations. A single OTU (the sebacinoid OTU549) was exclusively shared by the Italian population (O. patens) and the population of Gran Canaria (O. canariensis). Only the Italian population (O. patens) exhibited unique OTUs (15 OTUs: the tulasnelloid OTUs 60 and 77, the ceratobasidioid OTUs 84, 93, 101, 261, 271 and 317, the sebacinoid OTUs 448, 449, 468, 501, 508 and 510, and the serendipitoid OTU 358).

Among the OTUs shared by all populations, OTU6 (Tulasnella helicospora) was retrieved from all root samples (Figure 3), while the ceratobasidioid OTU83 was retrieved from 100% of O. canariensis roots samples from Gran Canaria and La Gomera and from 75% of O. patens samples.

The ceratobasidioid OTUs 230 and 273 were retrieved from 82.7% and 100% of root samples from Gran Canaria and La Gomera, respectively, while being absent in O. patens roots. On the contrary, the ceratobasidioid OTUs 84, 93, 101 and 261, the serendipitoid OTU448, the sebacinoid OTU501 and the tulasnelloid OTU60 were retrieved from 100% O. patens root samples while being absent in the other species.

OTU6 (T. helicospora) encompassed 89%, 69% and 81% of total OM fungal sequences from the Gran Canaria, La Gomera and Italian samples, respectively, while each of the other OTUs encompassed <7% of total OM fungal sequences, with the only exception of the ceratobasidioid OTU 83, which encompassed 18% of total OM sequences in La Gomera plants (Figure 4).

3.3. Amplification of Fungi From Soil

An average of 50,683 reads per soil sample were obtained from the Illumina sequencing. After the filtering and cleaning of the Illumina MiSeq reads obtained from 31 soil samples collected in seven areas within the Natural Park of Portofino (Liguria, Italy) at sites where O. patens was either present or not, high-quality sequences were clustered in 2664 non-singleton fungal OTUs (98% sequence identity), 134 of which were assigned to the “rhizoctonia” complex. In particular, 13, 35, 31 and 48 OTUs belonged to tulasnelloid, ceratobasidioid, sebacinoid and serendipitoid fungi, respectively.

Indicator species analysis pointed out that some OTUs were specifically associated to a single area, such as the serendipitoid OTU8207 (p = 0.015) for area 5 (Kulm), and the tulasnelloid OTU2063 and the ceratobasidioid OTU3669 (p = 0.003, p = 0.028, respectively) for area 6 (Madonna della Neve), the two areas where O. patens had been recently observed. While the tulasnelloid OTU2063 was closely related to OTU77 retrieved from O. patens roots, the other two OTUs from soil have not been found in any orchid species studied so far (Figure 1). Other OTUs were instead associated to areas currently lacking O. patens plants, namely the ceratobasidioid OTU 3223, and the serendipitoid OTUs 7687 and 8129 (p = 0.028, p = 0.004 and p = 0.015, respectively) for area 3 (Frate), and the serendipitoid OTUs 5522 (p = 0.011), 6707 (p = 0.010) 6838 (p = 0.001), 6897 (p = 0.010), 7079 (p = 0.001), 8205 (p = 0.005) and the sebacinoid OTUs 7899 (p = 0.008) and 8538 (p = 0.039) for area 7 (Sella Porcile). None of these OTUs clustered (with >80% identity) with OTUs retrieved from orchid roots (Supplementary Information Figures S1 and S2). Fifteen OTUs were specifically associated to two areas, only one to three, two to four and four to five areas.

Sixty-six percent, 40%, 50% and 45.4% of the tulasnelloid, ceratobasidioid, sebacinoid and serendipitoid OTUs retrieved in O. patens roots, respectively, were also found in soil. The latter tulasnelloid OTUs included OTU6 (T. helicospora). By contrast, several OTUs retrieved from soil were not closely related to any fungus associating with orchid roots (Figure 1; Supplementary Information Figures S1 and S2). Interestingly, Tulasnella helicospora (OTU1241, corresponding to plant OTU6) was not specifically associated to any area, being retrieved (although as a few reads) at several sampling sites. Its read numbers, however, varied across samples: indeed, at the two sites with O. patens plants (K4 at Kulm and MN3 at Madonna della Neve), this OTU accounted for 83% and 89%, respectively, of total OM fungal sequences, while in the remaining 29 sites, it accounted for <1% of total OM fungal sequences.

3.4. In Vitro Trials

In two months, germination was not observed on the asymbiotic medium M551 for any of the three taxa examined. Conversely, germination was observed for O. X fallax but also O. provincialis and O. patens, in the symbiotic system with the T. helicospora isolate, with a germination percentage of 56 ± 4.7% for the hybrid, and 47 ± 7.2% and 67.6 ± 3.9% for O. provincialis and O. patens, respectively.

4. Discussion

4.1. Specificity of Mycorrhizal Associations of the Two Orchis Sister Species

In this study, O. canariensis and O. patens were colonized by several ceratobasidioid, tulasnelloid, sebacinoid and serendipitoid fungi, some taxa being specifically associated to either species. These fungal groups, previously referred to as the “rhizoctonia” complex, are the main OM symbionts of photosynthetic grassland orchids, associating with both seedlings and adult plants [7,8,56]. The three tulasnelloid OTUs clustered with sequences from orchid roots (Figure 2). In particular, OTU77 clustered with sequences from Neottia cordata and Gymnadenia conopsea and Orchis spp. OTU60 clustered instead with sequences from Anacamptis morio, Orchis militaris and Orchis purpurea. OTU6 could instead be assigned to Tulasnella helicospora (UNITE SH1144039.08FU), and was included in a clade which also encompasses sequences from other Orchis spp., namely O. mascula, O. militaris and O. pauciflora, as well as from Cypripedium spp. and Corycium deflexum from Greece, Belgium, Sweden, China, California and South Africa, respectively. This cluster also includes a sequence from O. spitzelii from Italy (unpublished data) (Figure 2). However, the mycobiont spectra of both orchid species were dominated by Tulasnella helicospora, which occurred in 100% examined individuals and accounted for 69–89% of total OM fungal sequences, while the other fungal symbionts were less commonly observed (at least in terms of read numbers). Germination tests with the isolate obtained from an O. patens root indicated the effectiveness of this OM fungus at germinating seeds not only of its native host, but also of O. provincialis and their hybrid (O. X fallax) that occurs also in places where today O. patens is absent.

Tulasnelloid fungi are often the most frequently retrieved OM symbionts in both temperate and tropical orchid roots [1,8]. In particular, Tulasnella helicospora was found in the great majority of Orchis species studied up to date, both the generalist clade of anthropomorphic species—which includes, e.g., O. militaris, O. simia and O. anthropophora—and the clade of non-anthropomorphic species (as well as in other orchid genera, including phylogenetically distant orchids, e.g., Cypripedium spp.) [13,18]. However, in the orchid clade of non-anthropomorphic species of the Orchis genus, which comprises both narrowly distributed species (O. patens, O. canariensis, O. spitzelii) and more widespread species (such as O. provincialis and O. mascula), T. helicospora took on the role of the dominant, if not exclusive, fungal symbiont. These orchid species feature therefore apparent mycorrhizal generalism (the situation in which a plant species specializes on one or few fungal species that contribute unique resources, but also associates with other fungi that contribute functionally redundant resources), as observed in other orchid groups such as in the Cypripedium genus [18]. The fungal species replacements observed among the three examined orchid populations in our study are also consistent with such an interpretation, since apparent generalism should yield stabilizing selection for one or a few particular fungal species, with switching among other co-occurring fungi common [18]. The shared dominant symbiont (T. helicospora) in the two sister taxa O. patens and O. canariensis is therefore a synapomorphic feature, which provides further support to the idea of phylogenetic conservatism in mycorrhizal specificity in orchids [13].

4.2. Distribution of Tulasnella helicospora in Soil

Although putative OM fungi are rarely reported from soil in metabarcoding studies making use of generic fungal primers [57,58,59,60], when primers targeted to orchid-associated basidiomycetes [35] are used, these fungi are better represented [61]. OM fungi exhibit a highly patchy distribution in soil [15,16,31]. In our study, approximately 46% of fungi retrieved from Orchis patens roots (including T. helicospora) were also retrieved from soil samples collected in the Portofino Park (Liguria, Italy), including samples collected in areas where this orchid species has not been detected since 2014. T. helicospora features a wide geographical distribution, being reported from western and central Europe and South America [62]. Furthermore, the closely related O. provincialis, with which O. patens share the dominant symbiont [13], exhibits a larger distribution than the latter orchid species, also at a local scale (Portofino Park). Similar results were revealed in studies on the Australian genus Caladenia, where high specificity in mycorrhizal associations with sebacinoid fungi was found both in widespread and narrowly distributed orchid species, and the geographical range of the latter was far less extensive than for the OM fungi [10,63]. These observations suggest that the presence of the symbiotic fungus in the environment is not a limiting factor for the distribution of the host plant, that might be influenced by other biotic or abiotic factors. However, despite the global distribution of some OM fungal species, local scale availability of preferred strains might still constrain endemic orchid recruitment and distribution.

Within the Portofino Park, higher read numbers for the dominant symbiont (T. helicospora) were found in soil samples collected underneath O. patens roots than in samples collected at sites lacking adult plants of the same orchid. Although not statistically supported, this observation parallels the results of Kaur and co-authors [31], who found that the six Platanthera praeclara-specific OM fungal OTUs were nearly absent in bulk soil at locations unoccupied by the orchids but exhibited higher relative abundances in soil at orchid-occupied locations. These findings are also consistent with assessments of spatial variation in fungal abundance in soil using quantitative PCR (qPCR), which have indicated that OM fungal abundance declines rapidly with distance from the adult host plants [16,64]. The correspondence of orchid occurrence/abundance with that of its preferred fungal partners is generally considered as an evidence of OM fungal distribution limiting orchid distribution [15,16].

In this study, the analysis of soil samples yielded a number of OTUs assigned to the same families of OM fungi which also belong to, however are not closely related to, fungi found in orchids either at the study sites or elsewhere [see also 37]. Conversely, a number of OM fungi, including Tulasnella species, are exclusively found in orchid roots rather than in the surrounding soil or in background soil samples [31,37,61,65,66], and it has been hypothesized that these OM fungi could use their host plants as a “refuge” for survival and continued existence in the environment [67,68,69]. Additional soil samplings, as well as clarification of both species boundaries and intraspecific variation in tulasnelloid, ceratobasidioid, sebacinoid and serendipitoid fungi [70,71,72] are obviously needed to assess the possible existence, within the same family, of fungi adapted to the soil environment (free-living soil saprotrophs) as opposed to fungi associated to orchid roots/protocorms which are not very competitive in bulk soil and, for this reason, are rarely found outside the orchid rhizosphere. If the latter fungi are maintained by their orchid hosts, then the distribution of orchids would determine the distribution of their OM fungal associates, rather than the reverse. Under this scenario, the retrieval of OM fungi from soil would indicate the recent occurrence of orchid hosts (either adult plants or hypogeous, undetected protocorms), and the widespread geographic distribution of some of them (such as T. helicospora) might result from their broad host range.

5. Conclusions

In conclusion, the identification of fungal symbionts in adult plants of the two sister species O. patens and O. canariensis, as also compared with the closely related O. mascula and O. spitzelii, indicate high phylogenetic conservativism, which in this orchid clade has resulted in apparent mycorrhizal generalism. The main OM symbiont Tulasnella helicospora also induced the germination of seeds of O. patens and O. provincialis. Both the global distribution of T. helicospora and the larger distribution of O. provincialis suggest that the geographical pattern of the fungal symbiont does not limit the current spatial range of O. patens and O. canariensis at this scale. In fact, our findings apparently indicate that other factors, such as possibly soil and site properties (including pollinator availability) [73,74], actually constrain the distribution of these species. For instance, the establishment of orchid populations is affected by soil pH and moisture [74], which can be modified by habitat degradation and climate change. However, as fungal strains differ in their symbiotic compatibility by orchid life stages [75], research is needed to assess the functional importance of the various OM fungi across the orchid lifecycle. However, to fully understand the role of OM symbionts in the recruitment of endangered orchid species, the most dramatic piece of evidence missing is a clear comprehension of the ecology of the fungal species, including their environmental reservoirs and sensitivities to environmental factors [18].