# Coexistence of two sympatric cryptic bat species in French Guiana: insights from genetic, acoustic and ecological data

Ondine Filippi-Codaccioni, Marie-Pauline Beugin, Damien de Vienne, Elodie Portanier, David Fouchet, Cécile Kaerle, Lina Muselet, Guillaume Queney, Eric Petit, Corinne Regis, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

Ondine Filippi-Codaccioni, Marie-Pauline Beugin, Damien de Vienne, Elodie Portanier, David Fouchet, et al. Coexistence of two sympatric cryptic bat species in French Guiana: insights from genetic, acoustic and ecological data. BMC Evolutionary Biology, 2018, 18 (1), 10.1186/s12862-018-1289-8. doi:10.1186/s12862-018-1289-8

HAL Id: hal-01936895

https://hal.science/hal-01936895

Submitted on 22 Nov 2019

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH ARTICLE

Open Access

<figure><img src="image_3.png" /></figure>

CrossMark

# Coexistence of two sympatric cryptic bat species in French Guiana: insights from genetic, acoustic and ecological data

Ondine Filippi-Codaccioni<fnref n="1" />,<fnref n="2" />,³, Marie-Pauline Beugin<fnref n="1" />,<fnref n="2" />,⁴, Damien M. de Vienne<fnref n="1" />,<fnref n="2" />, Elodie Portanier<fnref n="1" />,⁵,⁶, David Fouchet<fnref n="1" />,<fnref n="2" />, Cecile Kaerle⁴, Lina Muselet⁴, Guillaume Queney⁴, Eric J. Petit⁷, Corinne Regis<fnref n="1" />,<fnref n="2" />, Jean-Baptiste Pons<fnref n="1" />,<fnref n="2" /> and Dominique Pontier<fnref n="1" />,<fnref n="2" />*

## Abstract

**Background:** The distinction between lineages of neotropical bats from the *Pteronotus parnellii* species complex has been previously made according to mitochondrial DNA, and especially morphology and acoustics, in order to separate them into two species. In these studies, either sample sizes were too low when genetic and acoustic or morphological data were gathered on the same individuals, or genetic and other data were collected on different individuals. In this study, we intensively sampled bats in 4 caves and combined all approaches in order to analyse genetic, morphologic, and acoustic divergence between these lineages that live in the same caves in French Guiana.

**Results:** A multiplex of 20 polymorphic microsatellite markers was developed using the 454-pyrosequencing technique to investigate for the first time the extent of reproductive isolation between the two lineages and the population genetic structure within lineages. We genotyped 748 individuals sampled between 2010 and 2015 at the 20 nuclear microsatellite loci and sequenced a portion of the cytochrome c oxidase I gene in a subset of these. Two distinct, non-overlapping haplogroups corresponding to cryptic species *P. alitonus* and *P. rubiginosus* were revealed, in accordance with previous findings. No spatial genetic structure between caves was detected for both species. Hybridization appeared to be quite limited (0.1–4%) using microsatellite markers whereas introgression was more common (7.5%) and asymmetric for mitochondrial DNA (mtDNA).

**Conclusions:** The extremely low rate of hybridization could be explained by differences in life cycle phenology between species as well as morphological and acoustical distinction between sexes in one or the other species. Taken together, these results add to our growing understanding of the nature of species boundaries in *Pteronotus parnellii*, but deserve more in-depth studies to understand the evolutionary processes underlying asymmetric mtDNA introgression in this group of cryptic species.

**Keywords:** *Pteronotus cf. parnellii*, Microsatellite, Mitochondrial DNA, Phylogeny, Echolocation, Population genetics, Asymmetric introgression, Neotropics

<figure><img src="image_4.png" /></figure>

BMC

© The Author(s). 2018 **Open Access** This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

<footnote>* Correspondence: dominique.pontier@univ-lyon1.fr</footnote>

<footnote n="1">¹University Lyon, Université Lyon 1, CNRS, Laboratoire de Biométrie et Biologie Evolutive UMR5558, F-69622 Villeurbanne, France</footnote>

<footnote n="2">²Université de Lyon, LabEx Ecofect, Nadine Cizaire, 92 rue Pasteur, CS 30122 69361 Lyon Cedex 07, France</footnote>

<footnote>Full list of author information is available at the end of the article</footnote>
## Background

Much natural diversity is morphologically hidden [1]. The detection of cryptic species, i.e., genetically divergent species previously classified as a single species due to morphological similarity [1, 2], has significantly increased for all major terrestrial and aquatic taxonomic groups and across biogeographical regions [2, 3] through the use of large-scale DNA sequencing approaches such as DNA barcoding [4, 5]. The discovery of this cryptic diversity has had profound implications for both evolutionary theory and future conservation decisions (see [6–8] for an example involving bats), especially in threatened ecosystems for which biodiversity has likely been underestimated. Morphologically similar species can indeed vary in geographic distribution and ecological requirements, and thus the conservation status among cryptic species belonging to the same species complex can differ (see [9] for an example in Hanuman langur *Semnopithecus entellus*).

Although these cryptic species can usually be identified due to mitochondrial DNA (mtDNA) differences, reproductive isolation is only scarcely investigated. There is indeed an implicit assumption that the molecular genetic divergence among these cryptic species can be taken as a surrogate for reproductive isolation [10–14]. While substantial mtDNA differences between similar or identical-looking individuals can reveal cryptic species, caution is needed as it does not tell us whether the two (or more) cryptic species are truly separate species (see [15] for an illustration of the danger of using only a limited amount of DNA to draw conclusions about evolutionary history). It may instead provide only a glimpse into the evolutionary past of animals that are now one population of interbreeding individuals. To confirm their reproductive isolation, it is necessary to use nuclear markers, which, unlike mtDNA, are inherited from both parents. This allows the question of whether interbreeding occurs or not between the cryptic species to be answered. Another important related question is then how truly cryptic species maintain species cohesion in an area of sympatry. Subtle isolating mechanisms, such as ecologic or recognition systems (e.g. [16–19]), may have evolved, permitting the co-existence of sympathically living cryptic species. An example of this is *Drosophila paulistorum* for which semi-species are morphologically similar, but have different courtship song patterns [20].

Cryptic species are fairly common in echolocating bats [21]. One of the first major discoveries was that the most widespread European bat species, the common pipistrelle (*Pipistrellus pipistrellus*) had [22] in fact been hiding two cryptic species: *P. pipistrellus* and the soprano pipistrelle (*P. pygmaeus*). Despite being morphologically

similar marked divergence in ecological requirements has been evidenced (e.g. [23]). The neotropical insectivorous bat *Pteronotus parnellii* (family Mormoopidae, subgenus *Phyllodia*) also comprises a particularly remarkable example of hidden diversity. Currently, nine cryptic species widespread in Middle America and the Caribbean [24–27] have been recognized within the *P. parnellii* complex. Phylogenetic analyses of mitochondrial genomes suggest that a split occurred in their maternal lineages ~1.1–2.8 millions years ago, while the entire complex would have shared a common ancestor ~2.5–6.1 millions years ago [24].

In French Guiana, the presence of two sympatric groups referred to as *Pteronotus* sp3 and sp4 has been proposed (sensu [24, 27]) using several tools: genetics – using mitochondrial molecular markers such as the cytochrome oxidase 1 (CO1), the cytochrome *b* (*Cyt b*), or the y-linked *Dby* genes [24, 25, 28] –, morphology [24, 28], and bioacoustics [28, 29]. Recently these two groups have been diagnosed as two distinct species, *P. rubiginosus* corresponding to *P*. sp4 [27] while *P*. sp3 has been named *P. Alitonus* [27]. However, these studies were limited by the number of samples when genetic and acoustic or morphological data were gathered on the same individuals, or because genetic and other data were collected on different individuals. Morphological studies failed to provide diagnostic characters to distinguish between members of this complex, but differences occur in their echolocation call frequency [28, 29]: *P. rubiginosus* uses a peak echolocation frequency of 53 kHz, *P. Alitonus* 59 kHz. It is evident that we still know very little about these newly discovered groups in terms of their preferred habitats, prey and roosting places, their breeding ecology or their population structure and potential for interbreeding - i.e., not enough information has been gathered as yet to formally describe them as new species.

