Insights into the Migration Routes and Historical Dispersion of Species Surviving the Messinian Crisis: The Case of Patella ulyssiponensis and Epizoic Rhodolith Lithophyllum hibernicum

The genus Patella (Patellogastropoda, Mollusca) is represented by a group of species exclusive to the Northeast Atlantic Ocean (including Macaronesian archipelagos) and Mediterranean Sea. The species Patella ulyssiponensis and Patella aspera are common in European waters, with the first inhabiting continental coast, and the second endemic to Macaronesian archipelagos. However, the acceptance of these two lineages as separate species is still highly debated. The red coralline species algae Lithophyllum hibernicum, distributed from Northeast Atlantic to the Mediterranean, is usually found as epilithic crusts or unattached forms (named rhodolith beds), although it also forms epizoic crusts on other organisms, e.g., shell surfaces. In order to study the historic dispersal and migration routes of the Patella ulyssiponensis-aspera complex, taxonomic, genetic and biogeographic approaches were employed based on haplotype network analyses and estimations for the most common recent ancestor (TMRCA), using Cytochrome Oxydase I. A synonymy for these two species is proposed, with the presence of a shared haplotype between the continental (P. ulyssiponensis) and insular (P. aspera) lineages, and with basis of morphological and nomenclatural data. We propose an evolutionary scenario for its dispersal based on a high haplotype diversity for the Mediterranean regions, indicating its possible survival during the Messinian Salinity Crisis (6–5.3 Mya), followed by a colonization of the Proto-Macaronesian archipelagos. The epizoic association of L. hibernicum on P. ulyssiponensis shell adult surface is recorded in this study, likewise the promotion of settlement conditions provided by these coralline algae to P. ulyssiponensis larvae, may explain the reach of P. ulyssiponensis distribution through rhodolith transportation.


Introduction
Patella Linnaeus, 1758 is a gastropod genus with 15 [1]. This group is characterized by "the coiling of the radula below the visceral mass and gonad" and "lack of pluricuspid cusp 1" [2]. Joseph Christiaens performed the first revision study of genus Patella, proposing the use of P. ulyssiponensis aspera for the Atlantic islands [3]. The last suggestion currently accepted was recommended by Weber and Hawkins, indicating P. ulyssiponensis for European coast shore, and P. aspera for Macaronesian archipelagos [4]. The P. ulyssiponensis-aspera

Data Collection
Samples were collected using spatulas (to remove the limpets from the rocks, Figure 2) and stored in thermal boxes with ice. In the laboratory, photographs of the organisms were taken using a Nikon Digital Sight D5-L1 camera and a stereomicroscope Nikon SMZ800 (Figures 3-7). Specimens were deposited at the Natural History Museum of the Iberian Peninsula (NatMIP-"Museu de História Natural da Península Ibérica"), Aquamuseu do Rio Minho, Vila Nova de Cerveira, Portugal.

Data Collection
Samples were collected using spatulas (to remove the limpets from the rocks, Figure 2) and stored in thermal boxes with ice. In the laboratory, photographs of the organisms were taken using a Nikon Digital Sight D5-L1 camera and a stereomicroscope Nikon SMZ800 (Figures 3-7). Specimens were deposited at the Natural History Museum of the Iberian Peninsula (NatMIP-"Museu de História Natural da Península Ibérica"), Aquamuseu do Rio Minho, Vila Nova de Cerveira, Portugal.

