http://orcid.org/0000-0001-5228-5368
Figures
Abstract
Heterocephalacria bachmannii is a lichenicolous fungus that takes as hosts numerous lichen species of the genus Cladonia. In the present study we analyze whether the geographical distance, the host species or the host secondary metabolites determine the genetic structure of this parasite. To address the question, populations mainly from the Southern Europe, Southern Finland and the Azores were sampled. The specimens were collected from 20 different host species representing ten chemotypes. Three loci, ITS rDNA, LSU rDNA and mtSSU, were sequenced. The genetic structure was assessed by AMOVA, redundance analyses and Bayesian clustering methods. The results indicated that the host species and the host secondary metabolites are the most influential factors over the genetic structure of this lichenicolous fungus. In addition, the genetic structure of H. bachmannii was compared with that of one of its hosts, Cladonia rangiformis. The population structure of parasite and host were discordant. The contents in phenolic compounds and fatty acids of C. rangiformis were quantified in order to test whether it had some influence on the genetic structure of the species. But no correlation was found with the genetic clusters of H. bachmannii.
Citation: Pino-Bodas R, Laakso I, Stenroos S (2017) Genetic variation and factors affecting the genetic structure of the lichenicolous fungus Heterocephalacria bachmannii (Filobasidiales, Basidiomycota). PLoS ONE 12(12): e0189603. https://doi.org/10.1371/journal.pone.0189603
Editor: Minou Nowrousian, Ruhr-Universitat Bochum, GERMANY
Received: October 9, 2017; Accepted: November 29, 2017; Published: December 18, 2017
Copyright: © 2017 Pino-Bodas et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: DNA sequences are available in the GenBank database, accession numbers are listed in S5 Table. All other relevant data are within the paper and its Supporting Information files.
Funding: This study has received funding from the European Union's Seventh Framework Programme for research, technological development and demonstration under grant agreement nº PIEF-GA-2013-625653, CLADOF. R. Pino-Bodas thanks the MINECO for a Juan de la Cierva-Incorporación nº 2015-23526 support.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The term “lichenicolous fungi” refers to a polyphyletic group of fungi specialized in living on lichens, whether as parasites, as commensals, or as saprophytes [1–2]. It is widely accepted that most of the species of lichenicolous fungi are highly specific, restricting themselves to a few host species that generally belong to a single genus [2–3]. Nevertheless, little is known about the factors that determine the host specificity and the genetic structure of the lichenicolous fungi. Very few researches using molecular tools have addressed these questions thus far. Molina et al. [4] studied the genetic variation of the lichenicolous basidiomyceteous fungus Marchandiomyces corallinus (Roberge) Diederich & D. Hawksw., that grows on host species of different lichen genera; they found that the genetic variation of this species was geographically structured. Werth et al. [5] studied the populations of Tremella lobariacearum Diederich & M.S. Christ. in Macaronesia on three host species of the genus Lobaria. In this case the authors found that the host species was the most important factor in the genetic structure of T. lobariacearum. Millanes et al. [6] investigated the genetic variation of Biatoropsis usnearum Räsänen, finding that it constitutes a complex of cryptic species, and that the host specificity seems to have influenced the speciation. Nadler [7] proposed that the host species can be the most determining factor for the genetic structure of the parasites. However, other possible factors, like the secondary metabolites of the host, could have a similar relevance in the lichenicolous fungi. Lichens synthesize a great variety of secondary metabolites, which are chemically variable and exclusive to them. These metabolites belong to aliphatic acids, lactones, quinones, dibenzofuranes, depsides, depsidones, terpenoids, xanthones, steroids, and carotenoids [8–9]. One of the biological functions attributed to secondary metabolites in the lichen hosts are the protection against the parasites [2, 10, 11], therefore they could play an important role in the genetic structure of lichenicolous fungi. It has been suggested that the different species of lichenicolus fungi are tolerant to only a limited number of these substances, and would select their hosts according to the substances they produce [12–13]. A different tolerance degree to several lichen compounds has been proved in fungi with a lichenicolous life style [14]; therefore, the secondary metabolites of the host could be a relevant factor of adaptive divergence for these organisms. To date, however, no research has been carried out in order to test whether any connection exists between the genetic variation of the lichenicolous fungi and the secondary metabolites of the host.
The genus Heterocephalacria (Filobasidiales, Basidiomycota) includes four species of gelatinous fungi, all of them mycoparasites, two of which are lichenicolous [15–16]. Heterocephalacria bachmannii (Diederich & M.S. Christ.) Millanes & Wedin is a widely distributed species, reported from Europe, Macaronesia, North America, and Asia [17–18]. This fungus induces the formation of gelatinous, red-brownish, more or less elongated galls that often cause deformations in the host thallus, but never kill it. Little is known of the life cycle of H. bachmannii–e.g. its mating system (homothallic or heterothallic) or dispersal mechanisms. This species parasitizes only species of the genus Cladonia [17] and it has been found so far on 39 species either on primary thalli or on podetia [18].
Most of the species of the genus Cladonia are terricolous, growing in areas of high light radiation and humidity [19]. The species that are hosts for H. bachmannii synthesize a great variety of secondary metabolites. These include aliphatic acids such as rangiformic and bourgeanic acids; depsides such as barbatic acid, depsidones such as fumarprotocetraric acid, the dibenzofuran usnic acid and terpenoid such as zeorin. In addition, a great number of species of Cladonia are chemically very variable [19–20], i.e. different specimens of the same species can synthesize different secondary metabolites. Though most of Cladonia species can reproduce sexually by means of apothecia, the majority of the species frequently lack apothecia, implying that asexual reproduction is predominant by dispersion of vegetative propagules, such as soredia, granules, squamules or thallus fragments [19]. Therefore, H. bachmannii could codisperse along with the host by these fragments carrying fungus galls; in this case we would expect a high congruence between the genetic structures of the parasite and its host. The system Cladonia–H. bachmannii is therefore suitable to study which of the agents determine the genetic structure of the lichenicolous fungi.
The aim of this study was to examine the genetic structure of H. bachmannii and the factors that determine it, assessing the importance of three potential factors: host species, host secondary metabolites, and geographical distribution.
Material and methods
Sampling and loci selection
The specimens were collected between July 2014 and July 2015 in areas belonging to three different biogeographical regions: the Azores Islands (Macaronesian region), Southern Europe (Mediterranean region) and Southern Finland (Hemiboreal region). On every locality, one to five specimens of H. bachmannii per host species were collected. The collected specimens were spaced at least 5 m from each other. To complete the sampling, specimens deposited in the herbaria H, LE, and MACB were selected, including eight specimens from America, and two from Asia. The new collections were deposited in H. In total, DNA sequences from 123 specimens of H. bachmannii were obtained (Table 1). The Cladonia host species were identified according to Ahti & Stenroos [20], by morphological and chemical study of the specimens.
The specimens from herbaria are indicated in brackets, N = number of specimens.
Three loci were selected for the population–based study of H. bachmannii: ITS rDNA, LSU rDNA, and mtSSU. Of these, ITS rDNA and LSU rDNA were selected according to the results of Millanes et al. [6], who recommend the use of these loci as barcodes in Tremellomycetes. The following loci were tested to determine which of them was the most informative at population level: l41 (with the primers L41F/L41R), IGS rDNA (with LR12R/5SRNA), atp6 (with ATP6-1/ATP6-2), rpb2 (with bRPB2-6F/bRPB2-11R), rpb1 (with RPB1-Af/RPB1-CR), ef1α (with 983F/1567). So far no sequences of H. bachmannii were obtained for any of these loci. Either the amplification failed or the obtained sequences corresponded to the lichen mycobiont. Therefore we decided to use mtSSU, the only additional locus that worked.
In addition, DNA was extracted from one of the host species, Cladonia rangiformis Hoffm., trying to amplify the same loci used for H. bachmannii, such as recomended by Vienne et al. [21]. However, mtSSU was uninformative for C. rangiformis, whereby we amplified IGS rDNA, a highly variable locus that we had previously tested in several species of the genus Cladonia [22].
DNA extraction, PCRs and sequencing
We cannot discard the coinfection of a single thallus of Cladonia by several strains of H. bachmannii. For this reason, we selected a single gall of each sample of H. bachmannii for the DNA extraction. The total DNA was extracted using E.Z.N.A. forensic DNA kit (Omega Bio-tek, Georgia, U.S.A.), according to the manufacturer's instructions. The DNA was eluted in the final step in 100 μl of elution buffer included in the kit. PCRs were carried out with Ready-to-Go-PCR Beads (GEHealthcare Life Sciences, UK). The volume of the reaction was 25 μl with 3 μl of extracted DNA. The primers used were: ITS1F/BasLSU3-3 [23–24] for ITS rDNA, BasLSU3-5/LR5 [24–25] for LSU rDNA, and MS1/MS2 [26] for mtSSU (primer sequences in S4 Table). The PCR programs for ITS rDNA and LSU rDNA were the same used by Millanes et al. [24], and for mtSSU: 95°C for 5 min; 35 cycles of 95°C for 30 s, 50°C for 30 s and 72°C for 1 min; with a final extension at 72°C for 10 min. The primers used to amplify C. rangiformis were: ITS1F/ITS4 [23, 26] for ITS rDNA, LROR/LR5 [25] for LSU rDNA, and IGSf/IGSr [27] for IGS rDNA. The PCR programs for ITS rDNA and IGS rDNA were described in Pino-Bodas et al. [28], and for LSU rDNA 95°C for 5 min; 30 cycles of 95°C for 30 s, 55°C for 30 s and 72°C for 1 min; with a final extension at 72°C for 10 min. PCR products were purified with E.Z.N.A. Ultra-Sep Gel Extraction Kit or ExoSAP-IT (USB Corporation, OH, USA). The sequencing was performed at Macrogen Europe service (www.macrogen.com), with the same primers as used for the PCR.
Secondary metabolites of the hosts
All the host specimens were studied by thin layer chromatography (TLC) according to the standardized procedures [29–30], with the solvent systems A and B. The secondary metabolites of C. rangiformis were additionally studied by high–performance liquid chromatography (HPLC) and ultra–performace liquid chromatography–mass spectrometry (UPLC–MS) according to [9, 31–32]. The fatty acids of C. rangiformis were studied by gas chromatography–mass spectrometry (GC–MS). Secondary metabolites extractions, HPLC, GC–MS and UPLC–MS protocols are included as Supporting Information (S1 File).
Phylogenetic analyses and haplotypes networks
The sequences were assembled in Sequencher 4.1.4 (GeneCodes, Ann Arbor, MI). The alignments were implemented in MAFFT [33] and BIOEDIT [34]. BLAST searches were done to verify the identity of the sequences. A few sequences corresponded to other Tremellomycetes genera and were deleted.
A phylogenetic analysis based on ITS rDNA and LSU rDNA was carried out to test whether H. bachmannii is monophyletic. To this end the following species were selected: two species of Hetereocephalacria, two species of Syzygospora, one species of Piskurozyma, three species of Filobasidium, three of Goffeauzyma (voucher specimens in Supporting Information, S1 Table). As outgroup Cystofilobasidium bisporidii, C. capitatum and C. ferigula were selected according to the results of Millanes et al. [26] and Weiss et al. [35]. The ambiguous regions were removed using Gblock [36] with the less stringent options. Each region and the combined dataset were analyzed by maximum likelihood (ML). The ML analyses were implemented using RAxML 7.0.3 [37] assuming the GTRGAMMA model. The node support was estimated with rapid bootstrap algorithm, using 1000 pseudoreplicates. Congruence between the loci was tested following Lutzoni et al. [38], manually checking the clades with at least 70% bootstrap support. No incongruity was detected.
Haplotype networks for each locus under statistical parsimony were constructed in TCS 1.21 [39], considering the gaps as a 5th character.
Genetic structure analyses
Summary statisticals including haplotype diversity and nucleotide diversity were calculated in DnaSP v.5 [40]. In order to assess the contribution of the different potential factors (host species, geographical origin and host secondary metabolites) to the overall genetic variation of H. bachmannii, analyses of molecular variance (AMOVA) were conducted in ARLEQUIN V3.5 [41]. In these analyses only the groups containing more than five specimens were considered. To examine the relative contribution of these three factors, redundancy analyses (RDA) and partial redundancy analyses (pRDA) were run in R (R Development core Team 2017), using the Vegan package [42]. For these analyses a binary matrix with the haplotypes of H. bachmannii was used as dependent matrix. Other three binary matrices containing data of host species, host secondary metabolites (classified in chemotypes), and the geographical origin (The Azores, Southern Europe, Southern Finland, America, Asia) were used as explanatory matrices. The variation explained by each variable group was estimated using adjusted R2 because the number of variables in each matrix was not the same. The statistical significance was assessed using a permutation–based ANOVA test with 2000 permutations.
A bayesian clustering algorithm to estimate genetically homogeneous groups was implemented in STRUCTURE 2.3.4 [43]. Although this method was implemented to infer the population structure using unlinked markers, several studies have proved that the SNPS from sequences data as independent loci are also suited [44–45]. The analyses were run assuming an admixture model, without consideration of locality or host species of the specimens and allele frequencies independently modelled. The analysis was performed with five runs per K value (number of clusters), the range of K was from 1 to 20 with 100.000 iterations discarded as burnin and 2000.000 iterations kept for each replicate. The optimum value of K was calculated with the Evanno method in STRUCTURE HARVESTER [46–47]. The results of different runs were combined using CLUMPP 1.1.2 [48] and the barplots were generated with DISTRUCT 1.1 [49]. To assess whether the specimens of H. bachmannii on different species hosts, on different chemotypes or from different geographical regions were randomly distributed across the clusters, Chi-square tests were performed (http://www.physics.csbsju.edu/stats/). Only the groups (host species, secondary metabolites and geographical regions) with a frequency ≥ 5 were included.
Host-parasite comparisons
To compare the genetic structure of the parasite and the host, we selected Cladonia rangiformis as host species for the following reasons: 1) the largest number of specimens of H. bachmannii was found on this host; 2) it has two different chemotypes, one with atranorin, rangiformic, and nor-rangiformic acids; another that, in addition to the above compounds, contains fumarprotocetraric and protocetraric acids, 3) the largest number of haplotypes of H. bachmannii found growing on this host.
The genetic variation for host and parasite was calculated in DnaSP v.5.
Three approaches were used in order to compare the genetic variation of H. bachmannii and C. rangiformis. Firstly, we performed partial Mantel test between pairwise [Fst/(1-Fst)] matrices of H. bachmannii and C. rangiformis correcting with geographic distances with 2000 random permutations using the VEGAN package in R. The pairwise geographical distance matrix was calculated with euclidean distances from the geographic coordinates using dist() function. Three different comparisons were tested: ITS–H. bachmannii/ITS–C. rangiformis, LSU–H. bachmannii/LSU–C. rangiformis and the combined dataset–H. bachmannii/combined dataset–C. rangiformis. Secondly, an analysis was run in STRUCTURE 2.3.4 with the same conditions specified above for the specimens of H. bachmannii and C. rangiformis. The congruence between the clusters of host and parasite defined by STRUCTURE was assessed by computing the proportion of individuals assigned together in the same cluster in both analyses. Contingency table analyses (http://www.physics.csbsju.edu/stats/) were used to test the association between the host chemotype and clusters generated by STRUCTURE. Variance analyses (ANOVA) in STATGRAPHICS 5.1 were carried out to study whether the clusters of H. bachmannii differed in the amount of phenolic compounds or fatty acids. Several ANOVAs were carried out, considering the total amount of substances (total amount of phenolic compounds, and total amount of fatty acids), the total amount of different subsets of fatty acids (saturated fatty acids, unsaturated acids and free fatty acids) and the amount of every phenolic substance separately (atranorin, fumarprotocetraric acid and protocetraric acid) and the amount of most abundant fatty acids separately (linoleic, oleic, palmitic and stearic acids). Bartlett test was used to check the variance homogeneity while Kolmogorov test checked the normality of the variables. All the variables had a homogeneous variance. Three variables were not normal (phenolic compounds total contents, fumarprotocetraric acid contents and protocetraric acid contents) and were analyzed by Kruskal-Wallis.
Finally, simple AMOVA analyses were made to study the genetic differentiation of H. bachmannii and C. rangiformis among geographic regions and host chemotypes.
Results
The 123 specimens of H. bachmannii were found on 20 host species (Table 1) with ten different chemotypes (Table 1). The chemotype containing fumarprotocetraric and protocetraric acids was the most abundant, found on 47 host specimens from ten species. HPLC and UPLC–MS analyses revealed 14 compounds for C. rangiformis (S1 File). The total contents of phenolic compounds varied from 19.1 to 61.5 mg/g of dry weight (d.w.). Atranorin was the major compound in all the populations except in one of them, where the contents of fumarprotocetraric acid exceeded that of atranorin. By means of GC–MS a total of 16 fatty acids methyl esters (FAME) derived from triglycerides and additionally two free fatty acids (FFA) were identified from C. rangiformis (S1 File).
In total 304, new sequences of H. bachmannii and 90 of C. rangiformis were generated (S5 Table). In the phylogenetic analysis, all the specimens of H. bachmanii formed one strongly supported clade (S1 Fig), into which another lichenicolous species, H. physciacearum was grouped. The haplotype and nucleotide diversity values are shown in Table 2.
For each of the three loci, the haplotypes were connected in a unique network with 95% of confidence (Fig 1). A total of 28 ITS rDNA haplotypes, 24 of LSU, and six of mtSSU were found. The distribution of haplotypes among host species revealed that six ITS rDNA haplotypes, nine LSU rDNA haplotypes, and two mtSSU haplotypes were present on multiple host species, while 22 ITS rDNA haplotypes, 19 LSU rDNA haplotypes and four mtSSU haplotypes restricted themselves to a unique host species. Five haplotypes of ITS rDNA, seven of LSU rDNA, and two of mtSSU were present in several geographical regions. Southern Europe showed the highest number of haplotypes: 19 in ITS rDNA, 14 in LSU rDNA, and three in mtSSU (Table 2). In the ITS rDNA haplotype network, the most abundant haplotype represented 32.1% of all the specimens, and it was found on several host species (Fig 1). A group of haplotypes restricted to C. rangiformis was separated by four mutational steps (Fig 1A). Another haplotype found on C. granulosa also separated from the rest by five mutational steps. In LSU rDNA network, the most abundant three haplotypes represented 44% of all the specimens, and they were found on several host species (Fig 1B). One haplotype restricted only to C. granulosa and another one only to C. macilenta were the most distant. In mtSSU haplotype network, the most abundant haplotype represented 58% of individuals, and it was present on 11 host species (Fig 1C).
Each circle represents a haplotype, the circle size is proportional to haplotype frecuency. Small circles represent haplotypes not observed in the data. The colors represent the different host species.
Genetic structure of Heterocephalacria bachmannii
The results of AMOVA analyses are shown in Table 3. Most of the varition was found among the chemotypes (33.3–73.1%) and among host species (29.0–68.5%) while the genetic differentiation among geographical regions explained less percentage of variation (7.0–40.3%).
Fig 2 shows the results of the redundancy analyses. The pRDAs show that the proportion of the genetic variation explained by the chemotypes, the host species or the geographical origin, when the other factors are kept under control, was very small. The greatest proportion of genetic variation was explained by the host species in conjunction with the chemotype (0.13 in ITS rDNA, 0.13 in LSU rDNA and 0.53 in mtSSU). The geographic origin on its own explained only a little genetic variation in ITS rDNA and LSU rDNA, while in mtSSU it explained a greater proportion than each of the other two factors. All the three variables, jointly taken, explained only a small proportion of the genetic variation in the three loci.
Fractions of genetic variance explained by host species, chemotype of the host and geographical origin.
The clustering algorithm implemented in STRUCTURE revealed that the optimal number of clusters was two (highest values of ΔK were obtained for K = 2). The 90.9% of individuals were assigned to one cluster with membership coefficients > 0.7. Ten speciemens could not be assigned to any of the clusters (membership coefficients < 0.7). The inviduals were nonrandomly distributed in the two genetic clusters (S2 Table). Different association with respect to the host species (Chi–square = 68.4, d.f. = 7, P–value ≤ 0.001), chemotypes (Chi–square = 78.5, d.f. = 6, P–value ≤ 0.001) and the geographical region (Chi–square = 14.4, d.f. = 4, P–value ≤ 0.01) was found. One cluster contained specimens on Cladonia rangiformis (n = 26), C. furcata (n = 1) and C. cervicornis (n = 1). All the specimens of this cluster were collected in Southern Europe and in the Azores. The other cluster was associated with the others host species, although five specimes were also associated with C. rangiformis as the host species (Fig 3). The specimens assembled in this cluster were collected in all geographic regions and chemotypes.
The vertical bars represent the individuals, the colors indicate the proportion of the genome assignable to each cluster. Populations are separated by vertical black lines. (A) Resulting clusters of H. bachmannii on all host species. The numbers represent the different hosts: 1. C. cervicornis, 2. C. crispata, 3. C. gracilis, 4. C. cornuta, 5. C. macroceras, 6. C. granulosa, 7. C. foliacea, 8. C. glauca, 9. C. ramulosa, 10. C. uncialis, 11. C. phyllophora, 12. C. pyxidata, 13. C. mitis, 14. C. stygia, 15. C. coniocraea, 16. C. furcata, 17. C. stereoclada, 18. C. squamosa, 19. C. macilenta and 20. C. rangiformis; (B) Resulting clusters of H. bachmannii collected on C. rangiformis; (C) Resulting clusters of C. rangiformis parasitized by H. bachmannii.
Comparison of Heterocephalacria bachmannii with Cladonia rangiformis
Table 4 shows the genetic variation of H. bachmannii (collected on C. rangiformis), as well as that of the host species itself. The results of partial Mantel test are shown in Table 5. Cladonia rangiformis and H. bachmannii pairwaise Fst values were not significantly correlated.
The STRUCTURE analysis carried out with H. bachmannii samples collected on C. rangiformis showed a ΔK peak at a K = 2, as did the STRUCTURE analysis of C. rangiformis (Fig 3). In both analyses all the specimens were assigned to one of the two clusters with values > 0.8, but the cluster composition was different. The clusters defined by STRUCTURE for H. bachmannii and C. rangiformis were less congruents, only 43.8% of specimens were assigned to the same cluster. The contingency table analyses indicated that there is not a significant association between the host chemotypes and the H. bachmannii clusters (S3 Table) generated in STRUCTURE (Chi–square = 0.406, d.f. = 1, P–value = 0.524). ANOVA and Kruskal–Wallis analyses did not find significant differences in the contents of the host secondary metabolites between the two clusters of H. bachmannii (Table 6).
The AMOVA analyses are shown in the Table 7. Significant genetic differenciation was found between H. bachmannii on different host chemotypes and geographical regions. However, the variance explained by the geographical origin (31.7% in ITS rDNA, 15.9% in LSU rDNA and 58.40% in mtSSU) was higher than the explained by the chemotypes (16.8% in ITS rDNA, 13.6 in LSU rDNA and 15.7% in mtSSU). The genetic variance of C. rangiformis from different geographical origins was only significantly different in LSU rDNA analysis, while no significant differences were found in any loci between specimens with different chemotypes.
Discussion
Owing to the high number of host specific lichenicolous fungi, it has been assumed that these fungi and their lichen hosts have coevolved [2]. Coevolution depends on the genetic variation and on the genetic structure of interacting species [50]. In a coevolution process, the parasite populations would be expected to have a structure determined by the hosts, and with a low gene flow among the host species [51]. Thus far, however, the investigations on the relations between lichenicolous fungi and lichens are still too poor. The genetic variation of H. bachmannii has been studied here for the first time in a broad geographic range. The data indicate that the genetic differentiation of H. bachmannii is related to the host secondary metabolites and the host species, while this species presents few geographic isolation. This study supports the presence of two genetic clusters of H. bachmannii specialized to different hosts.
Genetic variation and genetic structure of Heterocephalacria bachmannii
The number of haplotypes in ITS rDNA and LSU were greater in southern Europe than in the other regions, an expected result considering that the sampling and the number of host species were higher in this region than in the other ones (Table 1). Our results also show mtSSU was less variable than the nuclear loci, as to the number of haplotypes, haplotype diversity, and nucleotide diversity. This result is concordant with prior studies [52–56] conducted in other fungi, the authors considered that the mutational rate in the mitocondrial genome is lower than in the nuclear genome. Nevertheless, this result could be also attributed to the fact that the number of mtSSU sequences obtained was smaller than the number of sequences of the other two loci. To check whether this explanation was plausible, we calculated the number of haplotypes for the nuclear loci taking into account only the specimens for which a sequence of mtSSU was obtained (data not shown). However, the number of haplotypes (20 in ITS rDNA and 15 in LSU rDNA) continued being smaller in mtSSU.
According to our results, the host species and the host secondary metabolites are the most relevant factors in the genetic structure of H. bachmannii, while the populations of H. bachmannii from different geographical regions show slight genetic differentiation. It supports the presence of gene flow among the populations separate by long distances. Similar results have been found for other lichenicolous fungi of the genus Tremella [5, 57]; while the opposite occurs in Marchandiomyces corallinus [4]. According to Werth et al. [5] the generalist lichenicolous fungi might be geographically structured, while the more specialist species might be structured by the host species. Therefore, H. bachmannii is expected to be structured by its hosts, all of which belong to the genus Cladonia. In other basidiomyceteous non-lichenicolous fungi, populations not geographically structured have been frecuently found [58–62]. It is assumed that fungal spores have an ability for long distance dispersal. However in some cases, no correlation has been found between the genetic and geographic distances [59] and other hypotheses have been proposed. Kretzen et al. [59] proposed two hypotheses to explain the absence of geographic structure in ectomycorrhizal fungi of the genus Rhizopogon. The first holds that competence can prevent offspring from establishing close to parents. The second hypothesis asserts that the mating systems could strongly support the outcrossing, whereby the spores coming from other populations would have a higher probability of survival success. Owing to our limited knowledge about the reproductive biology of H. bachmanii, we can not be sure that either hypothesis explains our results. Therefore, further studies are required to determine dispersal mechanisms of this species.
The effect of the lichen secondary metabolites on the lichenicolous fungi had been previously studied. Lawrey [63] found that certain phenolic compounds inhibit the growth of these fungi. Other lichenicolous fungi can only colonize the host species if another fungus has previously degraded those compounds [13]. It is also proved that the thalli of parasitized lichens have a lower concentration of phenolic substances than the non-parasitized ones [64]. However, ours is the first study in which genetic divergence among the populations of a lichenicolous fungus has been proved to be associated with the host secondary metabolites. Werth et al. [5] demonstrated that the host species was the factor that explained most of the data variance in T. lobariacearum, suggesting that the host species could create a selective environment that only certain strains would have the ability to infect. This particular environment might be created by the secondary metabolites of the host species [5]. Therefore it was to be highly expected that the genetic structure of the parasite should be strongly influenced by the secondary metabolites of the host. Even though the AMOVA results indicated that secondary metabolites explained slightly more genetic variation of H. bachmannii than the host species, in the redundance analyses the conjuntion of both factors explained most of the genetic variation. The lack of resolution in our analyses may be due to the low intraspecific chemical variation of the hosts infected by H. bachmannii. Although several chemotypes are known for many of the host species studied [22], only one parasitized host species had several chemotypes (C. rangiformis), what makes difficult to separate both effects, “host species” and “host secondary metabolites”. We do not know if our sampling was biased towards one chemotype, or if H. bachmanii in some Cladonia species only parasitizes one chemotype. Despite this parasite is not a rare species [17–18], there are no data about its distribution on the chemotypes within the host species and the total range of secondary metabolites it tolerates.
It is also remarkable that the redundacy analyses showed high unexplained genetic variation in H. bachmannii. Therefore, other envarionmental factors could also be important in the genetic structure of this lichenicolous fungus.
Discordant genetic structure between Heterocephalacria bachmannii and Cladonia rangiformis
The results also point out that the genetic structures of the host and the parasite are different. The populations of H. bachmannii were genetically more variable in the Azores than in Southern Europe, while the populations of C. rangiformis were more variable in Southern Europe. In addition, the differentiation among populations from different geographical regions was higher for H. bachmannii than for C. rangiformis (Table 7), which reveals a lower gene flow among the populations of the parasite than among those of the host. The differences between the genetic structure of the parasite and their hosts are frequent [65–67], which could be due to differences in evolutionary rates, differences in dispersion rates or movement of H. bachmannii among the different potential hosts. Mantel test did not reveal any significant correlation between the genetic distances of the host and the parasite, which indicate the absence of codispersion. Therefore, the hypothesis that the galls of H. bahchmannii were dispersed together with thallus fragments of C. rangiformis is not supported by our data. The potential movement among different host species could explain why the geographical structure of H. bachmannii is more accentuated when only the specimens on C. rangiformis are considered. Other plausible explanation would be that the lineage of H. bachmannii (specialized on C. rangiformis) colonized once the Azores and subsequently has extended across the potential hosts in the islands.
The host chemotypes have less influence on the genetic structure of H. bachmannii on C. rangiformis (Table 7) than expected. These chemotypes differ in the presence or absence of fumarprotocetraric acid, a compound synthesized by a number of Cladonia species [19–20] that act as hosts for H. bachmannii. Most likely a large number of genotypes of H. bachmannii are tolerant to this substance. The genetic structure of the lichenicolous fungi is neither shaped by the amount of phenolic compounds nor by the concentration of fatty acids (Table 6), indicating that the two genetic clusters tolerate a wide range of these compounds.
No evidence of geographical isolation was found between the populations of the Southern Europe and those of the Azores in C. rangiformis. Cladonia species present different genetic structure patterns. While Printzen & Ekman [68] found that populations of C. subcervicornis separated by a few kilometers were genetically isolated, in other species the populations were geographically weakly structured [69–71]. The most extreme case was that of C. arbuscula and C. mitis, two species with a bipolar distribution, the specimens of which appear intermixed in the northern and the southern hemispheres [72]. It is assumed that the predominant reproduction in C. rangiformis is asexual by dispersion of thallus fragments, because the apothecia are rare. The propagule size can directly influence the population structure of the species. Since several studies show that thallus fragments are dispersed at close range by wind [73–74], one would expect that the species with this type of reproduction should have a marked geographic structure in its populations. However, this is not always the case; in species of Lobaria with a varying propagule size, no correlation has been found between the propagule size and the geographic structure [75]. In C. rangiformis, gene flow between the populations of the two geographic regions could have two explanations. First, long distance dispersion of thallus fragments by the wind or other vectors. In fact, some studies have proved that birds can act as dispersion vectors for lichen fragments [76–77]. Second, that although the sexual reproduction is less frequent, the ascospores (which easily disperse long distance by the wind) could have a higher establishment success than the thallus fragments. The success of the propagules in developing new thalli is highly depending on the substrate on which they fall, the substrate formed by bryophytes being the most suitable in C. mitis [78].
Supporting information
S1 File. Secondary metabolites of Cladonia rangiformis.
The analytical protocols, results, tables and references are decribed. Phenolic compounds were identified using HPLC and UPLC-MS, and fatty acids using GC-MS.
https://doi.org/10.1371/journal.pone.0189603.s001
(DOC)
S1 Fig. Maximum likelihood tree estimated from the concatenated dataset of ITS rDNA and LSU rDNA.
https://doi.org/10.1371/journal.pone.0189603.s002
(PDF)
S1 Table. Collection information for the specimens included in the phylogenetic analysis.
https://doi.org/10.1371/journal.pone.0189603.s003
(DOC)
S2 Table. Assignment of H. bachmannii specimens to the clusters inferred in STRUCTURE.
All the specimens had membership coefficients ≥ 0.7 and they were assigned to cluster 1 or 2 without uncertainty. Chemotypes, host species and geographical origin.
https://doi.org/10.1371/journal.pone.0189603.s004
(DOC)
S3 Table. Distribution of H. bachmannii specimens on C. rangiformis in the clusters generated by STRUCTURE: Chemotypes and geographical origin.
https://doi.org/10.1371/journal.pone.0189603.s005
(DOC)
S4 Table. List of primers used in PCR and sequencing reactions.
https://doi.org/10.1371/journal.pone.0189603.s006
(DOCX)
S5 Table. GenBank accession numbers of Heterocephalacria bachmannii and Cladonia rangiformis.
https://doi.org/10.1371/journal.pone.0189603.s007
(XLSX)
Acknowledgments
We thank Dr Tuulikki Laakso for helping with UPLC–MS, Prof Ana Rosa Burgaz, Prof Teuvo Ahti, Dr Mikhail P. Zhurbenko and Dr M. Kocakaya for providing specimens, and Dr J.V. Sandoval-Sierra for his helpful suggestions.
References
- 1. Hawksworth DL. Notes on British lichenicolous fungi. IV. Notes Roy Bot Gard Edinburgh 1982; 40: 375–397.
- 2. Lawrey JD, Diederich P. Lichenicolous fungi: interactions, evolution, and biodiversity. Bryologist 2003; 106: 80–120.
- 3.
Diederich P. Host-specificity and co-evolution in lichenicolous fungi. In The Fourth IAL Symposium, Progress and Problems in Lichenology at the Turn of the Millennium, Book of Abstracts 2000; 102: 523–527.
- 4. Molina MC, DePriest PT, Lawrey JD. Genetic variation in the widespread lichenicolous fungus Marchandiomyces corallinus. Mycologia 2005; 97: 454–463. pmid:16396353
- 5. Werth S, Millanes AM, Wedin M, Scheidegger C. Lichenicolous fungi show population subdivision by host species but do not share population history with their hosts. Fungal Biol. 2013; 117: 71–84. pmid:23332835
- 6. Millanes AM, Truong C, Westberg M, Diederich P, Wedin M. Host switching promotes diversity in host-specialized mycoparasitic fungi: uncoupled evolution in the Biatoropsis-Usnea system. Evolution 2014; 68: 1576–1593. pmid:24495034
- 7. Nadler SA. Microevolution and the genetic-structure of parasite populations. J Parasitol. 1995, 81: 395–403. pmid:7776124
- 8.
Culberson CF. Chemical and botanical guide to lichen products. Chapel Hill: University of North Carolina Press; 1969. 628 p.
- 9.
Huneck S, Yoshimura I. Identification of lichen substances. Berlin: Springer Berlin Heidelberg; 1996. 493 p.
- 10. Lawrey JD. Biological role of lichen substances. Bryologist 1986; 89: 111–122.
- 11. Stocker-Wörgötter E, Elix JA, Grube M. Secondary chemistry of lichen-forming fungi: chemosyndromic variation and DNA-analyses of cultures and chemotypes in the Ramalina farinacea complex. Bryologist 2004; 107: 152–162.
- 12. Hawksworth DL. The lichenicolous Hyphomycetes. Bull Brit Mus (Nat Hist), Bot. 1979; 6: 183–300.
- 13. Lawrey JD. Chemical interactions between two lichen-degrading fungi. J Chem Ecol. 1999; 26: 1821–1831.
- 14. Lawrey JD. Isolation, culture, and degradative behavior of the lichen parasite Hobsonia santessonii. Symbiosis 1997; 23: 107–116.
- 15. Liu XZ, Wang QM, Göker M, Groenewald M, Kachalkin AV, Lumbsch HT,et al. Towards an integrated phylogenetic classification of the Tremellomycetes. Stud Mycol. 2015; 81: 85–147. pmid:26955199
- 16. Ginns J. The genus Syzygospora (Heterobasidiomycetes: Syzygosporaceae). Mycologia 1986; 78: 619–636.
- 17. Diederich P. The lichenicolous heterobasidiomycetes. Bibl lichenol. 1996; 61: 1–198.
- 18. Zhurbenko M, Pino-Bodas R. Lichenicolous fungi growing on Cladonia, mainly from the Northern Hemisphere, with a worldwide key to the known species. Opusc Philolichenum 2017; 16: 188–266.
- 19. Ahti T. Cladoniaceae. Flora Neotropica Monograph 2000; 78: 1–363.
- 20.
Ahti T, Stenroos S. Cladoniaceae. In: Ahti T, Stenroos S, Moberg R, editors, Nordic Lichen Flora Vol. 5. Uppsala, Museum of Evolution, Uppsala University, Sweden; 2013. pp. 1–117.
- 21. Vienne DM, Refrégier G, López-Villavicencio M, Tellier A, Hood ME, Giraud T. Cospeciation vs host-shift speciation: methods for testing, evidence from natural associations and relation to coevolution. New Phytol. 2013; 198: 347–385. pmid:23437795
- 22. Pino-Bodas R, Martín MP, Burgaz AR, Lumbsch HT. Species delimitation in Cladonia (Ascomycota): a challenge to the DNA barcoding philosophy. Mol Ecol Resour. 2013; 13: 1058–1068. pmid:23437908
- 23. Gardes M, Bruns TD. ITS primers with enhanced specificity for basidiomycetes-application to the identification of mycorrhizae and rusts. Mol Ecol. 1993; 2: 113–118. pmid:8180733
- 24. Millanes AM, Diederich P, Ekman S, Wedin M. Phylogeny and character evolution in the jelly fungi (Tremellomycetes, Basidiomycota, Fungi). Mol Phylogenet Evol. 2011; 61: 12–28. pmid:21664282
- 25. Vilgalys R, Hester M. Rapid genetic identification and mapping of enzymatically amplified ribosomal DNA from several Cryptococcus species. J Bacteriol. 1990; 172: 4238–4246. pmid:2376561
- 26.
White TJ, Bruns T, Lee SB, Taylor JW. Amplification and direct sequencing of fungal ribosomal RNA genes for phylogenetics. In: Innis MA, Gelfand DH, White TJ, editors. PCR Protocols: A Guide to Methods and Applications SanDiego, Academic Press; 1990. pp. 315–322.
- 27. Wirtz N, Printzen C, Lumbsch HT. The delimitation of Antarctic and bipolar species of neuropogonoid Usnea (Ascomycota, Lecanorales): a cohesion approach of species recognition for the Usnea perpusilla complex. Mycol Res. 2008; 112: 472–484. pmid:18314319
- 28. Pino-Bodas R, Burgaz AR, Martín MP, Ahti T, Stenroos S, Wedin M, et al. The phenotypic features used for distinguishing species within the Cladonia furcata complex are highly homoplasious. Lichenologist 2015; 47: 287–303.
- 29.
White FJ, James PW. A new guide to microchemical techniques for the identification of lichen substances. British Lichen Society Bulletin 1985; 57 (supplement): 1–41.
- 30.
Orange A, James PW, White FJ. Microchemical methods for the identification of lichens. Twayne Publishers; 2001.
- 31. Huneck S, Feige GB, Schmidt J. Chemie von Cladonia furcata und Cladonia rangiformis. Herzogia 2004; 17: 51–58.
- 32. Yoshimura I, Kinoshita Y, Yamamoto Y, Huneck S, Yamada Y. Analysis of secondary metabolites from lichen by high performance liquid chromatography with a photodiode array detector. Phytochem Analysis 1994; 5: 197–205.
- 33. Katoh K, Standley DM. MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol Biol Evol. 2013; 30: 772–780. pmid:23329690
- 34.
Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. In Nucleic Acids Symposium Series 1999; 41: 95–98.
- 35.
Weiss M, Bauer R, Sampaio JP, Oberwinkler F. 12 Tremellomycetes and Related Groups. In McLaughlin DJ, Spatafora JW, editors. Systematics and Evolution. Springer Berlin Heidelberg. 2014; pp. 331–355.
- 36. Talavera G, Castresana J. Improvement of phylogenies after removing divergent and ambiguously aligned blocks from protein sequence alignments. Syst Biol. 2007; 56: 564–577. pmid:17654362
- 37. Stamatakis A. RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics 2006; 22: 2688–2690. pmid:16928733
- 38. Lutzoni F, Kauff F, Cox CJ, McLaughlin D, Celio G, Dentinger B, et al. Assembling the fungal tree of life: progress, classification, and evolution of subcellular traits. Am J Bot. 2004; 91: 1446–1480. pmid:21652303
- 39. Clement M, Posada DCKA, Crandall KA. TCS: a computer program to estimate gene genealogies. Mol Ecol. 2000; 9: 1657–1659. pmid:11050560
- 40. Librado P, Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics 2009; 25: 1451–1452. pmid:19346325
- 41. Excoffier L, Lischer HE. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour. 2010; 10: 564–567. pmid:21565059
- 42.
Oksanen J, Blanchet FG, Kindt R, Legendre P, Minchin PR, O’Hara RB, et al. vegan: community ecology package. R package version 2.0–3. Available at: http://CRAN.R-project.org/package=vegan. 2012.
- 43. Pritchard JK, Stephens M, Donnelly P. Inference of population structure using multilocus genotype data. Genetics 2000; 155: 945–959. pmid:10835412
- 44. Leavitt SD, Lumbsch HT, Stenroos S, Clair LLS. Pleistocene Speciation in North American Lichenized Fungi and the Impact of Alternative Species Circumscriptions and Rates of Molecular Evolution on Divergence Estimates. PLoS ONE 2013; 8: e85240. pmid:24386465
- 45. Falush D, Stephens M, Pritchard JK. Inference of population structure using multilocus genotype data: Linked loci and correlated allele frequencies. Genetics 2003; 164: 1567–1587. pmid:12930761
- 46. Earl DA, von Holdt BM. STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method. Conserv Genet Resour. 2012; 4: 359–361.
- 47. Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol Ecol. 2005; 14: 2611–2620. pmid:15969739
- 48. Jakobsson M. Rosenberg NA. CLUMPP: a cluster matching and permutation program for dealing with label switching and multimodality in analysis of population structure. Bioinformatics 2007; 23: 1801–1806. pmid:17485429
- 49. Rosenberg NA. distruct: a program for the graphical display of population structure: PROGRAM NOTE. Mol Ecol Notes 2003; 4: 137–138.
- 50. Dybdahl MF, Lively CM. The geography of coevolution: comparative population structures for a snail and itstrematode parasite. Evolution 1996; 50: 2264–2275. pmid:28565667
- 51.
Price PW. Evolutionary biology of parasites. Princenton: Princenton University Press, NJ; 1980. 256 p.
- 52. Fournier A, Fleer R, Yeh P, Mayaux JF. The primary structure of the 3-phosphoglycerate kinase (PGK) gene from Kluyveromyces lactis. Nucleic Acids Res. 1990; 18: 365. pmid:2326170
- 53. Saliola M, Schuster JR, Falcone C. The alcohol dehydrogenase system in the Kluyveromyces lactis. Yeast 1990; 6: 193–204. pmid:2190430
- 54. Clark-Walker GD. Contrasting mutation rates in mitochondrial and nuclear genes of yeast. Curr Genet. 1991; 20: 195–198. pmid:1657417
- 55. Cao Y, Zhang Y, Yu Z, Mi F, Liu C, Tang X, et al. Structure, gene flow, and recombination among geographic populations of a Russula virescens Ally from southwestern China. PLoS ONE 2013; 8: e73174. pmid:24069176
- 56. Zhan J, Pettway RE, McDonald BA. The global genetic structure of the wheat pathogen Mycosphaerella graminicola is characterized by high nuclear diversity, low mitochondrial diversity, regular recombination, and gene flow. Fungal Genet Biol. 2003; 38: 286–297. pmid:12684018
- 57. Lindgren H, Diederich P, Goward T, & Myllys L. The phylogenetic analysis of fungi associated with lichenized ascomycete genus Bryoria reveals new lineages in the Tremellales including a new species Tremella huuskonenii hyperparasitic on Phacopsis huuskonenii. Fungal Biol. 2015; 119: 844–856. pmid:26321732
- 58. Bergemann SE, Douhan GW, Garbelotto M, Miller SL. No evidence of population structure across three isolated subpopulations of Russula brevipes in an oak/pine woodland. New Phytol. 2006; 170: 177–184. pmid:16539614
- 59. Kretzer AM, Dunham S, Molina R, Spatafora JW. Patterns of vegetative growth and gene flow in Rhizopogon vinicolor and R. vesiculosus (Boletales, Basidiomycota). Mol Ecol. 2005; 14: 2259–2268. pmid:15969712
- 60. Kretzer AM, Molina R, Spatafora JW. Microsatellite markers for the ectomycorrhizal basidiomycete Rhizopogon vinicolor. Mol Ecol. 2000; 9: 1190–1191. pmid:10964248
- 61. Von Melzer C, Rothe GM. Genetic variability of two populations of the ectomycorrhizal fungus Xerocomus chrysenteron associated with European beech (Fagus sylvatica L.). Forest Genetics 2000; 7: 85–96.
- 62. Zhou Z, Miwa M, Hogetsu T. Polymorphism of simple sequence repeats reveals gene flow within and between ectomycorrhizal Suillus grevillei populations. New Phytol. 2001; 149: 339–348.
- 63. Lawrey JD. Chemical ecology of Hobsonia christiansenii, a lichenicolous hyphomyete. Am. J. Bot. 1993; 80:1109–1113.
- 64. Merinero S, Bidussi M, Gauslaa Y. Do lichen secondary compounds play a role in highly specific fungal parasitism?. Fungal Ecol. 2015; 14: 125–129.
- 65. Delmotte F, Bucheli E, Shykoff JA. Host and parasite population structure in a natural plant–pathogen system. Heredity 1999; 82: 300–308. pmid:10336705
- 66. Vercken E, Fontaine MC, Gladieux P, Hood ME, Jonot O, Giraud T. Glacial refugia in pathogens: European genetic structure of anther smut pathogens on Silene latifolia and Silene dioica. PLoS Pathog. 2010; 6: e1001229. pmid:21187901
- 67. James PM, Coltman DW, Murray BW, Hamelin RC, Sperling FA. Spatial genetic structure of a symbiotic beetle-fungal system: toward multi-taxa integrated landscape genetics. PLoS ONE 2011; 6: e25359. pmid:21991309
- 68. Printzen C, Ekman S. Local population subdivision in the lichen Cladonia subcervicornis as revealed by mitochondrial cytochrome oxidase subunit 1 intron sequences. Mycologia 2003; 95: 399–406. pmid:21156628
- 69. Robert S, Ravigne V, Zapater MF, Abadie C, Carlier J. Contrasting introduction scenarios among continents in the worldwide invasion of the banana fungal pathogen Mycosphaerella fijiensis. Mol Ecol. 2012; 21: 1098–1114. pmid:22256778
- 70. Yahr R, Vilgalys R, DePriest PT. Geographic variation in algal partners of Cladonia subtenuis (Cladoniaceae) highlights the dynamic nature of a lichen symbiosis. New Phytol. 2006; 171; 847–860. pmid:16918555
- 71. Park CH, Jeong G, Hong SG. Possible multiple introductions of Cladonia borealis to King George Island. Antarctic Sci. 2012; 24: 359–366.
- 72. Myllys L, Stenroos S, Thell A, Ahti T. Phylogeny of bipolar Cladonia arbuscula and Cladonia mitis (Lecanorales, Euascomycetes). Mol Phylogenet Evol. 2003; 27: 58–69. pmid:12679071
- 73.
Pyatt FB. Lichen propagules. In Ahmadjian V, Hale ME, editors. The Lichens. New York, Academic Press 1973; pp. 117–149.
- 74. Walser JC. Molecular evidence for limited dispersal of vegetative propagules in the epiphytic lichen Lobaria pulmonaria. Ame J Bot. 2004; 91: 1273–1276.
- 75. Werth S, Cheenacharoen S, Scheidegger C. Propagule size is not a good predictor for regional population subdivision or fine-scale spatial structure in lichenized fungi. Fungal Biol. 2014; 118: 126–138. pmid:24528636
- 76. Bailey RH, James PW. Birds and the dispersal of lichen propagules. Lichenologist 1979; 11: 105–106.
- 77. Galloway DJ, Aptroot A. Bipolar lichens: a review. Cryptogamic Botany 1995; 5: 184–191.
- 78. Roturier S, Bäcklund S, Sundén M, Bergsten U. Influence of ground substrate on establishment of reindeer lichen after artificial dispersal. Silva Fenn. 2007; 41: 269–280.