A cryptic species complex
Mayer and von Helversen [23
] showed that classification of bats based on morphological characteristics does not always correlate with mtDNA divergence. A similar pattern appears to be present here, where P. parnellii
on the mainland harbours at least four genetically distinct taxa. These taxa appear to be morphologically cryptic without obvious characters which differentiate them in the field, though statistically significant morphological differences, acoustic variation, and genetic divergence support their existence. The groups do not appear to match the distributions reported for sub-species [31
] (though see below). We have not included discrete character states in our morphological analysis (fur colour, banding pattern, etc.) as these are harder to quantify, but they may represent useful field characters for further investigation.
While any single uni-parentally inherited molecular marker has limitations (e.g. inability to assess hybridization), mtDNA can be a powerful tool for hypothesis generation in taxonomic research. In mammals, cytochrome b
has historically been employed for similar analyses (e.g. [7
]) though COI is favoured by DNA barcoding campaigns and has been employed extensively in Neotropical bats [8
] allowing for easy comparison of sequence divergence (in particular see Clare et al.
] for a review of COI sequence divergence in >9,000 individuals of 165 Neotropical bat species). Mitochondrial regions are frequently paired with one or more non-mitochondrial loci to test these hypotheses. Y-chromosome regions are ideal for comparison to mtDNA because they are fast evolving and non-recombining, but provide an exclusively paternal measure of gene flow. Particularly when male-biased gene flow is suspected (as is frequent in mammals; e.g. [60
]) y-chromosome DNA can provide a rigorous test of the patterns observed in mtDNA (see [9
]). Nuclear genes are also frequently used in phylogenetic analyses, but many are too slowly evolving for species level diagnosis, particularly when species are young. In addition, nuclear genes are bi-parentally inherited, raising problems of heterozygosity and recombination from multiple alleles. Nuclear genes can effectively support mitochondrial genes in phylogenetic reconstructions. In this case, the RAG2 region appears to provide limited but reliable support for a split between Groups (1+2) and (3+4) (Figure ) and thus supports a phylogenetic topology which unites these groups (Figure ). It should be noted that both the RAG2 and COI regions analyzed are protein-coding and thus subject to selection, however they have very different functions (immune system and electron transport chain respectively) resulting in different selection profiles. Additionally, neither region is known to be linked to any morphological trait, therefore they act as relatively independent markers.
In our analysis, all three genes provide evidence for a split between Groups (1+2) and Groups (3+4). The split between Groups 1 and 2 is not observed in the other two gene regions so we cannot evaluate the probability of male-biased gene flow but, as they are allopatric and acoustically distinct, the probability of hybridization is low. These two meet the criteria for the genetic species concept (GSC) advocated for mammals [7
]. The GSC evaluates species based on the Bateson-Dobzhansky-Muller model and permits small amounts of gene flow if the genetic groups are on independent evolutionary trajectories [7
]. Unlike the biological species concept, the GSC is applicable to allopatric populations and provides a framework for evaluating these groups. Finally, the y-chromosome supports the split between Groups 3 and 4. We were able to evaluate both mtDNA and the y-linked regions from these sympatric males and found little evidence of hybridization (n=1). This is lower than previously reported in other bat taxa, for example, Hoffman et al.
] proposed cryptic species in Uroderma bilobatum
and found two hybrids in 46 individuals and one potential F1, and in the European bats Myotis myotis
and M. blythii
introgression may be measured in 25% of individuals [61
]. The low level of hybridization measured here is acceptable under the GSC and is also generally permitted under a relaxed biological species concept.
Though the Dby
intron region has been amplified and sequenced in a wide variety of taxa [9
], some caution is required in interpreting the data. While the region is fast evolving it does not appear to evolve as fast as mitochondrial DNA. Additionally, like mtDNA, it may be subject to reduced variability from selective sweeps. An homologous region has been identified on the X-chromosome in some species though it is substantially divergent [62
]. We saw no evidence of co-amplification from the male X-chromosome and it is likely that it is either absent in P. parnellii
or that these primers preferentially bind to the target area, though no rigorous testing has been done to our knowledge. The family Mormoopidae contains two genera: Mormoops
, with three extant species, and Pteronotus
, with seven extant species [63
]. Of these, most are confined to Central America and the Antilles. P. parnellii
appears to be the oldest lineage in the genus, and perhaps the family [38
]. While all other species of Pteronotus
may have a Central American origin [38
] the origin of P. parnellii
is still unclear though our phylogeny suggests a South American origin (see below).
We employed fixed clock estimates of 2% and 5% to exploit the upper and lower calibration points commonly used for bat divergence time calculations in mammal mtDNA. More precise calibrations have been used in other bat taxa, but these are taxonomically specific and estimated for cytochrome b
. Given that COI is thought to evolve slightly slower than cytochrome b
], using exact calibrations calculated for cytochrome b
in bats would be inappropriate and we therefore employed minimum and maximum range bracketing. Previous calibrations have not been made for P. parnellii
or COI in bats, thus our divergence rate calculations are meant to approximate the upper and lower estimates of divergence time in these analyses and to minimize errors in the divergence time estimates if any specific calibration point was used. The similar estimates from the RAG2 region suggest these are appropriate boundaries.
Using fixed clock estimates (Figure ), Groups (1+2) and (3+4) diverged from a common ancestor 2.5- 6.1 MYBP during land bridge formation which ended ~3 MYBP [1
] and a similar estimate was recovered from the RAG2 region. Groups 1, 2, 3 and 4 all radiated during or after the rise of the Eastern Cordillera 2–3 MYBP [2
]. The most likely biogeographic scenario is that these groups have a South American origin as suggested by the phylogeny (a single invasion of Central America in the ancestor of Group 1) but invaded Central America during the great American interchange. The rise of the Eastern Cordillera effectively isolated the Central American population, giving rise to Group 1 in allopatry. However, our observations should be treated as preliminary. Additional samples from the Antilles region will be required to establish this definitively and to determine whether the invasion may have involved island hoping. Groups 2, 3 and 4 all exist within the Guyana Shield region which is one of the South American cratons [66
]. While no obvious geological event correlates to the divergence, specimens in Group 2 are associated with higher elevations within the Guyana Highlands [68
] though sampling here was minimal and the presence of the complex in Trinidad suggests more habitat variability. It is plausible that Group 2 has re-invaded South America either via the Antilles or a dispersal event over the Andes (additional sampling in Columbia, Venezuela, and the Antilles may resolve this). The divergence of Group 3 and 4 from a common ancestor occurred between 1.1 and 2.7 MYBP and both are associated with sympatric ranges at lower elevations within the Guyana shield region. While there is no contemporary pattern of geology or climatology that can provide a mechanism for this divergence, the region has experienced multiple climate shifts and continual forest change so divergence may have involved a past habitat restriction and vicariance event. Our sampling concentrates on the Central American and northern South American portions of the P. parnellii
range and includes individuals from ten countries. We have not sampled individuals in the southern extent of their range, particularly Bolivia, Peru, and Central Brazil, or populations in the Caribbean which includes distributions of four additional subspecies [31
]. Our intention was to sample areas encompassing all mainland subspecies however our results suggest that the distribution of these subspecies does not correspond with existing molecular divergences (see below). Additional sampling in these ranges will almost certainly uncover additional cryptic diversity and may help clarify the existence and ranges of these subspecies, their systematic status, and the origin and divergence of the complex.
Systematic considerations in the complex
The taxonomic status of P. parnellii
is exceedingly complex. Herd [31
] recognized nine subspecies; P. p. parnellii, P. p. gonavensis, P. p. portoricensis
and P.p. pusillus
in the Antilles and P. p. mexicanus, P. p. mesoamericanus, P. p. paraguanensis, P. p. fuscus
and P.p. rubiginosus
on the mainland. Of these, P. p. paraguanensis
was recently elevated to a species in Venezuela [64
]. The most important systematic question regarding the groups recognized in this study is whether any may be considered P. parnellii
sensu stricto. The type location for the taxon is Jamaica and is thought to encompass the subspecies P. p. parnellii
. While we have been unable to include a genetic sample from this location in our analysis, the Jamaican population is acoustically distinct. We have some limited morphological data from Jamaica (not included), which demonstrates that this population is comprised of much smaller individuals than those on the mainland. For example, the forearm measures in Jamaica are 53
mm±1.26 SD (range 43–55, n=58, S. Koenig pers. comm.) which is non-overlapping with any mainland population we have examined. Additionally, a comparison of our cranial measurements with those presented for subspecies in Smith [53
] indicate Group 1 best matches the measurements for P. p. mesoamericanus
with little overlap among our measurements and those for P. p. parnellii
. From this, we conclude that Group 1 corresponds with P. p. mesoamericanus
and an elevation of that name would be appropriate. The individuals in our remaining groups are even larger thus none corresponds morphologically to P. p
. Suggesting appropriate names for the remaining three groups is more difficult as their distribution does not correspond with subspecies distributions suggested by Herd [31
]. The remaining Antillean subspecies are even smaller than P. p. parnellii
and their distributions disjunct, therefore these names are unlikely to be appropriate; likewise the distribution of P. p. mexicanus
precludes a likely correspondence with Groups 2, 3 or 4. Of the remaining subspecies, P. p. paraguanensis
(hereafter referred to by its species status), P. p. fuscus,
and P. p. rubiginosus,
distributions from Herd [31
] do not correspond with our analysis. In Venezuela, Gutierrez and Molinari [64
] report that P. paraguanensis
is significantly smaller than either P. p. fuscus
or P. p rubiginosus
, with P. p. rubiginosus
the largest; however Gutierrez and Molinari [64
] agreed with Herd [31
] that only P. p. rubiginosus
is found east of the Rio Orinoco suggesting that these subspecies do not correspond with any specific group identified here unless the distributions are much larger than previously reported. Thus, while we are confident that none of our genetic groups can rightly be called P. parnellii
, we cannot, at this stage, suggest what names would be appropriate for Group 2, 3 and 4. In the interim, we suggest they be referred to as Pteronotus
species 2, 3 and 4 until more appropriate binomials can be established. Barring contradictory evidence, we further conclude that P. parnellii
be considered a taxon endemic to the Antilles. The most appropriate future analysis will be a molecular comparison between type material for these subspecies and the groups identified here with additional sampling in the Amazon and Antilles. See also Dávalos [69
] for a discussion of diversification in this family.
Divergence of morphological characters
Though forearm and cranium measurements are both correlated with overall size, we have treated them as independent in our statistical analysis as they have different applications. Cranial characters may only be measured accurately in extracted skulls (e.g. a museum collection), while forearm measurements can be easily obtained from live animals. As such, forearm length is frequently employed as a quick taxonomic field character. However, due to the high overlap in the forearm measurements of the groups, this would not be a reliable measure to separate these groups in the field. Our analysis indicates that the cranial characters advocated by Smith [53
] for the family can differentiate our groups with reasonable accuracy and may act as a useful tool for further investigation of subspecies throughout the entire range, though some overlap does exist thus molecular evidence is a more reliable definitive character. While the four groups can be largely discriminated by DFA (87% successful), there is a significant effect of latitude on the main principal component (PC1-size). Once we corrected for the latitudinal cline, DFA continued to show the groups separating out based on cranial morphological characters. Differences between the two sympatric groups, 3 and 4, became more apparent, while divergence between Groups 1 and 2 was not as obvious after correction (Figure ). This suggests that the morphological differences of Groups 3 and 4 to all other groups are not solely driven by the latitudinal cline, but this is less clear with Groups 1 and 2.
There was a strong effect of sex on both the size (PC1) and shape (PC2) of the skull suggesting some sexual dimorphism in all locations. The shape of the skull (PC2) did not show a latitudinal effect but did vary between Group 3 and all other groups, further suggesting that there are morphological differences among the groups that are not driven by latitude.
Acoustic divergences within and between groups
Our analysis suggests that there is significant acoustic variation in these groups, creating distinct calls which can differentiate among Groups 1, 2 and 3/4. There is a steady decrease in the frequency of the call with decreasing latitude which also corresponds to an increase in the size of many morphological components.
Interestingly, enough variation is present to further distinguish between some regions within Group 1. We cannot conclude whether there are acoustic differences between the sympatric Groups 3 and 4, as we are limited by relatively few recordings from the region and an inability to separate free-flying bats in field recordings. Further sampling will be required in combination with molecular analysis to determine whether these two groups share an identical echolocation call or whether an additional acoustic pattern exists in this region. If it does, it could prove a powerful tool in identifying the two groups in the field and this is a clear goal for future field studies (sound files available from the authors on request). In addition to echolocation calls, it would be extremely valuable to examine the role of non-echolocation vocalizations among these groups. The role of “social calls” in intra- and interspecific recognition has not been well documented; however the degree of variation in social call repertoire is extensive and may play an important role in mate recognition.
The acoustic variation observed here is not consistent with the pattern of harmonic hopping observed by Kingston and Rossiter [29
] in the R. philippinensis
but is similar to divergence patterns in the H. bicolor
]. Individuals are unlikely to recognize the calls of other groups if they come into auditory contact. Furthermore, the divergence between calls is small and, as in H. bicolor
, probably not sufficient to support ecological resource partitioning. As in H. bicolor
, social character displacement and selection for non-interference may have played a large role in the diversification of the P. parnellii
complex, particularly in northern South America where there are no obvious barriers to prevent contact between all three groups. If the Central American group has diverged in allopatry, it is likely that the variation in the constant frequency component has arisen by drift. It is not known what threshold of variation is used for social recognition in P. parnellii,
but estimates from Old World bats suggest that the constraints of high-duty cycle echolocation limit the possible variation more than in low-duty cycle species [33
]. Future research on inter- and intraspecific call divergence and the social aspects of call recognition and discrimination will prove interesting and may support an allopatric or sympatric origin of these groups.
Echolocation: a mode of reinforcement and speciation?
Theoretically, for speciation to occur despite substantial gene flow, reinforcement and a strong mate choice character are almost certainly required. In bats, non-morphological traits may be very important in speciation. In some potential cryptic species complexes, for instance Saccopteryx bilineata
(family Emballonuridae) [9
], olfactory cues may be particularly important. In bats that rely on echolocation, divergence between echolocation calls may reduce signal competition and influence prey selection [29
] but also change mate recognition [71
], creating an almost instantaneous form of pre-zygotic isolation. For reinforcement to directly influence speciation, a genetic association (linkage) between the ecological trait under selection and the mating signal is required [70
], a situation which should be disrupted by recombination [73
]. This problem is averted if the diverging trait is controlled by the same genetic loci which become fixed in the diverging populations (one allele model of Felsenstein [73
]). This is particularly effective when the characters governing both pre- and post-zygotic isolation are the same, so-called “magic traits” [74
], and direct selection has pleiotropic effects [72
]. Echolocation may meet both of these criteria [75
] and, in this context, changes in echolocation could rapidly lead to speciation.
Ecological character displacement (selection for increased ecological niche separation which indirectly influences mate choice [72
]) may also speed the process of differentiation. The divergence of echolocation call design may start as selection against signal interference between populations [29
] but may be taken over by selection for mate recognition. In this scenario, speciation can be swift in the presence of gene flow whether selection is initially acting on niche specialization or mate recognition [70
In our analysis, Groups 3 and 4 currently occupy sympatric distributions (Figure ). Elmer and Meyer [77
] outline four “gold standard” criteria for the hypothesis of sympatric speciation: 1) sympatric distributions, 2) a reciprocally monophyletic relationship between the taxa (sister species) with respect to others in the complex, 3) reproductive isolation and, 4) a setting where allopatric divergence is unlikely. Groups 3 and 4 appear to meet at least the first three criteria, but several lines of evidence could refute sympatric speciation. First, the relationship between individuals from the Brazilian Amazon and French Guiana is unknown. If a group from these areas is sister to Group 3 or Group 4 (raising the possibility of historical allopatric ranges), the hypothesis of sympatric speciation may be refuted. Additionally, pre-zygotic isolation should develop before post-zygotic isolation in cases of sympatric speciation [78
]. The presence of only one hybrid suggests that hybrids are rare but does not provide evidence for the cause. Searching for evidence of F2 individuals could illuminate this [78