Data Analyses
DNA samples were extracted from the foot of P. ulyssiponensis using an E.Z.N.A. Mollusc DNA Kit (Omega Bio-tek). Amplification of the mitochondrial gene cytochrome c oxidase I (COI-5P), was performed using the primer pairs LoboF1 and LoboR1 [26]. PCR    COI-5P sequences were edited and aligned using MEGA software [27] and verified for the presence of stop codons, insertions or deletions. Homology searches were performed in the BOLD database [28] and with BLAST [29] in the GenBank database [30]. Phylogenetic analyses for the Patella ulyssiponensis-aspera complex were carried using a dataset composed by the single sequence obtained in this study (OK465193 GenBank accession number) and 271 public sequences of Patella ulyssiponensis and Patella aspera (Table S1). The software JModelTest [31] was used to find the best-fit model for nucleotide substitutions, being this model Hasegawa-Kishino-Yano (HKY) with Gamma distributions. Maximum Likelihood (ML) tree was constructed with MEGA software using HKY model with Gamma distribution. Bootstraps were performed using 10,000 replicates. Bayesian Inference analyses were performed with Mr. Bayes software [32] using the General Time Reversible substitution model (GTR) with Gamma distribution, running for 35 million generations sampled every 1000 generations, using two chains with a temperature setting of 0.1. After confirming that the average standard deviation of the split frequencies was <0.01, the first 25% of the trees were discarded as burn-in and the remaining were used to construct a majorityrule consensus tree. Pairwise distances were calculated using the Kimura 2-parameter model. Haplotype genealogy analyses was performed using TCS software [33] with a 95% statistical parsimony connection limit using 226 P. ulyssiponensis and P. aspera sequences, with 46 sequences being removed from the original dataset due to lack of precise geographical location data. The final dataset, used only for haplotype genealogy analyses, was trimmed to a final size of 376 base pairs due to lack of resolution on the opening and closing segments of some sequences. Genetic differentiation among populations (European Atlantic continental, Mediterranean, Morocco, Azores and Madeira plus Canary sub-group) was calculated using an analysis of molecular variation (AMOVA) in GenAlEx [34]. F ST Haplotype (Hd) and nucleotide (π) diversity and Tajima's D and Fu's FS neutrality tests were performed using Arlequin 3.5 [35] for the populations mentioned above and for each Mediterranean population (Alboran, Balearic, Tyrrhenian, Ionian and Aegean).
The analyses between Atlantic Patella species were carried out by the construction of a Maximum Likelihood tree using a HKY + G substitution model, conducted with 10,000 bootstraps on MEGA software, from a dataset composed by 418 sequences of P. candei, P. lugubris, P. caerulea, P. pellucida, P. depressa, P. ferruginea, P. vulgata, P. rustica, P. ulyssiponensis and P. aspera. Public sequences used for tree construction available on Table S1.
The software BEAST v.1.7 [36] was used to estimate the most recent common ancestor (TMRCA) from a dataset consisting of one sequence of Eoacmaea profunda as the basal group, one sequence from one species for the genus Helcion, Cymbula and Scutellastra and one sequence of each species of Patella, except for Patella ulyssiponensis in which were used one sequence representing the continental basal haplotype, one sequence for the Azores archipelagos basal haplotype and one sequence from the most represented shared haplotype from Madeira and Canary archipelagos. The choice of these sequences was made in order to use the Yule speciation model obtaining a bifurcating pattern, chosen due to the usage of extant only lineages. The evolutionary model chosen was the HKY model with a proportion of invariable sites and four Gamma categories (HKY + I + G) chosen by JModelTest. The fossil calibration was performed by rooting the age of the clade containing the sister group of Patella using sequences of Scutellastra exusta, Helcion concolor and Cymbula safiana. This calibration was done through the oldest known fossil record of this paraphyletic group, "Patella" soyaensis Kase and Shigeta, 1996 from the Upper Cretaceous (75 Mya) of northern Japan, assigned to Scutellastra by Ridgway et al. [2]. A relaxed clock model was used with a lognormal distribution and a substitution rate of 1.2 ± 0.1% (s.d.) per million years, based on values for other Molluscs [37]. The analyses were performed by 20 million chains, sampling each 1000 generations with the initial 25% being discarded as burn-in. The software TRACER 1.5 [38] was used to check the effective size sample for each parameter (at least higher than 100), as well as for checking the confidence intervals for each TMRCA obtained from the 95% highest posterior densities interval. The software FigTree 1.4 was used to obtain the tree visualization.

Geological Maps
Ancient satellite Earth images (i.e., Upper Cretaceous, Lower Miocene, Messinian stage, and Lower Pliocene) were obtained by the GPlates software, version 2.2.0 [39]. Information of ocean currents old periods was obtained through Britannica website [40].

