Mitogenome and phylogenetic analyses support rapid diversification among species groups of small-eared shrews genus Cryptotis (Mammalia: Eulipotyphla: Soricidae).

The small-eared shrew genus Cryptotis is the third largest in the family Soricidae and occurs in North, Central, and northern South America. In Mexico and Central and South America, most species inhabit geographically isolated moist, montane habitats at middle and high elevations in a typical sky-island pattern. The 49 recognized species have been partitioned into as many as six species groups based on morphological and molecular phylogenetic studies. The relationships among these species groups are poorly resolved, and their evolutionary histories, including migration patterns and locomotor adaptations, remain unclear. Herein, we provide a new phylogeny incorporating complete mitochondrial genomes (mitogenomes) and supermatrix approach. We compared different evolutionary scenarios using approximately unbiased (AU), Kishino-Hasegawa (KH), and Shimodaira-Hasegawa (SH) statistical tests. The phylogenetic hypothesis based on mitogenomes revealed novel relationships supporting a basal position for the Cryptotis parvus-group in the genus, and a close relationship between C. gracilis and one clade of the C. thomasi-group. The former relationship is consistent with the least derived humerus morphology and northern distribution of the species. The latter relationship implies multiple migrations between Central and South America. The lack of fine resolution for the species group relationships may be due partly to the lack of taxon sampling. In contrast, multi-approach analyses suggest that the unresolved relationships may be a result of rapid diversification during the early stages of Cryptotis evolution.

INTRODUCTION

A robust phylogeny is a fundamental tool for testing evolutionary hypotheses, but one that remains unavailable for most animal taxa as few species have been sequenced and even fewer have more than a few short gene fragments available. Although some universally used gene fragments (e.g., cyt b and COI for mammals) provide genetic barcode-like information and are helpful for molecular identification of new and existing species, they often provide insufficient signals for phylogenetic reconstruction because of low mutation or high saturation rates.

The small-eared shrews (Mammalia: Eulipotyphla: Soricidae:Cryptotis) consist of 49 recognized species that are widely distributed across North America through Central America and the northern montane areas of South America (Woodman, 2019). In the last 20 years, systematic studies increased the number of species by 50%, including the discovery of at least 15 species in isolated mountain areas where they are vulnerable to habitat alteration and climate change. Cryptotis is now the third largest genus of the shrew family Soricidae after Crocidura (ca. 198 species) and Sorex (ca. 86 species) (Burgin & He, 2018). The recognized species of Cryptotis have been variously partitioned into three (Choate, 1970), four (Woodman & Timm, 1998), five (He et al., 2015), or six (Woodman, 2019) species groups based on morphological and poorly-resolved molecular phylogenetic studies (Supplementary Table S1). Species in two of these groups, i.e., Cryptotis goodwini-group andCryptotis goldmani-group, have evolved morphologies consistent with different levels of adaptation to a semi-fossorial lifestyle (Guevara, 2017; He et al., 2015; Woodman & Gaffney, 2014). At present, however, the relationships between and within these groups are incompletely resolved, and it remains unclear to what extent they evolved independently and convergently (He et al., 2015).

In the two most recent comprehensive phylogenetic studies of Cryptotis (Baird et al., 2018; He et al., 2015), support values among species groups were all below the level of significance (Bayesian posterior probabilities (PP)<0.9 and Maximum likelihood bootstrap support (BS)<50). In part, this likely reflects the fact that most species are represented by only one (in most cases mitochondrial cytb) or a few genes, making it difficult to resolve relationships except among closely related species. Those genetic relationships that have been revealed, however, have provided new and unexpected phylogenetic relationships that conflict with the traditional species-grouping hypotheses (He et al., 2015; Woodman, 2019). These results suggest that the genus has a highly complex evolutionary history.

Most species of Cryptotis inhabit humid montane habitats at mid to high elevations. Several species in Mexico and Central and South America occupy relatively small areas isolated by intervening lowlands, thus forming sky islands (He & Jiang, 2014; Heald, 1951; McCormack et al., 2009). The evolutionary history in this genus is likely affected by complex topography and periodic climatic fluctuations, resulting in migration and isolation. In the current study, we investigated whether complete mitochondrial genomes (mitogenomes) can be used to inform and stabilize the phylogeny of Cryptotis and our understanding of its evolutionary history. We also used a supermatrix approach to construct a new and more complete phylogeny of the genus.

MATERIALS AND METHODS

Sampling and mitogenome sequencing

