Introduction
Lichenized fungi form mutualistic relationships with photo-autotrophic organisms (photobionts), mainly green algae (Trebouxiophyceae and Ulvophyceae) and/or cyanobacteria. The lichen symbiosis has been highly successful within fungi, especially Ascomycota, with c. 19 400 currently accepted species (Lücking et al. Reference Lücking, Hodkinson and Leavitt2017) and an estimated diversity of more than 28 000 species (Lücking et al. Reference Lücking, Rivas Plata, Chaves, Umaña and Sipman2009; Leavitt et al. Reference Leavitt, Esslinger, Spribille, Divakar and Lumbsch2013). Additionally, lichens are commonly used to assess environmental disturbance, serving as bioindicators of air pollution, forest age and health, and climate change (Nimis et al. Reference Nimis, Scheidegger and Wolseley2002; Crespo et al. Reference Crespo, Divakar, Argüello, Gasca and Hawksworth2004; Giordani & Brunialti Reference Giordani, Brunialti, Upreti, Divakar, Shukla and Bajpai2015; Sujetovienė Reference Sujetovienė, Upreti, Divakar, Shukla and Bajpai2015; Sancho et al. Reference Sancho, Pintado and Green2019; Abas Reference Abas2021).
Recognizing phylogenetic relationships and delimiting species in lichens are crucial for ecological and conservation studies, assessing biotic diversity, and identifying factors driving diversification. They are also important for future investigations because phylogenetic differences may not be fully reflected in the phenotype.
Traditionally, to infer taxonomic boundaries in lichen-forming fungi, thin-layer chromatography (Culberson Reference Culberson1972), morphology and the expression of signature secondary metabolites, and isolation and identification of lichen substances have been used (Huneck & Yoshimura Reference Huneck and Yoshimura1996; Huneck Reference Huneck1999). However, these characters are highly variable and their homology has proved difficult to assess. For example, it has been shown that apothecial form and spore wall thickness have changed in parallel within Pertusaria. The Pertusaria-type ascus is plesiomorphic within the Pertusariaceae and thus cannot be used to circumscribe Pertusaria (Lumbsch & Schmitt Reference Lumbsch and Schmitt2001). However, with the use of molecular tools, it has been shown that phenotype-based taxonomy may not reflect natural groups, including cases in which morphologically distinct forms formerly recognized as distinct species are shown to represent a polymorphic species (see e.g. Boluda et al. Reference Boluda, Rico, Divakar, Nadyeina, Myllys, McMullin, Zamora, Scheidegger and Hawksworth2019), and other cases where several morphologically cryptic species are masked within a single nominal taxon (Del-Prado et al. Reference Del-Prado, Divakar, Lumbsch and Crespo2016).
Within lichen-forming Ascomycetes, Parmeliaceae (Lecanorales) constitutes one of the largest and best-studied families (Crespo et al. Reference Crespo, Lumbsch, Mattsson, Blanco, Divakar, Articus, Wiklund, Bawingan and Wedin2007, Reference Crespo, Divakar and Hawksworth2011; Thell et al. Reference Thell, Crespo, Divakar, Kearnefelt, Leavitt, Lumbsch and Seaward2012; Divakar et al. Reference Divakar, Crespo, Wedin, Leavitt, Hawksworth, Myllys, McCune, Randlane, Bjerke and Ohmura2015). This family is usually characterized morphologically by a certain type of ascoma ontogeny and the presence of an ascomatal structure known as a cupulate exciple (Henssen & Jahns Reference Henssen and Jahns1974). Parmeliaceae also comprises species which are frequently used in biomonitoring studies of atmospheric pollution, such as Parmelia sulcata, Flavoparmelia caperata, Parmotrema perlatum and Punctelia subrudecta (Hawksworth & Rose Reference Hawksworth and Rose1970; Crespo et al. Reference Crespo, Divakar, Argüello, Gasca and Hawksworth2004; De La Cruz et al. Reference De La Cruz, De La Cruz, Tolentino and Gioda2018). The application of phylogenetic analysis based on molecular (DNA) characters to delimit species allows us to determine a posteriori which types of phenotypic characters are good predictors of phylogenetic species and demonstrate how these characters evolve in this family and in lichenized fungi in general. These molecular data have led to the recognition of morphologically cryptic species, such as Parmelia serrana (Molina et al. Reference Molina, Crespo, Blanco, Lumbsch and Hawksworth2004), P. barrenoae (Divakar et al. Reference Divakar, Molina, Lumbsch and Crespo2005a), P. encryptata (Molina et al. Reference Molina, Divakar, Millanes, Sanchez, Del-Prado, Hawksworth and Crespo2011a) and Melanelixia californica (Divakar et al. Reference Divakar, Figueras, Hladun and Crespo2010), and conversely also to the union of species traditionally regarded as morphologically distinct (see e.g. Boluda et al. Reference Boluda, Rico, Divakar, Nadyeina, Myllys, McMullin, Zamora, Scheidegger and Hawksworth2019). Recently, coalescent-based species delimitation approaches have shown to be well suited to critically evaluate species boundaries in Parmeliaceae, as well as lichens in general (Leavitt et al. Reference Leavitt, Moreau, Lumbsch, Upreti, Divakar, Shukla and Bajpai2015). Furthermore, these methods can accurately display relationships when incomplete lineage sorting and gene tree heterogeneity hide phylogenetic relationships among species (Knowles & Carstens Reference Knowles and Carstens2007; Camargo et al. Reference Camargo, Morando, Avila and Sites2012). Commonly used methods to critically evaluate species boundaries include the poisson tree processes (PTP) model (Zhang et al. Reference Zhang, Kapli, Pavlidis and Stamatakis2013), the automatic barcode gap discovery (ABGD) (Puillandre et al. Reference Puillandre, Lambert, Brouillet and Achaz2012), the general mixed Yule coalescent model (GMYC) (Pons et al. Reference Pons, Barraclough, Gomez-Zurita, Cardoso, Duran, Hazell, Kamoun, Sumlin and Vogler2006; Monaghan et al. Reference Monaghan, Wild, Elliot, Fujisawa, Balke, Inward, Lees, Ranaivosolo, Eggleton and Barraclough2009), and SpedeSTEM (Ence & Carstens Reference Ence and Carstens2011).
Parmotrema (Massalongo Reference Massalongo1860) is one of the largest genera in the parmelioid group of the family Parmeliaceae. It includes c. 300 described species with an apparent centre of speciation in the Pacific Islands, tropical and subtropical regions of South America (Spielmann & Marcelli Reference Spielmann and Marcelli2020). The species of the genus are characterized by a pored epicortex, large thalli with broad lobes, a broad, naked marginal zone on the lower surface, and large, thick-walled, ellipsoid ascospores, sublageniform or filiform conidia (Elix Reference Elix1993), and (commonly) marginal cilia (Hale Reference Hale1974). Reproductive strategies vary among Parmotrema taxa. Sexual reproduction is restricted to characteristic fungal fruiting bodies (ascomata) producing ascospores. Ascospores are dispersed independently of the photosynthesizing partner (photobiont) and require reacquisition of the appropriate photobiont partner to re-establish the lichenized condition. Species within Parmotrema also commonly propagate asexually by means of vegetative diaspores, either isidia or soredia. These specialized vegetative reproductive propagules contain both fungal and algal symbionts, eliminating the need for independent acquisition of the appropriate photobiont.
In a molecular phylogeny of parmotremoid lichens (Blanco et al. Reference Blanco, Crespo, Divakar, Elix and Lumbsch2005), a single sample of both Parmotrema perlatum (Huds.) M. Choisy and Parmotrema crinitum (Ach.) M. Choisy, originally from Portugal, were included and these formed a well-supported monophyletic group, indicating that these specimens could be conspecific. Therefore, a critical evaluation of species boundaries is necessary.
Parmotrema crinitum is characterized by the presence of coralloid branched, apically ciliate isidia or often eciliate isidia (Divakar & Upreti Reference Divakar and Upreti2005), and the stictic acid complex in the medulla. According to Elix (Reference Elix1994), P. crinitum is a cosmopolitan species that is widespread in temperate, tropical regions and even in high humidity, sub-boreal forests (Elix Reference Elix1994; Kurokawa & Lai Reference Kurokawa and Lai2001). Many European countries have reported the presence of P. crinitum (Jablońska et al. Reference Jablońska, Oset and Kukwa2009) as have some Asian countries such as Japan (Yoshimura Reference Yoshimura1974), China (Wei Reference Wei1991) and Taiwan (Wang-Yang & Lai Reference Wang-Yang and Lai1976).
Parmotrema perlatum is a greenish grey foliose lichen, saxicolous or corticolous, loosely adnate, with rounded lobes. It can be recognized by its broad lobes, irregularly branched (7.0 cm), laterally overlapping, with frequent black cilia (0.20–)0.50–1.00 × 0.02(–0.10) mm, evenly distributed but less common in the lobe apices. Its soralia are marginal and linear, sometimes widely distributed or subcontinuous, concolorous with the thallus, and the medulla is white. The thallus contains atranorin, stictic acid, hypostictic acid, menegazziaic acid and norstictic acid (Spielmann & Marcelli Reference Spielmann and Marcelli2009). Parmotrema perlatum is a widespread species in temperate regions of the Northern and Southern Hemispheres: Asia (Hale Reference Hale1965; Kurokawa Reference Kurokawa1991; Kurokawa & Lai Reference Kurokawa and Lai2001), Europe (Hale Reference Hale1965), Africa (Hale Reference Hale1965; Swinscow & Krog Reference Swinscow and Krog1988), North America (Hale Reference Hale1965; Brodo et al. Reference Brodo, Sharnoff and Sharnoff2001; Nash & Elix Reference Nash, Elix, Nash, Ryan, Gries and Bungartz2002), Central America (Hale Reference Hale1965), South America (Hale Reference Hale1965), Australia (Kantvilas Reference Kantvilas2019) and New Zealand (Blanchon Reference Blanchon2013).
Morphologically, P. crinitum and P. perlatum can easily be separated based on the characteristics associated with their different dispersal strategies. Parmotrema crinitum has apically ciliate isidia or often eciliate isidia, and P. perlatum has marginal soralia. Parmotrema perlatum is one of the most common and widely utilized lichens in the Ayurvedic system of medicine and has been overexploited for uses in traditional medicine in India (Kumar & Upreti Reference Kumar and Upreti2001). In India, this species is currently considered threatened (Divakar & Upreti Reference Divakar and Upreti2005).
The aim of this study was to evaluate the phylogenetic relationships between P. crinitum and P. perlatum, and elucidate the possible monophyly between both species. Additionally, PTP, ABGD and genetic distance analyses were included to critically assess species boundaries. The ITS data matrix consisted of a total of 74 sequences of Parmotrema species, including 11 of P. crinitum and 35 of P. perlatum, plus two outgroup sequences of Crespoa carneopruinata. Morphological and chemical features of all samples were critically examined.
Material and Methods
Taxon sampling
Sequence data of the ITS locus were analyzed from 76 specimens, of which 11 sequences were of Parmotrema crinitum (two individuals were newly sequenced and nine sequences downloaded from GenBank) and 35 of Parmotrema perlatum (20 individuals were sequenced in this study and 15 sequences downloaded from GenBank). The new specimens were collected from distant geographical regions throughout the species distributions (Table 1). Additionally, 30 ITS sequences from 30 specimens, belonging to 14 different species (6 were newly sequenced and 24 sequences downloaded from GenBank), are included in this study to test the monophyly of both species within the genus. This total included Crespoa carneopruinata which was selected as outgroup since it has previously been shown to be the sister group of Parmotrema (Divakar et al. Reference Divakar, Crespo, Wedin, Leavitt, Hawksworth, Myllys, McCune, Randlane, Bjerke and Ohmura2015). Details of the material studied, including GenBank Accession numbers, are shown in Table 1. Species recognition was based mainly on morphological and chemical characters.
Parmotrema species used in this study with GenBank Accession numbers of the ITS sequences and voucher information for the specimens. Crespoa carneopruinata was included in the phylogenetic analyses as outgroup. * = newly generated sequences from specimens for which morphological and chemical characters are given in Table S1 (available online).
Chemistry and morphology
Thallus morphology of all new specimens of Parmotrema crinitum and P. perlatum included in the molecular analyses was studied using a Nikon SMZ-1000 stereomicroscope to identify morphological characteristics of the thallus: lobe shape, isidia, soralia, cilia and rhizines. This is because P. crinitum and P. perlatum are traditionally differentiated based on these features. Photographs were taken with a Nikon 105 mm f/2.8D AF Micro-Nikkor lens coupled to a Nikon D90 camera.
Spot tests were carried out on the medulla with usual chemical reagents such as aqueous potassium hydroxide (K), Steiner's stable paraphenylenediamine (PD) and aqueous calcium hypochlorite (C).
Thin-layer chromatography (TLC) was carried out following Orange et al. (Reference Orange, James and White2010). We used TLC solvent system C (200 ml toluene/30 ml acetic acid) (Elix & Ernst-Russell Reference Elix and Ernst-Russell1993), with concentrated acetone extracts at 50 °C spotted onto silica gel 60 F254 aluminum sheets (Merck, Darmstadt). The aluminum sheets were dried for 10 min in an acetic acid atmosphere to maximize resolution.
DNA extraction, PCR and sequencing
Total DNA was extracted from a single, clean (under a dissecting microscope) lichen lobe using the DNeasy Plant Mini Kit (Qiagen, Barcelona) with a slight modification to the manufacturer's instruction (Crespo et al. Reference Crespo, Blanco and Hawksworth2001). Genomic DNA (5–25 ng, quantified using a quantitative PCR machine) was used for PCR amplifications of the ITS region. Standard PCR amplifications were conducted in 25 μl reaction volumes containing 12.5 μl of Master Mix (50 units/ml of Taq DNA polymerase supplied in a proprietary reaction buffer (pH 8.5), 400 μM dATP, 400 μM dGTP, 400 μM dCTP, 400 μM dTTP, 3 mM MgCl2), and 1.5 μl of each primer at 10 μM, 4.5 μl of water (H2O) and 5 μl of DNA template. Fungal nuclear internal transcribed spacer (ITS) rDNA was amplified using the primer pair ITS1F (5ʹ [CTT GGT CAT TTA GAG GAA GTA A] 3ʹ) (Gardes & Bruns Reference Gardes and Bruns1993)/ LR1 (5ʹ [GGT TGG TTT CTT TTC CT] 3ʹ) (Vilgalys & Hester Reference Vilgalys and Hester1990). We also tried primer pair ITS1-LM (5ʹ [GAA CCT GCG GAA GGA TCA TT] 3ʹ) (Myllys et al. Reference Myllys, Lohtander and Tehler2001) and ITS2-KL (5ʹ [ATG CTT AAG TTC AGC GGG TA] 3ʹ) (Lohtander et al. Reference Lohtander, Källersjö and Tehler1998), but we had better results with primer pair ITS1F/LR1. For any failed samples, we tried a nested PCR using the primer pair ITS1F/LR1 for the first PCR, and 5 μl of this PCR product as DNA template for the second PCR using ITS1-LM/ITS2-KL. With this nested PCR, we aimed to amplify a smaller region than with the previous primers. However, it did not yield any additional products.
Amplification was run in an automatic thermocycler (Techne Progene, Jepson Bolton & Co., Watford, Herts, UK) using the following parameters: initial denaturation at 94 °C for 5 min followed by 35 cycles at 94 °C for 1 min, 54 °C (ITS1F/LR1) for 1 min, and 72 °C for 1 min 30 s, with a final extension at 72 °C for 10 min. Amplification products were visualized on 1% agarose gel stained with SYBR green DNA (Life Technologies Corporation, Grand Island, New York, USA), and subsequently purified using ExoSAP-IT (GE Healthcare, Chalfont St Giles, UK) according to the manufacturer's instructions. Sequencing was performed using BigDye Terminator reaction kit (ABI PRISM, Applied Biosystems, Waltham, Massachusetts, USA). Cycle sequencing reactions were performed with the same sets of primers used in the amplification step. Sequencing reactions were electrophoresed on a 3730 DNA Analyzer (Applied Biosystems) at the Unidad de Genómica (Parque Científico de Madrid).
Sequence alignments and phylogenetic analyses
Sequence fragments generated for this study were assembled and edited using the program SeqMan v. 7 (Lasergene R, DNASTAR, Madison, Wisconsin, USA). Sequence identity was confirmed using the megaBLAST search function in GenBank. For the alignment we prepared a matrix of 480 base pairs (bp). Then we used the program MAFFT v. 7 (Katoh & Standley Reference Katoh and Standley2013) with the parameters set to default. If sequences had different lengths, only the part shared by all the sequences was used, and after manually removing ambiguously aligned nucleotide positions, we kept 451 unambiguously aligned base pairs for the final matrix that we used as input for the phylogenetic reconstruction. The alignment was analyzed using maximum likelihood (ML) and a Bayesian Markov chain Monte Carlo (B/MCMC) approach. The ML analysis was performed using an online version of the program RAxML-HPC BlackBox v. 8.2.8 (Stamatakis Reference Stamatakis2006; Stamatakis et al. Reference Stamatakis, Hoover and Rougemont2008), implemented on the Cipres science gateway (https://www.phylo.org/portal2/home.action; Miller et al. Reference Miller, Pfeiffer and Schwartz2010). We used a GTRGAMMA model, which includes a parameter (Γ) for rate heterogeneity among sites but chose not to include a parameter to estimate the proportion of invariable sites (Stamatakis Reference Stamatakis2006; Stamatakis et al. Reference Stamatakis, Hoover and Rougemont2008). Support values were assessed using the ‘rapid bootstrapping’ option with 1000 replicates.
For the Bayesian reconstruction, we used MrBayes v. 3.2.1 (Ronquist & Huelsenbeck Reference Ronquist and Huelsenbeck2003). The analysis was performed assuming a discrete gamma distribution with six rate categories (GTR + G). The nucleotide-substitution model parameters were selected using the Akaike Information Criterion as implemented in jModelTest (Posada Reference Posada2008), molecular clock not assumed. A run with 10 million generations, starting with a random tree and employing 12 simultaneous chains, was executed. Trees were saved to a file every 200th generation. The first 2 million generations (i.e. 20 000 trees) were deleted as the ‘burn-in’ of the chains. We plotted the log-likelihood scores of sample points against generations using Tracer v. 1.5 (Rambaut & Drummond Reference Rambaut and Drummond2007) and determined that stationarity had been achieved when the log-likelihood values of the sample points reached an equilibrium value (Huelsenbeck & Ronquist Reference Huelsenbeck and Ronquist2001). The trees obtained before stationarity were discarded. Posterior probabilities (PPs) were obtained from the 50% majority-rule consensus of sampled trees after excluding the initial 20% as burn-in. Only clades that received bootstrap support ≥ 70% in the ML analyses and PP ≥ 0.95 were considered strongly supported. The phylogenetic tree was drawn using FigTree v. 1.2.3 (Rambaut Reference Rambaut2009) and modified with CorelDRAW v. 8.
Candidate species identification
In order to establish candidate species limits in the phylogenetic tree, three computational approaches were used:
1) Poisson tree processes (PTP): this method does not require an ultra-metric tree, as the transition point between intra- and inter-specific branching rates is identified directly using the number of nucleotide substitutions (Zhang et al. Reference Zhang, Kapli, Pavlidis and Stamatakis2013). PTP incorporates the number of substitutions in the model of speciation and assumes that the probability of a substitution giving rise to a speciation event follows a Poisson distribution. The branch lengths of the input tree are supposed to be generated by two independent classes of Poisson events, one corresponding to speciation and the other to coalescence. The ML phylogeny of the ITS dataset obtained with RAxML was used as the input trees to run PTP species delimitation analysis on the PTP web server (http://species.h-its.org/ptp/). We ran the PTP analysis for 100 000 MCMC generations, with a thinning value of 100 and a burn-in of 25%. Outgroup taxa were removed for species delimitation.
2) Automatic barcode gap discovery (ABGD): this is an automatic procedure that sorts the sequences into hypothetical species based on the barcode gap. This method automatically finds the distance where the barcode gap is located (Puillandre et al. Reference Puillandre, Lambert, Brouillet and Achaz2012). The ABGD method was carried out on the ITS dataset using the web interface at http://wwwabi.snv.jussieu.fr/public/abgd/abgdweb.html. Default parameters were chosen using Kimura 2-parameter (K2P) distances that correct for transition rate bias (relative to transversions) in the substitution process (Kimura Reference Kimura1980). The default for the minimum relative gap width was set to different values between 0.1 and 0.15.
3) Genetic distances: pairwise ML distances (given as the number of nucleotide substitutions per site) among the ITS rDNA sequences of Parmotrema crinitum and P. perlatum were calculated. Genetic distances were calculated with TREE-PUZZLE v. 5.2 (Schmidt et al. Reference Schmidt, Strimmer, Vingron and von Haeseler2002) using the GTR model of nucleotide substitution, assuming a discrete gamma distribution with six rate categories. The program generates an output file which consists of a triangular matrix with all pairwise distances between all the samples included. This matrix was visualized with Microsoft Excel 2010 and genetic distances between different specimens of P. crinitum and P. perlatum were manually identified. Candidate species were proposed based on the threshold of 0.016 substitutions per site (s/s) which separates intra- and interspecific distances in parmelioid lichens (Del-Prado et al. Reference Del-Prado, Cubas, Lumbsch, Divakar, Blanco, de Paz G, Molina and Crespo2010). Distance values in the matrix ≤ 0.016 s/s have been considered as the values between samples of a single species. By using the Excel filter, we separated values ≤ 0.016, providing for every specimen included in the analysis a group of specimens that share the values characterizing its species range.
Results
Morphological and chemical analyses
Morphological and chemical characters of the specimens (P. crinitum and P. perlatum) newly generated and included in the phylogenetic analysis are presented in Supplementary Material Table S1, available online. Parmotrema crinitum and P. perlatum could not be distinguished based on the colour of the thalli, presence of cilia (Fig. 1A & B) and colour of the lower surface (Fig. 1E & F). However, the reproductive structures allowed unequivocal separation into the two species; Parmotrema perlatum has soralia (Fig. 1C) and P. crinitum has isidia (Fig. 1D).
Main morphological characters of P. perlatum (A, C, E) and P. crinitum. (B, D, F), as indicated with arrows. A & B, black, simple cilia. C, marginal soralia concolorous with the thallus. D, isidia with apical cilia. E & F, black lower surface with brown naked zone peripherally. In colour online.
The TLC test showed that both species contained the stictic acid complex (stictic acid, menegazzic acid and hypostictic acid) and atranorin. The spot tests of P. crinitum were similar to those of P. perlatum (Fig. 2, Supplementary Material Table S1). However, the samples Parmotrema perlatum 34 and P. perlatum 35 lacked menegazzic acid (nos 10 and 12 in Fig. 2). The spot test was congruent with other samples of P. perlatum, and the morphology of both samples showed the same characteristics of P. perlatum, although the cilia were not as abundant as is usual in this species.
Thin-layer chromatography profile of Parmotrema species (included in the DNA analysis) in solvent system C. P is the control, Pleurosticta acetabulum. Lanes 1 & 7, P. dilatatum. Lane 2, P. flavotinctum. Lane 3, P. mellissii. Lane 4, P. paulense. Lanes 5, 6, 8, 11 & 19, P. crinitum. Lane 9, empty. Lanes 10 & 12, Parmotrema perlatum. Lanes 13, 14, 15, 16, 17, 18 & 20, P. perlatum; these have the same TLC result as P. crinitum but were distinguished from each other based on the morphology. a = norstictic acid; b = atranorin; c = protocetraric acid; d = unknown; e = hypostictic acid; f = stictic acid; g = menegazzic acid. In colour online.
Phylogenetic analysis
The ITS data matrix consisted of a total of 74 sequences of Parmotrema species, including 46 of P. crinitum and P. perlatum (11 sequences of P. crinitum and 35 sequences of P. perlatum), with two sequences of Crespoa carneopruinata were included as outgroup (Table 1). The final matrix used as input for the phylogenetic reconstruction contained 451 unambiguously aligned base pairs (bp). The tree reconstruction (Fig. 3) comprised the following species of Parmotrema, included to test the monophyly of P. crinitum and P. perlatum: Parmotrema clavuliferum, P. dilatatum, P. flavotinctum, P. grayanum, P. haitiense, P. hypoleucinum, P. mellissii, P. paulense, P. perforatum, P. pseudoreticulatum, P. reticulatum, P. robustum, P. subtinctorium and P. tinctorum.
Phylogenetic relationships between Parmotrema crinitum and P. perlatum. The green (thickened) line indicates the clade supported by maximum likelihood (ML) and Bayesian analyses; * = clades supported only by ML. ** = clades supported only by Bayesian analysis. Alphanumeric labels indicate clades and clusters. Species delimitation scenarios obtained from PTP, ABGD and genetic distances are indicated in columns to the right. In colour online.
The phylogenetic trees estimated from ML and Bayesian methods did not show any well-supported conflict; ML topologies are presented with bootstrap and B/MCMC analysis with posterior probability (ML bootstrap ≥ 70%; PP ≥ 0.95 in B/MCMC analysis) (Fig. 3). The LnL value was −1938.068210 for ML, and −2184.428 for the Bayesian analysis. The phylogenetic analysis showed that all samples of Parmotrema included in this analysis formed a monophyletic group. Within this genus the samples of Parmotrema crinitum and P. perlatum were grouped in one monophyletic group (clade A), with the exception of two samples of P. perlatum that were grouped in a separate clade B. While clade A comprised specimens of P. crinitum and P. perlatum collected from various geographical regions, clade B included only two specimens of P. perlatum collected from Tenerife (Canary Islands) and Lisbon (Portugal).
Although four well-supported monophyletic clusters can be recognized in clade A, the phylogenetic relationships among them remained unresolved (Fig. 3).
Cluster A1 (Fig. 3) contained 15 specimens of Parmotrema perlatum, most of them collected from Morocco, except four from the Canary Islands and two from Portugal and Turkey. Cluster A2 included eight specimens of Parmotrema crinitum from different geographical regions (Table 1). Cluster A3 grouped 18 samples of Parmotrema perlatum from the Canary Islands (Tenerife, Palma, Gomera), and cluster A4 included three samples of P. crinitum from the Canary Islands.
Identifying candidate species
We used the same RaxML tree obtained from the phylogenetic analysis to illustrate the delimitation of putative species recognized by the different approaches conducted with the ITS dataset. ABGD, PTP and genetic distance analyses applied to the ITS dataset detected two candidate species that corresponded to the well-supported clades A (including P. crinitum + P. perlatum) and B (including P. perlatum) obtained in the phylogenetic analysis (Fig. 3).
Discussion
Previous molecular phylogenetic studies have proposed the monophyly of Parmotrema crinitum and Parmotrema perlatum (Blanco et al. Reference Blanco, Crespo, Ree and Lumbsch2006; Crespo et al. Reference Crespo, Kauff, Divakar, del Prado, Pérez-Ortega, de Paz G, Ferencova, Blanco, Roca-Valiente and Núñez-Zapata2010). The primary goal of the present investigation was to evaluate the phylogenetic relationship between the two species using ML and Bayesian analyses. Additionally, we used different approaches to assess the species boundary, such as genetic distances based on the threshold of 0.015–0.017 s/s established by Del-Prado et al. (Reference Del-Prado, Cubas, Lumbsch, Divakar, Blanco, de Paz G, Molina and Crespo2010) to measure intra- and interspecific genetic distances in parmelioid lichens, and separate ranges of intra- and interspecific divergence. This threshold was established using both phylogenetically and morphologically well-delimited species (for details see Del-Prado et al. (Reference Del-Prado, Cubas, Lumbsch, Divakar, Blanco, de Paz G, Molina and Crespo2010)). We also used the poisson tree processes (PTP) model (Zhang et al. Reference Zhang, Kapli, Pavlidis and Stamatakis2013) and the automatic barcode gap discovery (ABGD; Puillandre et al. Reference Puillandre, Lambert, Brouillet and Achaz2012). Furthermore, we combined ecological, biogeographical, morphological and chemical data.
Most of the samples belonging to P. crinitum and P. perlatum were recovered in a well-supported monophyletic clade (clade A, Fig. 3), with the exception of two samples of P. perlatum, one from Tenerife (Canary Islands) and one from Lisbon (Portugal), that are separated in the monophyletic clade B, supported by ML and B/MCMC analysis. Morphologically, P. perlatum in clade A has more cilia compared to P. perlatum in clade B.
Thin-layer chromatography (TLC) revealed that the samples in clade B (nos 10 and 12 in Fig. 2) differ chemically from the samples of Parmotrema perlatum grouped in clade A. Several studies have shown that environmental factors, such as light, temperature, pH and culture media, can influence the secondary metabolism in lichens (BeGora & Fahselt Reference BeGora and Fahselt2001). However, we have samples from the same localities and same conditions as the samples in clade A. Furthermore, ABGD, PTP and genetic distance analyses supported clade B as a new sister group to clade A (Fig. 3), suggesting polyphyly, a common phenomenon in Parmeliaceae in general (Lumbsch & Leavitt Reference Lumbsch and Leavitt2011; Leavitt et al. Reference Leavitt, Divakar, Crespo and Lumbsch2016) and Parmotrema in particular (Divakar et al. Reference Divakar, Blanco, Hawksworth and Crespo2005b; Del-Prado et al. Reference Del-Prado, Divakar, Lumbsch and Crespo2016, Reference Del-Prado, Buaruang, Lumbsch, Crespo and Divakar2019; Widhelm et al. Reference Widhelm, Egan, Bertoletti, Asztalos, Kraichak, Leavitt and Lumbsch2016).
Clade A was split into 4 groups; however, the phylogenetic relationships among them remained unresolved (Fig. 3). Previous studies have shown that phylogenetic analyses alone are insufficient to explain phylogenetic relationships within Parmeliaceae. For example, based only on maximum parsimony and Bayesian analyses the sorediate P. sulcata was shown to belong to the same clade as the isidiate P. squarrosa (Molina et al. Reference Molina, Crespo, Blanco, Lumbsch and Hawksworth2004), leaving the authors unable to reach any conclusion regarding species boundaries. Similarly, phylogenetic analyses were insufficient to resolve genetic variability among Parmelia saxatilis specimens; whereas samples from distant geographical regions formed a monophyletic group, samples from neighbouring localities were separated (Crespo et al. Reference Crespo, Molina, Blanco, Schroeter, Sancho and Hawksworth2002; Molina et al. Reference Molina, Del-Prado, Divakar, Sánchez-Mata and Crespo2011b, Reference Molina, Divakar, Goward, Millanes, Lumbsch and Crespo2017). A study on Usnea perpusilla demonstrated that it was necessary to use a combined approach with molecular and morphological data to assess species boundaries in closely related and morphologically variable species (Wirtz et al. Reference Wirtz, Printzen and Lumbsch2008). For this reason, coalescent-based species delimitation analyses have been applied with the goal of explaining relationships among clades and delimiting species boundaries (Parnmen et al. Reference Parnmen, Rangsiruji, Mongkolsuk, Boonpragob, Nutakki and Lumbsch2012; Leavitt et al. Reference Leavitt, Esslinger, Spribille, Divakar and Lumbsch2013).
PTP, ABGD and genetic distance analyses supported clades A and B as putative distinct species. The values obtained were within the interspecific ranges (genetic distances as defined in Del-Prado et al. (Reference Del-Prado, Cubas, Lumbsch, Divakar, Blanco, de Paz G, Molina and Crespo2010)), which would support the two samples of P. perlatum in clade B as a separate species, different from clade A. Focusing on clade A, our species delimitation approaches (PTP, ABGD and genetic distance) supported P. crinitum and P. perlatum as one species. However, based only on one molecular marker (ITS) and the various phylogenetic analyses and species delimitation approaches used, can we consider Parmotrema perlatum conspecific with P. crinitum?
Similar cases have been reported in the Parmotrema perforatum species complex (Widhelm et al. Reference Widhelm, Egan, Bertoletti, Asztalos, Kraichak, Leavitt and Lumbsch2016), and authors have suggested that the phylogenetic relationships between sexual and asexual populations of this species group could be more complex than previously assumed. They also suggest that traditional tools based on reproductive mode and secondary metabolites (Culberson & Culberson Reference Culberson and Culberson1973) are no longer the key to identify species such as P. perforatum (Widhelm et al. Reference Widhelm, Egan, Bertoletti, Asztalos, Kraichak, Leavitt and Lumbsch2016).
A previous study on Umbilicaria (Ott et al. Reference Ott, Brinkmann, Wirtz and Lumbsch2004) focused on Umbilicaria kappenii and U. antarctica, which are distinguished only by their reproductive strategies. Umbilicaria antarctica propagates by thalloconidia and U. kappenii exhibits a variety of asexual propagules: soredia, adventive lobes and sorediate thallyles. To infer phylogenetic relationships between both species, the authors used molecular data from three loci. Results indicated that all samples morphologically referred to U. antarctica and U. kappenii form a monophyletic group and they proposed placing U. kappenii into synonymy with U. antarctica (Ott et al. Reference Ott, Brinkmann, Wirtz and Lumbsch2004). In Parmelia, Molina et al. (Reference Molina, Del-Prado, Divakar, Sánchez-Mata and Crespo2011b) rejected the previous hypothesis that P. sulcata and P. squarrosa form a monophyletic group (Molina et al. Reference Molina, Crespo, Blanco, Lumbsch and Hawksworth2004) and based on phylogenetic analyses and species delimitation approaches confirmed that P. sulcata is not conspecific with P. squarrosa. In addition, P. squarrosa is a reproductively isolated lineage and genetic distances clearly separate this from other Parmelia species (Del-Prado et al. Reference Del-Prado, Cubas, Lumbsch, Divakar, Blanco, de Paz G, Molina and Crespo2010). Within Usnea, studies have suggested that Usnea subfloridana was a secondary species derived from U. florida (Clerc Reference Clerc1984, Reference Clerc1987, Reference Clerc1997; Purvis et al. Reference Purvis, Coppins, Hawksworth, James and Moore1992). Usnea florida and U. subfloridana show many morphological similarities but they have different reproductive strategies. Usnea florida always displays many apothecia and produces no specialized asexual propagules. Usnea subfloridana has soralia, isidiomorphs and occasionally apothecia. Multilocus phylogenetic analyses based on sequences of the ITS, IGS and LSU regions of the nuclear ribosomal DNA and the gene coding for β-tubulin, Mcm7, RPB1 and RPB2 showed that specimens of the two morphospecies formed one monophyletic intermixed group, and not two groups corresponding to morphology (Articus et al. Reference Articus, Mattsson, Tibell, Grube and Wedin2002; Mark et al. Reference Mark, Cornejo, Keller, Flück and Scheidegger2016). This topology was further confirmed with a coalescent-based species delimitation approach (Mark et al. Reference Mark, Cornejo, Keller, Flück and Scheidegger2016). Authors have suggested that they could be conspecific but taxonomic conclusions must await further study. Moreover, a recent study using RADseq data suggests that closely related lichen species may need genome-wide data to test their species boundaries (Grewe et al. Reference Grewe, Lagostina, Wu, Printzen and Lumbsch2018).
Traditionally it was thought that asexually reproducing species in lichens, and in filamentous fungi in general, was an evolutionary dead end (Normark et al. Reference Normark, Judson and Moran2003). However, recent molecular studies have demonstrated that lineages with vegetative propagation can also present high genetic diversity (e.g. Parmelia sulcata (Molina et al. Reference Molina, Del-Prado, Divakar, Sánchez-Mata and Crespo2011b), Parmotrema reticulatum (Del-Prado et al. Reference Del-Prado, Divakar, Lumbsch and Crespo2016)). However, even in the absence of sexual reproduction, lichens can exchange photobionts and this process could provide opportunities for gene transfer (Piercey-Normore Reference Piercey-Normore2006; Nelsen & Gargas Reference Nelsen and Gargas2008).
While our study confirmed the monophyly of an intermixed clade of Parmotrema crinitum and P. perlatum, the taxonomic conclusion must await additional studies including more markers. The phylogenetic tree lacked ML or B/MCMC support for other widely accepted Parmotrema species, such as P. pseudoreticulatum and P. reticulatum. A more comprehensive taxon sampling and additional molecular markers will therefore be needed before making a formal taxonomic decision on the status of clade B.
Acknowledgements
This study was supported by the Spanish Ministerio de Ciencia e Innovacion (PID2019-105312GB-I00) and the Santander-Universidad Complutense de Madrid (PR87/19-22637 and G/6400100/3000).
Author ORCIDs
Ayoub Stelate, 0000-0001-8929-3046; David Alors, 0000-0001-7288-9521; Pradeep Kumar Divakar, 0000-0002-0300-0124; Ana Crespo, 0000-0002-5271-0157.
Supplementary Material
To view Supplementary Material for this article, please visit https://doi.org/10.1017/S0024282922000147.
Open access