Biogeography Patterns and Genetic Analyses
Trees obtained with Maximum Likelihood and Bayesian Inference shared a similar topology with a separation of two different lineages, one from the Northeast Atlantic and Mediterranean continental shelf and another from the Macaronesian archipelagos with some Mediterranean populations ( Figure 8 and Figure S1). On ML tree, bootstrap support values for the insular and continental clades are <70 (35 for the continental and 65 for the insular clades) ( Figure S1); nonetheless, Bayesian support values of these same clades are of 100 ( Figure S2). Within lineages, genetic distance was 0.5 ± 0.1% (s.e.) for the continental lineage and 0.7 ± 0.2% (s.e.) for the insular lineage. The distance between the two lineages was 3.5 ± 0.9% (s.e.). The specimen sequenced nested within the continental grade, being most similar to a specimen sampled in Asturias (North Spain). Interspecific distances range between 7% and 16%, with exception of the case between P. lugubris and P. candei which had an interspecific distance of 4% ( Figure 9). Haplotype genealogy revealed 85 different haplotypes, 15 unique for the Macaronesian region, 69 for the continental region and one shared between the Mediterranean, Madeira, Tenerife and Fuerteventura. Of those 69 haplotypes, five are exclusive to the Iberian Atlantic, six shared between Atlantic and Mediterranean, two exclusives to Morocco, plus one shared with Atlantic and Mediterranean populations, with the remaining being exclusive to the Mediterranean. Of the last ones, 36 are unique to Italy, four to Greece, two to Portugal, seven in the Spanish Mediterranean region, three for the Spanish region of Asturias and four for the French Mediterranean ( Figures 10 and 11, Table S2). According to AMOVA, 60% of the molecular variance was found between populations and 40% within populations.
Haplotype (Hd) and nucleotide (π) diversity and Tajima's D and Fu's FS neutrality tests for all regions are available in Table 1, and values for the Mediterranean regions are available in Table 2. F ST values show little differentiation between Atlantic and Mediterranean populations (F ST = 0.04) and high differentiation between Continental and Insular populations ( Table 3). Analyses of each Mediterranean population F ST values shows little differentiation, except for the Alboran Sea populations, which seem to have a moderate differentiation in structure comparing with the remaining Mediterranean regions (Table 4). Populations from the Balearic, Tyrrhenian and Ionian Seas seem to have gone through a recent demographic expansion, while the Alboran and Aegean Seas populations seem to be demographically stable ( Table 2). Estimated TMRCA values and corresponding confidence intervals are available on Table 5.    [46]. Type locality. "Ulyssiponem" [41], actual Lisbon, Portugal [3]. Diagnosis [5,84,119]: Top of shell (apex) located anteriorly to the middle ( Figure 3A); radiating ridges varying at single and triple ones ( Figure 3B). Tubercles arranged over the entire surface of the shell, in rows up to the apex ( Figure 3A). Black strips lines on surface from apex to basis. Border shell uniformly crenulate ( Figure 3A,B). Posterior margin (with lateral expansions) larger than anterior one ( Figure 3A,B). Internal shell region basally whitish with black bands (thinner in both ends), medially to upper yellowgreenish, to orangey apically ( Figure 3B). Radula with three unequal teeth: left tooth, with two protrusions (being one like a vestigial fourth tooth), and basally larger than central one. Yellowish to dark grey anterior cephalic tentacles basally inflated, anteriorly tapered ( Figure 4A,B). Body surface with a dark grey colour pattern ( Figure 3C). Pallial tentacles of mantle skirt whitish to translucent, organized in two series of distinct sizes ( Figure 3D), may contain adhering eggs ( Figure 4C). Uniformly yellow to orange foot with cream coloured edge ( Figures 4D and 6C).
Features in juvenile specimens. Juvenile individuals present slight distinctive morphological characteristics in relation to adult limpets: in individuals smaller than 3 mm, width of coloured and uncoloured strips are a uniform standard; in larger ones, dark ridges are broader, projecting somewhat beyond the shell perimeter. This change is related to the emergence and growth of the tentacles of the mantle (for details see [17]). These juvenile limpets may be found in pools, on red crustose algae ( Figure 5A,C), close to Corallina sp. and Mytilus sp. beds, as well as on adult shells surface. Young specimens have greyish feet, and as they grow, feed on red algae, then their feet start to take on an orange colour, typical of adults. Finally, mantle tentacles have a calcareous-type constitution.
Diagnosis [128,132]: Fruticose or spherical rhodoliths with up to 5 cm of diameter densely branched with cup-shaped tips up to 6 mm in diameter, or forming smooth crusts 1-5 mm thick, becoming irregular at the junction of the edges with contiguous thalli ( Figure 5A   Geographical distribution. Northeast Atlantic: England, Ireland, France and Iberian Peninsula; Mediterranean: Spain and France [128].
Ecological notes. Found in epilithic crust on the mid to low intertidal zone in rock pools, bedrock and on tidal canals, epizoic on mussels and other bivalves, such as Mytilus sp., and as rhodolith on intertidal or subtidal zones up to 15 m on the Atlantic Ocean. On the Mediterranean Sea, was only found subtidally as rhodolith and as epiphytic crust on other rhodoliths like Lithophyllum incrustans [128].

Association between P. ulyssiponensis and L. hibernicum
The association was verified by the presence of rhodolith L. hibernicum on shell surface of limpet P. ulyssiponensis populations ( Figure 6). The presence of this red calcareous alga was verified on all adult limpets, and one juvenile with 11.0 mm length from Belinho Beach (Esposende). It is also noted that juvenile limpets settle on the rhodolith surface [16], in a relationship similar to inquilinism or commensalism, and as the gastropods grow they feed on these calcareous algae, acquiring a more orange colouration. At the same time, rhodoliths begin to grow over the limpets, and may cover them completely behaving as epizoic. In speculative way, when limpets move to other environments, it is possible that transporting rhodoliths on your shell surface, may promote their dispersion.

Genetic Relationships of Two Patella ulyssiponensis Lineages
Analyses of sequence data from the P. ulyssiponensis-aspera complex revealed strong genetic differentiation between Atlantic continental and insular populations, however insular populations are remarkably similar to some Italian populations (Figure 8). AMOVA analyses showing a 60% molecular variance between populations indicates somewhat structured populations. With 15 unique haplotypes for the Macaronesian regions, nine of them being exclusive to the Azores archipelagos, three exclusives to Madeira, three exclusives to Tenerife and one haplotype shared between Madeira, Tenerife and Fuerteventura and Bonifacio, France (Figures 10 and 11, Table S2).
Insular haplotypes are divided in two groups, one formed between Azores archipelago and Tenerife, and the other formed between Tenerife, Fuerteventura and Madeira islands. The large number of exclusive haplotypes with high diversity present in Azores might indicate this archipelago as occurring during a possible first colonization, via Tenerife for continental P. ulyssiponensis. Madeira and Canary have smaller haplotype diversity and share a haplotype with the continent, indicating recent gene flow (or intermittent gene flow) or low times of genetic drift. Azores and the Canary/Madeira group present a starlike haplotype network with a high ratio of singletons, indicating population expansion from a small number of founders following a bottleneck [139]. Insular haplotypes are more related to Italian haplotypes from Lo Strangolato with one minimum nucleotide mutation ( Figure 10) between populations, although unsampled haplotypes may change the genealogy topology and solve the specific route of migration from continental to insular habitats. Sá-Pinto et al. suggested a differentiation between continental and insular lineages between 8.3 and 3.9 million years, with a basal haplotype for the Macaronesian Islands in Azores [140]. Nevertheless, our analysis with a larger dataset puts Tenerife as the basal haplotype, maintaining Azores the most differentiated group (Figure 10).
Although Mediterranean and Atlantic populations presented low differentiation indicating an ongoing genetic flow between these geographic areas, Mediterranean populations present a higher number of exclusive haplotypes, which might indicate the presence of barriers limiting gene flow, also suggested in [141], who pinpointed these barriers at the Almeria-Oran front and at South Italy, with the barrier at Almeria-Oran being supported by our own F ST calculations.
In the region of Agadir, Morocco (Atlantic) there are two exclusive haplotypes and one other haplotype shared between Portuguese, French Atlantic, Spanish populations, indicating a possible colonization from Iberian Peninsula (Figure 10). With the richest hap-lotype diversity on South Europe compared to northern Atlantic populations, the Western Mediterranean presents itself as strong candidate to this species ancestry. Although the differentiation between Atlantic continental and insular populations was found corroborating the analyses from Weber and Hawkins, and Sá-Pinto et al. [4,142], interspecific distances between Patella species range from 7 to 16% (with the exception of P. lugubris with P. candei) (Figure 9), and compared to the 3.5% between the continental and insular lineages of the P. ulyssiponensis-aspera complex, and the fact that Mediterranean populations are genetically very close to insular populations, even sharing a haplotype indicates a case of ongoing speciation process, rather than being two different species. This speciation process is also corroborated by a large negative value of Tajima's D and Fu's FS for the Macaronesian regions, indicating a recent expansion due to an excess of singletons or recent mutations. Although Tajima's D for Madeira and Canary is non-significative Fu's FS is, being this contradiction probably resolved with a higher sequence sampling for this region. Similarly, Tajima's D and Fu's FS show non-significative contradictory results for Morocco, which might indicate a recent bottleneck or may have resulted from small sample size and small collection range, being restrained to only the Agadir region.

Biogeographic Patterns
Nakano and Ozawa, based on the extant placement of Eoacmaea Nakano & Ozawa, 2007 at the base of Patellogastropoda, inferred by molecular data and worldwide distribution of its recent species, suggested that this group was widespread in the Tethys Sea from the Lower Jurassic to Middle Cretaceous, with the diversification of Patellogastropoda occurring between the Jurassic and Cretaceous periods [88]. The oldest known fossil of Patellidae is Patella costulata Münster, 1869, found in the St. Cassian Formation in northern Italy from the Triassic (210-230 Mya), however generic assignment was impossible due to poor preservation [2,143].
Due to the lack of fossil registers, the origin of the genus Patella was previously estimated, trough molecular divergence, to occur between the Upper Triassic and Upper Cretaceous, with a divergence from the rest of Patellidae estimated between 121-231 Mya [88] and between 69-169 Mya [144] with our estimations being similar to Koufopanou et al. results [144], ranging from 64.7-143 Mya (Figures 12 and 13A). Based on this time of origin and the distribution of its recent members, Patella spp. may have been the first Patellogastropoda to diversify in the Atlantic following its opening in the Jurassic, and may have been present in the northern Atlantic since the Mesozoic [88]. Following the closure of the Tethys Sea during the late Miocene, Patella may have been isolated from its sister group Scutellastra, suggesting a vicariant speciation [2] and explaining the mostly anti-tropical distribution of Patellidae ( Figure 13B). . Circle sizes are proportional to haplotype frequency for each locality. Grey circles represent haplotypes unique for that region, while coloured circles represent haplotypes shared between different sampling locations. Graphs constructed through rawgraphs.io.
Haplotype (Hd) and nucleotide (π) diversity and Tajima's D and Fu's FS neutrality tests for all regions are available in Table 1, and values for the Mediterranean regions are available in Table 2. FST values show little differentiation between Atlantic and Mediterranean populations (FST = 0.04) and high differentiation between Continental and Insular populations ( Table 3). Analyses of each Mediterranean population FST values shows little differentiation, except for the Alboran Sea populations, which seem to have a moderate differentiation in structure comparing with the remaining Mediterranean regions (Table 4). Populations from the Balearic, Tyrrhenian and Ionian Seas seem to have gone through a recent demographic expansion, while the Alboran and Aegean Seas populations seem to be demographically stable ( Table 2). Estimated TMRCA values and corresponding confidence intervals are available on Table 5. of ongoing speciation process, rather than being two different species. This speciation process is also corroborated by a large negative value of Tajima's D and Fu's FS for the Macaronesian regions, indicating a recent expansion due to an excess of singletons or recent mutations. Although Tajima's D for Madeira and Canary is non-significative Fu's FS is, being this contradiction probably resolved with a higher sequence sampling for this region. Similarly, Tajima's D and Fu's FS show non-significative contradictory results for Morocco, which might indicate a recent bottleneck or may have resulted from small sample size and small collection range, being restrained to only the Agadir region.

Biogeographic Patterns
Nakano and Ozawa, based on the extant placement of Eoacmaea Nakano & Ozawa, 2007 at the base of Patellogastropoda, inferred by molecular data and worldwide distribution of its recent species, suggested that this group was widespread in the Tethys Sea from the Lower Jurassic to Middle Cretaceous, with the diversification of Patellogastropoda occurring between the Jurassic and Cretaceous periods [88]. The oldest known fossil of Patellidae is Patella costulata Münster, 1869, found in the St. Cassian Formation in northern Italy from the Triassic (210-230 Mya), however generic assignment was impossible due to poor preservation [2,143].
Due to the lack of fossil registers, the origin of the genus Patella was previously estimated, trough molecular divergence, to occur between the Upper Triassic and Upper Cretaceous, with a divergence from the rest of Patellidae estimated between 121-231 Mya [88] and between 69-169 Mya [144] with our estimations being similar to Koufopanou et al. results [144], ranging from 64.7-143 Mya (Figures 12 and 13A). Based on this time of origin and the distribution of its recent members, Patella spp. may have been the first Patellogastropoda to diversify in the Atlantic following its opening in the Jurassic, and may have been present in the northern Atlantic since the Mesozoic [88]. Following the closure of the Tethys Sea during the late Miocene, Patella may have been isolated from its sister group Scutellastra, suggesting a vicariant speciation [2] and explaining the mostly anti-tropical distribution of Patellidae ( Figure 13B).    Figure 13C) and the Plio-Pleistocene glaciations. Crame noted a bipolar pattern among marine molluscs, suggesting three major phases of diversification: Jurassic-Cretaceous, Upper Palaeogene-Lower Neogene and Plio-Pleistocene-Recent [145].
The fossil record of the genus Patella is sparse, with only nine valid species according to WoRMS [1] (Table 6), although the full list of extinct species under Patella needs a revision due to the habit of assigning Patellidae taxa under Patella until the middle of the 19th century. From those nine extinct species, eight are from the Atlantic region, with three restricted to the Canary archipelagos and one from the Mediterranean, with extant dates ranging from the Oligocene to the Pleistocene. The disappearance of most of those species during the Miocene-Pliocene boundary might shed some light on the radiation of the extant species of Patella and its connection with the climactic change occurring within the Mediterranean and Macaronesian regions during this period.   Figure 13C) and the Plio-Pleistocene glaciations. Crame noted a bipolar pattern among marine molluscs, suggesting three major phases of diversification: Jurassic-Cretaceous, Upper Palaeogene-Lower Neogene and Plio-Pleistocene-Recent [145].
The fossil record of the genus Patella is sparse, with only nine valid species according to WoRMS [1] (Table 6), although the full list of extinct species under Patella needs a revision due to the habit of assigning Patellidae taxa under Patella until the middle of the 19th century. From those nine extinct species, eight are from the Atlantic region, with three restricted to the Canary archipelagos and one from the Mediterranean, with extant dates ranging from the Oligocene to the Pleistocene. The disappearance of most of those species during the Miocene-Pliocene boundary might shed some light on the radiation of the extant species of Patella and its connection with the climactic change occurring within the Mediterranean and Macaronesian regions during this period. Regarding the known fossils of P. ulyssiponensis, the oldest one was found in the Contrada Case Alte region in Italy, dating from the Lower Pleistocene (2.6-0.8 Mya) [86], with two more from the Upper Pleistocene from Santa Maria, Azores [48,50] and from the Canary islands Gran Canaria and Lanzarote [105].
Assuming the Jurassic-Cretaceous vicariance scenario as the origin of Patella, with the separation from its sister taxa (Scutellastra) through the breakup of Pangea and the Atlantic Ocean formation, followed by the closure of the Tethys seaway and successive closure of the connection between Indian and Atlantic Oceans during the Lower Miocene, the emergence of the Proto-Macaronesia archipelagos, the closure (before) and reopening (after this period) of the Mediterranean Sea during the Messinian Salinity Crisis and the closure of the passage on Isthmus of Panama (Table 7) as the radiation events for the genus Patella we propose the following scenarios for the origin and biogeography of Patella ulyssiponensis:

