3. Results

3.1. UCE sequencing

All combined, the Miseq runs produced on average 526,308 (SD ± 356,015) reads per specimen. The average number of loci per specimen after assembly was 7398 (SD ± 5099). After retaining only loci matching the UCE reference file and filtering for a 75% matrix completeness, the number of retrieved loci varied between 686 and 1860 across the five species complexes (Supplementary table S3). Resulting concatenated reads were on average 899,690 bp long; the shortest reads were obtained forN. goodeniana/succincta (i.e. 429,820 bp) and the longest forA. amieti/bicolor (i.e. 1,459,230 bp). In total, seven specimens did not pass the 90% missing data filter (Supplementary information S1) and were therefore excluded from downstream analyses. No significant departure form Hardy-Weinberg equilibrium was observed. The bi-allelic SNPs screening recovered on average 18,545 bi-allelic SNPs (SD ± 23,430) across all species complexes (Supplementary information S3).

3.2. Comparison of mitochondrial and nuclear analyses

3.2.1 Species complex 1: Andrena allosa/amieti/bicolor/montana

The phylogenetic trees performed on the COI and UCEs datasets provided very similar topologies for A. bicolor (Figure 1). Both trees showed distinct monophyletic clades with strong bootstrapping support (hereafter BS values) for the UCE tree. Only one specimen (i.e. “GBIFCH00135933”) sampled in Southern France was misplaced in the UCE tree compared to its position in the corresponding mitochondrial tree. For A. amieti , all specimens formed one monophyletic clade in the UCE tree, contrasting with the two paraphyletic clades in the COI tree. Furthermore, there was no apparent structuring in the UCE tree among mitochondrial lineages sampled in the alpine region. Specimens with large amounts of missing data (≥ 60%) exhibited longer branches in the UCE tree (e.g. “GBIFCH00131686”, “GBIFCH00136250”). Strong isolation by distance was found between the southern Italian and Alpine populations (R2 = 0.7221, p-value = 0.039; Supplemental information S4). These two individuals formed a distinct monophyletic clade sister to all alpine specimens in the UCE phylogenetic tree.
The DAPC with a prior knowledge on species assignments correctly reassigned membership for the majority of specimens (Supplementary information S5). The plotted cumulative variance of the eigenvalues suggested to retain the first eight principal components (conserving 70.3% of the total variance; Supplementary information S6). All specimens were correctly reassigned with a 100% membership probability for three taxa (i.e. A. allosa , Andrena sp3, A. montana ). For A. bicolor ML1, only one specimen (i.e. “GBIFCH00135933” from Southern France) revealed mixed membership probability, with 45.64% posterior probability (hereafter “pp”) to belong to ML1 and 54.36% to ML2. For A. bicolor ML2, one specimen (GBIFCH00117401) also showed mixed membership with a slight probability (i.e. 3.8%; Supplementary information S5) of belonging toA. bicolor ML1. The mixed membership probability for those two specimens are congruent with the PCA results where both specimens are found marginally away from the main A. bicolor ML2 aggregation (Supplementary information S7).
Between both lineages of A. amieti , genetic cluster assignments were much less supported and only the two specimens sampled in south Italy were assigned with a 100% probability to ML1. The lack of a clear separation between mitochondrial lineages for the alpine specimens is suggesting considerable levels of admixture.
The AMOVA (Table 1) depicted strong genetic difference (i.e. 43.43%; p-value ≤ 0.0001) between the two A. bicolor lineages but no significant difference between the two A. amieti lineages. The lowest, yet significant fixation index (Table 2) was obtained between both lineages of A. bicolor (Fst = 0.138). Nei’s genetic distance between both lineages within A. bicolor (i.e. 0.00231) was slightly higher that between A. allosa and A. amieti ML1 and ML2 (i.e. 0.00261 and 0.00299, respectively).
The GMYC analysis computed on all specimens identified nine clusters (Supplementary information S8): two clusters corresponding to the mitochondrial lineages found within A. bicolor , one cluster with with A. amieti and A. allosa merged together, and six clusters composed of only one specimen. The bGMYC analyses identified the same 9 clusters with a high probability (p = 0.95 - 1; Supplementary information S8). All other scenarios had very low posterior probabilities (pp = 0 - 0.05). The two parallel BPP analyses converged and were highly congruent (Supplementary information S9). Both runs depicted: (i) one tree model [((A. allosa, A. amietiML1+ML2), ((A. bicolor ML1, A. bicolor ML2), Andrena sp3 ))] with a posterior probability of ≥ 0.99; (ii) 5 delimited species (i.e. Andrena sp3 , A. bicolor ML1, A. bicolor ML2, A. allosa , A. amieti ML1+ML2), all with a posterior probability of 1; (iii) and a posterior probability of 1 for having 5 species present in the dataset. Finally, the DAPC analyses performed without a prior knowledge on species identifications identified K = 3 and K = 4 as best solution for A. amieti ML1+ML2 and A. bicolor ML1+ML2, respectively. Table 3 summarized the number of clusters found for each analysis.

