Punctelia borreri and P. subrudecta (Parmeliaceae) associate with a partially overlapping pool of Trebouxia gelatinosa lineages

Abstract

An increasing number of studies are describing the diversity of lichen phycobionts, which is leading to a better understanding of how lichen communities are assembled at different taxonomic, evolutionary and geographical scales. The present study explores the identity and genetic diversity of the microalgal partners of Punctelia borreri and P. subrudecta, two tropical and temperate parmelioid lichen fungi that often grow in temperate and Mediterranean forest ecosystems in Europe. Based on a specimen sampling distributed in two climatically divergent regions in the Iberian Peninsula, we found that these mycobionts are associated with Trebouxia gelatinosa, whose identity was also confirmed by an ultrastructural study of the pyrenoid. The bipartite network analysis indicated that each Punctelia species was associated with a different set of low frequency T. gelatinosa infraspecific lineages, whereas the two most abundant phycobiont lineages were shared between both mycobionts. Based on the current sampling, these two algal lineages occur exclusively in one of the two studied regions, which might point towards climate-driven, fine-tuned fungal-algal interactions. Finally, we documented visible symptoms of injury on the thalli in areas likely to have been impacted by air pollution.

Introduction

The family Parmeliaceae is amongst the most diverse lineages of lichenized ascomycete fungi, with c. 2800 species known worldwide (Thell et al. 2012; Jaklitsch et al. 2016; Stelate et al. 2022). These often form foliose or fruticose lichen thalli, and occur in a wide range of ecosystems including arctic-alpine (Goward & Myllys 2020; Xu et al. 2020), Antarctic (Øvstedal & Lewis-Smith 2001), tropical and temperate forests (Molina et al. 2017; Del-Prado et al. 2019) or even in fog-influenced desert regions such as the Namibian lichen fields, where some species thrive on quartz pebbles (Wirth 2010). 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. 2008; Nimis 2016; Nimis & Martellos 2022). 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. 2016; Widhelm et al. 2016; Lutsak et al. 2020) or apparently well-known species (Lendemer & Hodkinson 2010; Castellani et al. 2021; Stelate et al. 2022). 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. 2012), 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. 2019).

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. (2015) 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. 2021), 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. 2011; Garrido-Benavent et al. 2020). These works highlighted that switching to locally adapted phycobionts might have allowed lichenized fungi to widen their geographical distribution (Peksa & Škaloud 2011), although this does not seem to always be the case, as demonstrated in other lichenized fungal groups (Škvorová et al. 2022). Phycobiont switching has also been found to cause changes in reproduction type and phenotypic variability (Ertz et al. 2018). 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 1982). Punctelia currently hosts c. 45 species (Thell et al. 2012), and its distribution centre is probably subtropical South and North America and Africa (Krog 1982). Although the boundaries of some species in this genus have already been assessed phylogenetically (Crespo et al. 2004; Alors et al. 2016), there are limited data on the identity and range of associated phycobionts. The studies by Leavitt et al. (2015), Xu et al. (2020) and Kosecka et al. (2022) 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. (2020)) 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. 2017). 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 2003; but see Longán et al. 2000). 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. 2017). 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. (2000) and Smith et al. (2009). 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. (2018). 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 mm² 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 1993) and ITS4A (Larena et al. 1999); for the phycobiont, the primer pair was nr-SSU-1780 (Piercey-Normore & DePriest 2001) and ITS4T (Kroken & Taylor 2001). 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. 1997) 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. 2002; Katoh & Standley 2013), 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. (2016). 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. 2020). 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. 2015) 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 2006; Stamatakis et al. 2008), 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 2021) 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. 2012) 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. 2012), with two parallel, simultaneous four-chain runs executed over 5 × 10⁷ 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 2009) 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. (2020)). The TCS method implemented in PopART v. 1.7 (Templeton et al. 1992; Leigh & Bryant 2015) 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. 2008). Fungal haplotype specialization towards phycobiont nrITS haplotypes was estimated using the parameter dʹ (Blüthgen et al. 2006), 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 × 10⁶ 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. (2020), 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. (2020). 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. 2022; see also the gelatinosa-type pyrenoid of Friedl (1989)) 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 1986).

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. 2015; Xu et al. 2020; Moya et al. 2021; Kosecka et al. 2022). 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. 1980; Molins et al. 2018), 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. (2018). 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. 2014). 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. (2002).

The interaction between Trebouxia gelatinosa and P. subrudecta had already been reported by Piercey-Normore (2009) in a Canadian forest, and Friedl et al. (2000) 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. 2002) 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. 2011; Garrido-Benavent et al. 2020; Moya et al. 2021). 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. 2014).

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. 2017). Although genetically close, these two phycobiont strains might show contrasting physiological performances (Casano et al. 2011; Sadowsky et al. 2016). 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. 2006; Fernández-Mendoza et al. 2011; Rolshausen et al. 2018; Moya et al. 2021). 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. 2017; Ertz et al. 2018; Blázquez et al. 2021; Molins et al. 2021). 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. 2021). 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. 2011; Garrido-Benavent et al. 2020). 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. (2000) 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. 2000; Truong & Clerc 2003). 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 1983; Boonpragob & Nash 1990a, b; Ahroun et al. 2000). 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. 1997; Giordani et al. 2014). Future research must explore how those elements affect the different components of the holobionts in these foliose lichen communities (Barreno & Manzanera 2023).