Phylogenetic signal in the acoustic parameters of the advertisement calls of four clades of anurans
© Gingras et al.; licensee BioMed Central Ltd. 2013
Received: 18 December 2012
Accepted: 12 June 2013
Published: 1 July 2013
Anuran vocalizations, especially their advertisement calls, are largely species-specific and can be used to identify taxonomic affiliations. Because anurans are not vocal learners, their vocalizations are generally assumed to have a strong genetic component. This suggests that the degree of similarity between advertisement calls may be related to large-scale phylogenetic relationships. To test this hypothesis, advertisement calls from 90 species belonging to four large clades (Bufo, Hylinae, Leptodactylus, and Rana) were analyzed. Phylogenetic distances were estimated based on the DNA sequences of the 12S mitochondrial ribosomal RNA gene, and, for a subset of 49 species, on the rhodopsin gene. Mean values for five acoustic parameters (coefficient of variation of root-mean-square amplitude, dominant frequency, spectral flux, spectral irregularity, and spectral flatness) were computed for each species. We then tested for phylogenetic signal on the body-size-corrected residuals of these five parameters, using three statistical tests (Moran’s I, Mantel, and Blomberg’s K) and three models of genetic distance (pairwise distances, Abouheif’s proximities, and the variance-covariance matrix derived from the phylogenetic tree).
A significant phylogenetic signal was detected for most acoustic parameters on the 12S dataset, across statistical tests and genetic distance models, both for the entire sample of 90 species and within clades in several cases. A further analysis on a subset of 49 species using genetic distances derived from rhodopsin and from 12S broadly confirmed the results obtained on the larger sample, indicating that the phylogenetic signals observed in these acoustic parameters can be detected using a variety of genetic distance models derived either from a variable mitochondrial sequence or from a conserved nuclear gene.
We found a robust relationship, in a large number of species, between anuran phylogenetic relatedness and acoustic similarity in the advertisement calls in a taxon with no evidence for vocal learning, even after correcting for the effect of body size. This finding, covering a broad sample of species whose vocalizations are fairly diverse, indicates that the intense selection on certain call characteristics observed in many anurans does not eliminate all acoustic indicators of relatedness. Our approach could potentially be applied to other vocal taxa.
KeywordsAnurans Molecular phylogeny Phylogenetic signal Acoustics Advertisement calls Mitochondrial DNA 12 S Rhodopsin Bufo Hylinae Leptodactylus Rana
Although rich and varied, the vocalizations of anurans (frogs and toads) are relatively fixed and show no evidence of vocal learning [1, 2], and are thus presumably mostly genetically determined. Some basal level of general auditory stimulation may be necessary for the neural development of species-typical advertisement calls , but call structure itself appears to be highly stereotyped within a species. Moreover, anuran vocalizations, especially their advertisement calls, differ considerably across species and can often be used reliably to determine species [4, 5]. This leads straightforwardly to the hypothesis, inspired by Blair , that the degree of similarity between advertisement calls in anurans should be related to large-scale phylogenetic relationships, and that species that are evolutionarily distant would be expected, on average, to display vocalizations that are more dissimilar than species that are more closely related. However, a recent study by Tobias et al.  found only a weak phylogenetic signal in vocalizations of African clawed frogs, concluding that rapid evolution and frequent homoplasy can quickly erase acoustic indicators of phylogeny.
Detailed comparisons between differences in calling behavior and phylogenetic distances in vertebrates, involving up to 15 anuran species [8–10], 6 species of crested gibbons , and 11 deer species  can be found in the literature. However, these studies generally involved a small number of closely related species, with few exceptions such as a broader comparison involving 21 species of Bufonidae and Hylidae . There are, to our knowledge, no large-scale studies investigating the link between the acoustic similarity among advertisement calls and the phylogenetic distance on a large number of species representative of the broad geographic and taxonomic distribution of anurans. Moreover, previous studies typically relied upon specific acoustic features that were customized to the vocalizations of a particular clade, thereby restricting their applicability to a limited range of species.
The current study aimed to fill these lacunae, by comparing acoustic similarity and genetic distance for 90 species of anurans belonging to four clades with a wide geographic distribution: Bufo, Hylinae, Leptodactylus, and Rana. More precisely, we sought to find a set of low-level acoustic parameters applicable to distant clades that display a broad range of vocalization types (the expression “low-level acoustic parameters” refers to parameters that can be reliably extracted algorithmically in a variety of acoustic signals, such as spectral flux or dominant frequency). We then assessed the degree of autocorrelation between these acoustic parameters and phylogenetic dissimilarities to estimate the extent to which differences between call acoustics are linked to genetic divergence among clades.
This investigation is relevant to larger issues relating to the relative influence of various selective pressures on vocalization characteristics. If, in a taxon of non-vocal learners such as anurans, similarities in certain acoustic parameters were consistently related to phylogenetic distances, it would suggest that phylogenetic constraints are an important factor in determining those characteristics. On the other hand, when sexual selection is the main evolutionary force shaping vocalization parameters in anurans, the affected parameters might be expected to exhibit rapid, runaway-style selection, and the link between phylogenetic relatedness and acoustic similarity should be tenuous, except for closely-related species . Hence, our central aim is to understand whether some call parameters change slowly enough to retain a reliable signal of phylogeny across a broad range of clades and species.
To estimate phylogenetic distances, we chose the 12S rRNA region of the mitochondrial DNA (mtDNA), a region for which complete (or nearly complete) sequences are available for a large number of anurans and which has already been used to analyze relationships among hyloid frogs  and Malagasy reed frogs . Although some authors have questioned the suitability of mitochondrial DNA for phylogenetic inferences [16, 17], the extensive use of 12S in anuran phylogenetic studies [18, 19] means that sequences are available for numerous species, and suggests that it is appropriate for our purposes. However, a potential issue with using 12S mtDNA is that this gene is prone to rapidly accumulate mutations, thus possibly leading to saturation in nucleotide substitutions, which would potentially decrease the phylogenetic information contained in the dataset . Because genetic divergences in shallower relationships (such as intra-clade comparisons) could be expected to outweigh more distant relationships (such as inter-clade comparisons) in the case of a saturated gene, this could theoretically bias our analysis towards short-range relationships. We addressed this issue in two ways. First, we conducted a saturation test on the 12S mtDNA sequences to assess the degree of saturation present in our dataset. Second, we compared the results obtained using 12S mtDNA sequences to those obtained with the first exon of the rhodopsin gene, a nuclear gene, for a subset of 49 species for which both sequences were available. Nuclear protein-coding single-copy genes such as the rhodopsin gene have been shown to outperform mitochondrial sequences when analyzing deeper genetic divergences [21, 22]. Note that we do not expect 12S, rhodopsin, or other genes traditionally selected for phylogenetic analysis, to play any direct causal role in vocal production: these genes simply serve as proxies for overall phylogenetic distance. At present, the direct genetic determinants of acoustic dissimilarities in the vocalizations of any vertebrate species remain to be identified.
The acoustic properties of advertisement calls were quantified using mean values computed from a series of low-level acoustic parameters following the procedure previously developed for anurans in . The five parameters that this previous study showed to be most independent and informative were used: coefficient of variation of root-mean-square amplitude (CVA), dominant frequency (DF), spectral flux (SF), spectral irregularity (SI), and spectral flatness or tonality (TON). CVA refers to the standard deviation of the root mean square of the amplitude, divided by the true mean (note that, as a coefficient of variation, CVA is a unit-free measure and is thus independent of the absolute intensity of the sound or of the distance from the microphone), whereas DF (in Hz) represents the single frequency of maximal amplitude in the spectrum. SI is defined as the sum of the square of the differences in amplitude between adjoining partials . SF is a time-varying descriptor that corresponds to the Euclidean distance between two spectra. Lastly, TON is computed as the ratio between the geometric mean and the arithmetic mean of the power spectrum , with lower values representing more ‘spiky’ spectra, indicating the presence of strong partials (but not necessarily integer-multiple harmonics).
Crucially, these parameters can be measured from a very wide variety of sounds, unlike some traditional measures such as trill rate that may not be applicable to all species or vocalizations. Moreover, these parameters refer to well-characterized acoustic properties of the vocalizations that can be directly and automatically derived from recordings. These parameters correspond for the most part to spectral features, and global temporal patterning features are not considered in our analysis. However, CVA and SF capture local aspects of temporal variability, and are thus best characterized as spectro-temporal quantities.
The acoustic database for this study is derived entirely from published, commercially-available digital collections on CDs (listed in Additional file 1). The species in our sample were originally assigned to one of four genera on the basis of their identification by the recordists: Bufo, Hyla, Leptodactylus, and Rana. Several of the recordings originally assigned to Hyla species on these recordings have subsequently been reclassified to other genera, but remain within the Hylinae subfamily, a monophyletic clade. Here, we follow the nomenclature used by Pyron and Wiens , which is more conservative than Frost’s Amphibian Species of the World website , especially regarding Bufo and Rana which mostly retain the composition that they had prior to Frost et al.. However, given that the generic content of several anuran families remains in flux, we avoid referring specifically to “genera” or “subfamilies”, and will use the generic term “clades” to refer to our four taxonomic groupings for the remainder of the article. In any case, our analysis techniques use genetic distance as a proxy for phylogeny, and are thus robust to changes in systematic nomenclature.
Because body size imposes severe constraints on vocal signals, it is often closely related with the acoustical features of animal vocalizations. Indeed, an inverse relationship between body size and call frequency has been documented in numerous species of anurans (reviewed in [1, 27, 28]). In the case of the parameters selected for our analysis, both DF and TON were shown in a previous study to be inversely correlated with snout-vent length (SVL, a proxy for body size) in a sample of 136 species belonging to the four clades analyzed in the present study . Furthermore, body size is generally strongly autocorrelated with genetic distance (e.g., [30, 31]). It is thus necessary to partial out the relationship between body size and acoustic parameters prior to conducting a phylogenetic signal analysis. Here, we corrected for the effect of body size by first regressing the acoustic parameters on SVL using phylogenetic generalized least squares regression (PGLS) [32, 33], and then evaluating the autocorrelation between the regressed residuals and the genetic distances. PGLS regression was used instead of ordinary least squares regression because ignoring phylogeny in the size-correction procedure can lead to spurious results in the subsequent phylogenetic signal analysis [34, 35], and the PGLS method has been shown to be relatively robust to phylogenetic tree misspecification .
We tested for the presence of a phylogenetic signal in the size-corrected residuals of the acoustic parameters by applying three different statistical methods: Moran’s I test , Blomberg’s K, and the Mantel test . Both Moran’s and Mantel tests are general procedures for testing for spatial or genetic autocorrelation , whereas Blomberg’s K is a descriptive statistic based on a Brownian (random walk) model of trait evolution (BM), with a value of 1 corresponding to the degree of trait similarity expected under BM. A K of 0 indicates phylogenetic independence, whereas a K > 1 implies trait similarity greater than expected under BM . In contrast to Moran’s I and Blomberg’s K, which are suitable for univariate phenotypic traits, the Mantel test can be used to test for a correlation between a multivariate trait (corresponding for instance to a subset of the acoustic parameters examined here) and a distance matrix, and its power to detect a phylogenetic signal increases substantially with the number of traits . Here, trait distances were obtained by computing Euclidean distances on the standardized size-corrected residuals, and the best-fitting multivariate trait distance model was determined by a stepwise forward procedure (details given in the Methods section).
Although Moran’s I, Blomberg’s K, and the Mantel tests are closely related mathematically [41, 42], their performance depends heavily on the choice of genetic distance model . For this reason, we used three different models of genetic distance. The first model, based on the pairwise genetic distance between sequences, does not take into account phylogenetic relationships and can thus be seen as ‘phylogenetically naïve’. Distance matrices based on this model were used with both Moran’s and Mantel tests. The second model was based on Abouheif’s matrix of phylogenetic proximities [43, 44], which has been shown to be a powerful alternative to Blomberg’s K, especially for unresolved trees or when branch lengths are not accurate . Distances matrices based on this model were also tested with both Moran’s and Mantel tests. Finally, the genetic distance model used in Blomberg’s K is a variance-covariance matrix derived from the phylogenetic tree, where the covariance corresponds to the branch length from the root to the most recent common ancestor. Unlike the previous distance matrices that are not based on an explicit evolutionary model , the variance-covariance matrix is based on the BM model. All tests were conducted both on the entire sample of 90 species, and on each individual clade.
Testing for saturation in nucleotide substitutions in 12S mtDNA sequences
Testing for phylogenetic signal in acoustic parameters based on 12S mtDNA
Four acoustic parameters (DF, SF, SI, and TON), as well as SVL, were log-transformed to achieve normality (the values of the acoustic parameters computed for each recording, as well as the SVL values obtained from the literature, are listed in Additional file 3). All subsequent analyses were conducted on the log-transformed variables. All acoustic parameters (CVA, logDF, logSF, logSI, and logTON) were then phylogenetically regressed on logSVL using the PGLS method, to partial out the effect of body size . Separate regressions were conducted for each clade in the case of intra-clade analyses to account for the fact that the allometric relationship between body size and acoustic parameters may vary across clades .
Phylogenetic signal analysis (12S mtDNA) on size-corrected residuals
Moran’s I, pairwise distances
Moran’s I, Abouheif’s proximities
Blomberg’s K, variance-covariance matrix under BM
Mantel r, pairwise distances
Mantel r, Abouheif’s proximities
CVA, logDF, logSF
The intra-clade analyses were not as consistent across statistical methods, which may be due to the reduced power associated with smaller sample sizes [31, 42] and to the larger influence exerted by local discrepancies between different models of genetic distance. Nevertheless, we obtained several reliable results that did not depend on a particular statistical test or genetic distance model, such as the significant phylogenetic signals observed for CVA and logSF in Bufo, for CVA and logDF in Hylinae, and for logDF in Rana (the latter was detected by all methods except the Mantel test on pairwise distances). Importantly, Mantel tests on Abouheif’s proximities and the multivariate trait distances computed from CVA, logDF, and logSF yielded significant correlations for all four clades, indicating that a phylogenetic signal associated with this multivariate trait was detected on intra-clade genetic distances in all cases (note that the magnitude of the Mantel coefficient is often comparatively small even when statistically significant ).
Comparison between 12S mtDNA and rhodopsin
In order to assess the robustness of the phylogenetic signal detected in our acoustic parameters for these four clades and to validate our findings by comparing the results obtained on two genes displaying widely different mutation rates, we repeated our analysis on a subset of 49 species for which genetic data was available for both 12S mtDNA and rhodopsin (exon 1). Only three clades were represented: Bufo, Hylinae, and Rana. Given that rhodopsin is a very conserved gene, a number of closely related species harbored identical sequences for exon 1. The resolution of the phylogenetic tree derived from the rhodopsin sequence data was therefore relatively limited, rendering Blomberg’s K unsuitable for this analysis. Moreover, the small sample size led us to select methods that have been shown to be both powerful and relatively unaffected by the accuracy of the phylogenetic tree [40, 42]: Moran’s I (pairwise distances and Abouheif’s proximities), and the Mantel test based on Euclidean trait distances and Abouheif’s proximities.
Following the procedure outlined previously, all acoustic parameters (CVA, logDF, logSF, logSI, and logTON) were regressed on logSVL using the PGLS method, both for entire subset of 49 species and for each clade individually in the case of intra-clade analyses.
Comparison between 12S mtDNA and Rhodopsin on size-corrected residuals
Moran’s I, pairwise distances
Moran’s I, Abouheif’s proximities
Mantel r, Abouheif’s proximities
Results of intra-clade analyses were broadly consistent, with significance tests generally in agreement both between 12S and rhodopsin and between different models of genetic distance (pairwise distance or Abouheif’s proximities). Notably, a few significant intra-clade phylogenetic signals were detected by all tests on both genes, such as for CVA in Bufo and logDF in Hylinae (in line with the results observed with 12S on larger samples for both clades, see Table 1), as well as for logTON in Rana. Finally, Mantel tests on Abouheif’s proximities and the multivariate trait distances computed from logDF and logTON revealed significant correlations for Hylinae and Rana with both genes.
In this study, we tested for the presence of a phylogenetic signal in five acoustical features derived from recordings of anuran vocalizations. These acoustical parameters were size-regressed using the PGLS method to account for phylogeny . We first examined 90 species from four clades, using three different models of genetic distance based on 12S mtDNA sequences, after confirming that these sequences displayed little or no substitution saturation. A robust phylogenetic signal was detected in at least four acoustical parameters (CVA, logDF, logSF, and logTON) when considering all species. Intra-clade analyses were less consistent but nevertheless yielded multiple reliable results, such as a significant phylogenetic signal in CVA for Bufo and Hylinae, in logDF for Hylinae and Rana, and in logSF for Bufo. Notably, a multivariate trait computed from CVA, logDF, and logSF was significantly correlated with Abouheif’s matrix of phylogenetic proximities in all four clades.
We then investigated a subset of 49 species from three of these clades for which sequences were available for both rhodopsin and 12S mtDNA. Significant phylogenetic signals were observed over the entire subset for CVA, logDF, logSI, and logTON, for both 12S and rhodopsin. Significance tests on intra-clade analyses were generally in agreement between 12S and rhodopsin and between genetic distance models, and a reliable phylogenetic signal was found with both genes in CVA for Bufo, logDF for Hylinae and logTON for Rana.
Because the acoustic parameters used in this study were chosen on the basis of their ability to classify calls into one of the four clades studied here , the divergences between these parameters computed from advertisement calls of anurans belonging to different clades were expected to be greater than those computed from calls of anurans from the same clade. Indeed, among these five acoustic parameters, CVA, DF, and SF were previously found to be optimal in discriminating between the vocalizations of these clades , in line with the best-fitting multivariate trait distance obtained here. However, phylogenetic signals in some of these acoustic parameters were also observed at the intra-clade level in several cases, notably for CVA (Bufo and Hylinae) and logDF (Hylinae and Rana). These results do not follow automatically from our previous acoustic classification model, which only operated between clades, did not take into account body size and furthermore was generated without a priori knowledge of genetic distances . Thus, the significant intra-clade phylogenetic signals we found strongly suggest that phylogenetic proximity is, at least in these clades and for these acoustic features, a powerful and reliable predictor of the degree of acoustic similarity between advertisement calls of anurans.
A possible limitation of the approach followed in the present study is that it is based on only one or two recordings per species. Thus, the values for the acoustic parameters computed from the recordings represent a limited sample from a potentially broad range of values associated with our species. However, this is unlikely to significantly impact the general conclusions of our study, given that it is based on a large database and that it is concerned with detecting phylogenetic signals at the clade level. Additionally, the Mantel test on Abouheif’s proximities, in particular, has been shown to be relatively robust to sparse sampling and phenotypic variation within species . Nevertheless, follow-up studies may address this issue by sampling several recordings per species.
Polyploidy is fairly common in anurans, which can be problematic for phylogenetic analyses because polyploid taxa do not arise by ordinary cladogenesis . Furthermore, in the case of allopolyploidy, mitochondrial sequences, being inherited only from the maternal side, may not accurately reflect overall genetic distance. Because genetic distances were estimated from 12S mitochondrial DNA sequences in the present study, this could affect our analysis. However, only two species from our sample, Bufo viridis and Hyla chrysoscelis, are known polyploids (see , Additional file 1), suggesting that polyploidy had, at most, a limited impact on our findings.
Ryan  observed that characters regulated by behavior and physiology, such as call rate or amplitude modulation, are less conservative than characters that would require modifications of vocal morphology, such as some spectral characters. Supporting Ryan’s observation, and in line with our previous study , we found that the only parameters that correlated significantly with body size in our sample were DF and TON, spectral features that are likely to be largely determined by vocal morphology . However, after controlling for body size, we found that both spectral features and spectro-temporal features such as amplitude modulation (CVA in our analysis) and spectral flux exhibited reliable phylogenetic signals, suggesting that there is also a tendency toward phylogenetic conservation for acoustic characters whose link to morphology remains unclear (although Martin  related patterns of amplitude modulation in Bufo to the presence or absence of well-developed arytenoid valves). Furthermore, while body size may impose an indirect evolutionary constraint on acoustic parameters, and especially on spectral features such as DF [2, 5, 53], our results show that the trait similarity observed for these spectral features is not solely explained by body size similarity in related species. In that respect, Ryan  suggested that the evolution of call features such as DF might also be constrained by signaler-receiver interaction. For instance, modifications in the signal may be more likely to evolve in accordance with preexisting biases in the auditory system of the receiver . Thus, constraints at both production and perceptual levels may be operative. Finally, the values of Blomberg’s K were lower than 1 for all acoustic parameters and in all clades (see Table 1). Blomberg et al. observed that behavioral traits were generally more labile (i.e., characterized by K values below 1) than morphological or physiological traits (typically associated with higher K values). Our results thus indicate that the strength of the phylogenetic signal detected in our acoustic parameters is closer to that observed for behavioral traits than for morphological traits in Blomberg et al.’s study .
Given the central role played by advertisement calls in mate recognition and sexual selection in anurans, it is not surprising that bioacoustic features of these calls have been shown to be reliable taxonomical cues at the local species level. However, the rate at which these features evolve can vary between groups [13, 53, 55] and call evolution is not always closely related to phylogenetic distance. A recent study on African clawed frogs in fact reported very low phylogenetic signal in the underwater calls of this highly-derived clade . Our results here, studying airborne calls, clearly document the general tendency for more closely related species to exhibit more similar acoustic features, supporting the hypothesis that advertisement call acoustics are at least partly shaped by phylogenetic constraints. Of course, traditional research topics such as the potential influence of ecological constraints (e.g. the calling environment or the presence of sympatric species) [56–58], and other selective pressures, including especially sexual selection [14, 59, 60], on the evolution of anurans’ vocalizations should not be neglected. Nonetheless, our results suggest that phylogenetic constraints establish important and persistent ground rules, shaping the landscape within which adaptive calls evolution occurs. The clear phylogenetic signal uncovered in the present study, across a wide range of anurans with a broad geographical distribution, is consistent with this “phylogenetic constraint” hypothesis.
Because the analyses presented here were conducted using only a few basic parameters that are easily and automatically measured, and are applicable to a very wide range of sounds, the method described here may have broad biological relevance, applicable to many other clades. Our results indicate an robust relationship between acoustic similarity and genetic relatedness, even after partialling out the effect of body size, in a group of non-vocal learners whose vocalizations are fairly diverse, from the quasi-mechanical trills typical of many Bufo species to the whistle-like calls of Leptodactylus. It remains a topic for further research whether the acoustic parameters used to characterize advertisement calls in the present study would correlate with genetic distances for other anuran clades. Nevertheless, our results highlight the considerable potential of an approach based on versatile, low-level acoustic parameters, rather than handpicked, clade-specific characteristics. One key advantage of our approach is that it tends to be more resistant to over-fitting and does not require extensive parameterization or manual measurement of the acoustic properties of each vocalization, making these methods potentially useful for a wide range of vocalizations and animal species.
We found robust evidence for the presence of phylogenetic signals in several acoustic parameters derived from advertisement calls in a sample of 90 species representing four clades of anurans with a wide diversity of call acoustics. Moreover, these phylogenetic signals were detected both across the entire sample and, in several cases, within individual clades, using several measures of genetic distance. These results were generally validated on a subset of 49 species using genetic distances derived from both a highly variable mitochondrial region (12S mtDNA) and from a conserved protein-coding nuclear gene (rhodopsin). Because our methodology is based on general acoustic features found in most animal vocalizations, we suggest that the approach implemented here could fruitfully be applied to other vocal taxa.
Recordings and acoustic analysis
Recordings of anuran vocalizations were digitally copied as WAV files from commercially available CDs (see Additional file 1) comprising calls from North, South, and Central America, Europe, South Africa, Japan and Korea. Only advertisement calls representing single males were used. For each CD track, the longest continuous sequence containing only advertisement calls of a single male was extracted using the Praat software, version 5.1.44 . Two separate non-identical recordings (obtained from different CDs and thus presumably different individuals) were used for 41 species. Only one high-quality recording was available for the remaining 49 species.
The acoustic analysis was conducted using the MIR Toolbox 1.3.2 in MATLAB . Acoustic parameters were analyzed using a window of 40 ms with hop-size of 20 ms, according to the method described in .
Mean values for each recording were computed over all 40-ms frames for the following acoustic parameters: dominant frequency (DF), coefficient of variation of root-mean-square amplitude (CVA), spectral flux (SF), spectral flatness or tonality (TON), and spectral irregularity (SI). These five parameters were found to exhibit low multicollinearity in a set of 194 recordings that included the 131 recordings used in the current study . For the 41 species for which two recordings were available, the mean values were averaged over both recordings; the mean values computed from a single recording were used for the remaining 49 species. Detailed values for each recording, as well as the mean values for species for which two recordings were available, are given in Additional file 3.
Temperature data were available for some of the recordings and are provided in Additional file 3. These data were sparse and were not analyzed in the present study. In any case, temperature-induced effects on the acoustical properties of calls are presumably relatively small in comparison to taxonomical or ecological influences . Moreover, DF has rarely been shown to be affected by temperature .
Nucleotide sequence alignment
Our genetic dataset consisted of nucleotide sequences of approximately 930 bp from the 12S rRNA region of the mitochondrial genome (mtDNA) from 90 species belonging to four clades of anurans (Bufo, Hylinae, Leptodactylus, and Rana). As out-group, the same region of the mtDNA in the Siberian Salamander (Ranodon sibiricus, Family: Hynobiidae), California Tiger Salamander (Ambystoma californiens, Family: Ambystomatidae) and Western Spadefoot Toad (Pelobates cultripes, Family: Pelobatidae) was used.
In addition to the 12S mtDNA sequences, we analyzed a second dataset consisting of nucleotide sequences of 312 bp from exon 1 of the rhodopsin gene from 49 species belonging to three clades (Bufo, Hylinae and Rana). P. cultripes was used as out-group.
All 12S and rhodopsin sequences were obtained from GenBank, using search and extraction tools developed in Biopython , and the validity of obtained DNA sequences was checked using BlastSearch (National Center for Biotechnology Information). The nucleotide sequences were aligned in BioEdit software version 7.0 .
The Muscle (implemented in the software MEGA) , ClustalW2 , and Guidance methods  were all tested to obtain a reliable multiple sequence alignment. The ClustalW2 method was eventually used for the multiple sequence alignment on the 12S dataset with the following parameters: gap open penalties (GOP) = 10, gap extension penalties (GEP) = 0.2, gap distance = 5, and the UPGMA clustering method in BioEdit software version 7.0 . In the case of the rhodopsin gene, we used the Muscle Codon option, with GOP = −2.9, GEP = 0 and hydrophobicity multiplier = 1.2. Some sites were manually edited to maximize positional homology.
Phylogenetic tree reconstruction
The program jMODELTEST [69, 70] was used to identify the evolutionary models and other parameters. Based on the Akaike Information Criterion (AICc), the GTR+Γ+I model (general time-reversible)  with a proportion of invariant sites (28.6% invariant sites)  and rates at other sites varying according to a gamma distribution (gamma shape = 0.53, number of discrete gamma categories = 4) best described the 12S sequence data [73, 74]. However, the HKY model (Hasegawa-Kishino-Yano)  with gamma site heterogeneity model with 4 rate categories for the discrete approximation of the Gamma distribution rates among sites (HKY +Γ) (gamma shape = 0.31) best described the rhodopsin data.
The phylogenetic relationship among various anuran clades was reconstructed using the Monte-Carlo Markov Chain (MCMC) model implemented in BEAST version 1.6.2 . We used a Yule tree prior that assumes a constant (unknown) lineage birth rate for each branch in the tree. This model is suitable for trees describing the phylogenetic relationship among individuals from many different species. Wide uniform prior distributions were used as defaults. A relaxed clock model  with uncorrelated lognormal distribution was used to reconstruct the genetic phylogenies. The MCMC model was run for 30,000,000 generations (10,000,000 generations for rhodopsin) with the initial 1,000,000 steps discarded as burn-in. Trees and model parameters were sampled every 1000 steps thereafter. Effective Sample Sizes (ESS) for mean evolutionary rate, population size and posterior likelihood were found to be > 200 for all the models used. Subsequently, a single target tree with maximum clade credibility (MCC) and median node heights from a sample of trees produced by BEAST were constructed using TreeAnnotator v1.6.2 . The initial 1000 trees were discarded as burn-in. Abouheif proximities were computed from the phylogenetic trees using the function ‘proxTips’ in the package ‘adephylo’ in R .
Estimating pairwise genetic distances using maximum likelihood
Pairwise genetic distances for the 12S and rhodopsin nucleotide sequences were computed using the software MEGA . Both transition and transversion nucleotide substitutions were included in calculating the genetic distance. Gaps and missing data were treated with the partial deletion option, and the bootstrap method with 500 replications was employed.
Nucleotide substitution saturation test
The nucleotide substitution saturation test was conducted on the 12S mtDNA dataset using DAMBE [46, 47]. Because the test can only run on 32 taxonomic units and our sample included 90 species, 10,000 replications with random resampling of subsets of 4, 8, 16, and 32 species were conducted, following the methodology described in .
SVL values in mm for males, taken from the literature, were used as an estimate of male body size . Additional file 3 provides the median male SVL values for all 90 species included in this study, and Additional file 1 lists the sources used to obtain these data.
The phylogenetically corrected residuals of the acoustic parameters were computed using the ‘phyl_resid.R’ function in R . Following the methodology described in , the ‘phyl_resid.R’ function used a variance-covariance matrix representing phylogenetic relatedness under the Brownian trait evolution model, which was computed from the phylogenetic tree using the function ‘vcv.phylo’ in the package ‘ape’ in R . Multivariate trait Euclidean distances were computed on the standardized size-corrected residuals of the acoustic parameters using the ‘dist’ function from the package ‘stats’ in R (standardization was applied to give equal weight to all parameters).
The standard deviate of Moran’s autocorrelation coefficient  was computed using the ‘moran.test’ function from package ‘spdep’ in R . Statistical significance was assessed by Monte-Carlo simulations using the ‘moran.mc’ function from the package ‘spdep’. 10,000 permutations were conducted in each case. Geary’s C was computed in the same way using the functions ‘geary.test’ and ‘geary.mc’, respectively, from the package ‘spdep’. To confirm our analyses on Abouheif’s proximities conducted with the function ‘moran.test’, we used the function ‘abouheif.moran’ from the package ‘adephylo’ in R , and obtained nearly identical results (within rounding error).
Blomberg’s K was computed using the function ‘phylosig’ from the package ‘phytools’ in R . 10,000 permutations were conducted for significance tests.
Mantel tests were conducted using the ‘mantel’ function from the package ‘ecodist’ in R . 10,000 permutations were conducted for significance tests. To determine the best fit for multivariate trait Euclidean distances, a forward stepwise procedure was conducted in which an acoustic parameter was added to the multivariate trait distance model only if the Mantel correlation coefficient obtained on this distance model was significantly higher than the Mantel correlation coefficient obtained on a trait distance model that did not include this parameter. Because there is no formal test of significance to compare two Mantel correlation coefficients, we used confidence intervals as estimated by bootstrapping (the confidence intervals were also estimated with the ‘mantel’ function from the package ‘ecodist’). A given Mantel coefficient ‘A’ was considered as significantly higher than a coefficient ‘B’ when the value of ‘A’ was higher than 95% of the values of ‘B’ obtained on 10,000 bootstrapping iterations. Note that, to avoid over-fitting, the forward stepwise procedure was only conducted on the entire sample (90 species for the data presented in Table 1, and 49 species for the data presented in Table 2) and not on individual clades.
This research was supported by ERC Advanced Grant SOMACCA No. 230604 and a University of Vienna startup grant to W.T.F. We thank Walter Hödl, Carl Gerhardt, and four anonymous reviewers for very useful comments and suggestions, Jinook Oh for help with the Korean and Japanese frog CD notes, and Paolo Vatta for help with Biopython.
- Duellman W, Trueb L: Biology of amphibians. 1994, Baltimore: Johns Hopkins University PressGoogle Scholar
- Gerhardt HC: The evolution of vocalization in frogs and toads. Annu Rev Ecol Syst. 1994, 25: 293-324. 10.1146/annurev.es.25.110194.001453.View ArticleGoogle Scholar
- Dawson B, Ryan MJ: Early experience leads to changes in the advertisement calls of male Physalaemus pustulosus. Copeia. 2009, 2009: 221-226. 10.1643/CE-07-254.View ArticleGoogle Scholar
- Littlejohn MJ: Call differentiation in a complex of seven species of Crinia (Anura, Leptodactylidae). Evolution. 1959, 13: 452-268. 10.2307/2406128.View ArticleGoogle Scholar
- Ryan MJ: Constraints and patterns in the evolution of anuran acoustic communication. The evolution of the amphibian auditory system. Edited by: Fritzch B, Ryan MJ, Wilczynski W, Hetherington TE, Walkowiak W. 1988, New York: Wiley, 637-677.Google Scholar
- Blair WF: Mating call in the speciation of anuran amphibians. Am Nat. 1958, 92: 27-51. 10.1086/282007.View ArticleGoogle Scholar
- Tobias ML, Evans BJ, Kelley DB: Evolution of advertisement calls in African clawed frogs. Behaviour. 2011, 148: 519-549. 10.1163/000579511X569435.PubMed CentralPubMedView ArticleGoogle Scholar
- Wollenberg KC, Glaw F, Meyer A, Vences M: Molecular phylogeny of Malagasy reed frogs, Heterixalus, and the relative performance of bioacoustics and color-patterns for resolving their systematics. Mol Phylogenet Evol. 2007, 45: 14-22. 10.1016/j.ympev.2007.06.024.PubMedView ArticleGoogle Scholar
- Lehtinen RM, Wojtowicz EA, Hailey A: Male vocalizations, female discrimination and molecular phylogeny: multiple perspectives on the taxonomic status of a critically endangered Caribbean frog. J Zool. 2011, 283: 117-125. 10.1111/j.1469-7998.2010.00752.x.View ArticleGoogle Scholar
- Goicoechea N, De La Riva I, Padial JM: Recovering phylogenetic signal from frog mating calls. Zool Scr. 2010, 39: 141-154. 10.1111/j.1463-6409.2009.00413.x.View ArticleGoogle Scholar
- Thinh V, Hallam C, Roos C, Hammerschmidt K: Concordance between vocal and genetic diversity in crested gibbons. BMC Evol Biol. 2011, 11: 36-10.1186/1471-2148-11-36.PubMed CentralPubMedView ArticleGoogle Scholar
- Cap H, Deleporte P, Joachim J, Reby D: Male vocal behavior and phylogeny in deer. Cladistics. 2008, 24: 917-931. 10.1111/j.1096-0031.2008.00223.x.View ArticleGoogle Scholar
- Cocroft RB, Ryan MJ: Patterns of advertisement call evolution in toads and chorus frogs. Anim Behav. 1995, 49: 283-303. 10.1006/anbe.1995.0043.View ArticleGoogle Scholar
- Ryan MJ: Sexual selection and communication in frogs. Trends Ecol Evol. 1991, 6: 351-355. 10.1016/0169-5347(91)90224-L.PubMedView ArticleGoogle Scholar
- Darst CR, Cannatella DC: Novel relationships among hyloid frogs inferred from 12S and 16S mitochondrial DNA sequences. Mol Phylogenet Evol. 2004, 31: 462-475. 10.1016/j.ympev.2003.09.003.PubMedView ArticleGoogle Scholar
- Shaw KL: Conflict between nuclear and mitochondrial DNA phylogenies of a recent species radiation: What mtDNA reveals and conceals about modes of speciation in Hawaiian crickets. Proc Natl Acad Sci USA. 2002, 99: 16122-16127. 10.1073/pnas.242585899.PubMed CentralPubMedView ArticleGoogle Scholar
- Hurst GD, Jiggins FM: Problems with mitochondrial DNA as a marker in population, phylogeographic and phylogenetic studies: the effects of inherited symbionts. Proc Biol Sci. 2005, 272: 1525-1534. 10.1098/rspb.2005.3056.PubMed CentralPubMedView ArticleGoogle Scholar
- Frost DR, Grant T, Faivovich J, Bain RH, Haas A, Haddad CFB, de Sá RO, Channing A, Wilkinson M, Donnellan SC, Raxworthy CJ, Campbell JA, Blotto BL, Moler P, Drewes RC, Nussbaum RA, Lynch JD, Green DM, Wheeler WC: The amphibian tree of life. Bull Am Mus Nat Hist. 2006, 297: 1-370.View ArticleGoogle Scholar
- Pyron RA, Wiens JJ: A large-scale phylogeny of Amphibia including over 2800 species, and a revised classification of extant frogs, salamanders, and caecilians. Mol Phylogenet Evol. 2011, 61: 543-583. 10.1016/j.ympev.2011.06.012.PubMedView ArticleGoogle Scholar
- Hertwig S, De Sá RO, Haas A: Phylogenetic signal and the utility of 12S and 16S mtDNA in frog phylogeny. J Zoolog Syst Evol Res. 2008, 42: 2-18. 10.1111/j.1439-0469.2004.00225.x.View ArticleGoogle Scholar
- Springer MS, DeBry RW, Douady C, Amrine HM, Madsen O, de Jong WW, Stanhope MJ: Mitochondrial versus nuclear gene sequences in deep-level mammalian phylogeny reconstruction. Mol Biol Evol. 2001, 18: 132-143. 10.1093/oxfordjournals.molbev.a003787.PubMedView ArticleGoogle Scholar
- Hoegg S, Vences M, Brinkmann H, Meyer A: Phylogeny and comparative substitution rates of frogs inferred from sequences of three nuclear genes. Mol Biol Evol. 2004, 21: 1188-1200. 10.1093/molbev/msh081.PubMedView ArticleGoogle Scholar
- Gingras B, Fitch WT: A three-parameter model for classifying anurans into four genera based on advertisement calls. J Acoust Soc Am. 2013, 133: 547-559. 10.1121/1.4768878.PubMedView ArticleGoogle Scholar
- Jensen K, Andersen TH: Real-time beat estimation using feature extraction. Computer Music Modeling and Retrieval, Volume 2771. Edited by: Wiil UK. 2004, Berlin, Heidelberg: Springer Berlin Heidelberg, 13-22.View ArticleGoogle Scholar
- Dubnov S: Generalization of spectral flatness measure for non-Gaussian linear processes. IEEE Signal Process Lett. 2004, 11: 698-701. 10.1109/LSP.2004.831663.View ArticleGoogle Scholar
- Amphibian species of the world: an online reference. Version 5.5. http://research.amnh.org/vz/herpetology/amphibia/,
- Fitch WT, Hauser MD: Unpacking “honesty”: vertebrate vocal production and the evolution of acoustic signals. Acoustic Communication. Edited by: Simmons AM, Fay RF, Popper AN. 2002, New York: Springer, 65-137.Google Scholar
- Gerhardt HC, Huber F: Acoustic Communication in Insects and Anurans: Common Problems and Diverse Solutions. 2002, Chicago: University of Chicago PressGoogle Scholar
- Gingras B, Boeckle M, Herbst CT, Fitch WT: Call acoustics reflect body size across four clades of anurans. J Zool. 2013, 289: 143-150. 10.1111/j.1469-7998.2012.00973.x.View ArticleGoogle Scholar
- Freckleton RP, Jetz W: Space versus phylogeny: disentangling phylogenetic and spatial signals in comparative data. Proc Biol Sci. 2009, 276: 21-30. 10.1098/rspb.2008.0905.PubMed CentralPubMedView ArticleGoogle Scholar
- Blomberg SP, Garland T, Ives AR: Testing for phylogenetic signal in comparative data: behavioral traits are more labile. Evolution. 2003, 57: 717-745.PubMedView ArticleGoogle Scholar
- Felsenstein J: Phylogenies and the comparative method. Am Nat. 1985, 125: 1-15. 10.1086/284325.View ArticleGoogle Scholar
- Grafen A: The phylogenetic regression. Phil Trans R Soc B. 1989, 326: 119-157. 10.1098/rstb.1989.0106.PubMedView ArticleGoogle Scholar
- Revell LJ: Size-correction and principal components for interspecific comparative studies. Evolution. 2009, 63: 3258-3268. 10.1111/j.1558-5646.2009.00804.x.PubMedView ArticleGoogle Scholar
- Rohlf FJ: A comment on phylogenetic correction. Evolution. 2006, 60: 1509-1515.PubMedView ArticleGoogle Scholar
- Stone EA: Why the phylogenetic regression appears robust to tree misspecification. Syst Biol. 2011, 60: 245-260. 10.1093/sysbio/syq098.PubMed CentralPubMedView ArticleGoogle Scholar
- Moran PAP: Notes on continuous stochastic phenomena. Biometrika. 1950, 37: 17-23.PubMedView ArticleGoogle Scholar
- Mantel N: The detection of disease clustering and a generalized regression approach. Cancer Res. 1967, 27: 209-220.PubMedGoogle Scholar
- Borcard D, Legendre P: Is the Mantel correlogram powerful enough to be useful in ecological analysis? A simulation study. Ecology. 2012, 93: 1473-1481. 10.1890/11-1737.1.PubMedView ArticleGoogle Scholar
- Hardy OJ, Pavoine S: Assessing phylogenetic signal with measurement error: A comparison of Mantel tests, Blomberg et al.’s K, and phylogenetic distograms. Evolution. 2012, 66: 2614-2621. 10.1111/j.1558-5646.2012.01623.x.PubMedView ArticleGoogle Scholar
- Diniz-Filho JAF, Santos T, Rangel TF, Bini LM: A comparison of metrics for estimating phylogenetic signal under alternative evolutionary models. Genet Mol Biol. 2012, 35: 673-679. 10.1590/S1415-47572012005000053.PubMed CentralPubMedView ArticleGoogle Scholar
- Pavoine S, Ricotta C: Testing for phylogenetic signal in biological traits: The ubiquity of cross-product statistics. Evolution. 2013, 67: 828-840. 10.1111/j.1558-5646.2012.01823.x.PubMedView ArticleGoogle Scholar
- Abouheif E: A method for testing the assumption of phylogenetic independence in comparative data. Evol Ecol Res. 1999, 1: 895-909.Google Scholar
- Pavoine S, Ollier S, Pontier D, Chessel D: Testing for phylogenetic signal in phenotypic traits: new matrices of phylogenetic proximities. Theor Popul Biol. 2008, 73: 79-91. 10.1016/j.tpb.2007.10.001.PubMedView ArticleGoogle Scholar
- Münkemüller T, Lavergne S, Bzeznik B, Dray S, Jombart T, Schiffers K, Thuiller W: How to measure and test phylogenetic signal. Methods Ecol Evol. 2012, 3: 743-756. 10.1111/j.2041-210X.2012.00196.x.View ArticleGoogle Scholar
- Xia X, Xie Z, Salemi M, Chen L, Wang Y: An index of substitution saturation and its application. Mol Phylogenet Evol. 2003, 26: 1-7. 10.1016/S1055-7903(02)00326-3.PubMedView ArticleGoogle Scholar
- Xia X, Lemey P: Assessing substitution saturation with DAMBE. The Phylogenetic Handbook: A Practical Approach to DNA and Protein Phylogeny. Edited by: Lemey P, Salemi M, Vandamme A-M. 2009, Cambridge: Cambridge University Press, 615-630. 2View ArticleGoogle Scholar
- Legendre P, Fortin MJ: Spatial pattern and ecological analysis. Plant Ecol. 1989, 80: 107-138. 10.1007/BF00048036.View ArticleGoogle Scholar
- Geary RC: The contiguity ratio and statistical mapping. Inc Stat. 1954, 5: 115-127.Google Scholar
- Dutilleul P, Stockwell JD, Frigon D, Legendre P: The Mantel test versus Pearson’s correlation analysis: Assessment of the differences for biological and environmental studies. J Agric Biol Environ Stat. 2000, 5: 131-150. 10.2307/1400528.View ArticleGoogle Scholar
- Mable BK, Alexandrou MA, Taylor MI: Genome duplication in amphibians and fish: an extended synthesis. J Zool. 2011, 284: 151-182. 10.1111/j.1469-7998.2011.00829.x.View ArticleGoogle Scholar
- Martin WF: Evolution of vocalizations in the genus Bufo. Evolution in the Genus Bufo. Edited by: Blair WF. 1972, Austin, Texas: University of Texas, 279-309.Google Scholar
- Ryan MJ: Factors influencing the evolution of acoustic communication: biological constraints. Brain Behav Evol. 1986, 28: 70-82. 10.1159/000118693.PubMedView ArticleGoogle Scholar
- Ryan MJ, Fox JH, Wilczynski W, Rand AS: Sexual selection for sensory exploitation in the frog Physalaemus pustulosus. Nature. 1990, 343: 66-67. 10.1038/343066a0.PubMedView ArticleGoogle Scholar
- Richards CL: Has the evolution of complexity in the amphibian papilla influenced anuran speciation rates?. J Evol Biol. 2006, 19: 1222-1230. 10.1111/j.1420-9101.2006.01079.x.PubMedView ArticleGoogle Scholar
- Zimmerman BL: A comparison of structural features of calls of open and forest habitat frog species in the Central Amazon. Herpetologica. 1983, 39: 235-246.Google Scholar
- Kime NM, Turner WR, Ryan MJ: The transmission of advertisement calls in Central American frogs. Behav Ecol. 2000, 11: 71-83. 10.1093/beheco/11.1.71.View ArticleGoogle Scholar
- Amezquita A, Flechas SV, Lima AP, Gasser H, Hodl W: Acoustic interference and recognition space within a complex assemblage of dendrobatid frogs. Proc Natl Acad Sci USA. 2011, 108: 17058-17063. 10.1073/pnas.1104773108.PubMed CentralPubMedView ArticleGoogle Scholar
- Bee MA, Perrill SA, Owen PC: Male green frogs lower the pitch of acoustic signals in defense of territories: a possible dishonest signal of size?. Behav Ecol. 2000, 11: 169-177. 10.1093/beheco/11.2.169.View ArticleGoogle Scholar
- Wagner WE: Deceptive or honest signalling of fighting ability? A test of alternative hypotheses for the function of changes in call dominant frequency by male cricket frogs. Anim Behav. 1992, 44: 449-462.View ArticleGoogle Scholar
- Boersma P, Weenink D: Praat. doing phonetics by computer (Version 5.1.44), [Computer program]. [http://www.praat.org/]Google Scholar
- Lartillot O, Toiviainen P, Eerola T: A Matlab toolbox for music information retrieval. Data Analysis, Machine Learning and Applications. Edited by: Preisach C, Burkhardt H, Schmidt-Thieme L, Decker R. 2008, Berlin, Heidelberg: Springer Berlin Heidelberg, 261-268.View ArticleGoogle Scholar
- Hoskin CJ, James S, Grigg GC: Ecology and taxonomy-driven deviations in the frog call-body size relationship across the diverse Australian frog fauna. J Zool. 2009, 278: 36-41. 10.1111/j.1469-7998.2009.00550.x.View ArticleGoogle Scholar
- Cock PJA, Antao T, Chang JT, Chapman BA, Cox CJ, Dalke A, Friedberg I, Hamelryck T, Kauff F, Wilczynski B, de Hoon MJL: Biopython: freely available Python tools for computational molecular biology and bioinformatics. Bioinformatics. 2009, 25: 1422-1423. 10.1093/bioinformatics/btp163.PubMed CentralPubMedView ArticleGoogle Scholar
- Hall TA: BioEdit: a user-friendly biological sequence alignment editor for Windows 95/98/NT. Nucleic Acids Symp Ser. 1999, 41: 95-98.Google Scholar
- Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S: MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Mol Biol Evol. 2011, 28: 2731-2739. 10.1093/molbev/msr121.PubMed CentralPubMedView ArticleGoogle Scholar
- Thompson JD, Higgins DG, Gibson TJ: CLUSTAL W: Improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice. Nucleic Acids Res. 1994, 22: 4673-4680. 10.1093/nar/22.22.4673.PubMed CentralPubMedView ArticleGoogle Scholar
- Penn O, Privman E, Ashkenazy H, Landan G, Graur D, Pupko T: GUIDANCE: A web server for assessing alignment confidence scores. Nucleic Acids Res. 2010, 38: W23-W28. 10.1093/nar/gkq443.PubMed CentralPubMedView ArticleGoogle Scholar
- Guindon S, Gascuel O: A simple, fast, and accurate algorithm to estimate large phylogenies by maximum likelihood. Syst Biol. 2003, 52: 696-704. 10.1080/10635150390235520.PubMedView ArticleGoogle Scholar
- Posada D: jModelTest: phylogenetic model averaging. Mol Biol Evol. 2008, 25: 1253-1256. 10.1093/molbev/msn083.PubMedView ArticleGoogle Scholar
- Rodríguez F, Oliver JL, Marín A, Medina JR: The general stochastic model of nucleotide substitution. J Theor Biol. 1990, 142: 485-501. 10.1016/S0022-5193(05)80104-3.PubMedView ArticleGoogle Scholar
- Gu X, Fu YX, Li WH: Maximum likelihood estimation of the heterogeneity of substitution rate among nucleotide sites. Mol Biol Evol. 1995, 12: 546-557.PubMedGoogle Scholar
- Yang Z: Maximum-likelihood estimation of phylogeny from DNA sequences when substitution rates differ over sites. Mol Biol Evol. 1993, 10: 1396-1401.PubMedGoogle Scholar
- Yang Z: Maximum likelihood phylogenetic estimation from DNA sequences with variable rates over sites: approximate methods. J Mol Evol. 1994, 39: 306-314. 10.1007/BF00160154.PubMedView ArticleGoogle Scholar
- Hasegawa M, Kishino H, Yano T: Dating of the human-ape splitting by a molecular clock of mitochondrial DNA. J Mol Evol. 1985, 22: 160-174. 10.1007/BF02101694.PubMedView ArticleGoogle Scholar
- Drummond AJ, Rambaut A: BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evol Biol. 2007, 7: 214-10.1186/1471-2148-7-214.PubMed CentralPubMedView ArticleGoogle Scholar
- Drummond AJ, Ho SYW, Phillips MJ, Rambaut A: Relaxed phylogenetics and dating with confidence. PLoS Biol. 2006, 4: e88-10.1371/journal.pbio.0040088.PubMed CentralPubMedView ArticleGoogle Scholar
- Drummond AJ, Ho SYW, Rawlence N, Rambaut A: A rough guide to BEAST 1.4. 2007, University of Auckland: Auckland, New ZealandGoogle Scholar
- Jombart T, Balloux F, Dray S: adephylo: new tools for investigating the phylogenetic signal in biological traits. Bioinformatics. 2010, 26: 1907-1909. 10.1093/bioinformatics/btq292.PubMedView ArticleGoogle Scholar
- Paradis E, Claude J, Strimmer K: APE: Analyses of Phylogenetics and Evolution in R language. Bioinformatics. 2004, 20: 289-290. 10.1093/bioinformatics/btg412.PubMedView ArticleGoogle Scholar
- Bivand R: Implementing spatial data analysis software tools in R. Geogr Anal. 2006, 38: 23-40. 10.1111/j.0016-7363.2005.00672.x.View ArticleGoogle Scholar
- Revell LJ: phytools: an R package for phylogenetic comparative biology (and other things). Methods Ecol Evol. 2012, 3: 217-223. 10.1111/j.2041-210X.2011.00169.x.View ArticleGoogle Scholar
- Goslee SC, Urban DL: The ecodist package for dissimilarity-based analysis of ecological data. J Stat Softw. 2007, 22: 1-19.View ArticleGoogle Scholar
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.