Introduction
The family Parmeliaceae is amongst the most diverse lineages of lichenized ascomycete fungi, with c. 2800 species known worldwide (Thell et al. Reference Thell, Crespo, Divakar, Kärnefelt, Leavitt, Lumbsch and Seaward2012; Jaklitsch et al. Reference Jaklitsch, Baral, Lücking, Lumbsch and Frey2016; Stelate et al. Reference Stelate, Del-Prado, Alors, Tahiri, Divakar and Crespo2022). These often form foliose or fruticose lichen thalli, and occur in a wide range of ecosystems including arctic-alpine (Goward & Myllys Reference Goward and Myllys2020; Xu et al. Reference Xu, De Boer, Olafsdottir, Omarsdottir and Heidmarsson2020), Antarctic (Øvstedal & Lewis-Smith Reference Øvstedal and Lewis-Smith2001), tropical and temperate forests (Molina et al. Reference Molina, Divakar, Goward, Millanes, Lumbsch and Crespo2017; Del-Prado et al. Reference Del-Prado, Buaruang, Lumbsch, Crespo and Divakar2019) or even in fog-influenced desert regions such as the Namibian lichen fields, where some species thrive on quartz pebbles (Wirth Reference Wirth2010). In European temperate and Mediterranean forest ecosystems, epiphytic communities are usually dominated by parmelioid lichens forming relatively large thalli, such as members of the genera Flavoparmelia Hale, Parmelia Ach., Parmotrema A. Massal., Punctelia Krog, Pseudevernia Zopf or Usnea Dill. ex Adans. (Hawksworth et al. Reference Hawksworth, Blanco, Divakar, Ahti and Crespo2008; Nimis Reference Nimis2016; Nimis & Martellos Reference Nimis and Martellos2022). Although the identification of the vast majority of them is relatively straightforward, given the number of diagnostic macroscopic traits and secondary metabolite profiles, the use of DNA sequences in a phylogenetic framework has allowed species boundaries and regional occurrences to be (re-)assessed either in complex lineages (Alors et al. Reference Alors, Lumbsch, Divakar, Leavitt and Crespo2016; Widhelm et al. Reference Widhelm, Egan, Bertoletti, Asztalos, Kraichak, Leavitt and Lumbsch2016; Lutsak et al. Reference Lutsak, Fernández-Mendoza, Kirika, Wondafrash and Printzen2020) or apparently well-known species (Lendemer & Hodkinson Reference Lendemer and Hodkinson2010; Castellani et al. Reference Castellani, Bianchi, Coppi, Nascimbene and Benesperi2021; Stelate et al. Reference Stelate, Del-Prado, Alors, Tahiri, Divakar and Crespo2022). Based on the results of these and many other studies, it is clear that there is still room for further work based on single molecular markers, such as the nrITS fungal barcode (Schoch et al. Reference Schoch, Seifert, Huhndorf, Robert, Spouge, Levesque and Chen2012), to improve knowledge of the diversity of this fungal family on a regional scale, or in areas of special conservation concern (e.g. Singh et al. Reference Singh, Kukwa, Dal Grande, Łubek, Otte and Schmitt2019).
Phylogenetics based on single-locus datasets has also been essential for the identification of the microalgal partner associated with parmelioid lichens. A major work by Leavitt et al. (Reference Leavitt, Kraichak, Nelsen, Altermann, Divakar, Alors, Esslinger, Crespo and Lumbsch2015) showed that microalgal species in the genus Trebouxia De Puymaly (Chlorophyta) were the main phycobionts of parmelioid lichens across lineages and geographical regions. Further evidence of that specificity has been compiled later in studies focusing on particular groups (e.g. Parmelia saxatilis (L.) Ach. and P. sulcata Taylor species complexes; Moya et al. Reference Moya, Molins, Škaloud, Divakar, Chiva, Dumitru, Molina, Crespo and Barreno2021), or in phylogeographical assessments based on large specimen sampling. For example, diverse Trebouxia species and infraspecific lineages were found in association with Cetraria aculeata (Schreb.) Fr. and Pseudephebe M. Choisy spp. across an amphitropical distribution (Fernández-Mendoza et al. Reference Fernández-Mendoza, Domaschke, García, Jordan, Martín and Printzen2011; Garrido-Benavent et al. Reference Garrido-Benavent, Pérez-Ortega, de los Ríos and Fernández-Mendoza2020). These works highlighted that switching to locally adapted phycobionts might have allowed lichenized fungi to widen their geographical distribution (Peksa & Škaloud Reference Peksa and Škaloud2011), although this does not seem to always be the case, as demonstrated in other lichenized fungal groups (Škvorová et al. Reference Škvorová, Černajová, Steinová, Peksa, Moya and Škaloud2022). Phycobiont switching has also been found to cause changes in reproduction type and phenotypic variability (Ertz et al. Reference Ertz, Guzow-Krzemińska, Thor, Łubek and Kukwa2018). Due to the large species diversity in Parmeliaceae, there are still many information gaps regarding the fungal-algal relationships in most lineages.
One of those cases corresponds to Punctelia, a genus that was segregated from Parmelia based on a combination of several characters: pseudocyphellae ontogeny, conidia morphology, chemical profiles and geographical distribution (Krog Reference Krog1982). Punctelia currently hosts c. 45 species (Thell et al. Reference Thell, Crespo, Divakar, Kärnefelt, Leavitt, Lumbsch and Seaward2012), and its distribution centre is probably subtropical South and North America and Africa (Krog Reference Krog1982). Although the boundaries of some species in this genus have already been assessed phylogenetically (Crespo et al. Reference Crespo, Divakar, Argüello, Gasca and Hawksworth2004; Alors et al. Reference Alors, Lumbsch, Divakar, Leavitt and Crespo2016), there are limited data on the identity and range of associated phycobionts. The studies by Leavitt et al. (Reference Leavitt, Kraichak, Nelsen, Altermann, Divakar, Alors, Esslinger, Crespo and Lumbsch2015), Xu et al. (Reference Xu, De Boer, Olafsdottir, Omarsdottir and Heidmarsson2020) and Kosecka et al. (Reference Kosecka, Kukwa, Jabłońska, Flakus, Rodríguez-Flakus, Ptach and Guzow-Krzemińska2022) indicated that Punctelia rudecta (Ach.) Krog and P. hypoleucites (Nyl.) Krog are associated with several Trebouxia lineages included in ‘clade I’ (‘impressa’, sensu Muggia et al. (Reference Muggia, Nelsen, Kirika, Barreno, Beck, Lindgren, Lumbsch and Leavitt2020)) in North and South America, Africa, the Canary Islands, India and Japan. However, knowledge on phycobiont identity and interactions with other Punctelia, and especially lichenized species occurring in temperate and Mediterranean Europe, is virtually absent.
The present work examines phycobiont identity in Punctelia borreri (Sm.) Krog and P. subrudecta (Nyl.) Krog, two widespread species in Europe, and compares their range of associated phycobionts. To this aim, specimens from Atlantic and Mediterranean localities in the Iberian Peninsula were collected, and the nrITS region of both partners was amplified. Our hypothesis is that these two lichens would host a trebouxioid phycobiont and that it would probably belong to Trebouxia ‘clade I’, as previous works suggested. However, lineages of associated phycobionts might differ between Atlantic and Mediterranean areas, given the substantial differences in bioclimatic conditions (Rivas-Martínez et al. Reference Rivas-Martínez, Penas, del Río, Díaz González, Rivas-Sáenz and Loidi2017). Moreover, the ultrastructure of the algal cells is studied with transmission electron microscopy (TEM), and a phylogenetic study based on the nrITS of the mycobiont is conducted to ascertain their taxonomic identity, since thalli formed by the two Punctelia species have usually been difficult to distinguish morphologically (e.g. Truong & Clerc Reference Truong and Clerc2003; but see Longán et al. Reference Longán, Barbero and Gómez-Bolea2000). Finally, evidence of visible symptoms of thallus injury in these lichens, possibly due to air pollution, is illustrated and discussed.
Material and Methods
Fieldwork and morphological studies
In total, 41 specimens of Punctelia borreri and P. subrudecta were analyzed (Fig. 1). These were collected in 25 local populations, of which six and 19 were distributed in the northern and eastern Iberian Peninsula, respectively (Supplementary Material Table S1, available online). These two areas correspond with the Eurosiberian-Atlantic and Mediterranean biogeographical regions, respectively (Rivas-Martínez et al. Reference Rivas-Martínez, Penas, del Río, Díaz González, Rivas-Sáenz and Loidi2017). Two additional specimens from Richmond (London, UK) and Puglia (Italy) were also included for comparative purposes. Since performing thorough population genetics analyses was not the aim of the present study, up to five specimens per locality were sampled. Specimens growing on tree bark and rocks were collected. A preliminary taxonomic identification of each specimen followed Longán et al. (Reference Longán, Barbero and Gómez-Bolea2000) and Smith et al. (Reference Smith, Aptroot, Coppins, Fletcher, Gilbert, James and Wolseley2009). Thallus chemistry was assessed using a sodium hypochlorite solution (reagent C) on the thallus medulla; this generates a reddish stain in P. borreri and a pinkish one in P. subrudecta. The two species were further distinguished by the colour of the underside of the thalli: darker in P. borreri than in P. subrudecta. These characters were examined under a Leica M165C stereomicroscope. In addition, to document the visible symptoms of injury on thalli of a Punctelia sp., an image taken in 2011 in the coastal mountains of Sierra Calderona (Castellón, Valencia region, Spain) is provided (Fig. 2). The studied material is deposited in the Collection of Lichens VAL_Lich of the Department of Botany and Geology (University of Valencia, Spain).
Figure 1. Punctelia borreri and P. subrudecta thallus morphology and ultrastructure of the phycobiont. A, typical habitat of P. borreri on holm-oak (Quercus ilex subsp. rotundifolia) in the eastern Iberian Peninsula. B, P. borreri thallus when dry. C & D, P. subrudecta thalli with farinose soralia towards the centre and punctiform pseudochyphellae in young lobes. C, dry thallus. D, wet thallus. E–G, ultrastructural features by TEM of the microalga Trebouxia gelatinosa associated with P. subrudecta. E, cell ultrastructure showing a chloroplast (ch), pyrenoid (py) and cell wall (cw). F, detail of a pyrenoid showing the close association of pyrenoglobules (pg) to pyrenotubules (pt). G, fungal hyphae (hy) interacting with a phycobiont cell. Scales: E–G = 1 μm. In colour online.
Figure 2. Visible symptoms of injury in a Punctelia sp. thallus from Sierra Calderona (coastal mountains of Castellón, Valencia region, Spain) impacted by air pollutants. Note that the pink colouring and necrosis approximately start in the central areas of the thallus.
Transmission electron microscopy (TEM) examinations were performed on a portion of the P. subrudecta thallus (voucher utiel_a1) that was processed following Molins et al. (Reference Molins, Moya, García-Breijo, Reig-Armiñana and Barreno2018). Cross-sections of the phycobionts were observed with a Hitachi HT7800 (120 kV) electron microscope fitted with a high-speed CMOS camera (EMSIS XAROSA) at the SCSIE microscopy service of the University of Valencia.
Laboratory work and sequence processing
A small piece c. 1 mm2 was excised from a young thallus lobule in each specimen for DNA extraction. This procedure was carried out after careful examination of thalli under a stereomicroscope to avoid selecting areas infected by other microfungi or tiny clumps of epiphytic microalgae. Subsequently, total genomic DNA was isolated using the DNeasy Plant Mini Kit (Qiagen, Hilden, Germany) following the manufacturer's instructions. The nuclear internal transcribed spacer of the ribosomal DNA (nrITS) of the fungal and algal partners was selected for PCR amplification. In the mycobiont, the primer pair used was ITS1F (Gardes & Bruns Reference Gardes and Bruns1993) and ITS4A (Larena et al. Reference Larena, Salazar, González, Julián and Rubio1999); for the phycobiont, the primer pair was nr-SSU-1780 (Piercey-Normore & DePriest Reference Piercey-Normore and DePriest2001) and ITS4T (Kroken & Taylor Reference Kroken and Taylor2001). PCR amplifications were performed in a total volume of 15 μl, containing 7.5 μl of EmeraldAmp GT PCR Master Mix (Takara), 0.3 μl of each primer (10 μM), 0.6 μl of dimethyl sulfoxide (DMSO), 2 μl of DNA template and 4.3 μl of sterile Milli-Q water. Amplification by PCR of the fungal (algal) partners began with an initial denaturation step for 5(2) min at 95(94) °C, followed by 35(30) cycles of 30 s at 95(94) °C, 1 min (45 s) at 52(56) °C, and 1 min 30 s (1 min) at 72 °C, ending with a final elongation step at 72 °C for 10 min. PCR products were visualized on 2% agarose gels and purified using the NZYGelpure kit (NZYTech). Purified PCR products were sequenced with an ABI 3100 Genetic Analyzer using the ABI BigDye Terminator Cycle Sequencing Ready Reaction Kit (Applied Biosystems, Foster City, California). Raw electropherograms were manually checked, trimmed and assembled using SeqManII v. 5.07 (DNASTAR Inc.). GenBank Accession numbers are provided in Supplementary Material Table S1 (available online).
Sequence dataset construction and alignments
The mycobiont dataset was built using the newly obtained sequences of Punctelia borreri and P. subrudecta from the Iberian Peninsula, together with sequence data deposited in the GenBank database (https://www.ncbi.nlm.nih.gov/genbank/). Random sequences for each of the two Punctelia species were submitted to the BLAST online tool (Altschul et al. Reference Altschul, Madden, Schäffer, Zhang, Zhang, Miller and Lipman1997) to identify highly similar hits (i.e. sequence similarity ranging between 97–100%) in GenBank, irrespective of their taxonomic label. A total of 39 hits were downloaded and aligned with the new sequences using MAFFT v. 7.222 (Katoh et al. Reference Katoh, Misawa, Kuma and Miyata2002; Katoh & Standley Reference Katoh and Standley2013), as implemented in Geneious® v. 9.0.2, employing the FFT-NS-i ×1000 algorithm, the 200PAM/k = 2 scoring matrix, a gap open penalty of 1.5 and an offset value of 0.123. Flavopunctelia flaventior (Stirt.) Hale and F. soredica (Nyl.) Hale were also included in the alignment for rooting purposes in phylogenetic inference following Alors et al. (Reference Alors, Lumbsch, Divakar, Leavitt and Crespo2016). The resulting alignment was manually edited to trim ends of longer sequences that included part of the 18S-28S ribosomal subunits, and to replace gaps at the ends of shorter sequences with a code representing any base (‘N’). Partitions within the nrITS (i.e. spacers 1 and 2, and the 5.8S locus) were annotated.
The nrITS sequences of the two Punctelia-associated microalgae were also submitted to the BLAST online tool to identify GenBank hits with a high sequence similarity. The identity of the retrieved hits made it possible to infer the most likely phylogenetic location of the newly produced sequences within the most recently published Trebouxia phylogeny (Muggia et al. Reference Muggia, Nelsen, Kirika, Barreno, Beck, Lindgren, Lumbsch and Leavitt2020). Based on this work and the BLAST searches, 46 sequences were downloaded and aligned with our Punctelia spp. algal sequences using the MAFFT algorithm on the GUIDANCE2 web server (http://guidance.tau.ac.il/ver2/). First GUIDANCE2 calculated an overall reliability score for the nrITS alignment (Sela et al. Reference Sela, Ashkenazy, Katoh and Pupko2015) and then a separate alignment excluding columns with GUIDANCE2 scores < 95% was generated to minimize bias in phylogenetic reconstruction arising from ambiguously aligned regions. Further editing of the alignment was carried out as for the mycobiont (see above). No outgroup sequences were included in the algal dataset.
Phylogenetic analyses
Maximum likelihood (ML) phylogenetic trees were calculated for the mycobionts and phycobionts based on the previously inferred alignments. ML tree estimates were carried out with RAxML (Stamatakis Reference Stamatakis2006; Stamatakis et al. Reference Stamatakis, Hoover and Rougemont2008), using the GTRGAMMA nucleotide substitution model for nrITS1, nrITS2 and 5.8S. One thousand rapid bootstrap pseudoreplicates were conducted to evaluate nodal support. The resulting trees were drawn with the iTOL v. 5 web tool (Letunic & Bork Reference Letunic and Bork2021) and Adobe Illustrator CS5 was used for artwork. Tree nodes with bootstrap support (BS) values ≥ 70% were regarded as significantly supported. Additionally, the mycobiont alignment was used to infer a phylogeny under a Bayesian inference criterion. For this purpose, the program PartitionFinder v. 1.1.1 (Lanfear et al. Reference Lanfear, Calcott, Ho and Guindon2012) was first employed to select the following optimal partitioned evolutionary models: SYM + Γ (nrITS1 + 2) and JC (5.8S). The Bayesian phylogenetic reconstruction was then conducted with the program MrBayes v. 3.2.7 (Ronquist et al. Reference Ronquist, Teslenko, van der Mark, Ayres, Darling, Höhna, Larget, Liu, Suchard and Huelsenbeck2012), with two parallel, simultaneous four-chain runs executed over 5 × 107 generations starting with a random tree, and sampling after every 500th step. The first 25% of data was discarded as burn-in, and the 50% majority-rule consensus tree, and corresponding posterior probabilities (PP), were calculated from the remaining trees. Chain convergence was assessed by ensuring that values of the average standard deviation of split frequencies (ASDSF) dropped below 0.005, and values of the potential scale reduction factor (PSRF) approached 1.00. Tree nodes with PP ≥ 0.95 were regarded as significantly supported.
Inference of haplotype networks
The software DnaSP v. 5.10 (Librado & Rozas Reference Librado and Rozas2009) was used to infer myco- and phycobiont nrITS haplotypes considering sites with alignment gaps. Before conducting this step, however, alignments were trimmed to avoid short sequences with many ‘N's at their ends, and also columns with any ambiguous nucleotide. For the mycobiont alignment, this editing included the removal of the two outgroup sequences, a Punctelia specimen with the GenBank Accession MZ836018, which showed many ambiguous nucleotides, as well as 34 alignment columns. In the phycobiont, the alignment used for inferring haplotypes consisted of 44 sequences, including those newly generated (except for that of I1297 which was too short) as well as others downloaded from GenBank that corresponded with Trebouxia sp. OTU I05 (‘clade I’, sensu Muggia et al. (Reference Muggia, Nelsen, Kirika, Barreno, Beck, Lindgren, Lumbsch and Leavitt2020)). The TCS method implemented in PopART v. 1.7 (Templeton et al. Reference Templeton, Crandall and Sing1992; Leigh & Bryant Reference Leigh and Bryant2015) was then employed to infer genealogical relationships independently among mycobiont and phycobiont haplotypes under a statistical parsimony framework. Networks were artistically edited in Adobe Illustrator CS5; in the mycobiont network, haplotypes were colour-labelled according to different geographical sampling regions, and in the phycobiont networks the labelling followed the group of associated mycobionts according to GenBank metadata.
Interaction networks and specialization
Interaction networks between haplotypes of the two Punctelia species and phycobiont haplotypes were constructed using the function ‘plotweb’ in R package bipartite (Dormann et al. Reference Dormann, Gruber and Fründ2008). Fungal haplotype specialization towards phycobiont nrITS haplotypes was estimated using the parameter dʹ (Blüthgen et al. Reference Blüthgen, Menzel and Blüthgen2006), which can be calculated using the function ‘dfun’. This parameter ranges from 1 (high specialization) to 0 (no specialization) and offers a measure of how much a mycobiont lineage discriminates from a random selection of algal partners.
Results
Alignment statistics and phylogenetic analyses
A total of 41 mycobiont and 38 phycobiont nrITS sequences were generated de novo from Punctelia borreri and P. subrudecta specimens collected for the present study. No evidence of algal coexistence in a single lichen thallus was found. The fungal sequence alignment that included GenBank data totalled 82 sequences and 506 columns, of which 89 and 22 corresponded with variable and singleton sites, respectively. The algal sequence alignment obtained in GUIDANCE2 that included 84 sequences had a score of 0.986004 and was 752 bp in length. The number of variable and singleton sites was 186 and 34, respectively. ML analyses in RAxML estimated topologies with final log-likelihood values of −1300.4053 (mycobiont) and −3257.6684 (phycobiont). The MrBayes analysis executed with only the mycobiont alignment reached an ASDSF value of 0.005 after 3.93 × 106 generations, and Estimated Sample Sizes (EESs) were above 1900 for all parameters. Since the fungal RAxML and MrBayes topologies showed no supported conflicts, only the topology inferred with the former approach is presented in Fig. 3. The two Punctelia species were represented as two well-supported clades (BS ≥ 98%, PP = 1), and they differed in the amount of genetic structure within. The Bayesian approach generated support (PP ≥ 0.95) for many minor subclades within each Punctelia species compared to the RAxML analysis. Punctelia borreri was the species with the highest number of these subclades, whereas P. subrudecta included a large polytomy. One of the few well-supported (BS = 85%, PP = 1) minor subclades within P. borreri included Iberian Peninsula and Italian specimens. A subclade located at the base of P. borreri with BS = 73% included specimens from eastern Asia (South Korea, China and India). Interestingly, the P. borreri clade hosted several sequences already deposited in GenBank identified as P. subrudecta and P. perreticulata. GenBank sequences labelled as P. borreri were not found within the P. subrudecta clade in our study. The specimen collected in Richmond (London, UK) corresponded to P. subrudecta, while the Italian specimen fell within P. borreri (data not shown).
Figure 3. Evolutionary relationships among specimens of Punctelia borreri and P. subrudecta depicted by the RAxML phylogram constructed using an nrITS sequence dataset. Thickened branches indicate bootstrap support (BS, RAxML) ≥ 70% and/or posterior probabilities (PP, MrBayes) ≥ 0.95. The haplotype code is provided for each specimen except for those with sequences that were not used for haplotype inference; haplotype codes are the same as those shown in networks in Fig. 4. Data on the ecology and bioclimate of the newly sequenced specimens is provided as symbols to the right of the specimen's label. The deviant taxonomic identity in some specimens of P. borreri as it appears in GenBank is shown in quotation marks. GenBank Accession numbers or project extraction codes are provided and more information on the specimens can be found in Supplementary Material Table S1 (available online). Flavopunctelia flaventior and F. soredica are included in the phylogram for rooting purposes. In colour online.
The topology of the inferred algal phylogeny (Fig. 4) included all the newly generated phycobiont sequences within the Trebouxia sp. OTU I05 clade (BS = 96%). This algal lineage corresponds with Trebouxia gelatinosa Ahmadjian ex Archibald, according to Muggia et al. (Reference Muggia, Nelsen, Kirika, Barreno, Beck, Lindgren, Lumbsch and Leavitt2020), and it is included in Trebouxia clade I (‘impressa’ clade). In the estimated phylogeny, it was sister to the still undescribed Trebouxia sp. OTU I15. Clade structure within T. gelatinosa is almost negligible: a basal supported (BS = 90%) subclade, including five GenBank phycobiont sequences obtained from diverse lichens on both sides of the Atlantic, was sister to a large unsupported (BS = 34%) polytomic clade that hosted all Punctelia phycobionts. Interestingly, within this large clade, a supported clade (BS = 98%) was detected with some Punctelia phycobiont sequences obtained from the northern and eastern Iberian Peninsula (haps18–20) as well as South Korea.
Figure 4. Phylogram inferred with RAxML using nrITS sequence data depicting the phylogenetic relationships of phycobiont haplotypes (numbers within yellow boxes) in the analyzed Punctelia species. Codes within blue and green boxes correspond with species-level lineages in Trebouxia clade I as proposed by Muggia et al. (Reference Muggia, Nelsen, Kirika, Barreno, Beck, Lindgren, Lumbsch and Leavitt2020). Bootstrap support (BS) is indicated with different coloured branches. GenBank Accession numbers or project extraction codes are provided and more information on the specimens can be found in Supplementary Material Table S1 (available online). The sample with the code I1297 has no associated haplotype number because the sequence was too short and was therefore not included in the haplotype inference analysis.
Haplotype networks
The numbers of myco- and phycobiont haplotypes inferred after removing poor sequence data from the original datasets were 22 and 21, respectively. The mycobiont network showed two subnetworks (Fig. 5A), each corresponding to a Punctelia species, with two likely connections between them characterized by a high number of mutations. The P. borreri subnetwork had a rare and central Chinese haplotype (hap5) to which other more abundant haplotypes from Eurasia were connected. Haplotype sharing was not observed among geographically distant territories (e.g. eastern Asia and the Mediterranean); however, within each of these areas, haplotypes were shared between the Iberian Peninsula and Italy (haps1, 2) and among China, India and South Korea (haps4, 6). A close connection of the Iberian Peninsula hap12 to various Asian haplotypes was also found. The P. subrudecta subnetwork, however, differed from that of the other species by showing an abundant, central haplotype (hap15) distributed throughout Europe and connected to haplotypes restricted to the Iberian Peninsula, Germany, China and New Zealand. Moreover, the phycobiont haplotype network that was specifically estimated for Trebouxia gelatinosa and was colour-labelled according to the associated fungal partner, revealed several subnetworks separated by seven or more mutations (Fig. 5B). Some of these subnetworks included Punctelia phycobionts, except for hap5 which was shared between P. borreri and an unidentified species of Xanthoria (Fr.) Th. Fr. Haplotypes 6 and 18 were shared between both Punctelia mycobionts. Finally, more distantly related haplotypes (haps1, 2, 3, 21) associated with other parmelioid lichens, as well as members of other lichen fungal families (Physciaceae and Teloschistaceae).
Figure 5. Statistical parsimony networks for haplotypes of Punctelia borreri and P. subrudecta mycobionts (A) and phycobionts (B) based on nrITS sequence data. Haplotypes are coloured according to the geographical origin of samples (A, mycobiont network), or the identity of the associated fungal partner (B, phycobiont network). The sizes of the circles in the networks are proportional to the number of individuals bearing the haplotype; black-filled smaller circles indicate missing haplotypes and mutations are shown as hatch marks. Haplotype codes follow Figs 2 and 3.
Interaction networks and estimated dʹ parameter
The bipartite interactions plot showed that each Punctelia species was associated with a different set of T. gelatinosa lineages (haplotypes), and that the most abundant mycobiont haplotypes (hap15 in P. subrudecta, and hap1 in P. borreri) were also the ones that interacted with the highest number of low frequency (rare) phycobiont haplotypes (Fig. 6). Trebouxia gelatinosa haplotypes 6 and 18, which were the most abundant in the dataset, were shared between both Punctelia species. Based on the current sampling, the phycobiont hap6 also showed exclusive interactions with rare P. subrudecta haplotypes (hap16, 20, 21). Although in our dataset both P. borreri and P. subrudecta interacted exclusively with a number of phycobiont haplotypes, specialization (dʹ) was low in both species (Table 1); the highest specialization value (dʹ = 0.738) was shown by P. borreri hap12.
Figure 6. Interaction network structure between Punctelia spp. phycobiont (upper green boxes) and mycobiont (lower brown boxes) haplotypes. Link width is proportional to the number of specimens forming the association. The links of microalgal haplotypes associating with both Punctelia borreri and P. subrudecta are marked in deeper shades of grey. The bioclimate of localities for each microalgal haplotype is indicated by symbols: orange triangle (Mediterranean) and dark green circle (Atlantic). In colour online.
Table 1. Summary of the specialization parameter dʹ of the bipartite network in Punctelia borreri and P. subrudecta. This parameter ranges from 1 (high specialization) to 0 (no specialization) and offers a measure of how much a mycobiont lineage discriminates from a random selection of algal partners.
Phycobiont ultrastructure
The examination of the microalga of P. subrudecta (voucher utiel_a1) by TEM revealed an irregularly shaped, non-circular alga with a very thick cell wall (Fig. 1E). The cytoplasm was rather sparse compared to other green algae, with most of the cell volume occupied by a prominent multi-lobed chloroplast. Thylakoid membranes were grouped in scattered but evident grana. The chloroplast showed a single, large decolorans-type pyrenoid (Bordenave et al. Reference Bordenave, Muggia, Chiva, Leavitt, Carrasco and Barreno2022; see also the gelatinosa-type pyrenoid of Friedl (Reference Friedl1989)) with pyrenotubules traversing the matrix in a parallel arrangement, with numerous associated pyrenoglobuli (Fig. 1F). Mycobiont-phycobiont cell interactions were interpreted as simple wall-to-wall appositions and type-1 intraparietal haustoria (Fig. 1G; Honegger Reference Honegger1986).
Discussion
The results of the present study show that, at least in the studied localities of the Iberian Peninsula, Punctelia subrudecta and P. borreri are associated with Trebouxia microalgae. This finding is consistent with the results of previous works that indicated a high specific relationship between parmelioid lichen fungi and these eukaryote microalgae (e.g. Leavitt et al. Reference Leavitt, Kraichak, Nelsen, Altermann, Divakar, Alors, Esslinger, Crespo and Lumbsch2015; Xu et al. Reference Xu, De Boer, Olafsdottir, Omarsdottir and Heidmarsson2020; Moya et al. Reference Moya, Molins, Škaloud, Divakar, Chiva, Dumitru, Molina, Crespo and Barreno2021; Kosecka et al. Reference Kosecka, Kukwa, Jabłońska, Flakus, Rodríguez-Flakus, Ptach and Guzow-Krzemińska2022). The phylogenetic approach indicated that all genetic lineages could be assigned to Trebouxia gelatinosa. This microalga species also associates with Flavoparmelia caperata (L.) Hale and F. soredians (Nyl.) Hale (Slocum et al. Reference Slocum, Ahmadjian and Hildreth1980; Molins et al. Reference Molins, Moya, García-Breijo, Reig-Armiñana and Barreno2018), two foliose parmelioid lichens that often co-occur with Punctelia spp. as epiphytes on tree trunks. The pyrenoid ultrastructure revealed in the present study was very similar to those of T. gelatinosa analyzed by Molins et al. (Reference Molins, Moya, García-Breijo, Reig-Armiñana and Barreno2018). Interestingly, the phylogenetic study and the haplotype network also indicated a close genetic relationship between the Punctelia phycobiont haplotypes and those associated with lichenized fungi of the families Physciaceae and Teloschistaceae (Nyati et al. Reference Nyati, Scherrer, Werth and Honegger2014). In fact, several species in the genera Xanthoria and Physcia (Schreb.) Michx. are common in European temperate and Mediterranean forest ecosystems growing on tree trunks and twigs. These parmelioid, teloschistoid and physcioid lichens, which show different dispersion strategies, might be horizontally linked by using the same phycobiont species (or even genetic lineages), therefore forming a likely photobiont-mediated guild sensu Rikkinen et al. (Reference Rikkinen, Oksanen and Lohtander2002).
The interaction between Trebouxia gelatinosa and P. subrudecta had already been reported by Piercey-Normore (Reference Piercey-Normore2009) in a Canadian forest, and Friedl et al. (Reference Friedl, Besendahl, Pfeiffer and Bhattacharya2000) in an unspecified geographical locality (but probably European). That interaction, as well as the one linking the same microalgal species with P. borreri in South Korea and India, is confirmed by the metadata associated with Trebouxia sequences deposited in GenBank. This evidence suggests that the two lichens are highly selective (sensu Beck et al. Reference Beck, Kasalicky and Rambold2002) for this phycobiont taxon across a large geographical range in the Northern Hemisphere that includes various bioclimatic and biogeographical areas. These results agree with other observations in which identical phycobiont strains associated with a single mycobiont in distant geographical areas (Fernández-Mendoza et al. Reference Fernández-Mendoza, Domaschke, García, Jordan, Martín and Printzen2011; Garrido-Benavent et al. Reference Garrido-Benavent, Pérez-Ortega, de los Ríos and Fernández-Mendoza2020; Moya et al. Reference Moya, Molins, Škaloud, Divakar, Chiva, Dumitru, Molina, Crespo and Barreno2021). The fact that P. borreri and P. subrudecta usually reproduce by soredia might explain why these lichens retain the same microalgal species as their photosynthetic partner; however, genetic differentiation may still be occurring, for example, owing to intrathalline somatic mutations occurring during algal cell division (Dal Grande et al. Reference Dal Grande, Alors, Divakar, Bálint, Crespo and Schmitt2014).
The bipartite network provided crucial information to analyze mycobiont-phycobiont interactions at the geographical scale, even though it is very regionally focused. Based on the current sampling, the two most abundant T. gelatinosa haplotypes occurred exclusively in one of the two studied regions (hap6 in eastern and hap18 in northern Iberian Peninsula), which display radically different bioclimatic conditions (Rivas-Martínez et al. Reference Rivas-Martínez, Penas, del Río, Díaz González, Rivas-Sáenz and Loidi2017). Although genetically close, these two phycobiont strains might show contrasting physiological performances (Casano et al. Reference Casano, del Campo, García-Breijo, Reig-Armiñana, Gasulla, Del Hoyo, Guéra and Barreno2011; Sadowsky et al. Reference Sadowsky, Mettler-Altmann and Ott2016). Moreover, in each region, these phycobiont strains were shared by both Punctelia species. Our results, therefore, might point towards an association between the two studied Punctelia with locally adapted phycobionts, as already observed in other lichenized fungi (Blaha et al. Reference Blaha, Baloch and Grube2006; Fernández-Mendoza et al. Reference Fernández-Mendoza, Domaschke, García, Jordan, Martín and Printzen2011; Rolshausen et al. Reference Rolshausen, Dal Grande, Sadowska-Deś, Otte and Schmitt2018; Moya et al. Reference Moya, Molins, Škaloud, Divakar, Chiva, Dumitru, Molina, Crespo and Barreno2021). Local adaptation might also explain our finding of a different T. gelatinosa haplotype interacting with P. subrudecta in Richmond (London, UK). In addition, some algal nrITS haplotypes occurred in either the northern or eastern Iberian Peninsula, and this agrees with previous reports of a single Trebouxia species or even infraspecific lineages occurring in geographically and climatically very different regions, such as Macaronesia, the Iberian Peninsula, Central Europe and Scandinavia (e.g. Trebouxia lynnae; Singh et al. Reference Singh, Dal Grande, Divakar, Otte, Crespo and Schmitt2017; Ertz et al. Reference Ertz, Guzow-Krzemińska, Thor, Łubek and Kukwa2018; Blázquez et al. Reference Blázquez, Hernández-Moreno, Gasulla, Pérez-Vargas and Pérez-Ortega2021; Molins et al. Reference Molins, Moya, Muggia and Barreno2021). Ramalina species in the Macaronesian region also showed a pattern in which specimens associated indistinctly with different haplotypes of the same Trebouxia species (Blázquez et al. Reference Blázquez, Hernández-Moreno, Gasulla, Pérez-Vargas and Pérez-Ortega2021). Finally, the number of inferred phycobiont nrITS haplotypes was relatively high considering the limited number of specimens in our dataset. Although this limitation may bias our previous interpretations of fungal-algal interactions across the geographical scale, our results nevertheless fit with the general evidence that the Mediterranean region harbours a high diversity of phycobionts, probably due to the mosaic of existing macro- and microclimates (Fernández-Mendoza et al. Reference Fernández-Mendoza, Domaschke, García, Jordan, Martín and Printzen2011; Garrido-Benavent et al. Reference Garrido-Benavent, Pérez-Ortega, de los Ríos and Fernández-Mendoza2020). Future work will consider more Punctelia species and specimens to evaluate whether our interpretations are also supported in a wider European geographical context. Additionally, multilocus or even genome-scale datasets for the phycobionts are needed to refine species boundaries which might improve inferences from fungal-algal interaction patterns.
The phylogenetic study of the mycobionts based on nrITS data revealed a clear separation of the two Punctelia species, which aligns with previous findings by Longán et al. (Reference Longán, Barbero and Gómez-Bolea2000) on the discrimination of P. subrudecta and P. borreri solely on the basis of secondary metabolite profiling through thin-layer chromatography (TLC). The range of variation in other features, such as the lower cortex colour, thallus morphology and the size and morphology of pycnidia and conidia in the two Punctelia species often overlaps (Longán et al. Reference Longán, Barbero and Gómez-Bolea2000; Truong & Clerc Reference Truong and Clerc2003). With the currently available molecular data for many specimens, it would be worth checking whether these chemical, morphological and reproductive differences are really discriminatory and are supported phylogenetically.
Finally, the coloured reactions observed macroscopically on Punctelia thalli might be due to the impact of pollutants on thylakoid membranes of the trebouxioid algae, which release K+ ions that interact with lichen substances inducing these coloured stains (Sigal & Nash Reference Sigal and Nash1983; Boonpragob & Nash Reference Boonpragob and Nash1990a, Reference Boonpragob and Nashb; Ahroun et al. Reference Arhoun, Barreno and Ramis-Ramos2000). It is highly likely that the forest systems of the eastern Iberian Peninsula, which are exposed to the Mediterranean winds, are impacted by photo-oxidants and nitrogen deposition (Millán et al. Reference Millán, Salvador and Mantilla1997; Giordani et al. Reference Giordani, Calatayud, Stofer, Seidling, Granko and Fischer2014). Future research must explore how those elements affect the different components of the holobionts in these foliose lichen communities (Barreno & Manzanera Reference Barreno, Manzanera and Barba2023).
Acknowledgements
We thank J. G. Segarra-Moragues (UV, ES) for contributing to collections and also to two anonymous reviewers for their useful comments. This project was supported by Generalitat Valenciana (GVA, excellence in research: PROMETEO/2021/005 to EB) and the Ministerio de Ciencia e Innovación of the Spanish Government (PID2021-127087NB-I00 to IGB). We dedicate this article to Dr Pier Luigi Nimis in honour of his retirement.
Author ORCIDs
Isaac Garrido-Benavent, 0000-0002-5230-225X; María Reyes Mora-Rodríguez, 0000-0002-3044-2817; Salvador Chiva, 0000-0002-1615-6443; Simón Fos, 0000-0003-0976-5459; Eva Barreno, 0000-0002-2622-4550.
Competing Interests
The authors declare none.
Supplementary Material
The Supplementary Material for this article can be found at https://doi.org/10.1017/S0024282923000269.
Open access
