Biological Interaction as a Possible Ultimate Driver in the Local Extinction of Cedrus atlantica in the Iberian Peninsula

: The presence of Cedrus atlantica on the European continent, including, especially, the deter-mination of the time of its disappearance from the Iberian Peninsula, is one of the most controversial issues in recent decades regarding the successive extinction of conifers in the Western Mediterranean. This work propounds the possibility that C. atlantica and Pinus nigra could have co-habited in the past, mutually excluding each other in the areas with suitable conditions for both species, where, ultimately, the one that was the most competitive would have remained. The niche overlap in the two-dimensional ecological space was analyzed. In addition, the potential distribution of both species in the Western Mediterranean today and two past periods (Last Glacial Maximum and Mid-Holocene) was modeled to identify their common geographic area of distribution. The species showed very well differentiated niches and a distribution of their habitats virtually segregated by continents since the Mid-Holocene ( P. nigra in Europe and C. atlantica in Africa), which responds to differences in climatic afﬁnities. However, the contact of the bordering areas of their distributions in the Baetic mountain range suggests that C. atlantica could have maintained its presence in the Iberian Peninsula until recent times. P. nigra would have displace it in later stages due to its greater prevalence on the continent, so it would have had greater opportunities to occupy the available space.


Introduction
The presence of Cedrus atlantica (Endl.) Carrière on the European continent and the time of its disappearance from the Iberian Peninsula are issues that have raised controversy in recent decades regarding the successive extinction of conifers in the Western Mediterranean [1][2][3]. The presence of Cedrus in the Iberian Peninsula has been referenced from the Eocene to the Holocene [4] and widely in Europe since the Miocene [5]. However, the only fossil macro-remain that confirms its autochthonous origin corresponds to leaves of C. blombley found on the island of Rhodes dated to the Middle Pleistocene (500 ka) [6]. The absence of macro-remains (charred wood) of Cedrus on the European continent led to the statement that the presence of cedar pollen in the Iberian fossil record responded to its airborne uptake by winds from Africa [2,5,7,8]. However, the constancy of the pollen type in numerous deposits [1,3,7,[9][10][11][12], its proportion and the low dispersal power of grains [1,13,14]; the Cedrus-Abies assemblage in the records [10]; and the low correlation between the presence of pollen and the origin of the winds from Africa [15,16] question this possibility and support the hypothesis that Cedrus lived in southern Europe during the Holocene as a component of mixed conifers woods or as small copses [3].
The distribution dynamics of C. atlantica in Europe and Africa have been conditioned by climate, especially by the warmer temperatures in winter and the dryness of summer [17,18]. Thus, the warmest periods of the Quaternary, such as the last interglacial (200-140 ka) or the Holocene (11 ka-present), determined the migration of its populations towards the south and towards higher altitudes, in a progressive reduction in its distribution [17]. In the Iberian Peninsula, its last populations had to take refuge in the southern mountains of the Baetics as small, isolated nuclei [3,15]. In this context, it has been proposed that the definitive extinction of Cedrus in Europe during the Holocene could be motivated by human action [16,19,20], or even by competition with other species [3,16]. In this sense, an examination of the pollen diagrams of the available deposits on both sides of the Strait of Gibraltar (Baetic Mountain Range and Rif) revealed inverse patterns in the relative abundance of the Pinus and Cedrus pollen types; while, in the Baetics, Pinus dominates and Cedrus is scarce, in the Rif, Pinus appears meager where Cedrus is the dominant type [1,16,17,[21][22][23][24][25][26]. If the ecological affinities of the species are not mutually exclusive and in the absence of barriers to dispersal, this negative association between species suggests the possibility that species of both genera have interacted in the past, mutually excluding each other in areas where they shared a niche and where the one that was more competitive would have ultimately remained [27,28].
The objective of this study is to identify the interspecific interaction through a niche overlap analysis [29,30]. Pinus nigra Arnold was selected among the mountain conifers distributed in the Western Mediterranean due to its ecological affinity with C. atlantica [31][32][33][34]. Likewise, the distribution of both species in the Iberian Peninsula was explored and the location of the possible common area was identified, both currently and in two moments in the past, Last Glacial Maximum (LGM; 21 ka) and Holocene Climate Optimum (HCO; 6 ka). The variations in the distribution of both species would make it possible to understand their respective historical geographical dynamics.
The use of Species Distribution Models (SDMs) has been a complement to the traditional techniques used in Paleobiogeography that can help identify knowledge gaps derived from inherent geographic, temporal and taxonomic biases [35][36][37][38][39]. In fact, SDMs have been previously applied in the projections of different climatic conditions in the past with C. atlantica [14,40,41] and P. nigra [31,42,43].

Materials and Methods
To identify overlap in the potential geographic distribution of C. atlantica and P. nigra, Species Distribution Models (SDM) were used. These models were projected into the past (LGM, 21 ka; HCO, 6 ka). On the other hand, the relationship of the niches in the two-dimensional ecological space was revealed by a principal coordinate analysis (PCA), from which the equivalence and similarity of the niches of both species was evaluated.

Study Area
The study area was framed in the Western Mediterranean, where C. atlantica may have reached its maximum extension during the last glaciation (Würm, 21 ka) [5]. It covers the Iberian and Italic peninsulas, the Western Alps and the Dinarides, in Europe; the Atlas Mountains and Rif, in Africa; and Corsica, Sardinia, Sicily and the Balearic Islands.
For the purposes of the different analyses, the study area was defined by the coordinates lat 29 • N-46 • N; lon 12 • W-19 • E. For the niche overlap, we defined the surrounding environment of each species as the geographic space delimited by the extreme coordinates of the observations used ( Figure 1).

Species
Cedrus atlantica is currently distributed in mountainous areas of the Rif and Middle Atlas in Morocco and the Tell Atlas and Aures Mountains in Algeria [44], between 1300 and 2600 m of altitude, where annual rainfall ranges between 500 and 2000 mm and the minimum temperatures of the coldest month range from −1 to −8 • C [45]. For the analysis of this work, observations made directly in the field, in Morocco, were used; while in Algeria, the recognition of the species was conducted by photo interpretation (Figure 1). Diversity 2022, 14, x FOR PEER REVIEW 3 of 14 . Figure 1. Study area (outlined in red) showing the distribution of the observations for Pinus nigra (blue, n = 257) and Cedrus atlantica (red, n = 211). The surrounding environment of each species, marked by dotted lines, were applied in the niche overlap analysis.

Species
Cedrus atlantica is currently distributed in mountainous areas of the Rif and Middle Atlas in Morocco and the Tell Atlas and Aures Mountains in Algeria [44], between 1300 and 2600 m of altitude, where annual rainfall ranges between 500 and 2000 mm and the minimum temperatures of the coldest month range from −1 to −8 °C [45]. For the analysis of this work, observations made directly in the field, in Morocco, were used; while in Algeria, the recognition of the species was conducted by photo interpretation (Figure 1).
Pinus nigra is a species of medium and high mountain (1500-2500 m of altitude) of circum-Mediterranean distribution, which extends through southern Europe, from the Iberian Peninsula to Anatolia, and reaches the Caucasus, with occasional presence in North Africa, in the Rif and Tell Atlas. The subspecies P. nigra subsp. salzmannii (Dunal) Franco, the one contemplated in this study, is limited to the Western Mediterranean (Cévennes, Pyrenees, Iberian and Central Systems and Baetic mountain range in Europe; Rif and Tell Atlas in Africa). It inhabits rocky terrain where it presents greater resistance to drought than other species of medium-and high-mountain pine [18,32,46]. For the analysis of this work, 257 localities were randomly selected, at the working resolution (30 arc seconds), from the database of the Global Biodiversity Information Facility (GBIF; www.gbif.net). Redundant observations at working resolution were eliminated ( Figure  1).

Variables
The niche predictors were selected from the 19 bioclimatic variables of the 1960-1990 series for the present, from the WorldClim 1.4 portal, (www.worldclim.com, accessed on 1 February 2013). Preliminarily, a survey of the variables that best fit the distribution models of C. atlantica was carried out with MAXENT [47], implemented in the 'dismo' R package [48]. Those variables with the highest percentage of contribution and importance of permutation were selected alternatively when the model was adjusted [49], discarding the rest of the highly correlated variables (values higher than 0.7) in each choice (Table  A1). Pinus nigra is a species of medium and high mountain (1500-2500 m of altitude) of circum-Mediterranean distribution, which extends through southern Europe, from the Iberian Peninsula to Anatolia, and reaches the Caucasus, with occasional presence in North Africa, in the Rif and Tell Atlas. The subspecies P. nigra subsp. salzmannii (Dunal) Franco, the one contemplated in this study, is limited to the Western Mediterranean (Cévennes, Pyrenees, Iberian and Central Systems and Baetic mountain range in Europe; Rif and Tell Atlas in Africa). It inhabits rocky terrain where it presents greater resistance to drought than other species of medium-and high-mountain pine [18,32,46]. For the analysis of this work, 257 localities were randomly selected, at the working resolution (30 arc seconds), from the database of the Global Biodiversity Information Facility (GBIF; www.gbif.net). Redundant observations at working resolution were eliminated ( Figure 1).

Variables
The niche predictors were selected from the 19 bioclimatic variables of the 1960-1990 series for the present, from the WorldClim 1.4 portal, (www.worldclim.com, accessed on 1 February 2013). Preliminarily, a survey of the variables that best fit the distribution models of C. atlantica was carried out with MAXENT [47], implemented in the 'dismo' R package [48]. Those variables with the highest percentage of contribution and importance of permutation were selected alternatively when the model was adjusted [49], discarding the rest of the highly correlated variables (values higher than 0.7) in each choice (Table A1).
Finally, the following four predictor variables were selected: • Minimum temperature of the coldest month (bio6); • Temperature annual range (bio7); • Mean temperature of the wettest quarter (bio8); • Precipitation of the driest month (bio14).
To project the niche model in the geographical space, the climatic layers of WorldClim 1.4 were used, at a resolution of 30 arc seconds for current conditions and HCO (6 ka) and at a resolution of 2.5 arc minutes for LGM (21 ka). MIROC-ESM was chosen among the possible general circulation models [50], because it presented the best fit in its evaluation for the two moments in time in the past for C. atlantica (Table A2). In this evaluation, the available evidence of Cedrus atlantica in the Western Mediterranean from the fossil record for the two periods referred was used. In this review, the sites of the European Pollen Database [51] were considered, as well as compilations of works in the Iberian Peninsula [22] and in North Africa [14,17,21,[23][24][25][26].

Analysis
The niche of each species was modeled and projected in the geographical space using the MAXENT algorithm, included in the 'dismo' R' package [48], a machine learning-type algorithm that only requires presence records and a reference environment (background) to calibrate the models [52][53][54]. The potential habitat of each species was expressed through the binary transformation of suitability (suitable/unsuitable habitat for the species) applying, as threshold, the maximum of the sum (sensitivity + specificity). The overlap of the potential habitat of both species indicates the favorable areas for both species, thus where interaction between them is possible. There are numerous antecedents in the application of the SDM to the past and, especially, in the Iberian Peninsula [31,35,55].
The validation of the models was carried out by contrasting the value of the observed area under the curve (AUC) with the null model AUC [56]. The null model was built with the distribution of AUC values from 100 SDM iterations built with n randomly selected points (257 for P. nigra and 211 for C. atlantica) from the background. The null model AUC represents the 0.95 quantile in this distribution and the model is significantly valid (p < 0.05) if the observed AUC value is greater than the null model AUC.
On the other hand, the niche overlap in the two-dimensional environmental space was evaluated using the 'ecospat' R' package [57]. This is defined by the two main axes of a PCA built from the values of the four variables selected in the surrounding environment of C. atlantica and P. nigra. (Figure 2). This two-dimensional environmental space was gridded (100 × 100 cells) and the density surfaces for the observations were projected on this grid to represent the niches [29]. The relative positions of the niches of both species allowed us to infer possible interspecific relationships.  Niche overlap was measured with Schoener's D index (ranging from 0, no overlap, to 1, total overlap). The niche overlap was contrasted using the equivalence and similarity tests [29,30]. In each of them, the observed D-value is compared with a null distribution of D-values obtained as the result of the iteration of 100 superposition operations between pairs of sets of random observations. The equivalence test is aimed at identifying identical niches in the whole surrounding area. In this test, each of the two sets of observations of the null distribution is built by random selection among the whole surrounding area of the two species. The niches of the two species are admitted to be equivalent if the observed overlap value is significantly higher than the 0.95 quantile of the null distribution. That is, the two niches are more similar to each other than any two random sets of the whole surrounding area.
The similarity test assesses the similarity of one niche to another in relation to the surrounding environment in which it is found. The surrounding environment of each species is delimited by the extreme coordinates of the observations (Figure 1). In this test, the total of observations from which the pairs are drawn at random, for the null distribution, was constructed with the observations of one species and an extract from the surrounding environment of the opposite species (in a number equal to that of their observations). Thus, in this test, two asymmetric analyzes are carried out in which the observations of each species are compared with the surrounding environment of the species being confronted. Therefore, two species with environmental affinities would have higher overlap values than expected when comparing the niche of one species with the surrounding environment of the opposite. The niches of the two species are significantly similar if the D-value exceeds the 0.95 quantile of the null distribution [29,30].

Geographical Distribution
The four selected environmental variables defined the fundamental niche of both species with a high level of prediction of the optimal habitat for the current conditions. The validation of the models showed that the observed AUC values greatly exceeded those marked for a null distribution [56], with a significance of 95% (Table 1; Figure A1). According to the projections made in MAXENT, the potential habitats of C. atlantica and P. nigra are mostly geographically segregated from each other at present (Figure 3c). Currently, the optimal habitat for C. atlantica is distributed mainly in North Africa, from the High Atlas in the immediate vicinity of the Atlantic Ocean to the Tell Atlas and Aures Mountains, via the Rif and the Middle Atlas. In Europe, this species finds its optimal habitat in the Iberian Peninsula, Baetic Mountain Range, and in the Central System, extending to the Leon Mountains. The predicted distribution of P. nigra is framed in the European Mediterranean basin, from the Iberian and Central Systems and the Cantabrian Mountains, in the Iberian Peninsula; to the Apennines and western Alps, through the Pyrenees and Massif Central. Marginally, it is also represented in the Aures Mountains, Rif and Middle Atlas on the African continent. At present, C. atlantica and P. nigra would share their distribution in the Baetic Mountain Range and Central System in Europe and in the Middle Atlas and Aures Mountains in Africa (Figure 3c).
Compared with current conditions, the projection to the LGM (21 ka) shows the descent in altitude of the most favorable conditions for C. atlantica, where it would share the habitat with P. nigra (Figure 3a). This was revealed by a greater extension of the potential habitat of the species throughout the mountain systems of Africa and the prevalence of the plateaus in the Iberian Peninsula, as well as their presence in the mountain systems near the Tyrrhenian. In the Iberian Peninsula, the common areas would extend through the Iberian System and the inland areas close to the Cantabrian Mountains. In Africa, the common habitat would be linked to the High, Middle and Saharan Atlas and the Aures Mountains (Figure 3a). The warmer conditions of the HCO (6 ka) with respect to the LGM would have determined the ascent in altitude of the optimal conditions for both species (Figure 3b). This would have resulted in a drastic reduction in the potential area of C. atlantica and its isolation in a number of disjoint areas, largely coinciding with the current distribution of the species. The presence of the habitat of each of the species would have been virtually restricted to one continent, C. atlantica to Africa and P. nigra to Europe, with negligible presence on the opposite continent, where both species would have shared limits of habitat distribution marginally.

Niche Overlap
The analysis of niche overlap in the two-dimensional ecological space reflects the overlap with a D-value close to 0.37. According to the strict nature of the test, the equivalence test result was that both species niches were non-significantly equivalent (they did not match completely). Contrary to the equivalence test, both similarity tests indicated that both species niches were significantly similar to each other considering their different surrounding environment (Table 2; Figure A2). Table 2. Niche overlap of Pinus nigra and Cedrus atlantica and equivalence and similarity tests. Similarity implies two asymmetric analyses (niche→surrounding environment; 1 and 2 represent P. nigra and C. atlantica, respectively).
The arrangement of the density grids in the two-dimensional ecological space shows a relative displacement between the niches of both species along the second axis of the PCA. (Figure 4). This could be interpreted as a coincidence in the affinity of the species for the variables that determine the first axis of the PCA (precipitation of the driest month and mean temperature of the wettest quarter; Figure 2), while differences appear in the niche determined by the two remaining variables, minimum temperature of the coldest month, whose lowest values are better tolerated by C. atlantica (from −6 to −2 • C, compared to −3-2 • C, where P. nigra lives); and the temperature annual range, with C. atlantica better tolerating higher ranges in temperature (from 29 to 36 • C, compared to 23-26 • C in which P. nigra lives; Figures A3 and A4).

Discussion
According to our results, C. atlantica currently extends the distribution of its potential habitat in the Iberian Peninsula even though the species does not currently inhabit Europe. These results are consistent with other models made for C. atlantica [17,41]. During the LGM (21 ka), the extension of the habitat would occupy most of the center of the Iberian Peninsula, dominating the plateaus and mountainous areas of the Iberian System, the interior areas near the Cantabrian Mountains and most of the Baetic Mountain Range. The presence of Cedrus in the fossil records could support these projections to the past [11,58,59]. The increase in temperature during the HCO (6 ka) determined the increase in altitude of both species. Populations of C. atlantica were isolated in the mountain systems, which served as a refuge, from where the species spread to the present day. This is consistent with the dynamics of cedars in North Africa during the Holocene found in fossil records [1,15,17,21,23,24,26]. Therefore, the dynamics of the cedars in Europe have been linked to the climate since the last ice age, so their disappearance on the continent could have been due to a greater extent to a progressive reduction in their potential habitat [17,18].
This would not be incompatible with the existence of biotic relationships as a determining agent of local extinction in the redoubts present in the Baetic Mountain Range (Postigo-Mijarra et al., 2010). Bearing the possibility that C. atlantica and P. nigra interacted in the Iberian Peninsula in the past, occupying similar environmental spaces, they could mutually have excluded each other in mountainous areas where the more competitive species would have remained. The starting hypothesis arises from the fact that the distribution of each species is restricted to separate continents and that the Cedrus and Pinus spp. show some negative temporal association in the pollen profiles of the fossil records from both sides of the Mediterranean [1,16,17,[21][22][23][24][25][26]. However, both species show a similar ecological behavior (they form forests in the high Mediterranean mountains), which raises the question of whether and to what extent they share a niche. If so, the negative association of their presence (geographical and temporal) could be associated with biological interaction [27,28].
The results obtained indicate that the distribution of the potential habitat of C. atlantica and P. nigra overlapped both in Africa and in Europe during the Pleistocene, being able to maintain contact between their populations in their distribution limits. In addition, the niche overlap tests indicate that both species habitats are significantly similar to each other. In the geographical space, the results obtained indicate that the distribution of the potential habitat of C. atlantica and P. nigra overlapped both in Africa and in Europe during the Pleistocene, being able to maintain contact between their populations in their distribution limits. Warmer temperatures during the Holocene caused both species, C. atlantica and P. nigra, to migrate upwards [17]. In this situation, both species reduced their shared habitat and, even so, it is possible that their populations maintained contact. Thus, they could have been growing in the same stands as mixed coniferous forests; Cedrus could have formed part of pine forests as a companion species, or formed small separate copses in the same territory that P. nigra would have inhabited. Going up in altitude, the smaller distribution area of Cedrus in the Iberian Peninsula with respect to Pinus would make the species more sensitive to certain stochastic processes that would lead to habitat fragmentation and, in the extreme case, to its local disappearance, despite the existence of available habitat.
According to our results, the spatial segregation between species is much greater than that of its niches. Both species share their niches and both species have also shared habitat on the same continent since the Pleistocene. Nevertheless, their current distributions have drifted towards complete spatial segregation. It cannot be ruled out that this fact is due to the interaction between both species. In fact, the geographic segregation may have affected the results of the analysis by introducing a bias in the observations used to calibrate the distribution models. If there had been any observations of C. atlantica currently in Europe, the resulting model would likely have been different, with a much larger potential habitat distribution and greater overlap with the P. nigra niche.
In conclusion, C. atlantica and P. nigra show differentiated niches based on greater tolerance to extreme cold and continentality of Cedrus. Currently, the distribution of both species, segregated on different continents-P. nigra in Europe and C. atlantica in Africa-is virtually responding to dissociation from their optimal habitats. However, both species were able to share habitat in their distribution limits on the European continent until the end of the Pleistocene. The increase in temperatures of the Holocene Optimum Climate must have forced the rise of both species in altitude, causing the isolation of their populations in mountain refuges and the consequent loss of connectivity that has lasted until today. Although C. atlantica has continuously maintained its ideal habitat in the Baetic Mountains to date, it cannot be ruled out that the higher prevalence of P. nigra on mainland Europe has eventually displaced Cedrus.

Data Availability Statement:
The data used in this study are available upon request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.