In this study we have investigated colonies of *Pteronotus* sp. roosting in four caves in French Guiana during the period 2010 to 2015. We combined different approaches in order to investigate genetic, morphological and acoustic divergence between the two groups. For the first time we tested whether the identified groups in the caves interbreed using autosomal microsatellite markers and searching for hybrids in the sampled populations. Lastly, we explored the spatial genetic structure of *Pteronotus* sp., notably whether caves affect the genetic structure of the population of *Pteronotus* sp. For the DNA studies we used the 5' half of the mitochondrial CO1 and developed a panel of 20 specific microsatellites in one multiplex using the 454-pyrosequencing technique. The study's aim was to make comprehensive genetic comparisons that would reveal relationships between the two groups and to understand their microevolutionary
history - which might include hybridization and introgression patterns.

# Results

## Molecular analysis

### Species identification

The phylogenetic reconstruction (Fig. 1, Additional file 1: Table S1) showed that the 80 *Pteronotus* individuals sampled in the present study display two highly supported groups, corresponding to *P. alitonus* (hereafter *Pteronotus*

B) and *P. rubiginosus* (hereafter *Pteronotus A*) previously evidenced by [24, 27, 28].

Divergence of COI sequences was very small within the A and B groups (0.15 and 0.77% on average, respectively) but much larger between the groups (5.72%).

### Nuclear genetic analysis

We built a multiplex of 20 microsatellite markers (Additional file 1: Table S2). These 20 microsatellite markers were successfully amplified for 748 bat samples with a mean amplification rate of 98.7%. All markers were

<figure><img src="image_5.png" /></figure>
polymorphic, with the number of alleles ranging between 2 and 17.

Consistently with the phylogenetic reconstruction, the STRUCTURE analysis over 20 microsatellites combined with Evanno's method detected two sharply differentiated clusters (Additional file 1: Table S3). Among the 748 genotyped individuals, 325 belonged to *Pteronotus* A while 423 individuals belonged to *Pteronotus* B (see Table 1). The allelic diversity at the 20 microsatellite loci was comparable in both groups, but not identical (Table 2). The two groups were both present in all four caves sampled.

One locus showed significant signs of null allele frequency (PP08, *p* = 0.04) and four others appeared to be in Hardy-Weinberg disequilibrium for at least one of the two identified species (PP10 and PP16 for *Pteronotus* A, PP22 and PP24 for *Pteronotus* B, *p* < 10⁻⁴). We did not find any loci in linkage disequilibrium. Consequently, we ran STRUCTURE using 15 loci to redefine the two clusters, and considered these 15 markers in subsequent analyses. The reduction of the number of markers did not change the assignation of individuals to cluster A and B. Over these 15 microsatellite markers, the two species presented similar mean number of alleles per locus: 6.71 for *Pteronotus* A and 7.82 for *Pteronotus* B.

The two *Pteronotus* species were significantly differentiated with a *F*<sub>ST</sub> value of 0.139 [0.104–0.171]<sub>95%</sub>. Within species, only *Pteronotus* A males displayed significant genetic differentiation (global *F*<sub>ST</sub> values: 0.002, *p* = 0.007) and structure between caves (only significant between PA and SC caves; see pairwise *F*<sub>ST</sub> values, Additional file 1: Table S4). Additionally, we did not find any significant sex-biased dispersal (*p* > 0.05).

## Hybrid detection

Two thresholds, TP1 and TP2, were determined to identify hybrids based on the 15 microsatellite markers. They corresponded to the lowest *q*-value reached by simulated parents in each cluster. For each genotyped individual, either the mean *q*-value (conservative approach) or the

<table><caption>**Table 1** Number of *Pteronotus* A and B genotyped at each field session</caption><thead><tr><th rowspan="2">Year</th><th rowspan="2">Month</th><th colspan="3">Pteronotus A</th><th colspan="3">Pteronotus B</th><th rowspan="2">n(B)</th><th rowspan="2">n</th></tr><tr><th>F</th><th>M</th><th>NS</th><th>F</th><th>M</th><th>NS</th></tr></thead><tbody><tr><td>2010</td><td>August</td><td>0</td><td>1</td><td>0</td><td>1</td><td>3</td><td>4</td><td>0</td><td>7</td><td>8</td></tr><tr><td>2010</td><td>September</td><td>24</td><td>97</td><td>0</td><td>121</td><td>29</td><td>70</td><td>2</td><td>101</td><td>222</td></tr><tr><td>2011</td><td>July</td><td>7</td><td>63</td><td>0</td><td>70</td><td>26</td><td>37</td><td>0</td><td>63</td><td>133</td></tr><tr><td>2012</td><td>July</td><td>17</td><td>32</td><td>0</td><td>49</td><td>48</td><td>26</td><td>0</td><td>74</td><td>123</td></tr><tr><td>2015</td><td>October</td><td>24</td><td>54</td><td>6</td><td>84</td><td>71</td><td>105</td><td>2</td><td>178</td><td>262</td></tr><tr><td>Total</td><td></td><td>72</td><td>247</td><td>6</td><td>325</td><td>177</td><td>242</td><td>4</td><td>423</td><td>748</td></tr></tbody></table>

F Females, M Males, NS Not sexed, n(A) number of *Pteronotus* A, n(B) number of *Pteronotus* B, n total number of *Pteronotus*

<table><caption>**Table 2** Characteristics of the twenty microsatellite markers developed for *Pteronotus parnellii*</caption><thead><tr><th rowspan="2">Locus</th><th rowspan="2">Fst</th><th colspan="4">Pteronotus A</th><th colspan="4">Pteronotus B</th></tr><tr><th>Fis</th><th>Na</th><th>Ho</th><th>He</th><th>Fis</th><th>Na</th><th>Ho</th><th>He</th></tr></thead><tbody><tr><td>PP01</td><td>0.23</td><td>-0.095</td><td>4</td><td>0.465</td><td>0.424</td><td>0.049</td><td>7</td><td>0.730</td><td>0.768</td></tr><tr><td>PP02</td><td>0.17</td><td>-0.007</td><td>4</td><td>0.477</td><td>0.473</td><td>0.021</td><td>4</td><td>0.280</td><td>0.285</td></tr><tr><td>PP03</td><td>0.18</td><td>-0.070</td><td>3</td><td>0.465</td><td>0.434</td><td>-0.030</td><td>4</td><td>0.752</td><td>0.729</td></tr><tr><td>PP05</td><td>0.11</td><td>0.068</td><td>4</td><td>0.474</td><td>0.508</td><td>-0.006</td><td>5</td><td>0.711</td><td>0.706</td></tr><tr><td>PP08*</td><td>0.41</td><td>0.061</td><td>3</td><td>0.315</td><td>0.335</td><td>0.330</td><td>8</td><td>0.513</td><td>0.764</td></tr><tr><td>PP09</td><td>0.15</td><td>0.012</td><td>5</td><td>0.615</td><td>0.622</td><td>0.020</td><td>8</td><td>0.834</td><td>0.851</td></tr><tr><td>PP10<sup>#</sup></td><td>0.42</td><td>0.310</td><td>11</td><td>0.570</td><td>0.821</td><td>-0.002</td><td>11</td><td>0.839</td><td>0.837</td></tr><tr><td>PP12</td><td>0.09</td><td>0.004</td><td>3</td><td>0.108</td><td>0.108</td><td>0.062</td><td>4</td><td>0.374</td><td>0.398</td></tr><tr><td>PP13</td><td>0.13</td><td>NA</td><td>1</td><td>0.000</td><td>0.000</td><td>-0.038</td><td>2</td><td>0.260</td><td>0.250</td></tr><tr><td>PP14</td><td>0.16</td><td>0.044</td><td>10</td><td>0.674</td><td>0.704</td><td>-0.010</td><td>11</td><td>0.740</td><td>0.732</td></tr><tr><td>PP15</td><td>0.24</td><td>-0.043</td><td>9</td><td>0.834</td><td>0.798</td><td>0.011</td><td>8</td><td>0.707</td><td>0.714</td></tr><tr><td>PP16<sup>#</sup></td><td>0.19</td><td>0.110</td><td>4</td><td>0.482</td><td>0.539</td><td>0.003</td><td>8</td><td>0.861</td><td>0.862</td></tr><tr><td>PP17</td><td>0.15</td><td>-0.010</td><td>7</td><td>0.738</td><td>0.730</td><td>0.014</td><td>9</td><td>0.757</td><td>0.767</td></tr><tr><td>PP18</td><td>0.07</td><td>0.004</td><td>13</td><td>0.791</td><td>0.793</td><td>-0.022</td><td>17</td><td>0.856</td><td>0.837</td></tr><tr><td>PP21</td><td>0.09</td><td>0.030</td><td>4</td><td>0.517</td><td>0.532</td><td>0.013</td><td>8</td><td>0.770</td><td>0.779</td></tr><tr><td>PP22<sup>#</sup></td><td>0.08</td><td>-0.008</td><td>5</td><td>0.680</td><td>0.674</td><td>0.017</td><td>5</td><td>0.530</td><td>0.538</td></tr><tr><td>PP23</td><td>0.08</td><td>-0.024</td><td>4</td><td>0.683</td><td>0.666</td><td>-0.015</td><td>6</td><td>0.751</td><td>0.740</td></tr><tr><td>PP24<sup>#</sup></td><td>0.61</td><td>-0.043</td><td>6</td><td>0.441</td><td>0.422</td><td>0.211</td><td>4</td><td>0.186</td><td>0.236</td></tr><tr><td>PP26</td><td>0.06</td><td>0.020</td><td>5</td><td>0.397</td><td>0.404</td><td>0.004</td><td>6</td><td>0.621</td><td>0.622</td></tr><tr><td>PP27</td><td>0.01</td><td>0.009</td><td>4</td><td>0.443</td><td>0.446</td><td>-0.003</td><td>3</td><td>0.331</td><td>0.329</td></tr></tbody></table>