3.2.1. Species complex 2: Andrena barbareae/cineraria

Mitochondrial and nuclear phylogenies were discordant, with no clear separation between both species in the COI tree and two well-supported monophyletic clades corresponding to the two morphological species in the UCE tree (100% BS values; Figure 1). Results from the UCE phylogenetic tree were corroborated by the PCA in which both species were clearly separated (Supplementary information S7). The DAPC analyses correctly reassigned membership for all specimens with 100% probability (Supplementary information S5). The AMOVA revealed that 52.74% of the total observed variance could be explained by the species level. The GMYC and bGMYC provided similar results (Supplementary information S8). Both analyses grouped all specimens morphologically identified asA. cineraria in one clade. For A. barbareae , both analyses suggested the presence of three distinct clades. Both BPP runs were congruent and highly supported the presence of three species (pp = 1; Supplementary information S9). The runs however disagreed with respect to the phylogenetic relationships among the three species. The first run depicted three different possibilities for the species trees, with the most likely (pp = 1) tree being: [(A. barbareae , (A. cineraria , A. vaga)) ]. The second run depicted only one tree [((A. barbareae , A. cineraria) , A. vaga)],corresponding to the expected tree based on the phylogeny. In the first run this solution was supported at 75%. Congruent with the BPP analyses, the DAPC identified K = 3 (with outgroup) as the best solution (Supplementary information S10, Table 3). Clustering of all specimens corresponded to the morphological identifications.

3.2.2. Species complex 3: Andrena carantonica/trimmerana/rosae

Mitochondrial phylogenies recovered well-supported (BS 72-100%) monophyletic clades for A. carantonica and A. rosae , but not for A. trimmerana , which was composed of two clades forming a paraphyletic unit from which A. carantonica arose. One clade was composed of two specimens of A. trimmerana sampled in western Switzerland and was sister to the A. carantonica clade; support for this sister relationship was high (BS 93%; Figure 1). In contrast, all three species appeared as strongly supported monophyletic clades in the UCE tree (≥ 90% BS; Figure 1). Spring and summer generations ofA. rosae and of A. trimmerana were intermixed in both mitochondrial and UCE trees, supporting the view that A. stragulata and A. spinigera constitute the morphologically differentiated spring generation of A. rosae and A. trimmerana , respectively. Genetic distance between A. carantonica and A. trimmerana was relatively low (Nei’s D = 0.00061), although AMOVA and pairwise Fst depicted significant difference between both species (Table 1-2). The PCA with all three species showed no difference between A. carantonica and A. trimmerana , however when removing A. rosae from the analyses, both species were separated on the first two components (Supplementary information S7). The GMYC and bGMYC analyses failed to separate bothA. carantonica and A. trimmerana. In the DAPC both species were also not separated with K = 2 and K = 3 but were with K = 4 (Supplementary information S10). All three clustering scenarios had very closed BIC values. The BPP analyses however highly supported the presence of three distinct species, with the following tree topology [((A. carantonica , A. trimmerana ), A. rosae )].

3.2.3. Species complex 4: Andrena dorsata/propinqua

Strong mito-nuclear discordances were observed within this species complex. In mitochondrial trees, Swiss specimens formed two clusters corresponding to morphological identifications (Figure 1), but one specimen of A. propinqua (GBIFCH00133244) collected in southern France was sister to a well-supported clade containing all other specimens of A. dorsata and of A. propinqua (Figure 1). Our sampling also included one specimen of A. dorsata from this French site (GBIFCH00133243). Phylogenetic trees and PCAs based on UCEs recovered both species as separated clusters (Figure 1, Supplementary information S7); the French specimens were not particularly divergent. Both GMYC and bGMYC analyses, the BPP analyses and DAPC analysis successfully separated both species (Table 3, Supplementary information S8-S10). Taken together, these results indicate that A. dorsataand A. propinqua are valid species.

3.2.4. Species complex 5: Lasioglossum alpigenum/bavaricum/cupromicans

Comparison of mitochondrial and nuclear trees revealed mitochondrial-nuclear discordance for L. bavaricum and L. cupromicans (Figure 1). Both taxa were well delimited with highly supported monophyletic clades (95% bootstrap value) in the UCE tree but shared the same COI barcode. All L. alpigenum specimens clustered in a single monophyletic clade sister to both other taxa in both trees. Beside the GMYC analysis that over-clustered L. bavaricum andL. cupromicans (Supplementary information S8), all other analyses were congruent (Table 3) and supported the hypothesis of three distinct species as previously postulated based on morphology.

3.2.5. Species complex 6: Nomada goodeniana/succincta

COI and UCE trees depicted two well defined monophyletic clades (Figure 1). The bootstrap support for monophyly of N. goodeniana was low in the mitochondrial trees due to the presence of two slightly divergent specimens of N. goodeniana collected south of the Alps. In the nuclear trees, these two specimens clustered with other specimens ofN. goodeniana with high support values. The species delimitation tests also highly supported the hypothesis of two separated species (Table 3, Supplementary information S8-S10).