We sequenced 13 samples representing 11 recognized species of Cryptotis. We included two samples each of Cryptotis merriami and Cryptotis parvus as these species were previously discovered to contain two genetically distinct, but morphologically cryptic, lineages (He et al., 2015). Tissue samples were obtained from the Center for Conservation Genomics and the National Museum of Natural History, Smithsonian Institution, Washington, DC (loan no.: #2067019; Supplementary Table S2). We used a capture hybridization approach to obtain the complete mitogenomes (Chen et al., 2018). In brief, we extracted total DNA using a DNeasy Blood & Tissue Kit (Qiagen, USA) and sheared the DNA into small fragments to generate genomic DNA libraries. We generated biotin-linked homemade mitogenome probes using long-range polymerase chain reaction (PCR) amplicons. The DNA libraries and probes were incubated to capture mitochondrial libraries. The enriched libraries were amplified and sequenced using the Illumina high-throughput sequencing platform. We used FastQC v0.11.9 ( Andrews, 2010) and Trimmomatic v0.32 (Bolger et al., 2014) for quality control and data trimming, respectively, and mapped the reads to mitogenomes of Blarina brevicauda,Blarina hylophaga, Blarinella quadraticauda, Blarinella wardi, and Pantherina griselda using Geneious R11 v11.05 (Biomatters Ltd., New Zealand) (Ripma et al., 2014). Blarina is the sister genus to Cryptotis, and Blarinella+Pantherina is sister to Blarina+Cryptotis (He et al., 2018). We mapped the reads to each of the reference mitogenomes iteratively up to 25 times before generating the consensus sequences (Kearse et al., 2012). To confirm and improve our assemblies (i.e., reconciliation; Zimin et al., 2008), we repeated the mapping three times, and aligned the 15 consensus mitogenomes. Any missing data and incongruent positions were carefully checked by eye before generating the final consensus assemblies. Finally, we used the annotation transferring function in Geneious R11 to generate annotations for each mitogenome. The newly obtained mitogenomes are deposited in GenBank (accession Nos. MZ457409-MZ457421).

Mitogenome phylogeny and hypothesis testing

We used the 12S, 16S rRNA, and coding genes (except ND6, which is on the light chain of the mitogenome) to estimate maximum-likelihood (ML) and Bayesian gene trees (Duchêne et al., 2011). We included one mitogenome per species representing B. brevicauda, B. hylophaga, Bl. quadraticauda, Bl. wardi, and P. griselda as the outgroups for Cryptotis, and included six Sorex species and six Nectogalini shrews for comparison. We used Crocidura palawanensis andSuncus murinus, representing two ancient lineages in Crocidurinae, to root the tree (Supplementary Table S2; Hutterer et al., 2018). We grouped the data by gene and codon positions and used PartitionFinder v2.1.1 (Lanfear et al., 2017) to determine the best partitioning scheme under the GTR+G model based on a greedy algorithm, resulting in a nine-partition scheme (Supplementary Table S3). RAxML v8.2.12 (Stamatakis, 2014) was used to estimate the ML tree and CIPRES Science Gateway was used for implementation (Miller et al., 2015). We conducted rapid bootstrap analysis and searched for the best-scoring ML tree without the use of the BFGS searching algorithm (parameter: -f a --no-bfgs). BEAST v2.6 was used to estimate the Bayesian gene tree (Bouckaert et al., 2014). Mitogenome alignment was partitioned as mentioned above. We used a relaxed lognormal clock model, a Birth-Death model for the tree prior, and ran Markov chain Monte Carlo (MCMC) simulations for 50 million generations, with sampling every 5 000 generations. Analyses were conducted twice and Tracer v1.7 was used to examine the posterior distribution of each parameter in the log file to ensure that analyses reached a stationary state. The first 15% of MCMC samples were removed before the generation of the consensus tree. The RAxML and BEAST trees were identical except for two poorly supported nodes (see Results section). We also tested several alternative partitioning schemes and compared the results using the Shimodaira-Hasegawa (SH) test in RAxML (parameter: -f H) (Shimodaira & Hasegawa, 2001).

To examine potentially conflicting phylogenetic signals between genes, we calculated the partitioning Bremer support for each mitochondrial gene on each internal node of the best RAxML gene tree (Baker & Desalle, 1997) using PAUP v4.3.99.169.0 (Swofford, 2003) and a Tcl script (Göker et al., 2009). To test the rapid diversification hypothesis, we collapsed the two poorly supported nodes (i.e., C. mexicanus and C. goldmani; see Results section) on our best ML tree using TreeView v1.66 and generated all 15 possible dichotomic trees using the function “resolve_polytomy” in the ETE Toolkit v3.0 (Huerta-Cepas et al., 2016). We then calculated the site-wise log-likelihood supports using RAxML (-f G) and performed approximately unbiased (AU) (Shimodaira, 2002), Kishino-Hasegawa (KH) (Kishino & Hasegawa, 1989), and SH tests using CONSEL v0.20 (Shimodaira & Hasegawa, 2001).

Sequence data matrix and hypothesis testing

We downloaded two mitochondrial (cyt b, 16S rRNA) and two nuclear (ApoB and BRCA1) genes of the Cryptotis species and all outgroup species included in the mitogenome analyses from GenBank. We used these four genes and the mitogenomes to generate a gene matrix. The sequence manipulation is described in detail in the Supplementary Text. We obtained a sequence matrix of 136 samples and 14 955 bp in the alignment. The mitochondrial genes were partitioned as mentioned above and each nuclear gene was considered as one partition. To test alternative scenarios regarding the monophyly of known species groups, we constrained their monophyly, estimated the ML tree using RAxML, and performed AU, KH, and SH analyses as mentioned above.

RESULTS

Mitogenome phylogeny

The sequence alignments were uploaded to the GitHub repository (github.com/yinbinqiu/Cryptotis_phy). The mitogenome ML and Bayesian trees were congruent and generally well supported (BS>90, PP>0.95;Figure 1A). Relationships among all Sorex species were strongly supported, as were the relationships between Blarinellini and Blarinini and between Blarina and Cryptotis (i.e., BS=100, PP=1.0). Within Cryptotis, there was strong support for the monophyly of the C. nigrescens-group (C. nigrescens, C. mayensis, and C. merriami), C. parvus-group (C. parvus and C. tropicalis), and for a clade consisting of three species of the C. goodwini-group (C. lacertosus, C. mam, and C. oreoryctes), which has been recovered previously (Baird et al., 2018; He et al., 2015). Our analysis recovered two novel relationships supporting: (i) the C. parvus-group as one of the first branches in Cryptotis (BS=100, PP=1.0), and (ii) a close relationship between the C. gracilis and C. nigrescens-groups (BS=93, PP=1.0). These two relationships have not been observed in previous studies. He et al. (2015) and Baird et al. (2018) showed the C. mexicanus group (including C. magnus+C. phillipsii) as the first but weakly supported branch. Cryptotis gracilis was previously embedded in theC. goldmani-group but only supported by Bayesian analysis, which can be overestimated (Baird et al., 2018; He et al., 2015). We also observed several unresolved relationships: (i) the sister relationship of Soricini and Nectogalini+Anourosoricini in our outgroups was moderately supported, although this was not unexpected (see He et al., 2021); (ii) the relationships among Nectogalini water shrews were not resolved as observed by He et al. (2010); and (iii) within Cryptotis, the phylogenetic positions of two species, C. mexicanus and C. goldmani, were not resolved (BS<56, PP<0.94).

Figure 1

RAxML phylogenetic tree of Soricidae shrews estimated using mitogenomes (A) and concatenated alignment of mitochondrial and nuclear genes (B)

Branch lengths represent substitutions/site. A: Branch numbers refer to RAxML bootstrap support values (BS) and BEAST posterior probabilities (PP). Support values for a node with BS of 100 and PP of 1.0 are not shown. “-” indicates relationship was not recovered in BEAST analyses. All Cryptotis mitogenomes were obtained in the present study. B: Samples in a monophyletic clade representing the same species are collapsed. Branch numbers refer to BS. BS=100 is not shown and BS<50 is indicated as “-”. Species groups are based on Woodman (2019). Samples were assigned to eight clades best reflecting species groups (Supplementary Table S1). C. goodwini-group, C. mexicanus-group, and C. thomasi-group are paraphyletic, and C. gracilis is supported within a C. thomasi-group clade.

Because a suboptimal partitioning scheme may produce highly supported but incorrect nodes in a tree (Kainer & Lanfear, 2015), we first tested whether the phylogenetic relationships, especially the novel relationships, were due to a suboptimal partitioning scheme. We estimated the ML trees using two alternative partitioning schemes (Supplementary Table S3). All analyses resulted in similar topologies and did not affect the phylogenetic positions of either the C. parvus-group or C. gracilis (data not shown), suggesting that the novel relationships were not attributable to the partitioning scheme but may be the result of conflicting or poor phylogenetic signaling.

We next examined whether the unresolved relationships of C. mexicanus and C. goldmani may be due to conflicting phylogenetic signals or a lack of any signal based on partitioned Bremer values of all genes for each node (Supplementary Figure S1). None of the genes supported or rejected the position of C. goldmani (i.e., partitioned Bremer support (PBS)=0; Supplementary Table S4). The phylogenetic position of C. mexicanus was supported by seven genes, including the 12S rRNA and six coding genes (1≤PBS≤9), but was rejected by the other seven genes (−5≤PBS≤−1). We therefore split the mitogenome alignment into two sub-datasets based on PBS support for the position of C. mexicanus (i.e., PBS+ and PBS− alignments) and estimated the best ML trees individually. While the genes characterized by PBS+ recovered the same topology among Cryptotis species (Supplementary Figure S2A), the genes characterized by PBS− supported different phylogenetic relationships among C. gracilis, C. mexicanus, andC. goldmani with very low BS values (Supplementary Figure S2B). The PBS+ genes rejected this alternative tree based on SH analysis at a significance level of 0.05, whereas the PBS− genes could not significantly reject the best ML gene tree. Collectively, although we observed conflicting support over the phylogenetic position of C. mexicanus (but not C. goldmani) among the different genes, the genes causing the conflict did not strongly support an alternative phylogeny nor did they reject the best phylogenetic hypothesis. Thus, the poorly resolved relationships are unlikely due to strong conflicting signals.

We then asked whether the undetermined phylogenetic relationships were likely due to rapid diversification, and thus “hard polytomy”. We evaluated all alternative phylogenetic positions of C.mexicanus andC. goldmani using the AU, KH, and SH tests. Among the 14 alternative trees, only two supported C. mexicanus and C. goldmani on the basal branches after the C. parvus-group, and they were significantly worse than the best ML tree at the 0.05 level in all three tests (Supplementary Figure S3; Supplementary Table S5). Because the mitogenome data could not reject the alternative hypothesis, these results suggest that a rapid diversification scenario is plausible.

Multi-locus comprehensive phylogeny

Based on the mitogenome-nuclear gene concatenation tree (comprehensive tree hereafter), eight clades were recovered in Cryptotis (Supplementary Table S1), seven of which were well supported (Figure 1B; BS≥90); clade III, supporting the monophyly of the C. goldmani-group, was only weakly supported (BS=55). The comprehensive tree was congruent with the mitogenome gene tree in supporting both the basal position of theC. parvus-group (BS=78) and a close relationship between C. gracilis and the C. nigrescens-group (BS=90).

In addition to the non-monophyletic C. goodwini-group, two other hypothesized species groups, i.e.,C. mexicanus- and C. thomasi-groups, were each determined to be paraphyletic (Supplementary Table S1). We tested whether the monophyly of these two species groups could be rejected using the AU, KH, and SH tests. Unsurprisingly, the monophyly of the C. mexicanus-group could not be rejected as the BS values supporting the paraphyletic relationships were low (Supplementary Table S6). Congruent with the mitogenome results (Supplementary Figure S3), we also could not reject the sister relationship between the C. goldmani-group and C. goodwini-group, which were previously considered part of the same species group (He et al., 2015; Woodman & Timm, 1998). We also could not reject the grouping of C. gracilis with species from the C. thomasi-group (clade VIII) (Supplementary Table S6). Thus, based on the above analyses, the current six-species-group scenario is not violated (except for C. gracilis (see below)).

Fossorial morphology is becoming a common theme in the evolutionary trajectory of Cryptotis (He et al., 2015; Woodman & Wilken, 2019). Species in the widely separated clades III and VI (C. goldmani- andC. goodwini-groups, respectively) are characterized by enlarged forefeet and claws and a modified humerus. We tested whether these animals could instead be part of a monophyletic clade and thus support the hypothesis of a single evolutionary transition to fossoriality. Our analysis could not reject this hypothesis statistically (Supplementary Table S6); thus, whether there was a single or multiple trajectories toward greater fossoriality remains unresolved.

Species in the C. thomasi-group are mainly distributed in montane areas of northern South America but also in Costa Rica and Panama (Woodman & Timm, 2017). We identified a paraphyletic relationship for the C. thomasi-group, as reported previously (Zeballos et al., 2018), with monophyly statistically rejected (Supplementary Table S6). However, when constraining the monophyly of C. thomasi-group+C. gracilis (i.e., clades V and VIII), this hypothesis could not be rejected, even though the unconstrained comprehensive tree moderately supported a sister relationship between clades VII and VIII (Figure 1B; BS=76). Cryptotis gracilis is mainly distributed in the mountains of Panama and Costa Rica in southern Central America. Thus, our results suggest that either the C. thomasi-group species migrated to South America multiple times, or the ancestor of C. gracilis migrated in the reverse direction from South America to Central America. The latter is a plausible scenario given the distributions of two members of the C. thomasi group, i.e., C. endersi in Panama and C. monteverdensis in Costa Rica, (Pine et al., 2002; Woodman & Timm, 2017).

DISCUSSION

Our results revealed novel relationships supporting a basal position for the C. parvus-group. The humerus of C. parvus is the least derived among living Cryptotis (Woodman & Gaffney, 2014). Although this does not necessarily signify that C. parvus is the most primitive species of the genus, it is a plausible hypothesis.

The non-monophyletic relationships of the C. thomasi-group suggest their ancestors may have migrated to South America multiple times or that reverse migration to Central America also occurred. We have not yet included C. colombiana orC. brachyonyx, the two Colombian members of the C. nigrescens-group, in our taxon sampling. Combined with more robust results from analyses of the C. thomasi group, the C. nigrescens group provides important insights regarding patterns of mammalian migration between Central and South America. Exchanges between these two regions are likely to have occurred several times and could be more complicated than currently understood. Although divergence dating estimates could help clarify the timing of South American colonization (de Abreu-Jr et al., 2020), there are few Cryptotis fossils prior to the late Pleistocene. These are exclusively from the US, and their relationship with modern species is unclear (www.paleobiodb.org, last accessed 29 July 2021).

Despite the limited number of species sampled (11 of 49 recognized species), the use of mitogenomes undoubtably improved overall support of the phylogenetic relationships. Improving taxon sampling, especially for three of the eight major clades that were missed in the comprehensive gene tree, may better resolve these relationships. Although phylogenomic data are generally recommended, such data are more costly to obtain and may be hampered by the presence of non-orthologous sequences (Andermann et al., 2020). However, as the costs of next-generation sequencing (NGS) continue to decrease, using whole-genome shotgun sequencing to obtain mitochondrial sequences has become more economical (Gan et al., 2014). NGS also makes it possible to obtain complete mitogenomes from museum specimens up to 120 years old (de Abreu-Jr et al., 2020). This is highly recommended for Cryptotis as many species are only represented by old museum specimens (Woodman, 2019).

Several relationships, such as the positions ofC. mexicanus and C. goldmani (Figure 1A), could not be finely resolved, even with the availability of mitogenome data. This is due to insufficient rather than conflicting phylogenetic signals embedded within the data, as supported by partitioned Bremer analysis (Supplementary Figure S2) and SH tests (Supplementary Figure S3). Similarly, the relationships among species groups in the comprehensive tree were not well supported, nor could they be rejected statistically (Supplementary Table S6). The “hard polytomy”-like structure may be due to rapid diversification events. Many Cryptotis species inhabit high elevational habitats (Supplementary Figure S4) restricted to small montane areas and isolated by lowlands (Woodman, 2019; Zeballos et al., 2018). These sky islands can facilitate allopatric isolation and speciation (He & Jiang, 2014; McCormack et al., 2009). In addition, species can migrate and colonize new mountains during cool, humid periods (e.g., glacial periods), which could potentially result in rapid diversification, such as observed in the recent radiation of Crocidura shrews (Giarla & Esselstyn, 2015). The speciose Cryptotis group is a good model for understanding how geographic and climatic changes have shaped species diversity in the sky island mountains of Central and northern South America. In addition, their high elevational habitat and limited distribution means these animals are more vulnerable to the effects of anthropogenically induced global warming and habitat destruction. As such, more attention should be paid to the ecology and conservation of these enigmatic small montane mammals.

Supplementary data to this article can be found online.

Click here for additional data file.

COMPETING INTERESTS

The authors declare that they have no competing interests.

AUTHORS’ CONTRIBUTIONS

K.H. designed the study, analyzed the data, and wrote the manuscript; X.C. conducted the laboratory experiments; Y.B.Q. drew the figures, submitted the sequences to GenBank, and submitted the data to GitHub; Z.L. and J.E.M. collected samples; Z.L., W.Z.W., N.W., J.E.M., and X.P. revised the manuscript; X.P. supervised the work. All authors read and approved the final version of the manuscript.

We thank the curators and collection staff from the Smithsonian Institution, National Museum of Natural History, for access to tissue samples under their care. Any use of trade, product, or firm names is for descriptive purposes only and does not imply endorsement by the US government.