2.
Vicariant speciation event during the breakup of the supercontinent of Pangea, originating two lineages that eventually led to the radiation of the genus Patella on the Northeast Atlantic Ocean and Scutellastra, Helcion and Cymbula on Southwest Africa ( Figure 13B).

3.
Radiation of the recent species of the genus Patella followed by the closure of the Tethys Sea with Patella ulyssiponensis appearing somewhere between the Upper Oligocene and the Lower Miocene (6.3-33.2 Mya) ( Figure 13B). 4.
Since from this event period, three hypotheses are proposed:  Figure 13C) followed by a rapid population and genetic expansion in the Mediterranean after the Zanclean Flood. This hypothesis is mainly proposed due to the lack of fossil evidence and lack of haplotype sampling for some locations which could resolve the non-observed haplotypes in the haplotype network tree. Among the presented hypothesis we agree with hypothesis A, due to the presence of many exclusive haplotypes (including the most frequent one) in the Mediterranean, contrasting with the low haplotype diversity for the Northeast Atlantic, and the divergence date between continental and archipelago lineages occurring after the Zanclean Flood, according to TMRCA calculations

5.
Interruption of the gene flow between the populations from Azores and the populations from Madeira/Canary, possibly due to the global rise of sea level, consequence of the interglacial events during the Plio-Pleistocene (3.15-0.1 Mya), submerging the Proto-Madeira and Proto-Canary archipelagos. The disappearance of these islands may have cut the pathways of dispersal trough the entirety of the Macaronesian archipelagos ( Figure 13D).
As for the colonization of the Macaronesian archipelagos by P. ulyssiponensis, two possible pathways of dispersion are proposed:

1.
Transport by association with rhodoliths-The southeastern coasts of the Macaronesian archipelagos show abundant tempestite deposits that incorporate fossil rhodoliths, belonging to the Mid Miocene to Lower Pliocene [153], transported from eastern locations by storms [154]. Contrasting with the relative scarcity of rhodoliths in Madeira and Azores in present day [153], the Plio-Pleistocene might have rebalanced the circulation patterns, interrupting the gene flow between islands. In the present day, rhodolith debris is still found abundantly in Fuerteventura and seldomly found in Selvagens [153], which may explain the haplotypes shared between Madeira, Fuerteventura, Tenerife and Bonifacio (France).

2.
Larval dispersal-Although no proof of long-range larval dispersal for Patella ulyssiponensis has been demonstrated, this pathway of colonization cannot be discarded. Contrary to other species of Patella, larvae of Patella ulyssiponensis show an optimal development at higher temperatures (8-24 • C) and a higher larval longevity combined with a higher longevity in substratum absence [155]. All these aspects, associated with the presence of pathways via the Proto-Macaronesian archipelagos, may have an impact on long range larval dispersal during the warmer Upper Miocene ( Figure 13C).

Historical Context and Morphology: Relationships between Patella ulyssiponensis Gmelin, 1791 and Patella aspera Röding, 1798
Limpet gastropods of genus Patella were first erected by Johann Friedrich Gmelin in "Systema Naturae" [41], designating 237 species, including P. ulyssiponensis, that was probably named for its sampling site, "Ulyssiponem", actual Lisbon, Portugal. This species was erected based on drawing in "La conchyliologie" (pl. 2, Figure G2, [162]). After this, Peter Friedrich Röding established the species P. aspera likewise according to the same image of the "La conchyliologie" book [60], as well as using this same illustration and specific designation methodology, Jean-Baptiste de Lamarck named P. aspera [61], which nowadays is considered as a junior homonym of P. aspera Röding, 1798 [163]. Observing only this historical contextualization, it may be considered that both taxa belong to the same species, with the older authority that should be the one used for the specific designation, namely P. ulyssiponensis [45], which is in line with the "Principle of Priority"-Article 23 of International Code of Zoological Nomenclature [164,165].
Synonymy of P. ulyssiponensis and P. aspera was primarily suggested by Joseph Christiaens [3]. This author indicated that the 1780 s work referred inaccurately to Patagonia (South America) for limpet sampling (pl. 2, Fig. G2 [162]), which may have instead been collected in the Azores, Portugal [3]. Christiaens suggested the use of P. ulyssiponensis aspera for the subspecies of P. ulyssiponensis present on the Atlantic islands [3]. In the Checklist of British Marine Mollusca, the authors also considered P. aspera as a synonym of P. ulyssiponensis [166]. Sella, Robotti and Biglione also assumed that both species are the same in a study about genetic divergence among Mediterranean Patella spp. [167]. Posteriorly, Titselaar also considered P. aspera as a synonym of P. ulyssiponensis, and designated a neotype for a specimen of P. ulyssiponensis collected in Estoril, Portugal [46], which is 24 km away from the type locality (Lisbon). In a morphometric analysis of northwest Portuguese limpets, Cabral and Silva mentioned Röding and Lamarck's P. aspera as a synonym of P. ulyssiponensis [52]; Mauro et al. also considered both species as the same [168]. In relation to rhodolith-forming algae species, in 1897, Mikael Heggelund Foslie erected a variety for species Lithothamnion fasciculatum (Lamarck) Areschoug, 1852, identified as L. fasciculatum f. subtilis, with the remarkable characteristic of ramifications shorter than other varieties established in the same work (i.e., L. fasciculatum f. incrassata, f. gyrosa, and f. dilatata), and narrowed (1 to 1.5 mm width) compared to L. fasciculatum f. incrassata [137]. In 1906, the same author revised this identification and replaced the genus with Lithophyllum, re-naming for Lithophyllum hibernicum [120]. This last designation was, therefore, established and used in other subsequent studies (see Section 3.2 Taxonomic Account -Lithophyllum hibernicum).