Fst and Fis values correspond to Weir and Cockerham's estimates. The number of alleles is given in the Na column and the observed and expected heterozygosity are provided in the Ho and He columns, respectively. Loci with a star (*), a mu (μ) or a sharp (#) correspond to loci with significant null alleles frequency, loci in Hardy-Weinberg disequilibrium in species A or B, respectively. All these loci were excluded in population genetics analysis

lower bound (relaxed approach) of the credibility interval returned by STRUCTURE was compared to those thresholds to determine whether the individual was hybrid. The simulation procedure led us to consider 0.90 and 0.76 as values for TP1 and TP2. Accordingly, one putative hybrid was detected with the “conservative” approach while 30 putative hybrids were detected with the “relaxed” approach. A further simulation study revealed that the probability to identify parental individuals as hybrids, either with the conservative or relaxed approach, was null with these thresholds used in a configuration where the number of hybrids is lower than the number of parental forms. No hybrid was detected using the computer program NEWHYBRIDS.

Among individuals for which we have sequenced mitochondrial DNA (n = 80), all 42 individuals assigned to *Pteronotus* B with nuclear DNA had *Pteronotus* B mitochondrial sequences. In contrast, 6 individuals out of the 38 assigned to *Pteronotus* A with nuclear DNA were found to have *Pteronotus* B mitochondrial sequences (see Fig. 1).
## Morphology

Only 23 of the 30 individuals identified as putative nuclear hybrids with the relaxed approach have morphometric data. Excluding these 23 individuals to be conservative as well as the 6 mtDNA hybrids (see above), *Pteronotus A* (*n* = 166) had a mean forearm length of 64.37 (95% CI [64.13;64.56]) vs 62.06 (95% CI [61.90;62.22]) for *Pteronotus B* (*n* = 300) (Fig. 2a). This difference between both groups was found significant even after correcting for the potential confounding effect of sex (F = 293.17, df = (1462), *p* < 10⁻¹⁶).

The global comparison between sexes showed non-homogeneous forearm length in at least one group (F = 8.38, df = (2462), *p* = 2.7 × 10⁻⁴). The post-hoc analysis revealed that the effect of sex was significantly different between the two groups (F = 5.90, df = (1462), *p* = 1.5 × 10⁻²) and was only found significant in group B (*t* = 4.07, df = 298, *p* = 6 × 10⁻⁵, forearms being longer in females; in group A: *t* = -0.38, df = 164, *p* = 0.73).

Due to the absence of F1 and to a small number of morphological data for putative hybrids (Additional file 1: Figure S5), we could not compare hybrid and parental morphological features.

## Echolocation

Among the 262 (73 *Pteronotus* A, 171 *Pteronotus* B, 13 putative nuclear hybrids identified with the relaxed approach and 5 mtDNA hybrids) microphone recorded

individuals, 146 (of which 15 are hybrids) were precisely measured for the Frequency of Maximal Energy (FME) with BatSound (FME dataset). Excluding the 15 hybrid individuals of the FME dataset and the four outliers from group B (pure individuals), *Pteronotus A* (*n* = 60) emitted at a call frequency situated around 53 kHz (53.73 kHz, 95% CI [53.51;53.96]) while *Pteronotus B* (*n* = 67) emitted at a call frequency situated around 59 kHz (59.64 kHz, 95% CI [59.49; 59.78]) (Fig. 2b). This difference between both groups was found significant even after correcting for the potential confounding effect of sex (F = 2080.8, df = (1116), *p* < 10⁻¹⁶).

The global comparison between sexes showed non-homogeneous echolocation frequency in at least one group (F = 6.32, df = (2116), *p* = 2.4 × 10⁻³). The post-hoc analysis revealed that the effect of sex was significantly different between the two groups (F = 10.68, df = (1116), *p* = 1.4 × 10⁻³) and was only found significant in group A (*t* = 2.92, df = 53, *p* = 5.2 × 10⁻³, frequencies being higher in females; in group B: *t* = -1.37, df = 63, *p* = 0.18).

Again, due to the absence of F1 and to a small number of acoustic data for putative hybrids (see Additional file 1: Figure S5), we could not compare hybrid and parental acoustic features (see Additional file 1: Figure S5).

## Reproductive periods

Regarding *Pteronotus A*, only two females (over 24) were found to be pregnant in July with one in MA and one in PA, and one lactating female was found in October (over 24). No juveniles were found during the different field periods. Concerning *Pteronotus B*, pregnant (27 over 74) and lactating (15 over 74) females were found in July and post-lactating (33 over 103) females between late August and October in the MA and PA caves (Additional file 1: Table S6). Ten juveniles were captured in September in the PA cave.

## Discussion

In this study, molecular, acoustic, and morphological differentiation has been analyzed in *Pteronotus* sp. collected from four caves in French Guiana. Using COI - a mtDNA marker commonly employed for species barcoding - our phylogenetic reconstruction of *Pteronotus* reveals two distinct, non-overlapping haplogroups corresponding to cryptic species *P. alitonus* and *P. rubiginosus*, in accordance with previous findings (see [24, 25, 27, 28]). The smallest divergence between *P. alitonus* and *P. rubiginosus* is still almost 3 times larger (2.7) than the largest within group divergence – a typical value for interspecific differences in mammals [5, 30]. We found good agreement between mitochondrial and nuclear microsatellite markers, i.e. the two mtDNA lineages encompass two distinct nDNA genetic clusters. The *F*<sub>st</sub> value of 0.139 between *P*.

<figure><img src="image_6.png" /></figure>
*alitonus* and *P. rubiginosus* indicates a strong population genetic differentiation in the 15 microsatellite markers used.

## Reproductive season

*P. alitonus* and *P. rubiginosus* appear to have an overlapping reproductive season. However the very low proportions of pregnant and lactating *P. rubiginosus* females compared to *P. alitonus*, together with the absence of juveniles, suggest a small shift in the reproductive season between the two species. Alternatively, we cannot exclude the possibility that *P. rubiginosus* uses breeding sites other than caves (e.g., buildings, underside of bridges), potentially reducing the competition with *P. alitonus* for food resources. The identification of a colony with pregnant (9 over 24 females) and lactating (14 over 24 females) females of *P. rubiginosus* in the roof of a village house in June 2017 may confort this second hypothesis. It was also suggested that *P. alitonus* prefers to forage in more dense forest than *P. rubiginosus* [27]. Such spatial separation of foraging habitats may reduce competition between the two species thereby explaining their co-existence in the same area.

## Phenotypic differentiation

When comparing forearm length between the two cryptic species, we found a broad overlap between measurements showing that the species distinctions did not precisely match the genetic differences (see also [24, 28, 31]). We also showed that *P. rubiginosus* is slightly larger than *P. alitonus*, as previously reported by de Thoisy et al. [28] in the same area, in López-Baucells et al. [31] in a diversity of habitats in the Central Brazilian Amazon. Interestingly, Pavan et al. [27] showed that this phenotypic difference is enhanced in areas where the two species occur in symmetry – a result which is expected in species undergoing character displacement [27, 32]. Furthermore, in our study we observed that females have a slightly larger forearm length than males in *P. alitonus*. Sexual dimorphism had already been demonstrated by Clare et al. [24] in skull size, with males having larger ones. Such subtle size differences between males and females suggest that *P. alitonus* may vary in sexual selection intensity over its distribution area – as it was observed for *Rhinolophus ferrumequinum* ([33]; see also [34]). However, such dimorphism may be closely linked with habitat types and thus with natural selection rather than sexual selection. For example, the existence of habitat-specific sexual dimorphism has been reported in *Anolis* lizards [35]. Thus, females of *P. alitonus* may have evolved larger forearms than males because it may confer reproductive advantages [34, 36] depending on their environmental conditions (e.g., food availability, competitors).

