Spatial analysis
Occurrence records were obtained from 287 observations from the Global
Biodiversity Information Facility (211 occurrences validated by
morphology, http://www.gbif.org/) and data from herbarium
specimens (76 occurrences). For modelling purposes, the dataset was
reduced to keep only georeferenced data.
We used environmental niche modelling (ENM) approaches to reconstruct
the potential distribution of the four main clades uncovered in the
phylogenomic analyses (see Results) under current and past climatic
conditions using the maximum entropy algorithm implemented in the R
package ‘Maxent’ (Phillips, 2021). Nineteen bioclimatic variables were
extracted from the Bioclim dataset, provided by WorldClim 1.4 in a
GIS-based raster format (2.5 min resolution). The correlations between
environmental variables were determined with a Pearson’s correlation
matrix and subsequent realization of a dendrogram cluster for its
visualization (Supplementary Data Figure S1). We selected a different
set of uncorrelated variables for each geographical region with a high
percentage of contribution (PC): bio4 (temperature seasonality), bio8
(mean temperature of wettest quarter), bio9 (mean temperature of driest
quarter), bio15 (precipitation seasonality) and bio16 (precipitation of
wettest quarter) for the circum-Mediterranean region; bio3
(isothermality), bio6 (min temperature of coldest month), bio8 (mean
temperature of wettest quarter), bio14 (precipitation of driest month),
bio15 (precipitation seasonality) and bio16 (precipitation of wettest
quarter) for the Eastern Mediterranean region; bio3 (isothermality),
bio4 (temperature seasonality), bio16 (precipitation of wettest quarter)
and bio18 (precipitation of warmest quarter) for the Macaronesian
region. ENM analyses were carried out under current climatic conditions,
and projected to climatic conditions of the Mid Holocene (MH, about
6,000 years ago), the Last Glacial Maximum (LGM, about 22,000 years ago;
Braconnot et al., 2007) and the Last Interglacial (LIG,
~120,000 - 140,000 years BP; Otto-bliesner et
al., 2006) using the palaeoclimatic ‘Community Climate System Model’
(CCSM; Gent et al., 2011). Layers were cropped to represent the
distribution range of each phylogenetic group (i.e., DC1, DC2 and DC3
for D. communis , and D. orientalis ; see Results) to
maximize the reliability of the results and discard false occurrences
using the package ‘raster’ (version 3.5-15; R v.4.0.5). We used the
Schoener’s D and Hellinger’s I indices as implemented in ENMtools
v.1.0.4 to evaluate niche overlap (Warren et al., 2008, 2010).
Equivalence and similarity tests with 1,000 replicates were carried out
to assess if the overlap between ENMs is higher than expected under
randomized ENMs.
A multivariate ordination analysis (principal component analysis, PCA)
was carried out using uncorrelated bioclimatic variables obtained from
WorldClim using the packages ‘ade4’, ‘factoextra’, ‘magrittr’, ‘dismo’
and ‘HH’ in R v.4.0.5. The correlation analysis was performed with a
Pearson correlation matrix and the subsequent visualization and
selection of variables using a dendrogram cluster. The selection of
uncorrelated variables with the highest contribution to the PCA were 10:
bio1 (annual mean temperature), bio2 (mean diurnal range), bio3
(isothermality), bio7 (temperature annual range), bio8 (mean temperature
of wettest quarter), bio9 (mean temperature of driest quarter), bio10
(mean temperature of warmest quarter), bio12 (annual precipitation),
bio15 (precipitation seasonality) and bio19 (precipitation of coldest
quarter).