Limpets and Rhodoliths
The existence of epizoic algae on the hard shell of molluscs is a small, sparse and poorly understood research topic; thus, the ecological importance of this relationship is sometimes overlooked [169]. Calcareous red algae are able to grow on stone shores (epilithic forms), as well as on biological structures, i.e., epiphytic types growing on other algae, and epizoic carbonate algae on alive and dead animals: on shell fragments, gastropods, bivalves, barnacles [170], sponges, corals, or sea urchins [171,172]. In relation to genus Lithophyllum, associations between Lithophyllum incrustans Philippi, 1837 and a snail shell [173], or between Lithophyllum hibernicum Foslie, 1906 and the bivalve Mytilus galloprovincialis Lamarck, 1819 [174] have already been observed ( Figure 2B,C). The association between P. ulyssiponensis and L. hibernicum is described for the first time. An ecological relationship between the two species, benefiting the limpet at a certain stage of its life cycle (settlement) and even as camouflage to protect it from predation should be further investigated.
Haplotype's network of this limpet indicates that the distribution and historical dispersion originated from Mediterranean Sea, until arrival to European Atlantic and Macaronesian islands (Figures 10 and 11). For L. hibernicum, the actual distribution ranges from the Mediterranean Sea to European Atlantic, but there are no records for the Atlantic islands. Nonetheless, fossil registers may help to understand the rhodolith-limpet relationship. For example, Johnson et al. identified many Pleistocene stratigraphic profiles containing rhodoliths (including Lithophyllum sp.) at the Cape Verde Islands [175]. They also verified the presence of Patella sp. (and other molluscs) in same strata of rhodoliths [175], similarly to the first work performed in the same study area by Charles Darwin, that found many shells fossils, including limpet specimens, with rhodoliths formation [176].