As previously reported for French Guiana and in the Central Amazon [24, 28, 29, 31], the two cryptic species correspond to two entirely distinct phonic types – one displaying frequencies of maximum energy around 53 kHz (*P. rubiginosus*), and the other around 59 kHz (*P. alitonus*). Jiang et al. [37] and Lin et al. [38] suggested that variations of 5–7 kHz do not impact the ability to detect prey and thus, should not affect resource use. The two species occupy the same caves, which are a limiting resource in the studied area, suggesting that they coexist without major ecological competition either for roost caves or for prey, as already discussed for other cryptic species such as *Rhinolophus mehelyi* and *R. euryale* [39]. We thus hypothesize that the different acoustic calls may have evolved to facilitate intraspecific and interspecific communication/recognition, rather than to facilitate resource partitioning. Though unlikely [40], a neutral explanation for this observation cannot entirely be ruled out.

Furthermore, *P. rubiginosus* - which have longer forearms - use lower peak frequency of echolocation calls (53 kHz) than *P. alitonus* (59 kHz), in accordance with the mass-signal frequency allometry (e.g. [41, 42]). However, this size-dependent effect in call frequencies cannot explain the inter-sex variation observed in *Pteronotus* sp. While *P. alitonus* females have slightly larger forearm length than males both genders emit at the same call frequency. On the contrary, although *P. rubiginosus* bats do not exhibit dimorphism in forearm length between the sexes, we did find that females emit at higher frequencies than males. By contrast with *P. alitonus*, sexual selection may act more on echolocation call frequencies in males than body size in *P. rubiginosus*. Such variation in call frequencies between the sexes has not been reported in previous studies on *Pteronotus* species [24, 28]; however, the number of individuals sampled was limited and no information about the sampling period was given. Although we have no idea if males were active or not during the sampling periods for both species (handling animals giving rise to fluctuating testis size we had no reliable indicator), we suggest that certain bat species - at least *P. rubiginosus* - can change their call frequency during the mating period. For example, Grilliot et al. [43] have observed in *Eptesicus fuscus* that bats can modulate their call according to functional context, making them monomorphic for activities such as foraging, but dimorphic during mating activity. Acoustics may thus play an important role in sex recognition in *Pteronotus* and/or indicate some aspects of male condition or quality at least in *P. rubiginosus* (as highlighted in *R. mehelyi* [44]). Such differences in calls may limit mating between species in *Pteronotus*. Some studies have suggested that changes in echolocation frequency (for instance *Rhinolophus philippinensis*) are associated with assortative mating, and ultimately reproductive isolation
and speciation, regardless of external barriers to gene flow [18].

## Hybridization

We found evidence of limited hybridization between the two species using the 15 microsatellite markers: only 1 individual (conservative approach) and 30 individuals (relaxed approach) out of 748 showed signs of hybridization. On the long term however, such events leave traces that are best observed in the mitochondrial genome: 6 individuals out of 80 contain COI sequence of the other species. The most striking finding is that in all cases mtDNA introgression has occurred asymmetrically, from *P. alitonus* to *P. rubiginosus* (6 out of 38 *P. rubiginosus* have mtDNA from *P. alitonus* while all 42 *P. alitonus* contain the mtDNA of *P. alitonus*). Our proportion of mtDNA hybrids in *P. rubiginosus* (0.158) is higher than that previously reported by Clare et al. [24] who identified - in the lowlands of the Guyana Shield - only one hybrid out of 61 specimens (35 *P. alitonus* and 26 *P. rubiginosus*), but similarly, this hybrid recovered the mitochondrial DNA of *P. alitonus* (giving a proportion of hybrids of 0.04). Similar results have been reported in bats between the two sibling bat species *Myotis myotis* and *Myotis blythii* [45], and *Rhinolophus sinicus septentrionalis* and *R. s. sinicus* [46]. Many different scenarios have been suggested to explain biased hybridization in mtDNA such as sex-biased dispersal, asymmetry in mate choice, differential production of offspring, demographic dynamics of local and colonizing species, differential selection, or a combination of these effects [47]. Differentiating between the possible causes is complex but some of them can be tentatively discussed to explain this unexpected asymmetry of mtDNA lineages between *P. alitonus* and *P. rubiginosus*.

Based on our analysis of genetic data at the scale of the colony, sex-biased dispersal is an unlikely explanation for the observed asymmetry since there is no evidence that males and females - of both species - disperse at different distances.

The observed unidirectional introgression of mtDNA could alternatively reflect a propensity of mating between *P. rubiginosus* males and *P. alitonus* females but not the opposite (*P. rubiginosus* females and *P. alitonus* males). Such asymmetry in mating between the sexes has been observed in insects [48] and in birds [49] but to our knowledge very few data regarding mating behavior are available for bats (see, e.g., [50–52]) and no data exist in particular for our two species. For example Bogdanowicz et al. [53] proposed that swarming behavior at swarming sites, where high number of bats belonging to several species meet, is an important factor that could explain hybridization among the three bat species they studied. Considering that both species coexist in the same habitats, form mixed maternity colonies in the caves, that their

reproductive season seems to overlap, that they differ in call frequency, as well as the extremely low rate of hybridization on microsatellite markers, we may consider asymmetry in mating between the sexes as anecdotal. Phenotypic differences between *P. alitonus* and *P. rubiginosus*, both in size and regarding the presence or absence of sexual dimorphism (see also [24, 27, 28]), could help species in the recognition of conspecifics and participate in the limitation of hybridization. As already mentioned, recognition between sexual partners can also be achieved through acoustic cues such as echolocation. This acoustic dimorphism may also help avoiding hybridization by erroneous recognition of conspecifics.

