Sequence variation and phylogenetic analyses
From 160 wedge-tailed sabrewings we obtained mtDNA sequences that contained 81 polymorphic sites (38 in the ATPase and 43 in the control region) of the 1,407 bp of the genes analyzed. No insertions or deletions were present; therefore, variants were identified based solely on nucleotide substitutions. Sixty-three haplotypes were identified for the 22 wedge-tailed sabrewing populations. Haplotype and nucleotide diversity are summarized in Table . Haplotype proportion relative to the number of samples per population was more than 50% in most cases, indicating high haplotype diversity. No differences between populations were evident. Nucleotide diversity was low, indicating little variation between sequences from the same population.
Population genetic variability of Campylopterus curvipennis
The consensus tree obtained from Bayesian inference clustered the haplotypes into two main well-supported clades (posterior probabilities of 1.0), corresponding to disjunct western (curvipennis, SMO) and eastern (pampa, YUC) groups on either side of the Isthmus of Tehuantepec (Figure ). Haplotypes of the excellens group from the Tuxtlas region were clustered in a well-supported clade nested within the poorly resolved SMO clade.
Figure 1 Phylogenetic relationships among mtDNA haplotypes based on Bayesian inference. Values above branches denote posterior probabilities, and numbers at the tip of the branches indicate distinct haplotypes. Outgroups were Campylopterus rufus, C. hemileucurus (more ...)
Assuming a constant molecular clock and rates of 2 and 5% divergence per My, for the ATPase coding region BEAST estimated that the SMO and TUX clades diverged from the YUC clade 1.47 (0.35-3.42) and 0.52 (0.13-1.21) Mya, respectively. These results suggest that the split between the Sierra Madre Oriental and Yucatan Peninsula clades may have occurred during the mid-Pleistocene. With respect to SMO, for the TUX clade the constant clock TMRCA estimated divergence times of 614,000 (0.12-1.55 Mya) and 202,000 (0.05-0.47 Mya) years ago for 2 and 5% substitution/My, respectively. This implies that the split between the TUX and SMO clades may have occurred more recently in the late Pleistocene.
Phylogeographic and genetic structure
The haplotype network showed a strong phylogeographic structure among three groups (SMO, TUX, YUC), but not among populations or areas within the SMO (nSMO, cSMO, sSMO) (Figure ). Of the 63 haplotypes obtained, 44 were private to SMO populations. The most frequent haplotype was shared by samples from all of the SMO populations, and the rest were low frequency haplotypes (Figure ). Seven haplotypes were private to the Tuxtlas region, and are separated from the SMO haplotypes by five mutational steps. The remaining haplotypes (12) were exclusively found in populations from the Yucatan Peninsula, and formed a separate network (Figure ).
Figure 2 Sampling localities, geographic distribution, and genealogical relationships of mtDNA haplotypes. A. Pie charts indicate the frequency of occurrence of haplotypes in each population, and colors correspond to those shown in the haplotype network below. (more ...)
The AMOVA results revealed significant genetic differentiation at every hierarchical level. When grouped by geographic areas (nSMO, cSMO, sSMO, TUX, and YUC; see Table and field procedures for the definitions) most of the variation (85%) was explained by differences among areas, whereas variation among populations within areas (0.92%) and variation within populations (13.94%) only explained a small percentage of the total variation (Table ). The same pattern was observed when grouped by subspecies (curvipennis
, SMO; excellens
, TUX; and pampa
, YUC) with the highest percentage of the variation being explained by differences among groups (90.41%) and the lowest by differences among populations within groups (0.55%; Table ). ΦCT
values for both analyses were very high (Table ) indicating high levels of genetic differentiation among areas and groups. Pairwise comparisons of FST
among sampling localities ranged from 0.003 for Aqm/Ama to 0.95 for UG/Gar and were significant in most cases between sampling localities from the TUX region and the other localities, and for the YUC region and the other localities [Additional file 1
]. Comparisons between sites located in the same geographic area (SMO, TUX, YUC) were not significant [Additional file 1
Analysis of molecular variance (AMOVA) for the control region and ATPase 6-8 mtDNA genes of Campylopterus curvipennis
Across all sampling localities, the number of alleles per locus varied from 4 to 20, and within localities across all loci this number varied from 13 to 62 (mean alleles per locus, 1.3-6.2) (Table ). Observed heterozygosity values did not consistently deviate from H-W equilibrium. Only four localities were not in H-W equilibrium after Bonferroni corrections at locus CACU13-2 and one locality at locus CACU13-7 (Table ) probably due to the presence of null alleles. No significant linkage disequilibrium was detected in any of the population-loci comparisons after Bonferroni corrections.
Significant genetic subdivision was detected among sampling localities (global RST
estimate ± SE, 0.085 ± 0.0021, P
< 0.0001), but not among sampling localities within the fragmented areas of the SMO (global RST
estimate ± SE, -0.020 ± 0.0021, P
> 0.05). Pairwise FST
were quantitatively similar but we only report RST
values because the stepwise mutation model implemented in this estimate is more appropriate for microsatellites [Additional file 1
]. Pairwise RST
values among sampling localities ranged from 0.001 for Xil/Ord to 0.42 for Risc/Nov. As with the pairwise FST
values for the mitochondrial genes, estimates of RST
were significant for almost all comparisons between localities from TUX versus
the rest of the localities and comparisons between localities from YUC versus
the rest after Bonferroni corrections. Values from comparisons between sampling localities within SMO were not statistically significant. Comparisons between samples from Bec and Gar (from the YUC region) versus
the rest were not significant in most cases after Bonferroni corrections probably due to sample size, affecting statistical power.
The phenogram constructed with pairwise RST
, recovered sampling localities from the YUC group clustering together in a basal position followed by the TUX and the SMO groups [Additional file 2
]. Within the SMO group, samples from the northern limit of the distribution (Tamaulipas) are clustered together, as are samples from San Luis Potosi (central distribution). However, localities from the southern distribution (Veracruz and Puebla) are scattered throughout the tree [Additional file 2
Results of the STRUCTURE analysis corroborated population substructure. Using all genotypes together the log likelihood was highest for K = 6; however, when ΔK is calculated the break in the slope of the distribution of L(K) was at K = 2 (Figure ). One cluster includes populations from the Yucatan Peninsula (YUC), and the other includes samples from populations from the TUX and SMO groups (Figure ). However, in the clustering pattern at K = 3, samples were assigned with high probabilities to three clusters, corresponding to the YUC, TUX and SMO groups, which is congruent with the sequence mtDNA data results. In further analyses of these data from which the YUC samples were first excluded, STRUCTURE assigned samples to two groups (TUX and SMO); in the analysis where the TUX samples were excluded two clusters were detected for the SMO samples with no evidence of genetic clustering at this level. These results suggest a hierarchical clustering pattern where YUC samples are highly divergent and samples from TUX are sub-structured in a cluster that contains samples from SMO.
Figure 3 Spatial genetic structure determined by a Bayesian assignment analysis using STRUCTURE. A. Posterior assignment probabilities of 160 individuals of Campylopterus curvipennis to an optimal number of K = 2 and K = 3. Each individual is represented by a (more ...)
Among-group comparisons of gene flow (M) from the MIGRATE analysis ranged from 0.49 (TUX to YUC) to 1.98 (TUX to SMO), however, all comparisons were less than or not significantly greater than 1.0 (Table ). None of the comparisons among genetic groups indicate contemporary gene flow, suggesting that the Isthmus of Tehuantepec is preventing gene flow between populations on either side of the isthmus, and that the more recently divergent populations from the Tuxtlas region and SMO also have interrupted gene flow. Estimates of M among populations were only significantly greater than 1.0 between those within the SMO (M values ranged between 1.4 and 2.3) and the direction of gene flow was asymmetric in most cases. The direction was mainly northwards from sSMO to cSMO, and only in few cases did gene flow reached populations on the northern limit of the distribution (Ciel, GF). Estimates of gene flow between populations of different genetic groups were not significantly greater than 1.0.
Estimates of M (mutational corrected migration) from the MIGRATE analysis of microsatellites among genetic groups
Our simulations to test whether genetic groups might have originated in the face of gene flow were consistent across replicates, and produced confident posterior probability peaks for the parameters estimated. When testing for migration following the split between the SMO and TUX groups, migration rates were higher than 1 in the SMO-TUX direction (m = 1.416, 90% highest posterior density (HPD), 0.004-4.21) and close to 1 in the opposite direction TUX-SMO (m = 0.91, 90% HPD, 0.007-4.47), suggesting that divergence took place in the presence of gene flow from the SMO to the Tuxtlas region. In contrast, low migration rates following the split were estimated for populations east and west of the Isthmus of Tehuantepec (west: m = 0.52, 90% HPD, 0.003-1.27; east: m = 0.36, 90% HPD, 0.002-1.23), indicating that the population split occurred in the absence of gene flow, i.e. this geographic barrier was not permeable to migrants. Assuming a mutation rate for microsatellites of 2.96 × 10-3, estimates of effective population size indicate that the ancestral population was significantly larger (3,924 individuals, 90% HPD, 2,598-19,005) than the populations after divergence. The TUX group had the lowest population size (103.88 individuals 90% HPD, 38.01-259.3) followed by YUC (175.67 individuals 90% HPD, 96.28-282.93) and SMO (326.01 individuals 90% HPD, 203.54-474.66).
Morphological and acoustical variation
The MANOVA showed significant lineage differences in male morphology (Wilks' Lambda, F6,222
= 33.77, P
= 0.0001), and these differences were significant for each morphological variable (one-way ANOVAs: wing chord F2,113
= 16.71, P
= 0.0001; tail length F2,113
= 48.45, P
= 0.0001; culmen F2,113
= 57.58, P
= 0.0001, SMO n
= 106, TUX n
= 6, YUC n
= 16; Figure ). The relationship between culmen length and wing chord (as a measure of body size; r
= 0.38, P
< 0.0001) was further examined for possible allometric effects using the residuals of this relationship, and the differences between lineages were maintained (F2,113
= 48.33, P
< 0.0001). For females the MANOVA also showed significant lineage differences (Wilks' Lambda, F6,42
= 3.89, P
< 0.005), and the results of the univariate ANOVAs were all significant (wing chord F2,23
= 10.17, P
< 0.001; tail length F2,23
= 11.94, P
< 0.0005; culmen F2,23
= 4.65, P
< 0.05, SMO n
= 21, TUX = 11, YUC n
= 5; Figure ). Overall these results are congruent with previous subspecific designation based on bill and body size [31
], indicating that individuals from the excellens
group (TUX) are significantly larger than those of the curvipennis
(SMO) and pampa
(YUC) groups, and that pampa
individuals have shorter bills.
Figure 4 Morphological traits taken for males and females of C. curvipennis lineages. Boxplots show the 10th, 25th, 50th (median), 75th, and 90th percentiles. Values above the 90th and below the 10th percentile are plotted as open circles. Illustrations correspond (more ...)
A total of 344 syllable types were detected across populations. Songs were very versatile and no successive syllable repetitions were detected except for vocalizations recorded in the YUC region from the C. c. pampa lineage, where individuals commonly repeated one of the syllables three times in succession (Figure ). Based on vocal similarity measures (syllable type sharing), individuals from each sampling locality were clustered accordingly (Figure ). In addition, the analysis showed individuals from the YUC lineage group together in a basal position followed by individuals from the TUX and SMO lineages (Figure ). On average populations from SMO shared a lower proportion of syllable types with the TUX and YUC lineages (0.062 ± 0.018), than within any other locality from the SMO (0.145 ± 0.037). Despite the great syllable diversity observed across populations and regions, there were some syllables shared by all recorded individuals, whereas other were shared only by members from populations from the SMO. Regarding acoustic measurements of the common syllable, the first three PCs accounted for 53% of variation (PC1 = 24.4%, PC2 = 18.2%, PC3 = 10.4%). One-way ANOVAs yielded significant group differences in the first PC (PC1, F2,74 = 107.05, P < 0.0001), but not in the second and third (PC2, F2,74 = 2.13, P > 0.05; PC3, F2,74 = 1.01, P > 0.05 SMO n = 66, TUX = 6, YUC n = 5). PC1 was mainly explained by minimum and peak frequency of note 1, duration of note 3, and frequency range, peak frequency and duration of note 4.
Figure 5 Dendrogram generated by cluster analysis of a presence/absence matrix of syllable types from recordings of individuals. Colored lines correspond to individuals recorded at different geographic locations: green, purple and blue correspond to sites located (more ...)
Comparisons between morphological, acoustic, habitat-related, and genetic distances
Mantel tests showed a strong positive correlation between genetic (pairwise FST) and morphological distances for males (r = 0.42, P < 0.05), even when controlled for geographic distance (r = 0.38, P < 0.05). However, the relationships between pairwise RST and morphology were not significant. In the case of females none of the relationships were significant. Regarding acoustic distances, Mantel test showed a strong positive correlation between genetic and song sharing distance (FST: r = 0.72, P < 0.05; RST: r = 0.71, P < 0.05), and also for the genetic and common syllable distance (FST: r = 0.72, P < 0.01; RST: r = 0.72, P < 0.01). These relationships were maintained when geographic distance was accounted for (song sharing: FST, r = 0.50, P = 0.05; RST, r = 0.51, P < 0.05; common syllable: FST, r = 0.44, P < 0.05, RST, r = 0.44, P = 0.01). These analyses indicate that more genetically divergent lineages shared fewer syllable types, and differed acoustically in the common syllable, independently of distance. Lastly, we did not find a significant relationship between morphological and acoustic distances versus habitat-related (climate and topography) distances for either males or females.
The role of drift and selection: a coalescent test
Coalescent simulations [30
] of morphology and mtDNA sequence data from three populations with large sample sizes of sequences corresponding to the three genetic groups in which morphological characters are fixed (Tux, Ciel, Nov) were used to differentiate the roles of selection and genetic drift in morphology evolution. If morphological characters, assumed to be encoded by nuclear genes, have sorted significantly faster than mtDNA haplotypes, then the hypothesis of divergent selection rather than that of drift is supported. The observed s
value for reconstructed trees was equal to 2 when considering three populations for morphology characters. The upper 95% CL for time since population divergence was 3Ne
generations assuming a dichotomous branching model of divergence and 3.3Ne
generations assuming a simultaneous model. When branch lengths were reduced four times to mimic nuclear genes, an s
value of 2 occurred in one third of the simulated trees for dichotomous branching, and in one fifth of the simulated trees for simultaneous branching. There is a high probability (P
> 0.05) that nuclear genes would be fixed in these populations under neutrality, suggesting that drift cannot be rejected as a possibility for morphological divergence [Additional file 3
The coalescent simulations of the song and mtDNA sequence data from four populations with large sample sizes of sequences in which song characters are fixed (Ciel, Ord, Tux and Nov) were used to differentiate the roles of selection and genetic drift in song evolution. The observed s
value for reconstructed trees was equal to 8 when considering four populations for vocal characters. The upper 95% CL for time since population divergence was 0.48Ne
generations assuming a dichotomous branching model of divergence and 0.68Ne
generations assuming a simultaneous model. When branch lengths were reduced four times, an s
value of 3 occurred in none of simulated trees, suggesting that there is a very low probability (P
< 0.0001) that nuclear genes would be fixed in different populations under neutrality [Additional file 3
]. In this case the hypothesis of divergent selection rather than that of drift is supported.