Conclusions
According to morphological and genetic analyses, both P. ulyssiponensis and P. aspera are the same species, existing in continental and insular shore environments, influenced by the Northeast Atlantic Ocean and Mediterranean Sea, and probably surviving the Messinian Salinity Crisis event. The current distribution and ecological release may be associated with the Atlantic Ocean tides after the opening of the Strait of Gibraltar (end of Messinian Crisis). In this context, among the presented hypothesis, we agree with the first hypothesis, due to the many exclusive haplotypes (including the most frequent one) for the Mediterranean, contrasting with the low haplotype diversity for the Northeast Atlantic, and their divergence date between continental and archipelago lineages occurring after the Zanclean Flood, according to TMRCA calculations. Therefore, it is very possible that limpets survived the Messinian event, following susceptible dispersions on Atlantic islands.
Rhodolith-forming algae L. hibernicum may influence the settlement and dispersal of limpet larvae, but on the other hand possesses an epizoic relationship on the adult shell surface of P. ulyssiponensis. This probable association can currently be considered as a profitable case of co-evolution of both species.
The fossil record of Patellogastropoda is poorly known, due to the rarity of their fossils, mostly because of their occurrence on wave exposed rocky shores in which conditions for preservation are limited. This scarce knowledge is also compounded by a lack of shell microstructure examination, causing genera and family misplacements. Consequently, the divergence time estimated from genetic distances and calibrations based on a small number of fossil records can result in its underestimation or overestimation. Concerning the origin and lineages divergence times of P. ulyssiponensis, an expansion of sampling locations might provide some insights and solve some patterns, especially by having a higher sampling effort in the Iberian Atlantic, Macaronesia archipelagos and Moroccan coast, and by the inclusion of sequences from the African Mediterranean, Black Sea and Northern Atlantic Regions, combined with a compilation of datasets based on multiple genes.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/hydrobiology1010003/s1, Figure S1: Maximum Likelihood tree obtained for COI-5P sequences of Patella ulyssiponensis, using HKY + G substitution model. Value of nodes corresponds to bootstrap support. Figure S2: Bayesian consensus tree of COI-5P for Patella ulyssiponensis. Node values correspond to Bayesian posterior probabilities. Table S1: List of COI-5P sequences for Patellogastropoda with GenBank accession numbers.