Finally, the asymmetry of introgression could also arise if the mtDNA of *P. alitonus* has some selective advantage that promotes its introgression in an alternative environment (a scenario proposed to explain the spread of mtDNA in some hare populations, [54] or in a foreign genetic background [55]. The mitochondrial genome plays a central role in cellular energy provision and can affect different life-history traits like life span and fertility, and/or behavior such as activity and exploration (e.g. [56]). An increasing number of studies have indicated that mtDNA appears to be under selection (e.g. [57–59]) but arguments in favor of this hypothesis are, as yet, far from being established for *Pteronotus*. That asymmetric hybridization of mtDNA occurs at a large geographical scale (Guyana Shield) but at different frequencies, can however give some support to this hypothesis.

## Genetic structuring of populations

For the first time, these two cryptic species have been studied using autosomal microsatellites. These new markers also allowed us to shed light on the genetic patterning of the populations. In our dataset, the level of genetic differentiation is an order of magnitude higher between species than within species ($F_{ST}$ = 0.13 versus $F_{ST}$ = 0.002, respectively), suggesting the existence of high gene flow between caves for the two species. Thus the cave does not appear to be a structuring unit as the populations of *P. alitonus* and *P. rubiginosus* seem to form two unique populations covering distances up to 80 km (from the Nouragues region to the Kaw Mountains). One possible explanation for this lack of genetic structure is that the distance between caves is lower than the dispersal distance of both species. Thus, the structuring unit may be an ensemble of caves rather than each cave.

## Conclusions

In agreement with previous studies, our genetic, morphological, and acoustic results justify the existing classification
of *Pteronotus* into two species - *P. alitonus* and *P. rubiginosus* - in French Guiana. Much work remains to be done to increase our knowledge of the evolutionary mechanisms that generate the process of speciation in the *Pteronotus* complex and to identify the causes of incongruence between mitochondrial and nuclear data. Furthermore, it remains to determine whether call differences have resulted from food competition or from intraspecific recognition in the context of sexual selection. Investigation of foraging behavior and diet of each species is crucial to analyze whether call frequency differences between the two species influence habitat or prey preference. We also need to identify the reproductive period of both species as well as their reproductive sites, and in particular whether *P. alitonus* and *P. rubiginosus* are both able to modulate their echolocation call during the mating period. The presence of a sexual dimorphism in echolocation in *P. rubiginosus* but not in *P. alitonus* could indeed reflect a slight shift in their breeding phenology, thereby limiting current genetic interaction between the two species. Future studies should also look deeper into the relationship between morphological, call frequencies and habitat use by sympatric and allopatric populations of *Pteronotus* species.

# Methods

## Study area and bat capture

*Pteronotus* bats were sampled from 4 caves in the tropical rainforests of French Guiana (Fig. 3): Mathilde (MA), Scierie (SC), Parfums (PA) and Montagne des gouffres (MG). There were four sampling sessions: late August/September 2010, July 2011, July 2012, October 2015, except for MG which was only sampled in October 2015. Distances between caves vary from 12 to 15 km between MA and SC to 80 km between MA or SC and PA. MG is 55 km from PA and 25 km from SC and MA.

Bats were caught in a Two-Bank HarpTrap (AUSTBAT Research Equipment, Victoria, Australia; catching surface of 4.0 m²) as they left the caves for foraging. We placed them individually in cotton bags until sampling. Adults were distinguished from juveniles by trans-illumination of the cartilaginous epiphyseal plates in the phalanges [60]. We estimated the reproductive status of each bat by examining the development of testes and caudae epipidymides for males, and the development of mammary glands and of nipples for females. Pregnant females were identified through palpation. Lactation was confirmed by gentle squeezing the mammary glands and nipples. Non-juvenile individuals were classified as inactive or active males, and inactive, pregnant or lactating females. We measured forearm length (with a dial caliper) and body weight (with an electronic scale) to the nearest 0.05 mm and 0.1 g, respectively. Tissue sample for DNA analysis was collected from the wing membrane (patagium) using a 3-mm diameter biopsy punch (Kai Industries, Gifu, Japan) and

preserved in 70% ethanol solution until DNA extraction. Bats were released at the place of capture after sampling.

## Molecular analysis

### Identification of microsatellites via high-throughput sequencing

Ten individuals were randomly chosen in order to build a pool of DNA. This pool was then used to design microsatellite markers based on the GS-FLX® method [61]. This method consists of firstly, a fragmentation of the genomic DNA. The method then proceeds by an enrichment of microsatellite sequences through the addition of primers with common STR patterns (TG, AAC, AAG, ACAT, TC, AAG, ACG, ACTC). Finally, the enriched DNA is amplified with a High Fidelity Taq polymerase. The bands thus obtained were quantified in order to have a minimum of 5 ng of genomic DNA. The resulting sequences were analyzed with QDD software [62].

### Multiplex building and amplification

Thirty microsatellite markers were chosen among the set of 3487 markers identified via GS-FLX®. A Simplex Step was performed in order to verify their proper amplification. More steps of primer design and primer concentration adjustment followed in order to build a multiplex of 20 microsatellite markers. During this step of multiplex building, tetranucleotides were privileged because we expected them to have a lower mutation rate - although this is uncertain given the high interspecies variability [63, 64] - and because of their higher legibility in subsequent computer analyses.

### Sequencing and genotyping

A total of 748 *Pteronotus* specimens out of the 1349 captured were genotyped at 20 microsatellite markers (see Additional file 2: Table S7). Total genomic DNA was extracted using purification column kits (Nucleospin 96 Tissue, Macherey-Nagel) following the manufacturer's instructions and in the presence of positive and negative controls. PCR reactions were carried out in 96-well microplates with three negatives and positive amplification controls to verify lack of contamination, in a total volume of 10 µl (2 ng/µl of DNA, 1.58 µl of a mix of primers at a concentration between 0.06 and 0.5 µM in the final PCR, 5 µl of 2X Mastermix). The samples were first denatured at 95 °C for 5 min. Then, 30 cycles followed (denaturation at 95 °C for 30s; hybridization at 58 °C for 90s; elongation at 72 °C for 30s) and a final elongation step at 60 °C for 30 min. PCR products were resolved on a ABI PRISM 3130XL capillary sequencer (Applied Biosystems) with formamide (denaturing conditions) and an internal size marker (600 liz; Applied Biosystems) in one migration. The electrophorograms were analyzed using GENEMAPPER 4.1
<figure><img src="image_7.png" /><figcaption>**Fig. 3** Location of the sampled sites in French Guiana. MA: Mathilde cave, SC: Scierie cave, PA: Parfums cave, MG: Montagne des gouffres cave</figcaption></figure>
(Applied Biosystems/Life Technologies) independently twice by different operators. Results were then compared, and ambiguous loci were set to missing data.

We also amplified the cytochrome c oxydase I (COI) mitochondrial marker (500 bp) for 80 individuals randomly chosen using the primers Bat-COI-01-F (5'-TGAG CAGGAATAGTAGGCAC-3') and Bat-COI-03-R (5'-CG GCAGGGTCAAAGAATGTG-3'). Sequencing reactions were performed with the corresponding amplifying primers using a BigDye Terminator Cycle Sequencing Kit v.2.0 (Applied Biosystems, USA) and run on an ABI 3730 automated sequencer (Applied Biosystems, USA).

## Species identification and hybridization detection

In order to identify the species to which the sampled specimens belonged, we built a tree based on the COI sequences of these 80 individuals and added 28 samples obtained in previous studies [24, 25, 28] for which COI sequences were retrieved from Genebank (see Additional file 1: Table S1 for the detailed list of GenBank sequences used). *Pteronotus gymnonotus* and *P. personatus* were used as outgroups.

COI sequences for the 108 specimens (80 + 28) were aligned using MACSE [65] with default parameters. The alignment was trimmed with trimAl v1.2rev59 [66] using the *-gappyout* option. Phylogenetic reconstruction was achieved using the IQTree Maximum Likelihood method [67] after automatic selection of the best fitting model (X + F + G4) and 1000 ultrafast bootstrap replicates were performed to get support values associated with the branches.

The genetic distance between COI sequences was also computed for all pairs of sequences in order to compare within and between COI divergence. This was done in R language [68] using the *dist.dna* function of the "ape" package v.5 [69] with the "raw" method (simple computation of the proportion of sites that differ between each pair of sequence).

The group identification of all individuals was also carried out using the microsatellites we developed specifically for *Pteronotus*, along with the detection of hybrids, using the Bayesian software STRUCTURE 2.3.4 [70, 71]. We used the admixture model with correlated allele frequencies between populations without any prior information on the putative group affiliation of individuals in all analyses. The program was run with a Monte-Carlo Markov chain length of 1,000,000 after a burn-in of 300,000 iterations. For the identification of groups, we first validated the number of clusters that best described the data following the method of Evanno et al. [72] implemented in STRUCTURE HARVESTER online web 0.6.94 [73]. For this purpose, we ran STRUCTURE for values of *K*, the number of clusters, ranging from one to four, twenty independent runs

being carried out for each value of *K*. The outputs from the 20 runs corresponding to the optimal value of *K* were then compiled using CLUMPP [74] to access consensus individual probabilities of assignment (*q*-values), and their 90% credibility interval, to the different clusters identified.

For the detection of hybrids, we selected the 80 individuals presenting the highest *q*-values in each cluster and used them to simulate individuals belonging to different hybrid classes (300 simulated individuals from parental classes, F1, F2, first generation backcrosses) using the function *hybridize* from the package *adegenet* [75] implemented in R [76]. These individuals were then pooled together and analysed using STRUCTURE to obtain their *q*-values for the *K* clusters identified in the real population. Based on these *q*-values, we defined two thresholds (two thresholds because *K* = 2, see Results section) corresponding to the lowest *q*-value reached by a simulated parent in each cluster (TP1 and TP2). These thresholds were then used to categorize the sampled individuals either as a parent or a hybrid either on the basis on their mean *q*-value (a conservative approach that detects individuals that are clearly hybrids) or the lower bound of their 90% credibility interval (a relaxed approach that detects all individuals presenting cues of hybridization as hybrids). We ran additional simulations in order to assess whether hybrids detected with the relaxed approach may be misclassified parental forms. To do so, we simulated again hybrid classes (300 individuals from each parental population, 10 F1, 15 backcrosses to each population) using the function *hybridize* and analysed the simulated individuals with STRUCTURE. We repeated the simulations thirty times and counted the number of parental forms detected as hybrids with the methods previously described. For this additional study, we simulated a few hybrids (10 F1, 15 backcrosses to each parental form) to consider a scenario consistent with the proportion of hybrids detected in the sampled population. We used the relaxed approach to select unambiguously parental groups for the morphologic and acoustic analyses. The detection of hybrids was also performed using the computer program NEWHYBRIDS 1.1 [77] This analysis was carried out using Jeffreys' prior with a burn-in period of 50,000 and the MCMC chain of 100,000 iterations.

## Genetic variability and population genetics of *Pteronotus* sp.

Once the group of each individual was determined, microsatellite loci were examined for null (non-amplifying) alleles in each species using the program Micro-checker v.2.2.3 [78]. The significance of null allele frequencies was assessed using a binomial exact test following De Mêeus et al. [79]. We also tested for
deviation from Hardy-Weinberg equilibrium and linkage disequilibrium for each pair of loci with FSTAT 2.9.3.2 [80]. P-values were adjusted using a Bonferroni correction. All subsequent analyses were performed without the loci for which we detected significant signs of null allele frequencies or deviance from Hardy-Weinberg equilibrium.

Allelic richness and Weir and Cockerham's *F*<sub>ST</sub> estimates for each group were assessed using FSTAT 2.9.3.2 [80]. In order to test the role of the caves in the structuration of the bat populations, we also calculated pairwise *F*<sub>ST</sub> between each pair of caves. Exact G-tests of population differentiation [81] were performed to test for differentiation significance. Finally, sex-biased dispersal tests were achieved based on the *F*<sub>ST</sub> because it is considered to be the most powerful statistic available for this test [82]. A thousand randomizations were made for each species and unilateral tests were performed setting either male or female as the dispersing sex.

## Morphology

Only adults were measured (presented in Additional file 2: Table S7). We performed a two-way ANOVA to determine whether there were differences in forearm (FA) length between sexes and between groups. Because sexual dimorphism may be quite different between species [24, 28], we tested the effect of sex by comparing the nested models Sex\*Group and Group using the classical Fisher's test of linear models. Rejecting the null hypothesis would then reveal the presence of a sexual dimorphism in at least one species. In that case, post-hoc *t*-tests comparing males and females were performed in each group considered independently and the hypothesis that sexual dimorphism is identical in both groups was tested by comparing models Sex\*Group and Sex +Group. Recorded body weights were not analysed because of the high variability according to time (before or after being fed) and period of capture (parturition state). Statistical analyses were performed with R version 3.3.3, and significant data were determined as having *p*-values less than 0.05.

## Echolocation frequency

Echolocation calls were recorded in 2015 from 262 *Pteronotus* bats held at 30 cm from a ZoomH2 microphone (Zoom Corporation, Japan) linked to a Pettersson D240X ultrasonic detector (Pettersson Elektronik AB, Uppsala, Sweden) (presented in Additional file 2: Table S7). Signals were analyzed with Bat Sound Pro 3.4 ultrasound analysis software (Pettersson Elektronik AB, Uppsala, Sweden) based on spectrograms with a Hanning window at a sampling frequency of 44.1 kHz and a Fast Fourier transformation (FFT) size of 512. From two to five calls were analyzed for each bat, and for each recorded call, the harmonic containing most energy (Frequency of Maximal

Energy: FME) was identified from the power spectrum and measurements taken from the constant frequency (CF) component of the call (see [28, 83]). Among genotyped individuals, a two-way ANOVA was applied to evaluate differences in frequency between groups and sexes. We followed the same procedure as described above for forearm length. For the echolocation variable, five outlying points were identified in *Pteronotus* B (with very low echolocation value compared to the rest of the sample in this group). These outliers were removed from the analysis.

## Reproductive periods

In order to have an inkling as to the reproductive periods of the two groups of *Pteronotus* identified, we analysed the reproductive status for each individual as well as their age (juvenile or adult) during the different sampling periods, without any statistical test due to the low number of individuals sampled for each period.

## Additional files

**Additional file 1: Table S1.** Table. List of external COI sequences used in this study. The sequence in bold corresponds to type-locality of *P. rubiginosus*, caught in MatoGrosso, Brazil compared to *P. sp4* in De Thoisy et al. (2014). **Table S2.** Characteristics of microsatellite markers. **Figure S3.** Graphical representation of assignation probabilities for *Pteronotus* A and *Pteronotus* B sampled. Each individual (x axis) is represented by a vertical bar divided in two parts according to its assignation probability (y axis) in each of the two clusters. **Figure S4.** Graphical representation of the difference between the Conservative and the Relaxed approach for the detection of hybrids. **Table S5.** Pairwise *F*<sub>ST</sub> values between caves for each species and both sexes. **Table S6.** Life cycle of *Pteronotus* A and B. PRE = Pregnant, LAC = Lactating females, PLAC = Post-lactating females, NS = No status (non pregnant, non-lactating, non post-lactating females, and non breeding males). n(A): number of *Pteronotus* A (number of adult females in parentheses); n(B): number of *Pteronotus* B (number of adult females in parentheses); n: total number of *Pteronotus*. (DOCX 319 kb)

**Additional file 2: Table S7.** Data used in the study, including sex, forearm length, FME, genotype at 20 microsatellite markers for 748 *Pteronotus* sp. Data in bold are individuals identified as hybrids using microsatellite (N = 30) or mitochondrial data (N = 6). (XLSX 253 kb)

### Abbreviations

*C*: Degrees Celsius; CO1: cytochrome oxidase 1; CytB: cytochrome b; DNA: Deoxyribonucleic acid; FME: Frequency of Maximal Energy; mtDNA: Mitochondrial DNA; PCR: Polymerase chain reaction

### Acknowledgements

We are particularly grateful to all volunteers and field workers who have helped us during the field sessions: V. Alt, M. Bely, G. Chagneau, M. Dorfiac, S. Dufour, C. Gizardin, G. Leblanc, M. Leuchtmann, E. Loufti, A. Le Guen, and L. Trebucq. We thank A. Pavé for his scientific support in French Guiana. We thank M. Barataud, R. Buckland, and two anonymous referees for their helpful comments.

### Funding

Financial support was received from the CNRS (to DP) and by the LABEX ECOFECT (ANR-11-LABX-0048) of Université de Lyon, within the program "Investissements d'Avenir" (ANR-11-IDEX-0007) operated by the French National Research Agency (ANR). The funding agencies of this study play no role in the design, collection, analysis, and interpretation of data and writing of this manuscript.
**Availability of data and materials**

Data (including genotypes, FME and forearm length of all the individuals analyzed during the current study) are included within the article and its additional files.

**Authors' contributions**

JBP and DP conceived the study and supervised the field work. MPB, DF, EP1, EP2, OFC and DP carried out the statistical and genetic analyses, and DDV the phylogenetic analyses. LM, CK and GQ developed the microsatellite markers, and CR genotyped the individuals. OFC, MPB, DF, DDV and DP wrote the manuscript. All authors read, edited and approved the final version of the manuscript.

**Ethics approval and consent to participate**

Authorization for bat capture in French Guiana was provided by the Ministry of Ecology, Environment, and Sustainable development over the period 2015–2020 (approval no. C692660703 from the Departmental Direction of Population Protection (DDPP, Rhône, France). All methods (capture and animal handling) were approved by the MNHN, the SFEPM and the DEAL-Guyane.

**Consent for publication**

Not applicable

**Competing interests**

The authors declare that they have no competing interests.

**Publisher's Note**

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Author details**

¹University Lyon, Université Lyon 1, CNRS, Laboratoire de Biométrie et Biologie Evolutive UMR5558, F-69622 Villeurbanne, France. ²Université de Lyon, LabEx Ecofect, Nadine Cizaire, 92 rue Pasteur, CS 30122 69361 Lyon Cedex 07, France. ³SEISE 26 bis Barrouil, 33720 Illats, France. ⁴ANTAGENE, Animal Genomics Laboratory, 6 allée du Levant, 69890 La Tour de Salvagny (Lyon), France. ⁵Université de Lyon, VetAgro Sup - Campus Vétérinaire de Lyon, 1 Avenue Bourgelat, BP 83, F-69280 Marcy l'Etoile, France. ⁶Office National de la Chasse et de la Faune Sauvage - Unité Faune de Montagne, 147 Route de Lodève, Les Portes du Soleil, F-34990 Juvignac, France. ⁷UMR ESE, Ecology and Ecosystem Health, INRA, Agrocampus Ouest, 65 rue de Saint-Brieuc, 35042 Rennes Cedex, France.

Received: 22 June 2018 Accepted: 2 November 2018

Published online: 20 November 2018

**References**

1. Bickford D, Lohman DJ, Sodhi NS, Ng PKL, Meier R, Winker K, et al. Cryptic species as a window on diversity and conservation. Trends Ecol Evol. 2007; 22:148–55.

2. Pfenninger M, Schwenk K. Cryptic animal species are homogeneously distributed among taxa and biogeographical regions. BMC Evol Biol. 2007;7:1–7.

3. Beheregaray LB, Taccone A. Cryptic biodiversity in a changing world. J Biol. 2007;6(4):9.

4. Hebert PDN. Cywinska a, ball SL, deWaard JR. biological identifications through DNA barcodes. Proc R Soc B Biol Sci. 2003;270:313–21. Available from: https://doi.org/10.1098/rspb.2002.2218.

5. Hebert PDN, Ratnasingham S, de Waard JR. Barcoding animal life: cytochrome c oxidase subunit 1 divergences among closely related species. Proc R Soc B Biol Sci. 2003;270:S96–9. Available from: https://doi.org/10.1098/rbl.2003.0025.

6. Ashrafi S, Rutishauser M, Ecker K, Obrist MK, Arlettaz R, Bontadina F. Habitat selection of three cryptic *Plecotus* bat species in the European Alps reveals contrasting implications for conservation. Biodivers Conserv. 2013;22:2751–66.

7. Rutishauser MD, Bontadina F, Braunisch V, Ashrafi S, Arlettaz R. The challenge posed by newly discovered cryptic species: disentangling the environmental niches of long-eared bats. Divers Distrib. 2012;18:1107–19.

8. Witt JDS, Threloff DL, Hebert PDN. DNA barcoding reveals extraordinary cryptic diversity in an amphipod genus: implications for desert spring conservation. Mol Ecol. 2006;15:3073–82.

9. Karanth KP. Molecular systematics and conservation of the langurs and leaf monkeys of South Asia. J Genet. 2010;89:393–9.

10. Ballard JWO, Whitlock MC. The incomplete natural history of mitochondria. Mol Ecol. 2004;13:729–44.

11. Coyne JA, Orr HA. Patterns of speciation in *Drosophila*. Evolution (N Y), 1989; 43:362–81.

12. Funk DJ, Omland KE. Species-level Paraphyly and Polyphyly: frequency, causes, and consequences, with insights from animal mitochondrial DNA. Annu Rev Ecol Evol Syst. 2003;34:397–423. Available from: https://doi.org/10.1146/annurev.ecolsys.34.011802.132421.

13. Karanth KP. Species complex, species concepts and characterization of cryptic diversity: vignettes from Indian systems. Curr Sci. 2017;112:1320–4.

14. Paterson ID, Mangan R, Downie DA, Coetzee JA, Hill MP, Burke AM, et al. Two in one: cryptic species discovered in biological control agent populations using molecular data and crossbreeding experiments. Ecol Evol. 2016;6:6139–50.

15. Miller W, Schuster SC, Welch AJ, Ratan A, Bedoya-Reina OC, Zhao F, et al. Polar and brown bear genomes reveal ancient admixture and demographic footprints of past climate change. Proc Natl Acad Sci. 2012;109:E2382–90. Available from: https://doi.org/10.1073/pnas.1210506109.

16. Höbel G, Gerhardt HC. Reproductive character displacement in the acoustic communication system of green tree frogs (*Hyla cinerea*). Evolution. 2003;57: 894–904.

17. Honda-Sumi E. Difference in calling song of three field crickets of the genus *Teleogryllus*: the role in premating isolation. Anim Behav. 2005;69:881–9.

18. Kingston T, Rossiter SJ. Harmonic-hopping in Wallacea's bats. Nature. 2004; 429:9–12.

19. Puechmaille SJ, Gouilh MA, Piyapan P, Yokubol M, Me KM, Bates PJ, et al. The evolution of sensory divergence in the context of limited gene flow in the bumblebee bat. Nat Commun. 2011;2:573.

20. Ritchie MG, Gleason JM. Rapid evolution of courtship song pattern in *Drosophila willistoni* sibling species. J Evol Biol. 1995;8:463–79.

21. Jones GJ. Acoustic signals and speciation: the roles of natural and sexual selection in the evolution of cryptic species. Adv Study Behav. 1997;26:317–54.

22. Barlow KE, Jones G. Roosts, echolocation calls and wing morphology of two phonic types of *Pipistrellus pipistrellus*. Zeitschrift für Saugetierkunde. 1999; 64:257–68.

23. Davidson-Watts I, Walls S, Jones G. Differential habitat selection by Pipistrellus pipistrellus and Pipistrellus pygmaeus identifies distinct conservation needs for cryptic species of echolocating bats. Biol Conserv. 2006;133:118–27.

24. Clare EL, Adams AM, Maya-Simões AZ, Eger JL, Hebert PDN, Fenton MB. Diversification and reproductive isolation: cryptic species in the only New World high-duty cycle bat, *Pteronotus parnellii*. BMC Evol Biol. 2013;13:26.

25. López-Wilchis R, Flores-Romero M, Guevara-Chumacero LM, Serrato-Díaz A, Díaz-Larrea J, Salgado-Mejia F, et al. Evolutionary scenarios associated with the *Pteronotus parnellii* cryptic species-complex (Chiroptera: Mormoopidae). Acta Chiropterologica. 2016;18:91–116. Available from: https://doi.org/10.3161/15081109ACC2016.18.1.004.

26. Pavan AC, Marroig G. Integrating multiple evidences in taxonomy: species diversity and phylogeny of mustached bats (Mormoopidae: *Pteronotus*). Mol Phylogenet Evol. 2016;103:184–98 Available from: https://doi.org/10.1016/j. ympev.2016.07.011.

27. Pavan AC, Bobrowiec PED. R PA. Geographic variation in a south American clade of mormoopid bats, *Pteronotus* (Phyllogdia), with description of a new species. J Mammal. 2018;99:624–45.

28. De Toisy B, Pavan AC, Delaval M, Lavergne A, Luglia T, Pineau K, et al. Cryptic diversity in common mustached bats *Pteronotus cf. parnellii* (Mormoopidae) in French Guiana and Brazilian Amapa. Acta Chiropterologica. 2014;16:1–13 Available from: https://doi.org/10.3161/ 150811014X683228.

29. Barataud M, Giosa S, Leblanc F, Rufray V, Disca T, Tillon L, et al. Identification et écologie acoustique des chiroptères de Guyane française. Le Rhinolophe. 2013;19:103–45.

30. Johns GC, Avise JC. A comparative summary of genetic distances in the vertebrates from the mitochondrial cytochrome b gene. Mol Biol Evol. 1998;15:1481–90.

31. López-Baucells A, Torrent L, Rocha R, Pavan AC, Bobrowiec PED, Meyer CFJ. Geographical variation in the high-duty cycle echolocation of the cryptic common mustached bat *Pteronotus cf. rubiginosus* (Mormoopidae). Bioacoustics. 2018;27:341–57.
32. Pfennig K. Character displacement: ecological and reproductive responses to a common evolutionary problem. Q Rev Biol. 2009;84:253–76.

33. Wu H, Jiang T, Huang X, Lin H, Wang H, Wang L, et al. A test of Rensch's rule in greater horseshoe bat (*Rhinolophus ferrumequinum*) with female-biased sexual size dimorphism. PLoS One. 2014;9:e86085.

34. Wu H, Jiang T, Huang X, Feng J. Patterns of sexual size dimorphism in horseshoe bats : testing Rensch 's rule and potential causes. Sci Rep. 2018; 1–13 Available from: https://doi.org/10.1038/s41598-018-21077-7.

35. Butler MA, Schoener TW, Losos JB. The relationship between sexual size dimorphism and habitat use in greater Antillean *Anolis* lizards. Evolution. 2000;54:259–72.

36. Myers P. Sexual dimorphism in size of Vespertilionid bats. Am Nat. 1978;112: 701–11.

37. Jiang T, You Y, Liu S, Lu G, Wang L, Wu H, et al. Factors affecting geographic variation in echolocation calls of the endemic *Myotis davidii* in China. Ethology. 2013;119:881–90.

38. Lin A, Jiang T, Kanwal JS, Lu G, Luo J, Wei X, et al. Geographical variation in echolocation vocalizations of the Himalayan leaf-nosed bat: contribution of morphological variation and cultural drift. Oikos. 2015;124:364–71.

39. Russo D, Almenar D, Aihartza J, Goiti U, Salsamendi E, Garin I. Habitat selection in sympatric *Rhinolophus mehelyi* and *R. euryale* (Mammalia: Chiroptera). J Zool. 2005;266:327–32.

40. Zhang J. Neutral theory and phenotypic evolution. Mol Biol Evol. 2018;35: 1327–31. Available from: https://doi.org/10.1093/molbev/msy065.

41. Jones G. Scaling of echolocation call parameters in bats. J Exp Biol. 1999; 202:3359–67.

42. Thiagavel J, Santana SE, Ratcliffe JM. Body size predicts echolocation call peak frequency better than gape height in Vespertilionid bats. Sci Rep. 2017;7:1–6 Available from: https://doi.org/10.1038/s41598-017-00959-2.

43. Grilliot ME, Burnett SC, Mendonça MT. Sexual dimorphism in big Brown bat (*Eptesicus fuscus*) ultrasonic vocalizations is context dependent. J Mammal. 2009;90:203–9. Available from: https://doi.org/10.1644/07-MAMM-A-161.1.

44. Puechmaille SJ, Borissov IM, Zsebok S, Allegrini B, Hizem M, Kuenzel S, et al. Female mate choice can drive the evolution of high frequency echolocation in bats: a case study with Rhinolophus mehelyi. PLoS One. 2014;9(7):e103452.

45. Berthier P, Excoffier L, Ruedi M. Recurrent replacement of mtDNA and cryptic hybridization between two sibling bat species *Myotis myotis* and *Myotis blythii*. Proc R Soc B Biol Sci. 2006;273:3101–23. Available from: https://doi.org/10.1098/rspb.2006.3680.

46. Mao X, He G, Zhang J, Rossiter SJ, Zhang S. Lineage divergence and historical gene flow in the Chinese horseshoe bat (*Rhinolophus sinicus*). PLoS One. 2013;8(2):e56786.

47. Toews DPL, Brelsford A. The biogeography of mitochondrial and nuclear discordance in animals. Mol Ecol. 2012;21:3907–30.

48. Shaw KL, Lugo E. Mating asymmetry and the direction of evolution in the Hawaiian cricket genus Laupala. Mol Ecol. 2001;10:751–9.

49. Stein AC, Uy JA. Unidirectional introgression of a sexually selected trait across an avian hybrid zone: a role for female choice? Evolution. 2006;60: 1476–85.

50. Angell RL, Butlin RK, Altringham JD. Sexual Segregation and Flexible Mating Patterns in Temperate Bats. PLoS One. 2013;8: https://doi.org/10.1371/journal.pone.0054194.

51. Bogdanovicz W, Piksa KTA. Genetic structure in three species of whiskered bats (genus *Myotis*) during swarming. 2012;93:799–807.

52. Watt E, Fenton MB. DNA fingerprinting provides evidence of discriminate suckling and non-random mating in little brown bats *Myotis lucifugus*. Mol Ecol. 1995;4:261–4.

53. Bogdanowicz W, Piksa K, Tereba A. Hybridization hotspots at bat swarming sites. PLoS One. 2012;7(12):e53334.

54. Melo-Ferreira J, Boursot P, Suchentrunk F, Ferrand N, Alves PC. Invasion from the cold past: extensive introgression of mountain hare (*Lepus timidus*) mitochondrial DNA into three other hare species in northern Iberia. Mol Ecol. 2005;14:2459–64.

55. Loren H. Rieseberg. Evolution: replacing genes and traits through hybridization. Curr Biol. 2009;19:119–22.

56. Ballard JWO, Melvin RG, Miller JT, Katewa SD. Sex differences in survival and mitochondrial bioenergetics during aging in Drosophila. Aging Cell. 2007;6: 699–708.

57. Bazin E, Glemin S, Galtier N. Population size does not influence mitochondrial genetic diversity in animals. Science. 2006;312:570–2.

58. Grant WS, Spies IB, Canino MF. Biogeographic evidence for selection on mitochondrial DNA in North Pacific walleye Pollock *Theragra chalcogramma*. J Hered. 2006;97:571–80.

59. Stewart JB, Freyer C, Elson JL, Wredenberg A, Cansu Z, Trifunovic A, et al. Strong purifying selection in transmission of mammalian mitochondrial DNA. PLoS Biol. 2008;6:e0063–71.

60. Anthony ELP. Age determination in bats. Ecological and behavioral methods for the study of bats (T. H. Kunz, ed.). Washington, D.C.: Smithsonian Institution Press; 1988. Pp. 47–58. https://jhupbooks.press.jhu.edu/content/ecological-and-behavioral-methods-study-bats.

61. Malausa T, Gilles A, Meglècz E, Blanquart H, Duthoy S, Costedoa C, et al. High-throughput microsatellite isolation through 454 GS-FLX titanium pyrosequencing of enriched DNA libraries. Mol Ecol Resour. 2011;11:638–44.

62. Meglècz E, Costedoa C, Dubut V, Gilles A, Malausa T, Pech N, et al. QDD: a user-friendly program to select microsatellite markers and design primers from large sequencing projects. Bioinformatics. 2009;26:403–4.

63. Ellegren H. Microsatellites: simple sequences with complex evolution. Nat Rev Genet. 2004;5:435–45.

64. Vigouroux Y, Jaqueth JS, Matsuoka Y, Smith OS, Beavis WD, Smith JSC, et al. Rate and pattern of mutation at microsatellite loci in maize. Mol Biol Evol. 2002;19:1251–60.

65. Ranwez V, Harispe S, Delsuc F, Douzery EJP. MACSE: multiple alignment of coding SEquences accounting for frameshifts and stop codons. PLoS One. 2011;6(e9):e22594.

66. Capella-Gutiérrez S, Silla-Martinez JM, Gabaldón T. trimAl: A tool for automated alignment trimming in large-scale phylogenetic analyses. Bioinformatics. 2009;25:1972–3.

67. Nguyen LT, Schmidt HA, Von Haeseler A, Minh BQ. IQ-TREE: a fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Mol Biol Evol. 2015;32:268–74.

68. R Core Team. R: a language and environment for statistical computing. 2015.

69. Paradis E, Claude J, Strimmer K. APE: analyses of phylogenetics and evolution in R language. Bioinformatics. 2004;20:289–90.

70. Pritchard JK, Stephens M, Donnelly P. Inference of population structure using multilocus genotype data. Genetics. 2000;155:945–59.

71. Falush D, Stephens M, Pritchard JK. Inference of population structure using multilocus genotype data: linked loci and correlated allele frequencies. Genetics. 2003;164:1567–87.

72. Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol Ecol. 2005;14:2611–20.

73. Earl DA, von Holdt BM. STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method. Conserv Genet Resour. 2012;4:359–61.

74. Jakobsson M, Rosenberg NA. CLUMPP: a cluster matching and permutation program for dealing with label switching and multimodality in analysis of population structure. Bioinformatics. 2007;23:1801–6.

75. Jombart T. Adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics. 2008;24:1403–5.

76. R Core Team. R: a language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. 2016. http://www.R-project.org/.

77. Anderson EC, Thompson EA. A model-based method for identifying species hybrids using multilocus genetic data. Genetics. 2002;160:1217–29.

78. Van Oosterhout C, Hutchinson WF, Wills DPM, Shipley P. MICRO-CHECKER: software for identifying and correcting genotyping errors in microsatellite data. Mol Ecol Notes. 2004;4:535–8.

79. de Meeüs T, Béati L, Delaye C, Aeschlimann A, Renaud F. Sex-biased genetic structure in the vector of Lyme disease, Ixodes Ricinus. Evolution. 2002;56:1802–7. Available from: https://doi.org/10.1111/j.0014-3820.2002.tb00194.x.

80. Goudet J. Computer Note J Hered. 1995;86:485–6.

81. Goudet J, Raymond M, De Meeüs T, Rousset F. Testing differentiation in diploid populations. Genetics. 1996;144:1933–40.

82. Goudet J, Perrin N, Waser P. Tests for sex-biased dispersal using bi-parentally inherited genetic markers. Mol Ecol. 2002;11:1103–14.

83. Schnitzler HU, Kalko EKV. In: Kunz TH, Racey PA, editors. How Echolocating bats search and find food. in Bats: phylogeny, morphology, echolocation, and conservation biology, vol. 365. Washington, D.C.: Smithsonian Institution Press; 1998. p. 183–96.