Jairo Arroyave a, *, Dalia Angélica De La Cruz-Fernández a, b
a Universidad Nacional Autónoma de México, Instituto de Biología, Tercer Circuito s/n, Ciudad Universitaria, Coyoacán, 04510 Ciudad de México, Mexico
b Universidad Nacional Autónoma de México, Posgrado en Ciencias del Mar y Limnología, Circuito Exterior s/n, Ciudad Universitaria, Coyoacán, 04510 Ciudad de México, Mexico
*Corresponding author: firstname.lastname@example.org (J. Arroyave)
Received: 18 August 2020; accepted: 15 January 2021
Four of the 7 species of Rhamdia present in Mexico stand out for being microendemic and also troglobitic, that is, for being restricted to their type-locality caves and for exhibiting a distinctive phenotype characterized by ocular reduction/loss and body depigmentation. Diagnosis and recognition of Mexican troglobitic forms as distinct species, however, appears to be primarily based on their regressive troglomorphic phenotype and highly localized geographic distributions. To test the adequacy of its current taxonomy, we investigated patterns of genetic and phenotypic variation in Mexican troglobitic Rhamdia in a phylogenetic context. Our results indicate that external morphology does not allow for unambiguous differential diagnoses and robust distinction among troglobitic species. Similarly, beyond typical regressive troglomorphic traits, troglobitic species do not differ greatly in external morphology from their most closely related congener, the epigean species Rhamdia laticauda. From a phylogenetic perspective, continued recognition of troglobitic species implies a deep and generalized paraphyly in R. laticauda. Despite the evidence presented herein, we refrain from making nomenclatural decisions until we can unambiguously ascertain that our findings are indeed explained by phylogeographic structure in R. laticauda, instead of by recent divergence and subsequent speciation of cave-dwelling lineages from this widespread epigean species.
Keywords: Rhamdia reddelli; Rhamdia zongolicensis; Rhamdia macuspanensis; Rhamdia laluchensis; Troglomorphism; Cave fishes; DNA barcoding; Geometric morphometrics
Evidencia genética y morfológica pone en duda la validez de las especies mexicanas troglobias del bagre neotropical del género Rhamdia (Siluriformes: Heptapteridae)
Cuatro de las 7 especies de Rhamdia presentes en México se destacan por ser microendémicas y también troglobias, es decir, por estar restringidas a sus cuevas/localidades tipo y exhibir un fenotipo distintivo caracterizado por la reducción/pérdida ocular y despigmentación corporal. El diagnóstico y reconocimiento de formas troglobias como especies distintas, parece estar basado principalmente en su fenotipo troglomorfo regresivo y distribuciones geográficas altamente restringidas. Con el fin de evaluar la idoneidad de su taxonomía actual, investigamos los patrones de variación genética y fenotípica de las especies troglobias de Rhamdia en México en un contexto filogenético. Nuestros resultados indican que la morfología externa no permite diagnosticar diferencial e inequívocamente ni distinguir de manera robusta a las especies troglobias. Similarmente, más allá de los rasgos troglomórficos regresivos típicos, las especies troglobias no difieren substancialmente en morfología externa de su congénere más cercanamente relacionado, la especie epigea Rhamdia laticauda. Desde una perspectiva filogenética, el reconocimiento continuado de las especies troglobias implica una parafilia profunda y generalizada en R. laticauda. A pesar de la evidencia aquí presentada, nos abstenemos de tomar decisiones nomenclaturales hasta que podamos establecer de forma definitiva que nuestros hallazgos, efectivamente, se explican solo por la estructura filogeográfica en R. laticauda, en lugar de por divergencia reciente y posterior especiación de linajes cavernícolas, a partir de esta especie de superficie ampliamente distribuida.
Palabras clave: Rhamdia reddelli; Rhamdia zongolicensis; Rhamdia macuspanensis; Rhamdia laluchensis; Troglomorfismo; Peces de cueva; Códigos de barras de la vida; Morfometría geométrica
© 2021 Universidad Nacional Autónoma de México, Instituto de Biología. Este es un artículo Open Access bajo la licencia CC BY-NC-ND
Troglobitic organisms are those restricted to hypogean (i.e., subterranean) habitats throughout their life cycle, and as a result of this particular ecology and distribution they tend to evolve an exceptional phenotype characterized by the reduction —or even complete loss— of body pigmentation and ocular structures (Romero & Paulson, 2001; Wilkens & Strecker, 2017). The troglobitic phenotype is therefore considered by many an adaptation to cave life, as well as a classic example of regressive and convergent evolution (Porter & Crandall, 2003). Evidence from developmental biology and molecular evolution studies, however, suggests that both selection and genetic drift have played a role in the evolution of troglomorphism (Rétaux & Casane, 2013). Mexico is not only a megadiverse country in a broad sense, but it also contains one of the most diverse subterranean faunas of the world (Reddell, 1981), which includes troglobitic species of both invertebrates and vertebrates. In Mexico, invertebrate troglobites comprise species from a number of arthropod groups such as crustaceans, insects, arachnids, diplopods and chilopods, while vertebrate troglobites are taxonomically restricted to teleost fishes (Nicholas, 1962). Although troglomorphism has evolved multiple times in the tree of life, fishes are one of the groups of organisms with the largest number of taxa and lineages independently adapted to cave life (Helfman et al., 2009; Proudlove, 2006). Of the 200+ troglobitic species of fish worldwide, 10 are found in Mexico (Proudlove, 2006), 9 of which are endemic, and 1 of which, Astyanax mexicanus (or A. jordani, depending on the author), is arguably the most popular model in studies of evolutionary cave adaptation mechanisms (Gross, 2012; Miller, 2005). Notably, Mexican-endemic cave fishes comprise representatives from phylogenetically distant lineages, and while the majority can be considered primary freshwater fishes, they also include descendants from brackish/marine ancestors, such as the eleotrid Caecieleotris morrisi (Gobiiformes: Eleotridae) and the viviparous brotula Typhlias pearsei (Ophidiiformes: Dinematichthydae) (Møller et al., 2016; Walsh & Chakrabarty, 2016). In contrast to A. mexicanus, very little is known about most of Mexico’s troglobitic ichthyofauna, even on fundamental aspects such as their taxonomy and basic ecology.
Notably, 4 of the 10 troglobitic fish species present in Mexico are members of the genus Rhamdia Bleeker 1858 (Siluriformes: Heptapteridae), a catfish widely distributed in the Neotropical region, from Mexico to Argentina (Hernández et al., 2015; Perdices et al., 2002). Rhamdia is a taxonomically complex group that lacks synapomorphies and is therefore of questionable monophyly (Bockmann, 1998; DoNascimiento et al., 2004; Garavello & Shibatta, 2016). In the latest revision of the genus, Silfvergrip (1996) consolidated its diversity to only 11 species (of the 100+ considered as valid at the time), arguing that, by increasing the geographical coverage of comparative material, phenotypic discontinuities that initially allowed to differentiate between putative species faded into a spectrum of intraspecific morphological variation. Subsequent authors, however, in view of patterns of morphological, ecological, and genetic variation, did not adopt several of Silfvergrip’s nomenclatural proposals (Hernández et al., 2015; Perdices et al., 2002; Weber et al., 2003; Weber & Wilkens, 1998). Currently, Rhamdia comprises 27 valid species (Fricke et al., 2020), 6 of which are troglobitic and endemic to their respective country of occurrence: the Brazilian Rhamdia enfurnada (Bichuette & Trajano, 2005), the Venezuelan Rhamdia guasarensis (DoNascimiento et al., 2004), and the Mexican Rhamdia reddelli (Miller, 1984), Rhamdia zongolicensis (Wilkens, 1993), Rhamdia macuspanensis (Weber & Wilkens, 1998), and Rhamdia laluchensis (Weber et al., 2003). Four of the 7 species of Rhamdia present in Mexico are in fact troglobitic, with Rhamdia laticauda, Rhamdia guatemalensis, and Rhamdia parryi being the only epigean, or surface forms (Miller, 2005). Mexican troglobitic Rhamdia species occur in karstic systems in the southern portion of the country, each being restricted to its own type-locality cave (Fig. 1). Although troglobitic forms of Rhamdia have been known to exist in Mexico since the 1930s (Hubbs, 1936, 1938), the discovery and description of the 4 currently recognized cave-dwelling species ultimately resulted from specimens collected throughout the 1970s and 1980s by American and Italian explorers during speleological —not biological— surveys (Mosier, 1984; Robertson, 1983). Despite reports of troglobitic populations from other caves (Mosier, 1984), these were never investigated from a taxonomic perspective. Various authors believe that the troglobitic forms of Rhamdia in Mexico derive from the epigean species R. laticauda (Greenfield et al., 1982; Perdices et al., 2002; Silfvergrip, 1996; Wilkens, 1993), which includes a hypogean population in Belize described under the category of subspecies (R. laticauda typhla) (Greenfield et al., 1982). The evidence in support for the hypothesis that R. laticauda populations gave rise to the troglobitic Rhamdia species in Mexico, however, is presently scarce and insufficient.
The taxonomic history of Rhamdia as it relates to its Mexican troglobitic species has not been altogether straightforward, as evidenced by a series of past synonymies and the ensuing taxonomic ambiguities. Silfvergrip (1996) was the first to cast doubt on the validity of the troglobitic species described at the time, and synonymized R. zongolicensis and R. reddelli with R. laticauda arguing a lack of diagnostic characters beyond the typical reductions/losses associated with troglomorphism. Subsequently, in the only study to date that has investigated phylogenetic relationships and species boundaries in Middle American Rhamdia, Perdices et al. (2002) concluded that R. reddelli —the only troglobitic species included in their analysis— ought to be synonymized with R. laticauda in order to maintain a taxonomy congruent with phylogeny. More recently, Miller (2005) regarded R. zongolicensis a junior synonym of R. reddelli, while acknowledging the validity of R. macuspanensis and R. laluchensis. Interestingly, in his description of R. zongolicensis, Wilkens himself admitted that morphological differences between R. zongolicensis and the previously described R. reddelli were “minute because of convergent evolution” (Wilkens, 1993, p. 375).
It is thus self-evident that one of the main problems with the existing taxonomy of Rhamdia lies in the insufficiency of the evidence offered to justify the recognition of troglobitic species. Although the presence of morphological discontinuities is traditionally used in species descriptions as indirect evidence of divergence and reproductive isolation, the main morphological discontinuities offered for the differentiation between epigean and hypogean forms of Rhamdia are precisely those related to the regressive troglobitic phenotype, which has resulted in a lack of consensus on the exact diversity of the group and the abovementioned taxonomic ambiguities. That the presence of characters related to cave life is sufficient evidence for the distinction of species is questionable (Miller, 2005), even more so when there appears to be significant intrapopulation variation in those characters, manifested in the coexistence of normal (i.e., eyed and pigmented), troglomorphic, and phenotypically intermediate individuals (Mosier, 1984; Sbordoni et al., 1986; Wilkens, 2001). Indeed, the same type of clinal variation has been documented in A. mexicanus, and none of its troglobitic populations is considered a distinct species (Keene et al., 2015; Miller, 2005). Furthermore, current taxonomic designations do not consider molecular evidence in a phylogeographic/phylogenetic context, which for more than a decade has become a standard methodological approach of integrative taxonomic studies (Dayrat, 2005; Padial et al., 2010).
Given the manifest need for testing the adequacy of its current taxonomy, this study investigates comprehensively, for the first time, patterns of genetic and phenotypic variation in Mexican troglobitic Rhamdia in a phylogenetic context. We leverage novel collections of comparative material, including fresh specimens and tissue samples from all 4 valid cave species and other cave and surface populations/species, to shed light on the hypothesis that troglobitic forms derive from populations of the epigean species R. laticauda that opportunistically colonized underground environments in different regions, but never completely diverged/speciated from ancestral epigean populations.
Materials and methods
Voucher specimens and associated tissue samples were obtained through dedicated collecting trips to the type-locality caves of the 4 troglobitic species in southern Mexico. Whenever possible, epigean habitats —such as rivers and streams— adjacent to or in the vicinity of the type-locality caves were sampled for surface Rhamdia populations/species (Fig. 1). Hypogean specimens were primarily caught by means of minnow traps and dip nets, whereas surface forms were collected using electrofishing, cast nets, and hooks and lines. Some of the sampled cave systems were considerably deep (>50 m and up to 200 m) and not easily accessible. In those cases, sampling was accomplished with the assistance of professional speleologists well trained in vertical caving techniques (Minton & Droms, 2019). Fishes were collected and euthanized prior to preservation in accordance with recommended guidelines for the use of fishes in research (Nickum et al., 2004) and under collecting permits SGPA/DGVS/04259/17 and SGPA/DGVS/05375/19 issued by the Mexican Ministry of Environment and Natural Resources (Secretaría de Medio Ambiente y Recursos Naturales; Semarnat). Stress and suffering were minimized by limiting handling prior to euthanasia, which was accomplished by means of the anesthetic tricaine mesylate (MS-222). Tissue samples (fin clips) were taken in the field, immediately preserved in 95% ethanol, and later cryopreserved at -80 °C. Voucher specimens were fixed in a 10% solution of formalin and subsequently transferred to 70% ethanol for long-term storage. All newly collected specimens were curated, cataloged, and deposited in the Colección Nacional de Peces (CNPE) of the Universidad Nacional Autónoma de México (UNAM). Overall, we collected a total of 129 Rhamdia specimens, including representatives of the 4 Mexican troglobitic species (R. reddelli: n = 8; R. zongolicensis: n = 9; R. macuspanensis: n = 5; R. laluchensis: n = 27), of 4 undescribed hypogean populations from the Sierra de Zongolica, Veracruz, herein treated as R. cf. zongolicensis (n = 43), and of the primarily epigean R. laticauda (n = 22) and R. guatemalensis (n = 15) (Table 1).
Novel comparative DNA sequence data —mtDNA markers cytochrome c oxidase subunit I (COI) and cytochrome b (CYB)— were generated from 27 samples, as follows: R. reddelli (n = 2), R. zongolicensis (n = 2), R. cf. zongolicensis (n = 6), R. macuspanensis (n = 2), R. laluchensis (n = 2), R. laticauda (n = 8), and R. guatemalensis (n = 8). Total genomic DNA was extracted from fresh tissue samples using the Qiagen DNeasy Tissue Extraction Kit following the manufacturer’s protocol. Amplification and sequencing of COI and CYB were carried out using the primer pairs LCO1490/HCO2198 (Folmer et al., 1994) and CytbL14841/CytbH15915 (Irwin et al., 1991; Kocher et al., 1989), respectively. DNA sequencing was carried out at Laboratorio de Secuenciación Genómica de la Biodiversidad y de la Salud (Instituto de Biología, UNAM), in-house Sanger sequencing facilities. Contig assemblage and sequence editing was performed using Geneious Prime 2020.0.4 (https://www.geneious.com). In addition to the genetic data generated from specimens collected in the field during the course of this study, we mined available COI and CYB data from Rhamdia generated in previous studies from public repositories such as The Barcode of Life Data System (BOLD) (Ratnasingham & Hebert, 2007) and GenBank. Overall, we generated/mined COI and CYB sequence data for a total of 109 specimens from 14 Rhamdia species, including all 7 distributed in Mexico. GenBank/BOLD accession numbers for both novel and preexisting sequences and their corresponding voucher specimen data are referenced in table 2.
The monophyletic status, phylogenetic placement, and relationships of Mexican troglobitic Rhamdia species —with respect to each other and to epigean congeners throughout the geographic range of the genus— was investigated via phylogenetic analysis. Specifically, we used 3 different datasets based on the abovementioned comparative DNA sequence data —a COI matrix (645 bp × 53 ingroup terminals), a CYB matrix (1110 bp × 84 ingroup terminals), and a COI +CYB concatenated matrix (1755 bp × 28 ingroup terminals)— to infer phylogenetic relationships in Rhamdia, emphasizing representation of Middle American diversity and Mexican troglobitic forms (Table 2). DNA sequence data from each gene/marker were independently aligned via multiple sequence alignment using the software MUSCLE (Edgar, 2004) under default parameters. The concatenated matrix was taxonomically limited to terminals (voucher specimens) for which both COI and CYB data were available, assembled by linking together individual alignments using the software 2matrix (Salinas & Little, 2014). Best-fit models of nucleotide substitution were statistically selected using the software jModelTest2 (Darriba et al., 2012) based on the Akaike Information Criterion and the following settings for the computation of likelihood scores: number of substitution schemes = 3, base frequencies = +F, rate variation = +I and +G (with nCat = 4), base tree for likelihood calculations = ML optimized, and base tree search = best. Individual gene alignments (i.e., COI and CYB datasets) were each analyzed partitioned by codon position, while the concatenated alignment was analyzed partitioned by both gene and codon position. Maximum likelihood (ML) inference of phylogeny for each dataset was carried out with the software RAxML-HPC BlackBox (v. 8.2.12) (Stamatakis, 2006) as implemented through the CIPRES Science Gateway (v. 3.3) (http://www.phylo.org). Nodal support was assessed using bootstrap percentage values calculated via RAxML’s rapid bootstrapping heuristics (Stamatakis et al., 2008). All trees were rooted at Pimelodella chagresi (COI: MG937103; CYB : AY036748). To further investigate species boundaries in Middle American Rhamdia —with emphasis on troglobitic forms— we analyzed comparative COI sequence data in a phenetic context (i.e., overall genetic similarity) under a DNA barcoding threshold-based approach (Hebert et al., 2003, 2004; Ward, 2009). To this end, we first computed a matrix of Kimura 2-parameter (K2P) genetic distances from the abovementioned COI matrix using the R (www.r-project.org) package “ape” (Paradis et al., 2004; Popescu et al., 2012). We used K2P-corrected genetic distances to allow broader comparisons as this metric represents the standard in DNA barcoding literature (Čandek & Kuntner, 2015). The resulting corrected distance matrix was used as input to generate a dendrogram using the hierarchical clustering algorithm Unweighted Pair Group Method with Arithmetic mean (UPGMA) (Sokal & Michener, 1958) as implemented in the R package phangorn (Schliep, 2011).
Variation in the degree of troglomorphism (i.e., body depigmentation and eye reduction) in hypogean populations was documented and assessed qualitatively by imaging external morphology dorsally, based on a sample of 4 ethanol-preserved specimens per locality, representing the spectrum of variation at each hypogean locality. Localities consisted of all 4 type-locality caves plus 2 caves from the Sierra de Zongolica, Veracruz, harboring undescribed troglobitic populations (i.e., Sótano de Popocatl and Sumidero de Cotzalostoc). Variation in external morphology and overall body shape among cave-dwelling species (including R. cf. zongolicensis from Sótano de Popocatl and Sumidero de Cotzalostoc) and available closely related epigean species (R. laticauda and R. guatemalensis) was investigated quantitatively under a multifarious approach comprising meristics, traditional morphometrics (TM), and landmark-based geometric morphometrics (GM). Variation in overall body shape was assessed on both lateral and dorsal profiles. For meristic and morphometric analyses, novel collections were complemented with 6 voucher specimens of troglobitic species from historical collections, namely R. reddelli (CNPE-IBUNAM 2440; n = 3), R. macuspanensis (CNPE-IBUNAM 22287; n = 1), and R. laluchensis (CNPE-IBUNAM 18455; n = 2). On the other hand, undescribed hypogean populations from “El Naranjal” and “Piedras Blancas” caves (Sierra de Zongolica, Veracruz) were not considered for morphological analyses, due to their small sample sizes (n = 3 and n = 2, respectively).
List of Rhamdia specimens (and associated data) collected during the course of this study and included in genetic and morphological analyses. Mexican states abbreviated using ISO 3166 codes.
|R. reddelli||Cave||CNPE-IBUNAM 23807||JA854-861||8||Cueva Nacimiento Río San Antonio||18.48||-96.64||Cañada San Antonio||OAX||Papaloapan|
|R. zongolicensis||Cave||CNPE-IBUNAM 23800||JA785-791||7||Cueva Ostoc||18.61||-96.89||Zongolica||VER||Papaloapan|
|CNPE-IBUNAM 23806||JA838-839||2||Cueva Ostoc||18.61||-96.89||Zongolica||VER||Papaloapan|
|R. cf. zongolicensis||Cave||CNPE-IBUNAM 23799||JA763-784||22||Sumidero de Cotzalostoc||18.70||-96.92||Totolacatla||VER||Papaloapan|
|CNPE-IBUNAM 23802||JA794-801||8||Sótano de Popocatl||18.76||-96.95||Coetzala||VER||Papaloapan|
|CNPE-IBUNAM 23817||DADF2018.1-2018.7||7||Sótano de Popocatl||18.76||-96.95||Coetzala||VER||Papaloapan|
|CNPE-IBUNAM 23816||DADF2015||1||Sótano de Popocatl||18.76||-96.95||Coetzala||VER||Papaloapan|
|CNPE-IBUNAM 23801||JA792-793||2||Cueva Piedras Blancas||18.67||-96.90||Zongolica||VER||Papaloapan|
|CNPE-IBUNAM 23803||JA802-804||3||Cueva del Naranjal||18.78||-96.96||Naranjal||VER||Papaloapan|
|R. macuspanensis||Cave||CNPE-IBUNAM uncat.||JA742||1||Grutas de Agua Blanca||17.62||-92.47||Macuspana||TAB||Usumacinta|
|CNPE-IBUNAM 23814||JA932-935||4||Grutas de Agua Blanca||17.62||-92.47||Macuspana||TAB||Usumacinta|
|R. laluchensis||Cave||CNPE-IBUNAM uncat.||JA743-747||5||Sótano de La Lucha||17.06||-93.90||La Lucha||CHP||Grijalva|
|CNPE-IBUNAM 23811||JA898-919||22||Sótano de La Lucha||17.06||-93.90||La Lucha||CHP||Grijalva|
|R. laticauda||Surface||CNPE-IBUNAM 23804||JA805-809||5||Río Atoyac||18.91||-96.78||Atoyac||VER||Papaloapan|
|CNPE-IBUNAM 23805||JA831-834||4||Río Hoyotempa||18.60||-96.89||Cuahuixtlahuac||VER||Papaloapan|
|CNPE-IBUNAM 23810||JA886-896||11||Nacimiento Río Tonto (Huixtla)||18.59||-96.88||Aticpac||VER||Papaloapan|
|CNPE-IBUNAM 23809||JA864||1||Río San Antonio||18.46||-96.63||Cañada San Antonio||OAX||Papaloapan|
|CNPE-IBUNAM 23813||JA928||1||Río Agua Blanca||17.07||-93.87||La Lucha||CHP||Grijalva|
|R. guatemalensis||Surface||CNPE-IBUNAM 23808||JA862-863||2||Río San Antonio||18.46||-96.63||Cañada San Antonio||OAX||Papaloapan|
|CNPE-IBUNAM 23812||JA897, 920-927, 929||10||Río Agua Blanca||17.07||-93.87||La Lucha||CHP||Grijalva|
|Cave||CNPE-IBUNAM 23815||JA943-945||3||Grutas de Coconá||17.56||-92.93||Teapa||TAB||Grijalva|
Rhamdia samples (catalog and voucher) used for genetic analyses and their corresponding GenBank/BOLD accession numbers. Acronyms of ichthyological collections follow Sabaj (2020). Superscripts following accession numbers correspond to the sequence source as follows: 1) this study, 2) Arroyave et al. (2020), 3) Bermingham et al. (Unpublished), 4) BOLD (Unpublished), 5) de Carvalho et al. (2011), 6) Usso et al. (2019), 7) García-Dávila et al. (2015), 8) Hernández et al. (2015), 9) Perdices et al. (2002), 10) Pereira et al. (2013), 11) Ríos et al. (2017), and 12) Valdez-Moreno et al. (2009). Countries and provinces/states abbreviated using ISO 3166 codes.
|Rhamdia reddelli||CNPE-IBUNAM 23807||JA854||MX||OAX||Papaloapan||A||18.48||-96.64||MW4760781||MW4675851|
|Rhamdia reddelli||CNPE-IBUNAM 23807||JA856||MX||OAX||Papaloapan||A||18.48||-96.64||–||MW4675861|
|Rhamdia reddelli||MNCN 2781||MNCN-2781||MX||OAX||Papaloapan||A||18.48||-96.64||–||AY0366979|
|Rhamdia zongolicensis||CNPE-IBUNAM 23800||JA786||MX||VER||Papaloapan||A||18.61||-96.89||MW4760791||MW4675871|
|Rhamdia zongolicensis||CNPE-IBUNAM 23800||JA787||MX||VER||Papaloapan||A||18.61||-96.89||–||MW4675881|
|Rhamdia cf. zongolicensis||CNPE-IBUNAM 23799||JA764||MX||VER||Papaloapan||A||18.70||-96.92||MW4760751||MW4675761|
|Rhamdia cf. zongolicensis||CNPE-IBUNAM 23799||JA765||MX||VER||Papaloapan||A||18.70||-96.92||–||MW4675771|
|Rhamdia cf. zongolicensis||CNPE-IBUNAM 23801||JA792||MX||VER||Papaloapan||A||18.67||-96.90||–||MW4675791|
|Rhamdia cf. zongolicensis||CNPE-IBUNAM 23801||JA793||MX||VER||Papaloapan||A||18.67||-96.90||–||MW4675801|
|Rhamdia cf. zongolicensis||CNPE-IBUNAM 23802||JA794||MX||VER||Papaloapan||A||18.76||-96.95||–||MW4675811|
|Rhamdia cf. zongolicensis||CNPE-IBUNAM 23803||JA802||MX||VER||Papaloapan||A||18.78||-96.96||–||MW4675781|
|Rhamdia macuspanensis||CNPE-IBUNAM 23814||JA932||MX||TAB||Usumacinta||A||17.62||-92.47||–||MW4675821|
|Rhamdia macuspanensis||CNPE-IBUNAM 23814||JA934||MX||TAB||Usumacinta||A||17.62||-92.47||MW4760761||MW4675831|
|Rhamdia macuspanensis||CNPE-IBUNAM 23814||JA935||MX||TAB||Usumacinta||A||17.62||-92.47||MW4760771||MW4675841|
|Rhamdia laluchensis||CNPE-IBUNAM 23811||JA900||MX||CHP||Grijalva||A||17.06||-93.90||MW4760671||MW4675701|
|Rhamdia laluchensis||CNPE-IBUNAM 23811||JA901||MX||CHP||Grijalva||A||17.06||-93.90||–||MW4675711|
|Rhamdia laticauda||CNPE-IBUNAM 22978||A05||MX||CHP||Lacantún||A||16.10||-91.01||MW4760691||–|
|Table 2. Continued|
|Rhamdia laticauda||CNPE-IBUNAM 22978||F04||MX||CHP||Lacantún||A||16.10||-91.01||MW4760701||–|
|Rhamdia laticauda||CNPE-IBUNAM 22978||G04||MX||CHP||Lacantún||A||16.10||-91.01||MW4760711||–|
|Rhamdia laticauda||CNPE-IBUNAM 22978||H04||MX||CHP||Lacantún||A||16.10||-91.01||MW4760721||–|
|Rhamdia laticauda||CNPE-IBUNAM 23011||D04||MX||CHP||Lacantún||A||16.11||-90.94||MW4760731||MW4675731|
|Rhamdia laticauda||CNPE-IBUNAM 23804||JA805||MX||VER||Papaloapan||A||18.91||-96.78||–||MW4675741|
|Rhamdia laticauda||CNPE-IBUNAM 23805||JA831||MX||VER||Papaloapan||A||18.60||-96.89||MW4760681||MW4675721|
|Rhamdia laticauda||STRI 11528||STRI-11528||PA||1||Róbalo||A||9.07||-82.29||MG9371773||–|
|Rhamdia laticauda||CNPE-IBUNAM 23810||JA887||MX||VER||Papaloapan||A||18.59||-96.88||MW4760741||MW4675751|
|Rhamdia laticauda||STRI 16400||STRI-16400||PA||2||Antón||P||8.60||-80.14||MG9371803||–|
|Rhamdia laticauda||MNCN 10266-GU||MNCN10266-GU||GT||CQ||Lempa||P||14.53||-89.42||–||AY0367179|
|Rhamdia laticauda||MNCN 3499-MEX||MNCN3499-MEX||MX||VER||Catemaco||A||18.62||-96.87||–||AY0366969|
|Rhamdia laticauda||STRI 2686||STRI-2686||PA||2||Indio||A||8.65||-80.12||MG9371753||–|
|Rhamdia laticauda||MNCN 46-IST||MNCN46-IST||PA||4||Chiriqui Viejo||P||8.42||-82.68||–||AY0367289|
|Rhamdia laticauda||STRI 6873||STRI-6873||PA||1||Róbalo||A||9.04||-82.29||MG9371743||–|
|Rhamdia laticauda||MNCN 653-MEX||MNCN653-MEX||MX||OAX||Papaloapan||A||17.68||-97.09||–||AY0366989|
|Rhamdia laticauda||MNCN GU-3||MNCNGU-3||MX||CHP||Usumascinta||A||16.10||-91.77||–||AY0367049|
|Rhamdia laticauda||STRI 10084||STRI-10084||PA||2||Cocle del Norte||A||8.70||-80.45||MG9371783||AY0367329|
|Table 2. Continued|
|Rhamdia laticauda||STRI 11527||STRI-11527||PA||1||Róbalo||A||9.07||-82.29||MG9371763||AY0367319|
|Rhamdia laticauda||STRI 13659||STRI-13659||NI||AS||Escondido||A||12.01||-84.67||MG4962173||AY0367119|
|Rhamdia laticauda||STRI 13805||STRI-13805||NI||CO||San Juan||A||12.00||-85.17||–||AY0367139|
|Rhamdia laticauda||STRI 14536||STRI-14536||NI||SJ||San Juan||A||11.52||-84.83||–||AY0367129|
|Rhamdia laticauda||STRI 2161||STRI-2161||CR||A||San Juan||A||10.91||-85.21||MG4962143||AY0367149|
|Rhamdia laticauda||STRI 2162||STRI-2162||CR||A||San Juan||A||10.91||-85.21||–||AY0367159|
|Rhamdia laticauda||STRI 218||STRI-218||CR||L||Sixaloa||A||9.60||-82.80||MG4962133||AY0367279|
|Rhamdia laticauda||STRI 219||STRI-219||CR||L||Sixaloa||A||9.60||-82.80||–||AY0367169|
|Rhamdia laticauda||STRI 2686||STRI-2686||PA||2||Indio||A||8.65||-80.12||–||AY0367339|
|Rhamdia laticauda||STRI 6345||STRI-6345||PA||4||Chiriqui Viejo||P||8.67||-82.85||–||AY0367309|
|Rhamdia laticauda||STRI 662||STRI-662||PA||4||Chiriqui Viejo||P||8.76||-82.83||MG9371813||AY0367299|
|Rhamdia laticauda||STRI 722||STRI-722||PA||2||Antón||P||8.42||-80.25||MG9371793||AY0367349|
|Rhamdia laticauda||STRI 7821||STRI-7821||GT||AV||Usumacinta||A||15.83||-90.33||–||AY0367059|
|Rhamdia laticauda||STRI 7822||STRI-7822||GT||AV||Usumacinta||A||15.83||-90.33||–||AY0367069|
|Rhamdia laticauda||STRI 8099||STRI-8099||GT||PE||Belize||A||16.96||-89.36||–||AY0367079|
|Rhamdia laticauda||STRI 8100||STRI-8100||GT||PE||Belize||A||16.96||-89.36||–||AY0367089|
|Rhamdia laticauda||STRI 8199||STRI-8199||GT||IZ||Polochic||A||15.54||-88.90||–||AY0367099|
|Table 2. Continued|
|Rhamdia laticauda||STRI 8284||STRI-8284||GT||AV||Polochic||A||15.32||-89.86||–||AY0367109|
|Rhamdia laticauda||STRI 8345||STRI-8345||HN||CP||Ulua||A||14.65||-88.88||–||AY0367219|
|Rhamdia laticauda||STRI 8393||STRI-8393||HN||CP||Copán||A||14.86||-89.12||–||AY0367209|
|Rhamdia laticauda||STRI 8585||STRI-8585||HN||YO||Aguan||A||15.48||-86.67||–||AY0367229|
|Rhamdia laticauda||STRI 8729||STRI-8729||HN||OL||Patuca||A||14.40||-85.93||–||AY0367249|
|Rhamdia laticauda||STRI 8794||STRI-8794||HN||OL||Patuca||A||14.58||-86.29||–||AY0367239|
|Rhamdia laticauda||STRI 8965||STRI-8965||HN||FM||Choluteca||P||14.01||-87.00||–||AY0367259|
|Rhamdia laticauda||STRI 8966||STRI-8966||HN||FM||Choluteca||P||14.01||-87.00||–||AY0367269|
|Rhamdia nicaraguensis||STRI 2122||STRI-2122||CR||G||Higuerón||P||10.34||-85.08||MG4962153||AY0367189|
|Rhamdia nicaraguensis||STRI 2121||STRI-2121||CR||G||Higuerón||P||10.34||-85.08||MG4962163||AY0367199|
|Rhamdia parryi||MNCNGU 178||MNCNGU-178||GT||ES||Acomé||P||14.32||-90.98||–||AY0367029|
|Rhamdia parryi||STRI 7723||STRI-7723||GT||SU||Nahualate||P||14.50||-91.41||–||AY0367009|
|Rhamdia parryi||STRI 7724||STRI-7724||GT||SU||Nahualate||P||14.50||-91.41||–||AY0367019|
|Rhamdia parryi||STRI 7776||STRI-7776||GT||ES||María Linda||P||14.20||-90.71||–||AY0367039|
|Rhamdia saijaensis||IMCN-TE00018||IMCN-TE00018||CO||–||San Juan||P||–||–||–||KM4890748|
|Rhamdia saijaensis||IMCN-TE00019||IMCN-TE00019||CO||–||San Juan||P||–||–||–||KM4890758|
|Rhamdia saijaensis||IMCN-TE00020||IMCN-TE00020||CO||–||San Juan||P||–||–||–||KM4890768|
|Table 2. Continued|
|Rhamdia saijaensis||IMCN-TE00021||IMCN-TE00021||CO||–||San Juan||P||–||–||–||KM4890778|
|Rhamdia saijaensis||IMCN-TE00022||IMCN-TE00022||CO||–||San Juan||P||–||–||–||KM4890788|
|Rhamdia cinerascens||STRI 12103||STRI-12103||EC||–||Daule||P||–||–||–||AY0367359|
|Rhamdia cinerascens||STRI 12104||STRI-12104||EC||–||Daule||P||–||–||–||AY0367369|
|Rhamdia guatemalensis||CIUA 8494||CIUA-8494||CO||CAL||Cauca||P||5.06||-75.71||UDEA184-184||–|
|Rhamdia guatemalensis||CIUA 8669||CIUA-8669||CO||ANT||Atrato||P||6.40||-76.71||UDEA187-184||–|
|Rhamdia guatemalensis||CNPE-IBUNAM 23233||JA012||MX||YUC||YP aquifer||A||20.58||-89.61||MW0025601||MW0104461|
|Rhamdia guatemalensis||CNPE-IBUNAM 23316||JA323||MX||ROO||YP aquifer||A||20.28||-87.48||MW4760621||MW4675651|
|Rhamdia guatemalensis||CNPE-IBUNAM 23808||JA862||MX||OAX||Papaloapan||A||18.46||-96.63||MW4760631||MW4675661|
|Rhamdia guatemalensis||CNPE-IBUNAM 23812||JA897||MX||CHP||Grijalva||A||17.07||-93.87||MW4760641||MW4675671|
|Rhamdia guatemalensis||CNPE-IBUNAM 23812||JA920||MX||CHP||Grijalva||A||17.07||-93.87||MW4760651||MW4675681|
|Rhamdia guatemalensis||CNPE-IBUNAM 23815||JA943||MX||TAB||Grijalva||A||17.56||-92.93||MW4760661||MW4675691|
|Rhamdia guatemalensis||ECOCH 5785||MX528||GT||AV||Cahabón||A||15.78||-90.34||EU75195912||–|
|Rhamdia guatemalensis||ECOCH 7593||MXV-771||BZ||2||–||A||17.28||-88.66||MXV775-154||–|
|Rhamdia guatemalensis||STRI 1207||STRI-1207||CR||P||Cañas||P||10.35||-85.17||MG4962183||AY0366719|
|Rhamdia guatemalensis||STRI 13519||STRI-13519||NI||CI||Estero Real||P||12.95||-86.85||MG4962203||–|
|Rhamdia guatemalensis||STRI 1669||STRI-1669||PA||KY||Mandinga||A||9.47||-79.09||MG9371883||AY0366859|
|Table 2. Continued|
|Rhamdia guatemalensis||STRI 2165||STRI-2165||CR||A||San Juan||A||10.91||-85.21||MG4962213||AY0366609|
|Rhamdia guatemalensis||STRI 1569||STRI-1569||CO||CHO||Atrato||A||5.62||-76.74||–||AY0366899|
|Rhamdia guatemalensis||STRI 816||STRI-816||CO||BOL||Magdalena||A||9.36||-74.72||–||AY0366929|
|Rhamdia guatemalensis||STRI 7831||STRI-7831||GT||AV||Usumacinta||A||15.83||-90.33||–||AY0366289|
|Rhamdia guatemalensis||STRI VZ-1||STRI-VZ-1||VE||V||Lake Maracaibo||A||8.84||-71.98||–||AY0366949|
|Rhamdia quelen||MZUEL 14359-4477||MZUEL14359-4477||AR||E||Lower Paraná||A||–||–||KU8456876||–|
|Rhamdia quelen||LBP 29359||LBP29359||BR||SP||Upper Paraná||A||-22.64||-44.62||GU70211910||–|
|Rhamdia quelen||MCP 44783||MCP44783||BR||MG||São Francisco||A||-19.47||-43.87||HM9060115||–|
|Rhamdia quelen||IIAP 470||IIAP-470||PE||–||Amazon||A||–||–||KT9524557||–|
|Rhamdia quelen||ROM uncat.||TZGAA016-06||GY||BA||Waini||A||7.50||-59.55||TZGAA016-064||–|
|Rhamdia quelen||STRI 517||STRI-517||PE||MDD||Amazon||A||-11.80||-71.44||–||AY0367419|
|Rhamdia quelen||STRI 2308||STRI-2308||AR||W||Paraná||A||-28.53||-58.91||–||AY0367429|
|Rhamdia voulezi||UEL 16503||BRH013-16||BR||PR||Iguaçu||A||-25.62||-52.62||BRH013-164||–|
|Rhamdia voulezi||UEL 16503||BRH014-16||BR||PR||Iguaçu||A||-25.62||-52.62||BRH014-164||–|
|Rhamdia voulezi||UEL 16503||BRH015-16||BR||PR||Iguaçu||A||-25.62||-52.62||BRH015-164||–|
|Rhamdia branneri||UEL 16503||BRH001-16||BR||PR||Iguaçu||A||-25.62||-52.62||BRH001-164||–|
|Rhamdia branneri||UEL 16504||BRH002-16||BR||PR||Iguaçu||A||-25.62||-52.62||BRH002-164||–|
A total of 6 meristic and 38 mensural (linear) traits were recorded from a sample of 122 specimens: 116 of the 129 newly collected plus the abovementioned 6 from historical collections (Table 3; Fig. 2). Eight of the 129 newly collected specimens were unavailable for meristics and TM analyses because they were cryopreserved (-80 °C) immediately after photographic documentation and tissue sampling, and are being kept ultrafrozen for future genomic studies. Selection of meristic and linear traits was informed by relevant previous taxonomic work (Hernández et al., 2015; Silfvergrip, 1996). Linear traits were measured with a handheld digital caliper to the nearest 0.01 mm, standardized as a proportion of either standard length (SL) or head length (HL), and log-transformed for use in TM analyses. Exploratory analysis of linear data was conducted via Principal Component Analysis (PCA), seeking to visualize general patterns of morphological variation among species (including 2 undescribed cave-dwelling populations provisionally identified as R. cf. zongolicensis) and to identify morphometric variables that contribute to the distinction of groups within the multivariate space. To test for significant differences in overall body shape among species/populations, a Multivariate Analysis of Variance (MANOVA) was conducted using the resulting first 3 principal components (PCs) as dependent variables, and species/populations as independent (predictor) variables. Lastly, post hoc pairwise comparisons (Bonferroni corrected) were conducted as follow-up tests of differentiation (Hotelling’s T2 test) between pairs of species/populations. Both PCA and MANOVA analyses (including post hoc comparisons) were carried out with the software PAST (v. 4.04) (Hammer et al., 2001).
For GM analyses, we employed homologous landmarks of the types 1 and 2, as well as homologous semilandmarks. Type-1 landmarks are those clearly defined by a structure or insertion (e.g., pectoral-fin insertion), whereas those of the type 2 are defined more ambiguously, often as inflection points in curved structures (e.g., posterior margin of the operculum). On the other hand, semilandmarks refer to points along a contour used to represent homologous curves or surfaces (Bookstein, 1997; Gunz et al., 2005). To reduce the criterion of bending energy of the thin-plate spline surface, semilandmarks are allowed to slide along pre-determined vectors (sliders). After sliding, landmarks and semilandmarks can be treated the same way in subsequent statistical analyses (Gunz et al., 2005). The sample for lateral GM (lGM) analysis consisted of 135 specimens: the 129 specimens listed in table 1 plus the abovementioned 6 from historical collections. Specimens were photographed post-mortem but prior to preservation (except for the 6 ethanol-preserved from historical collections) in a standardized lateral view indicating 15 homologous landmarks (LMs): 14 type 1 (LMs 1-14; mostly located at the insertion of structures) and 1 type 2 (LM 15; posteriormost point of the operculum). Additionally, 7 equidistant semilandmarks (LMs 16-22) were used to capture variation in the dorsal cephalic profile, from the anteriormost tip of the head to the insertion of the dorsal fin (Fig. 3a). For dorsal GM (dGM) analysis, we used a total of 8 landmarks —7 type 1 (LMs 1-7; mostly insertion points) and 1 type 2 (LM 8; posteriormost point of occipital process)— from a smaller subset consisting of 114 specimens but with representation of all sampled localities/populations (Fig. 3b). Just as with meristic and linear trait data, landmark selection was informed by consideration of taxonomically important characters reported in previous studies (Hernández et al., 2015; Silfvergrip, 1996). Homologous landmarks and semilandmarks were digitized with the computer programs tpsDig2 (v.2.31) and tpsUtil (v.1.70), respectively (Rohlf, 2015). The sliders for each semilandmark were defined as tangent vectors to the outline in the position of the point. This tangent was defined as the chord between the previous and next points to the semilandmark along the contour. Shapes were aligned using generalized Procrustes analysis as implemented in the software tpsRelw (v. 1.67) (Rohlf, 2015), removing non-shape variation due to specimen size, location, and orientation. Overall patterns of shape variation were visualized with a PCA paired with thin-plate spline analysis (Bookstein, 1989) to generate warp grids for displaying shape deformations at extremes of principal component (PC) axes relative to mean shape. Relative warps (RWs) analysis was carried out in the software tpsRelw (v. 1.67) (Rohlf, 2015). We tested for measurement error (Bailey & Byrnes, 1990) due to digitizing by means of a one-way Procrustes ANOVA based on 2 alternative digitization schemes for a sample of 30 individuals (1,000 randomized residual permutations) with the function ‘procD.lm’ of the R package geomorph (Adams & Otárola-Castillo, 2013). The resultant estimates of repeatability and measurement error were 98.75% and 1.25%, respectively, and the “individual” main effect was highly significant (p < 0.001). These results imply that the variation among individuals greatly exceeds the measurement error due to digitizing, confirming that variation in residuals was not due to digitizing error but to biological variation. To further investigate group structure and differentiation in multivariate body shape data, a Canonical Variate Analysis (CVA) was performed on both lateral- and dorsal-view Procrustes shape coordinates. By finding shape features that maximize the discrimination among pre-defined groups (relative to the variation within groups), CVA is often used for assessment of group separation and thus to support the identification of distinct species. The statistical significance of pairwise differences in mean shapes was assessed by means of permutation tests based on both Procrustes and Mahalanobis distances and 10,000 permutations per test. CVA and permutation tests were conducted with MorphoJ (v. 1.07) (Klingenberg, 2011).
Meristic and linear morphometric variables used to assess external morphological variation.
|Body depth at dorsal fin||BD|
|Snout to posterior tip of occipital process||S-OP|
|Snout to insertion of dorsal fin||S-DF|
|Snout to insertion of adipose fin||S-AdF|
|Snout to insertion of pectoral fin||S-PcF|
|Snout to insertion of pelvic fin||S-PvF|
|Snout to anal cavity||S-A|
|Snout to insertion of anal fin||S-AF|
|Pectoral fin spine length||PcFSL|
|Dorsal fin spine length||DFSL|
|Inter-maxillary barbel distance||IMBD|
|Body width at pectoral fins||BW|
|Inter-pelvic fin distance||IPvFD|
|Dorsal-fin base length||DFL|
|Adipose-fin base length||AdFL|
|Anal-fin base length||AFL|
|Caudal-fin base length||CFL|
|Insertion of dorsal fin to insertion of adipose fin||DF-ADF|
|Insertion pectoral fin to insertion of pelvic fin||PcF-PvF|
|Insertion pelvic fin to insertion of anal cavity||PvF-A|
|Anal cavity to insertion of anal fin||A-AF|
|Posterior insertion of anal fin to origin of caudal fin||pAF-CF|
|Posterior insertion of adipose fin and caudal fin||pAdF-CF|
|Insertion of dorsal fin to operculum end||DF-O|
|Insertion of dorsal fin to insertion of pectoral fin||DF-PcF|
|Table 3. Continued|
|Insertion of dorsal fin to insertion of pelvic fin||DF-PvF|
|Insertion of adipose fin to insertion of pectoral fin||AdF-PcF|
|Insertion of adipose fin to insertion of pelvic fin||AdF-PvF|
|Insertion of adipose fin to insertion of anal fin||AdF-AF|
|Posterior insertion of adipose fin to posterior insertion of anal fin||pAdF-pAF|
|Caudal-fin upper lobe rays||UCFr|
|Caudal-fin lower lobe rays||LCFr|
For all 3 datasets, the best-fit substitution model was the general time-reversible (GTR) with rate heterogeneity among sites modeled by the gamma distribution (Γ) and a proportion of invariant sites (I). The phylogenies resulting from analysis of the COI, CYB, and COI +CYB datasets are presented in figure 4. In all phylogenies, Mexican species of the “R. laticauda-group” (Weber & Wilkens, 1998; Wilkens, 2001) are resolved as forming a well-supported clade —demarcated by a gray dotted line— that also includes the Costa Rican/Nicaraguan Rhamdia nicaraguensis. Previously undocumented populations from non-type-locality caves in the Sierra de Zongolica, Veracruz (R. cf. zongolicensis), are also resolved as nested within the abovementioned clade, closely related to —but not reciprocally monophyletic with— R. zongolicensis. Regardless of dataset, the phylogenetic placement of troglobitic species (including R. cf. zongolicensis) renders R. laticauda non-monophyletic. The monophyly of R. laticauda is similarly challenged by the placement of the epigean R. nicaraguensis and R. parryi. On the other hand, the monophyly of troglobitic species is not contested by the phylogeny resulting from analysis of the dataset in which multiple individuals of each of them were sampled (i.e., CYB) (Fig. 4b). All phylogenies reveal a modest level of lineage differentiation —as implied by branch lengths— in the troglobitic R. laluchensis and R. macuspanensis. Conversely, the degree of lineage divergence between the troglobitic R. reddelli and R. zongolicensis (including R. cf. zongolicensis) and R. laticauda samples from Mexico (particularly from the Papaloapan basin in Veracruz) appears to be negligible.
The resulting COI UPGMA dendrogram displaying clustering of samples according to the (K2P-corrected) genetic distances between them is presented in figure 5. The greatest divergences are those between R. guatemalensis samples and the remaining sampled Rhamdia species (up to 12.5%). Genetic distances among species of the “R. laticauda-group” do not exceed 4% for any pairwise comparison, with the majority being under 2.5%. The most COI -divergent of the troglobitic species is R. macuspanensis, at ~4% from the remaining troglobitic forms, which do not differ among them in more than ~ 2%. Barcodes of R. reddelli and R. zongolicensis (including R. cf. zongolicensis from the Cotzalostoc cave) are effectively identical, and barely divergent (< 1%) from those of R. laticauda samples from the Papaloapan basin.
A broad spectrum of variation in the degree of troglomorphism (i.e., eye reduction and depigmentation) was observed in most troglobitic species, including populations from non-type-locality caves in the Sierra de Zongolica (Fig. 6). In the case of R. zongolicensis, none of the sampled individuals displayed complete eye loss and depigmentation (Fig. 6b). Notably, our sample of R. macuspanensis did not exhibit varying degrees of troglomorphism, but instead only individuals either fully or non-troglomorphic; that is, syntopy of cave and surface forms (Figs. 6c, 7). Taxonomic designation of the non-troglomorphic R. macuspanensis specimen (CNPE-IBUNAM 23814, JA935) was based on multiple sources of evidence such as morphometric and meristic traits, genetics (COI barcodes), and geographic distribution (co-occurrence in the type-locality cave).
The distribution of meristic counts in troglobitic forms (all 4 described species plus 2 undescribed populations/species) and the epigean species R. laticauda and R. guatemalensis is presented in table 4. The PCA on linear data resulted in 37 PCs that explain 100% of the variation, with only the first 3 individually accounting for more than 10%, and together explaining 52.6% of the overall variance. A graphical representation of these results is presented in figure 8. The first PC is primarily associated with variation in barbel length, being the maxillary barbels those exhibiting the greatest variation. Despite the considerable overlap among species/populations along this PC, it can be observed that samples from the epigean R. guatemalensis and the hypogean R. macuspanensis cluster at the upper half (implying longer barbels), whereas those from the epigean R. laticauda and the hypogean R. cf. zongolicensis are mostly confined to the lower half (implying shorter barbels). Samples from the remaining species (R. reddelli, R. zongolicensis, and R. laluchensis) are highly overlapping and cluster at intermediate values along this gradient of variation. The second PC is primarily associated with variation in body depth. With the exception of R. macuspanensis, which displays a relatively robust body more akin to that of epigean species, samples from troglobitic species/populations cluster primarily towards the negative side of this axis, meaning that they have considerably slenderer bodies. Lastly, the third PC is mostly associated with variation in the length of the interdorsal space. The observed pattern implies that R. guatemalensis is the most divergent in this respect, having the shortest interdorsal space with respect to R. laticauda and the troglobitic forms, which together display almost complete overlap and are therefore indistinguishable as far as this character is concerned. Despite the apparent lack of clear differentiation among species/populations implied by the general patterns of morphometric variation uncovered via PCA (except for R. guatemalensis and R. macuspanensis) (Fig. 8), MANOVA results indicate that in fact there are statistically significant differences among groups, a finding that holds across most pairwise comparisons (Tables 5, 6). Notably, R. laticauda was found to be significantly different in overall body shape from all troglobitic species/populations. In contrast, no significant differences were found between R. reddelli and R. zongolicensis.
The PCA of lateral-view body shape in the sampled specimens indicates that 40 principal components/RWs explain 100% of the variation, with the first 2 explaining 41.92% of the overall shape variance. The largest component of shape variation (RW1), accounting for 25.76% of the total variance, is associated primarily with changes in body depth (at the level of pelvic fin) and, to a lesser extent, with head length (Fig. 9a). Except for R. macuspanensis, samples from troglobitic species/populations cluster towards the positive side of RW1 (i.e., slenderer, longer-headed bodies), whereas samples from the epigean species R. guatemalensis and R. laticauda scatter on the negative side (i.e., deeper, shorter-headed bodies). The second-largest component of shape variation (RW2) explains 16.66% of total variance and is associated mainly with changes in interdorsal space and, to a smaller degree, with head depth (Fig. 9a). Therefore, individuals located on the positive side of RW2 display a longer interdorsal space coupled with a shorter adipose fin, as well as a slenderer head compared to those on the negative side of this RW. Troglobitic species display considerable overlap on this axis, with most of the clustering at intermediate values (near zero), and with R. laluchensis displaying the widest range of variation (possibly due to having the largest sample size). Although R. laticauda displays the largest values along this axis (i.e., largest interdorsal space and slenderest head), these partially overlap with those of troglobitic forms, particularly at the lower half of the range. On the other hand, the PCA of dorsal-view body shape resulted in 12 RWs that explain 100% of the variation, with the first 2 explaining 66.07% of the overall variance. In this case, the largest component of variation (RW1) explains most of the total variance (52.34%) and is associated with variation in head size and shape and in the relative position of the pectoral and pelvic fins (Fig. 9b). Specifically, extreme values along this axis describe either individuals with longer and sharper heads and paired fins closer together (on the negative side) or with shorter and blunter heads and paired fins more apart (on the positive side). With the exception of R. macuspanensis, which displays a broad range of values that extend on both negative and positive sides of this axis, samples from described troglobitic species cluster on the negative side and therefore have comparatively larger heads and paired fins closer together. On the other hand, undescribed troglobitic populations from the Sierra de Zongolica (Popocatl and Cotzalostoc) display intermediate values along this axis, while samples from R. laticauda occur mostly on the positive side. Variation in the second-largest component of dorsal-view shape variance (RW2; 13.71%) is mainly explained by the shape of the trunk, particularly as it relates to the position of the pelvic fins, with fishes on the positive side of the axis displaying posteriorly narrower bodies (Fig. 9b). Except for R. guatemalensis, whose individuals cluster narrowly towards the central region, sampled populations display a broad and highly overlapping range of variation along this axis, implying that this shape component is not particularly useful to morphologically differentiate species. The results from CVA on both lateral and dorsal profiles are basically equivalent to the abovementioned general patterns of morphological variation implied by PCA (Fig. 10). Lateral-shape CVA indicates that 7 canonical variates (CV) explain 100% of the variation, the first 2 of which account for 63.15% of the overall shape variance (Fig. 10a). The largest component of lateral shape variation (CV1; 35.51% of total variance) is primarily associated with changes in body depth and head length. Lower values on this axis indicate deeper bodies and shorter heads (such as in R. guatemalensis), while higher values specify slenderer bodies with larger heads (such as in most troglobitic species) (Fig. 10a). The second-largest component of lateral shape variation (CV2; 27.64% of total variance) is associated mainly with changes in interdorsal space. For the most part, troglobitic species occur at intermediate values along this axis, between R. guatemalensis (at the positive extreme) and R. laticauda (at the negative extreme) (Fig. 10a). On the other hand, CVA of dorsal shape data resulted in 7 CVs that explain 100% of the variation, of which the first 2 account for 80.59% of the overall shape variance (Fig. 10b). The largest component of dorsal shape variation (CV1; 67.59% of total variance) is mainly related to head size and body width at pelvic-fin origin. Most troglobitic species and populations fall on the positive side of this axis (i.e., have longer heads and narrower bodies) with considerable overlap, while the epigean species and the troglobitic R. macuspanensis are mostly restricted to the negative side (i.e., have smaller heads and thicker bodies) (Fig. 10b). The second-largest component of dorsal shape variation (CV2; 13% of total variance) is related to changes in head and trunk width. The substantial overlap among species along this axis, however, entails no discernible pattern of differentiation based on this shape attribute (Fig. 10b). Overall, with few exceptions, all pairwise comparisons of lateral shape resulted in statistically significant differences (Table 7).
Summary of the variation in meristic trait data registered for all 4 troglobitic species, 2 undescribed cave-dwelling populations (R. cf. zongolicensis1 = Sumidero de Cotzalostoc; R. cf. zongolicensis2 = Sótano de Popocatl), and the epigean R. laticauda and R. guatemalensis. Total number of individuals in parenthesis. Meristic formula: spine counts in Roman numerals; (soft) ray counts in Arabic. Abbreviations: DFr = dorsal-fin rays; PvFr = pelvic-fin rays; PcFr = pectoral-fin rays; AFr = anal-fin rays; UCFr = upper caudal-fin rays; LCFr = lower caudal-fin rays.
|R. reddelli (10)||I 6 (10)||6 (10)||I 8 (2), I 9 (6), I 10 (2)||11 (2), 12 (1), 13 (3), 14 (2), 15 (1), 16 (1)||8 (10)||9 (2), 10 (8)|
|R. zongolicensis (8)||I 5 (1), I 6 (7)||6 (8)||I 8 (1), I 9 (7)||11 (2), 12 (2), 13 (3), 14 (1)||7 (1), 8 (7)||10 (7), 11(1)|
|R. macuspanensis (5)||I 6 (5)||6 (5)||I 8 (3), I 10 (2)||10 (1), 11 (2), 13 (1), 14 (1)||8 (5)||9 (1), 10 (4)|
|R. laluchensis (25)||I 6 (25)||6 (25)||I 7 (1), I 8 (1), I 9 (10), I 10 (11), I 11 (2)||11 (6), 12(13), 13(4), 14(2)||7 (2), 8(22), 9(1)||9 (6), 10 (18), 11 (1)|
|R. cf. zongolicensis1 (22)||I 6 (22)||4 (1), 6 (21)||I 8 (10), I 9 (11), I 10 (1)||11 (8), 12 (11), 13 (2), 14 (1)||8 (21), 9 (1)||9 (3), 10 (18), 11 (1)|
|R. cf. zongolicensis2 (15)||I 6 (15)||6 (15)||I 8 (4), I 9 (9), I 10 (2)||12 (2), 13 (4), 14 (5), 15 (3), 16 (1)||7 (3), 8 (12)||8 (2), 9 (7), 10 (6)|
|R. laticauda (22)||I 6 (20), I 7 (2)||6 (22)||I 7 (1), I 8 (9), I 9(12)||9 (1), 10 (3), 11 (10), 12 (6), 13 (2)||6 (1), 7 (4), 8 (16)||7 (2), 10 (20)|
|R. guatemalensis (15)||I 6 (15)||6 (15)||I 7 (1), I 8 (12), I 9 (2)||10 (1), 11 (5), 12 (8), 13 (1)||8 (14), 9 (1)||10 (14), 11(1)|
Results of MANOVA on linear trait data.
|Test statistic||df1, df2||F||p-value|
|Wilks’ lambda||0.06074||21, 322.2||25.4||< 0.0001|
|Pillai’s trace||1.707||21, 342||21.5||< 0.0001|
P-values for all possible pairwise morphometric comparisons between the Rhamdia species/populations investigated in this study (Hotelling’s T2 test; α = 0.05). Values with asterisk indicate no significant differences (i.e., p-value > α). R. cf. zongolicensis1 = Sumidero de Cotzalostoc; R. cf. zongolicensis2 = Sótano de Popocatl.
|R. guatemalensis||R. laticauda||R. laluchensis||R. macuspanensis||R. reddelli||R. zongolicensis||R. cf. zongolicensis1|
|R. laticauda||< 0.0001|
|R. laluchensis||< 0.0001||< 0.0000|
|R. macuspanensis||0.0421||0.0002||< 0.0000|
|R. reddelli||< 0.0001||< 0.0001||4.1082*||0.0023|
|R. zongolicensis||< 0.0001||< 0.0001||0.0057||0.0014||2.2947*|
|R. cf. zongolicensis1||< 0.0001||< 0.0001||0.0009||0.0001||0.0009||0.0001|
|R. cf. zongolicensis2||< 0.0001||< 0.0001||0.0011||< 0.0001||0.029||0.0638*||0.0668*|
In contrast, pairwise comparisons of dorsal shape resulted in a smaller number of significantly different species pairs (Table 8). Notably, these results imply that R. zongolicensis does not differ significantly in lateral profile from either R. reddelli or R. macuspanensis, and that R. macuspanensis does not differ significantly in dorsal profile from any of the remaining species/populations investigated.
Although our focus was the pattern of relationships of Mexican troglobitic Rhamdia in the context of the Middle American radiation of the genus, our results broadly conform with previous molecule-based hypotheses of relationships inferred for non-troglobitic Rhamdia across the entire geographic range of the genus (Hernández et al., 2015; Perdices et al., 2002) (Fig. 4). The phylogenetic evidence we present here implies that Mexican troglobitic species are part of a clade that also includes the epigean species R. laticauda, R. parryi, and R. nicaraguensis. The inferred phylogenetic placement of troglobitic forms, as well as that of R. parryi and R. nicaraguensis, however, renders R. laticauda irreconcilably non-monophyletic. Under this phylogenetic framework, continued recognition of R. laticauda as currently delimited is problematic and undesirable. Notwithstanding, our phylogenetic results (CYB phylogeny) also fail to refute the monophyly of all recognized troglobitic species (with respect to R. laticauda), while revealing a modest level of lineage differentiation in R. laluchensis and R. macuspanensis. Granted, such a degree of lineage distinction does not necessarily entail interspecific differentiation, for it is also consistent with intraspecific geographic genetic structuring under the hypothesis that troglobitic species are in fact cave-adapted populations of R. laticauda. Comparatively larger mtDNA divergences between epigean and hypogean populations, have been reported for the Mexican tetra A. mexicanus (Ornelas-García et al., 2008), and yet both cave and surface forms are presently regarded as a single species by many authors (Gross, 2012; Hausdorf et al., 2011; Miller, 2005; Strecker et al., 2012). Thus, our phylogenetic findings offer support to the hypothesis that Mexican troglobitic species of Rhamdia, as well as the epigean R. parryi and R. nicaraguensis, effectively represent lineages of the more widespread species R. laticauda.
Whether these findings warrant a taxonomic overhaul, however, deserves further scrutiny. Reciprocal monophyly is a property that lineages acquire as they separate and diverge from one another, and therefore it is widely considered as evidence in support for the process of speciation (de Queiroz, 2007). Although a number of authors have argued for the existence and recognition of paraphyletic species (e.g., Crisp & Chandler, 1996; Kizirian & Donnelly, 2004), paraphyly of “ancestral” species, while certainly a possibility (particularly in speciation by peripheral isolation), would only be transient under the assumption that species ultimately evolve into monophyletic lineages (de Queiroz, 2007). Accordingly, our phylogenetic results imply that Mexican troglobitic forms of Rhamdia represent either highly localized, “exclusive” cave-dwelling populations of the widespread and primarily “non-exclusive” epigean species R. laticauda (“exclusiveness” sensu Wiens & Penkrot, 2002) or recently diverged species at very early stages of their speciation from R. laticauda. We recognize that our genetic data might be insufficient to unambiguously distinguish between these 2 contrasting hypotheses. Therefore, until these alternative scenarios can be properly tested, caution should be exercised when drawing taxonomic conclusions based solely on the phylogenetic evidence presented herein. Our intention is not to prioritize taxonomic concerns —such as the monophyly of binominals— over the recognition of diversity, but to shed light on the evolution and systematics of the genus Rhamdia by testing the species status ascribed to long-recognized Mexican troglobitic populations using multiple lines of evidence, including phylogenetic.
Notably, this is not the first time that taxonomic/nomenclatural issues involving R. laticauda and closely related species have been raised. In his revision of the genus, Silfvergrip (1996) concluded that the troglobitic species described at the time —R. zongolicensis and R. reddelli— as well as the epigean R. parryi, were junior synonyms of R. laticauda. Although he considered R. nicaraguensis valid and distinct from R. laticauda, the main character he proposed to distinguish them (i.e., shorter adipose fin) had been previously shown to display considerable overlap between the 2 species (Bussing, 1987; Silfvergrip, 1996, p. 92). Similarly, in the only previous study that sampled a Mexican troglobitic species of Rhamdia for phylogenetic analysis, the authors concluded that “preserving the species status of R. reddelli would suggest that more than 20 additional species should be recognized in the R. laticauda clade alone if a taxonomic goal is to fairly represent evolutionary history” (Perdices et al., 2002, p. 182). It should be noted that, based on their resultant phylogeny, the same conclusion would apply to R. parryi and R. nicaraguensis as it applies to R. reddelli (Perdices et al., 2002, p. 177, Fig. 2).
According to our results, among the troglobitic species, R. reddelli and R. zongolicensis are by far the most closely related phylogenetically. This conforms with the taxonomic assessment of Miller (2005) —who regarded them as synonyms— and with the expectations from geography. The type-locality caves of these 2 species are in relative proximity to one another (~ 30 km) and are part of the same hydrological basin (Papaloapan). Along this line of reasoning, our finding that R. laticauda samples from the Papaloapan are the most closely related to R. reddelli and R. zongolicensis (including R. cf. zongolicensis) is not only unsurprising but also supportive of the hypothesis that ancestral populations of R. laticauda from the region colonized hypogean habitats resulting in the evolution of the cave-adapted troglobitic forms present in the Sierra de Zongolica (Veracruz) and neighboring karstic systems (Cueva del Nacimiento del Río San Antonio, Oaxaca). The timing and other specifics about this colonization, however, have yet to be determined. Our findings also imply that R. laluchensis appears to be the earliest divergent of the Mexican troglobitic species, but details on the history of colonization of the Sótano de La Lucha (Chiapas) —as well as of the Grutas de Agua Blanca (Tabasco) in the case of R. macuspanensis— are presently unknown. Notwithstanding, our phylogenetic results suggest that the known diversity of troglomorphic Rhamdia in southern Mexico can be explained by at least 3 cave colonization events by the epigean species R. laticauda leading to the establishing of troglomorphic populations in caves of Chiapas (R. laluchensis), Tabasco (R. macuspanensis), and Veracruz/Oaxaca (R. zongolicensis and R. reddelli), respectively. A larger geographic sampling of epigean populations will be required to shed further light on the apparently complex history of cave colonization by Rhamdia in southern Mexico.
Traditional DNA barcoding species identification relies on the use of thresholds set to differentiate between intraspecific variation and interspecific divergence (Hebert et al., 2003, 2004). In fishes, the proposed ~ 3% COI sequence divergence heuristic threshold for conspecifics (Ward, 2009) appears to hold across a broad range of teleost lineages (Arroyave et al., 2019; Decru et al., 2016; Lara et al., 2010; Lowenstein et al., 2011; Pereira et al., 2011, 2013; Valdez-Moreno et al., 2009). Our phenetic analysis of COI sequence variation in Rhamdia uncovered unexpectedly low divergences (mostly < 3%) among species of the “R. laticauda-group” inclusive of R. nicaraguensis, effectively impeding unambiguous identification and distinction of most troglobitic species via DNA barcoding. Even the largest divergences among samples from this multispecies group (~ 4%) pale in comparison with the divergences between R. laticauda and demonstrably different species such as R. guatemalensis (~ 12%) or the cis-Andean R. quelen (~ 10%). The uncovered patterns of DNA barcode variation in Rhamdia are therefore consistent with the hypothesis that currently recognized troglobitic species effectively correspond to cave-adapted lineages/populations of the more widespread, epigean R. laticauda. Because of the limitations of the method, however, in and of themselves these results should not be interpreted as sufficient evidence of conspecificity among troglobitic forms.
P-values of permutation tests (α = 0.0001) of pairwise differences in mean body shape (lateral profile) based on Procrustes (standard typeface) and Mahalanobis (bold) distances. Values with asterisk indicate no significant differences (i.e., p-value > α). R. cf. zongolicensis1 = Sumidero de Cotzalostoc; R. cf. zongolicensis2 = Sótano de Popocatl.
|R. guatemalensis||R. laticauda||R. laluchensis||R. macuspanensis||R. reddelli||R. zongolicensis||R. cf. zongolicensis1||R. cf. zongolicensis2|
|R. guatemalensis||–||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001|
|R. laticauda||< 0.0001||–||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001|
|R. laluchensis||< 0.0001||< 0.0001||–||< 0.0001||0.0001||< 0.0001||< 0.0001||0.0011*|
|R. macuspanensis||< 0.0001||< 0.0001||< 0.0001||–||0.0002*||< 0.0001||< 0.0001|
|R. reddelli||< 0.0001||< 0.0001||< 0.0001||< 0.0001||–||0.036*||< 0.0001||< 0.0001|
|R. zongolicensis||< 0.0001||< 0.0001||< 0.0001||0.0002*||0.0049*||–||< 0.0001||< 0.0001|
|R. cf. zongolicensis1||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||–||0.0002*|
|R. cf. zongolicensis2||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||< 0.0001||–|
P-values of permutation tests (α = 0.0001) of pairwise differences in mean body shape (dorsal profile) based on Procrustes (standard typeface) and Mahalanobis (bold) distances. Values with asterisk indicate no significant differences (i.e., p-value > α). R. cf. zongolicensis1 = Sumidero de Cotzalostoc; R. cf. zongolicensis2 = Sótano de Popocatl.
|R. guatemalensis||R. laticauda||R. laluchensis||R. macuspanensis||R. reddelli||R. zongolicensis||R. cf. zongolicensis1||R. cf. zongolicensis2|
|R. guatemalensis||–||0.001*||< 0.0001||0.1707*||< 0.0001||0.0002*||0.0501*||0.1051*|
|R. laticauda||0.0001||–||< 0.0001||0.1155*||< 0.0001||< 0.0001||< 0.0001||< 0.0001|
|R. laluchensis||< 0.0001||< 0.0001||–||0.0031*||< 0.0001||0.0002*||< 0.0001||0.0018*|
|R. reddelli||< 0.0001||< 0.0001||0.0001||0.0011*||–||0.0011*||< 0.0001||< 0.0001|
|R. zongolicensis||< 0.0001||< 0.0001||0.0001||0.0007*||0.0001||–||< 0.0001||0.0022*|
|R. cf. zongolicensis1||< 0.0001||< 0.0001||< 0.0001||0.006*||< 0.0001||< 0.0001||–||0.0465*|
|R. cf. zongolicensis2||0.0001||< 0.0001||0.0031*||0.0025*||< 0.0001||0.0083*||0.0144*||–|
In support for the recognition of R. reddelli and R. zongolicensis, Wilkens argued that “considerable genetic divergences exist between them” and therefore “morphological similarity is the result of convergent evolution and the cave catfishes should be treated as separate species” (Wilkens, 2001, p. 260). Our results from both phylogenetic and phenetic analyses of comparative mtDNA data, however, imply that the degree and pattern of genetic divergences among Mexican troglobitic populations is considerably lower than the expected for samples belonging to different species (“contra” Wilkens, 2001). Remarkably, some of the evidence brought up by Wilkens to support this claim consisted of unpublished data: “the evolution of separate gene pools is proven by mtDNA-analysis (unpublished data)” (Wilkens, 2001, p. 260). Our own analyses of mtDNA data strongly dispute this statement. The remaining genetic evidence cited by Wilkens to endorse the species status of R. reddelli and R. zongolicensis came from previous work in which laboratory crosses between these species and the epigean R. laticauda resulted in F1s consisting exclusively of individuals from the heterogametic sex (in this case females), conforming to Haldane’s Rule (Wilkens, 1993). Hybrid dysfunction under Haldane’s Rule, however, applies to early stages of speciation (Coyne & Orr, 1989; Haldane, 1922). Furthermore, most research on Haldane’s rule has been conducted on Drosophila, thus making generalizations and extensions to other taxa, especially poorly studied ones and with heterogametic females, suspect at best (Johnson, 2008). While we acknowledge that our results are also consistent with patterns expected at early stages of the speciation process, we believe that current available evidence appears insufficient to assert absolute divergence, reproductive isolation, and lineage differentiation worthy of species status in any form of troglobitic Rhamdia from southern Mexico.
Despite previous reports about the coexistence of individuals with varying degrees of troglomorphism in populations of Mexican troglobitic Rhamdia, the implications and generality of this pattern had not been critically investigated (Mosier, 1984; Sbordoni et al., 1986; Wilkens, 2001). Our qualitative assessment of intraspecific variation in ocular development and degree of pigmentation not only corroborates the aforementioned reports, but also reveals that this phenomenon is more widespread than previously documented, effectively occurring in all cave populations examined (Fig. 6), and even manifested through previously unreported syntopy of epigean and troglobitic forms in R. macuspanensis (Fig. 7). Regarding the abovementioned syntopy of forms, we tentatively attribute the absence of intermediate phenotypes in R. macuspanensis to sampling error (effectively the smallest sample size at n = 5).
Intermediate troglobitic phenotypes have been documented both in the wild and in F1 laboratory crosses between surface and cave forms of the Mexican tetra A. mexicanus, and their existence attributed to either recent divergence from surface populations (i.e., not enough time since hypogean colonization for fixation of troglomorphic traits) or present-day hybridization between cave and surface fish (Wilkens, 2016). All populations of troglobitic Rhamdia investigated in this study live in pools inside caves that are supplied, either perennially or intermittently, by surface water (e.g., streams, surface runoff) and/or groundwater from karstic springs. Therefore, it is reasonable to assume that individuals from the epigean R. laticauda have the potential to “hybridize” with cave-dwelling fish. Remarkably, complete replacement of a cave-dwelling troglomorphic population by its epigean ancestor has been previously documented in Rhamdia, namely R. quelen from the Cumaca cave in Trinidad, West Indies (Romero et al., 2002). Originally described as Caecorhamdia urichi Norman 1926, this troglomorphic population was subsequently synonymized with R. quelen precisely because of its lack of diagnostic features other than the absence of eyes and depigmentation (Mees, 1974; Silfvergrip, 1996). Varying degree of troglomorphism had been reported for this population since the 1950s and, by the early 2000s, the troglobitic phenotype had been completely substituted with the surface form, possibly prompted by changes in precipitation regimes and the ability of epigean individuals to outcompete troglobitic ones (Romero et al., 2002). Without entirely discounting the possibility of intraspecific phenotypic plasticity (Bilandžija et al., 2020), we hypothesize that the observed clinal variation in troglomorphism in Mexican cave-dwelling species/populations of Rhamdia (Figs. 6, 7) is the result of crossings between epigean and hypogean forms rather than recent divergence followed by only partial fixation of troglomorphism. Further research into this subject (e.g., QTL analysis of genome-wide data), however, is required to properly test this hypothesis and better understand the genetic basis of regressive phenotypic loss in Mexican cave-dwelling Rhamdia.
Of the 6 meristic variables investigated, dorsal-fin and pelvic-fin ray counts were the ones that exhibited the least amount of variation across all species/populations evaluated (6 rays in almost all samples). In contrast, anal-fin ray counts displayed the widest range of variation (9-16), although with most values concentrated at 11-13 rays, and a slight tendency of troglobitic forms to display a comparatively higher ray count. The remaining variables, while displaying more intermediate ranges, still exhibited a tendency to concentrate most records on a single value (i.e., a clearly defined mode). Overall, none of these meristic traits display a non-overlapping distribution with potential to diagnose or distinguish among troglobitic species, nor between these and the epigean R. laticauda. These findings are consistent with a recurrent issue in catfish systematics, that is, the fact that meristic features, while not necessarily uninformative, are plagued by considerable overlap when a significant number of taxa is compared (Silfvergrip, 1996; p. 120).
Although exploratory analysis of body shape variation as implied by linear data (Fig. 8) did not seem to unambiguously differentiate among troglobitic species, as well as between these and the epigean R. laticauda, MANOVA results did support the existence of significant differences among most of the species/populations investigated (Tables 5, 6). Similarly, GM analyses initially suggested a lack of absolute and clear-cut differentiation in overall body shape among troglobitic forms and between them and R. laticauda (Figs. 9, 10), but subsequently demonstrated the existence of statistically significant differences, at least in the lateral-view component of body shape (Table 7). Whether such differences constitute proof of species-level differentiation (vs. intraspecific variation), however, requires further scrutiny. On one hand, statistically significant differences in body shape between troglobitic species/populations and R. laticauda are ultimately phenetic and therefore cannot be unequivocally attributed to historical processes such as speciation. On the other hand, such differentiation appears to be mainly related with phenotypically plastic features associated with life in hypogean habitats (e.g., body depth, head length).
Despite being a character system explaining a considerable proportion of the observed variation in linear trait data (PC1 in figure 8), barbels are not consistently longer in troglobitic forms of Mexican Rhamdia with respect to epigean species. This finding is certainly contrary to our expectations, since barbel elongation has long being considered a typical constructive adaptation in troglobitic fishes (Romero & Green, 2005; Soares & Niemiller, 2020). Without discounting other possible explanations, this counterintuitive finding could be attributed to the relative recency of the divergence between R. laticauda and hypogean species/populations (none of the cave-dwelling forms has on average shorter barbels than R. laticauda). Furthermore, the presence of comparatively longer barbels in the epigean R. guatemalensis is unsurprising when considering its strong nyctophilia (Arroyave et al., 2020; Wilkens, 2001).
Among the shape attributes that appear to offer some support for the morphological distinction between R. laticauda and currently recognized troglobitic species are body depth (PC2, RW1, and CV1 in figures 8a, 9a, 10a, respectively) and head size and shape coupled with the relative position of the pectoral and pelvic fins (RW1 and CV1 in figures 9b, 10b, respectively). While variation in body depth seems potentially useful to discriminate troglobitic species (except R. macuspanensis) from the epigean R. laticauda, this shape attribute is potentially problematic because cave fishes tend to have slenderer bodies than their epigean counterparts due to the fact that caves are often nutrient limited (Langecker, 2000; Venarsky et al., 2014). Notably, our findings that troglobitic species tend to have longer and sharper heads relative to R. laticauda are in agreement with each and every one of the original diagnoses, which also list the presence of a weak and short occipital process as diagnostic (Miller, 1984; Wilkens, 1993; Weber & Wilkens, 1998; Weber et al., 2003). The fact that troglobitic forms have longer heads begs the question of whether head enlargement in troglobitic Rhamdia, just like barbel elongation, constitutes an instance of constructive troglomorphism and convergent adaptation to hypogean habitats, and therefore of questionable taxonomic utility. While typical constructive troglomorphic traits involve improved chemosensory capacities through elongated sensory structures (Jones & Culver, 1989), head enlargement might be driven solely by functional necessity, to accommodate hypertrophied sensory receptors (olfactory, tactile, lateral-line, etc.) and brain centers (Poulson, 1963). Thus, even if only developmentally and mechanically correlated with traits under positive selection in hypogean habitats (enhanced non-visual sensory systems), head enlargement would be an indirect byproduct of selection, and so, effectively a convergent adaptation.
This study shows that the validity of currently recognized troglobitic species of Mexican Rhamdia is unequivocally supported by neither genetic nor morphological evidence. Meristic and morphometric variation appear insufficient to unambiguously and robustly diagnose troglobitic species and distinguish them from each other. Furthermore, beyond typical regressive troglomorphic traits, troglobitic species/populations do not differ greatly in external morphology from their most closely related congener, the epigean species R. laticauda. Although statistically significant differences in body shape exist between troglobitic species/populations and R. laticauda, such differences might not be the result of historical processes such as speciation but of convergent adaptation in phenotypically plastic traits instead. Based on the phylogenetic evidence presented here, continued recognition of all troglobitic species, as well as of the epigean R. nicaraguensis and R. parryi, implies a deep and generalized paraphyly in R. laticauda. Although the most parsimonious solution to this paraphyly is for the species involved to be synonymized with R. laticauda, we must admit that our geographic sampling of R. nicaraguensis and R. parryi is limited and based on vouchers unavailable for direct examination. Moreover, and more relevant to our research question, we recognize that the phylogenetic pattern uncovered herein is also consistent with a very recent divergence and early speciation of cave-dwelling populations from R. laticauda. Distinguishing between recent speciation of cave-dwelling lineages and intraspecific genetic structuring in R. laticauda (shaped by ongoing —albeit restricted— gene flow between epigean and hypogean populations), would therefore be the ultimate test for the validity of the current taxonomy. Discrimination of population and species boundaries in this system, however, will likely require additional comparative genetic data —ideally with genome-wide representation (e.g., genome-wide SNP data)— and state-of-the-art analytical approaches (e.g., coalescent-based species delimitation methods). Accordingly, in spite of the morphological and molecular evidence presented herein —which unarguably casts doubt on the validity of Mexican troglobitic species of Rhamdia— we refrain from making nomenclatural decisions until the abovementioned alternative evolutionary scenarios can be properly tested.
This research was financially supported by the Universidad Nacional Autónoma de México (UNAM) through a “Programa de Apoyo a Proyectos de Investigación e Innovación Tecnológica (PAPIIT)” grant (IA202119) to JA. We are especially grateful to speleologists Hugo Salgado and Victor Cruz for their invaluable assistance in the challenging endeavor of sampling fishes from deep caves. We would also like to thank Pedro Pablo Cruz Cano and Manuel Fuentes Ayohua from Zongolica (Veracruz), and Benjamín Hernandez from La Lucha (Chiapas), for their hospitality and instrumental support during our collecting trips to these regions. Likewise, we thank our colleagues and friends Christian de Jesús Narcia Rico (UNICACH) and Carlos Garita Alvarado (UNAM) for their assistance and support in the field during our fish-collecting campaigns. Special thanks to Andrea Jiménez Marín (Laboratorio Nacional de Biodiversidad, LaNaBio) for her assistance during the generation of comparative DNA sequence data. This work is dedicated to the memory of Adán Gómez González (RIP), Mexican ichthyologist and dear friend who inspired and motivated JA to pursue research on troglomorphic Rhamdia.
Adams, D. C., & Otárola-Castillo, E. (2013). geomorph: an R package for the collection and analysis of geometric morphometric shape data. Methods in Ecology and Evolution, 4, 393–399. https://doi.org/10.1111/2041-210X.12035
Arroyave, J., Martinez, C. M., Martínez-Oriol, F. H., Sosa, E., & Alter, S. E. (2020). Regional-scale aquifer hydrogeology as a driver of phylogeographic structure in the Neotropical catfish Rhamdia guatemalensis (Siluriformes: Heptapteridae) from cenotes of the Yucatán Peninsula, Mexico. Freshwater Biology, 2020, 1–17. https://doi.org/10.1111/fwb.13641
Arroyave, J., Martinez, C. M., & Stiassny, M. L. J. (2019). DNA barcoding uncovers extensive cryptic diversity in the African long-fin tetra Bryconalestes longipinnis (Alestidae: Characiformes). Journal of Fish Biology, 95, 379–392. https://doi.org/10.1111/jfb.13987
Bailey, R. C., & Byrnes, J. (1990). A new, old method for assessing measurement error in both univariate and multivariate morphometric studies. Systematic Biology, 39, 124–130. https://doi.org/10.2307/2992450
Bichuette, M. E., & Trajano, E. (2005). A new cave species of Rhamdia (Siluriformes: Heptapteridae) from Serra do Ramalho, northeastern Brazil, with notes on ecology and behavior. Neotropical Ichthyology, 3, 587–595. https://doi.org/10.1590/S1679-62252005000400016
Bilandžija, H., Hollifield, B., Steck, M., Meng, G., Ng, M., Koch, A. D. et al. (2020). Phenotypic plasticity as a mechanism of cave colonization and adaptation. eLife, 9, e51830. https://doi.org/10.7554/elife.51830
Bockmann, F. A. (1998). Análise filogenética da família Heptapteridae (Teleostei, Ostariophysi, Siluriformes) e redefinição de seus gêneros (PhD. Thesis). Universidade de São Paulo, São Paulo, Brasil.
Bookstein, F. L. (1989). Principal warps: thin-plate splines and the decomposition of deformations. IEEE Transactions on Pattern Analysis and Machine Intelligence, 11, 567–585. https://doi.org/10.1109/34.24792
Bookstein, F. L. (1997). Landmark methods for forms without landmarks: morphometrics of group differences in outline shape. Medical Image Analysis, 1, 225–243. https://doi.org/10.1016/S1361-8415(97)85012-8
Bussing, W. A. (1987). Peces de las aguas continentales de Costa Rica. Editorial Universidad de Costa Rica.
Čandek, K., & Kuntner, M. (2015). DNA barcoding gap: reliable species identification over morphological and geographical scales. Molecular Ecology Resources, 15, 268–277. https://doi.org/10.1111/1755-0998.12304
Coyne, J. A., & Orr, H. A. (1989). Patterns of speciation in Drosophila. Evolution, 43, 362–381. https://doi.org/10.1111/j.1558-5646.1989.tb04233.x
Crisp, M. D., & Chandler, G. T. (1996). Paraphyletic species. Telopea, 6, 813–844. http://dx.doi.org/10.7751/telopea19963037
Darriba, D., Taboada, G. L., Doallo, R., & Posada, D. (2012). jModelTest2: more models, new heuristics and parallel computing. Nature Methods, 9, 772–772. https://doi.org/10.1038/nmeth.2109
Dayrat, B. (2005). Towards integrative taxonomy. Biological Journal of the Linnean Society, 85, 407–415. https://doi.org/10.1111/j.1095-8312.2005.00503.x
de Carvalho, D. C., Oliveira, D. A., Pompeu, P. S., Leal, C. G., Oliveira, C. et al. (2011). Deep barcode divergence in Brazilian freshwater fishes: the case of the São Francisco River basin. Mitochondrial DNA, 22 (Suppl.), 80–86. https://doi.org/10.3109/19401736.2011.588214
de Queiroz, K. (2007). Species concepts and species delimitation. Systematic Biology, 56, 879–886. https://doi.org/10.1080/10635150701701083
Decru, E., Moelants, T., Gelas, K. D., Vreven, E., Verheyen, E., & Snoeks, J. (2016). Taxonomic challenges in freshwater fishes: a mismatch between morphology and DNA barcoding in fish of the north-eastern part of the Congo basin. Molecular Ecology Resources, 16, 342–352. https://doi.org/10.1111/1755-0998.12445
Do Nascimiento, C., Provenzano, F., & Lundberg, J. G. (2004). Rhamdia guasarensis (Siluriformes: Heptapteridae), a new species of cave catfish from the Sierra de Perijá, northwestern Venezuela. Proceedings of the Biological Society of Washington, 117, 564–574. https://doi.org/10.1590/S1679-62252005000400016
Edgar, R. C. (2004). MUSCLE: a multiple sequence alignment method with reduced time and space complexity. BMC Bioinformatics, 5, 113. http://dx.doi.org/10.1186/1471-2105-5-113
Folmer, O., Hoeh, W. R., Black, M. B., & Vrijenhoek, R. C. (1994). DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Molecular Marine Biology and Biotechnology, 3, 294–299.
Fricke, R., Eschmeyer, W. N., & Van der Laan, R. (2020). Eschmeyer’s catalog of fishes: genera, species, references. Accessed Jul 18, 2020. http://researcharchive.calacademy.org/research/ichthyology/catalog/fishcatmain.asp
Garavello, J. C., & Shibatta, O. A. (2016). Reappraisal of Rhamdia branneri Haseman, 1911 and R. voulezi Haseman, 1911 (Siluriformes: Heptapteridae) from the rio Iguaçu with notes on their morphometry and karyotype. Neotropical Ichthyology, 14. https://doi.org/10.1590/1982-0224-20140111
García-Dávila, C., Castro-Ruiz, D., Renno, J. F., Chota-Macuyama, W., Carvajal-Vallejos, F. M., Sanchez, H. et al. (2015). Using barcoding of larvae for investigating the breeding seasons of Pimelodid catfishes from the Marañon, Napo and Ucayali rivers in the Peruvian Amazon. Journal of Applied Ichthyology, 31, 40–51. https://doi.org/10.1111/jai.12987
Greenfield, D. W., Greenfield, T. A., & Woods, R. L. (1982). A new subspecies of cave-dwelling pimelodid catfish, Rhamdia laticauda typhla from Belize, Central America. Brenesia, 19/20, 563–576.
Gross, J. B. (2012). The complex origin of Astyanax cavefish. BMC Evolutionary Biology, 12, 1–12. https://doi.org/10.1186/1471-2148-12-105
Gunz, P., Mitteroecker, P., & Bookstein, F. L. (2005). Semi-
landmarks in three dimensions. In D. E. Slice (Eds.), Modern morphometrics in physical Anthropology. Developments in Primatology: progress and prospects. Boston, MA: Springer. https://doi.org/10.1007/0-387-27614-9_3
Haldane, J. B. S. (1922). Sex ratio and unisexual sterility in hybrid animals. Journal of Genetics, 12, 101–109. https://doi.org/10.1007/BF02983075
Hammer, Ø., Harper, D. A., & Ryan, P. D. (2001). PAST: Paleontological statistics software package for education and data analysis. Palaeontologia electronica, 4, 9.
Hausdorf, B., Wilkens, H., & Strecker, U. (2011). Population genetic patterns revealed by microsatellite data challenge the mitochondrial DNA based taxonomy of Astyanax in Mexico (Characidae, Teleostei). Molecular Phylogenetics and Evolution, 60, 89–97. https://doi.org/10.1016/j.ympev.2011.03.009
Hebert, P. D. N., Ratnasingham, S., & deWaard, J. R. (2003). Barcoding animal life: cytochrome c oxidase subunit 1 divergences among closely related species. Proceedings of the Royal Society B: Biological Sciences, 270, S96–S99. https://doi.org/10.1098/rsbl.2003.0025
Hebert, P. D. N., Stoeckle, M. Y., Zemlak, T. S., & Francis, C. M. (2004). Identification of birds through DNA barcodes. Plos Biology, 2, e312. https://doi.org/10.1371/journal.pbio.0020312
Helfman, G., Collette, B. B., Facey, D. E., & Bowen, B. W. (2009). The diversity of fishes: biology, evolution, and ecology. Hoboken, New Jersey: Wiley-Blackwell.
Hernández, C. L., Ortega-Lara, A., Sánchez-Garcés, G. C., & Alford, M. H. (2015). Genetic and morphometric evidence for the recognition of several recently synonymized species of trans-Andean Rhamdia (Pisces: Siluriformes: Heptapteridae). Copeia, 103, 563–579. https://doi.org/10.1643/CI-14-145
Hubbs, C. L. (1936). Fishes of the Yucatan peninsula. In A. S. Pearse, E. P. Creaser, & F. G. Hall (Eds.), The cenotes of Yucatán: a zoological and hydrographic survey (pp. 157–287). Washington D.C.: Carnegie Institution of Washington Publications.
Hubbs, C. L. (1938). Fishes from the caves of Yucatán. Publications of the Carnegie Institution of Washington, 491, 261–295.
Irwin, D. M., Kocher, T. D., & Wilson, A. C. (1991). Evolution of the cytochrome b gene of mammals. Journal of Molecular Evolution, 32, 128–144. https://doi.org/10.1007/bf02515385
Johnson, N. A. (2008). Haldane’s rule: the heterogametic sex. Nature Education, 1, 1–6.
Jones, R., & Culver, D. C. (1989). Evidence for selection on sensory structures in a cave population of Gammarus minus (Amphipoda). Evolution, 43, 688–693. https://doi.org/10.2307/2409074
Keene, A., Yoshizawa, M., & McGaugh, S. E. (2015). Biology and evolution of the Mexican cavefish. Cambridge, Massachusetts: Academic Press. https://doi.org/10.1016/C2014-0-01426-8
Kizirian, D., & Donnelly, M. A. (2004). The criterion of reciprocal monophyly and classification of nested diversity at the species level. Molecular Phylogenetics and Evolution, 32, 1072–1076. https://doi.org/10.1016/j.ympev.2004.05.001
Klingenberg, C. P. (2011). MorphoJ: an integrated software package for geometric morphometrics. Molecular Ecology Resources, 11, 353–357. https://doi.org/10.1111/j.1755-0998.2010.02924.x
Kocher, T. D., Thomas, W. K., Meyer, A., Edwards, S. V., Pääbo, S., Villablanca, F. X. et al. (1989). Dynamics of mitochondrial DNA evolution in animals: amplification and sequencing with conserved primers. Proceedings of the National Academy of Sciences, 86, 6196–6200. https://doi.org/10.1073/pnas.86.16.6196
Langecker, T. G. (2000). The effects of continuous darkness on cave ecology and cavernicolous evolution. In H. Wilkens, D. C. Culver, & W. F. Humphreys (Eds.), Subterranean ecosystems (pp. 135–157). Amsterdam: Elsevier.
Lara, A., León, J. L. P., Rodríguez, R., Casane, D., Côté, G., Bernatchez, L. et al. (2010). DNA barcoding of Cuban freshwater fishes: evidence for cryptic species and taxonomic conflicts. Molecular Ecology Resources, 10, 421–430. https://doi.org/10.1111/j.1755-0998.2009.02785.x
Lowenstein, J. H., Osmundson, T. W., Becker, S., Hanner, R., & Stiassny, M. L. J. (2011). Incorporating DNA barcodes into a multi-year inventory of the fishes of the hyperdiverse Lower Congo River, with a multi-gene performance assessment of the genus Labeo as a case study. Mitochondrial DNA, 22, 52–70. https://doi.org/10.3109/19401736.2010.537748
Mees, G. F. (1974). The Auchenipteridae and Pimelodidae of Suriname (Pisces, Nematognathi). Zoologische Verhandelingen, 132, 1–246.
Miller, R. R. (1984). Rhamdia reddelli new species, the first blind pimelodid catfish from Middle America, with a key to the Mexican species. Transactions of the San Diego Society of Natural History, 20, 135–144.
Miller, R. R. (2005). Freshwater fishes of Mexico (with the collaboration of WL Minkley and SM Norris). Chicago: The University of Chicago Press.
Minton, M., & Droms, Y. (2019). Exploration of caves —vertical caving techniques. Chapter 49. In W. B. White, D. C. Culver, & T. Pipan (Eds.), Encyclopedia of caves (Third edition). (pp. 420–425). Cambridge, Massachusetts: Academic Press.
Møller, P. R., Knudsen, S. W., Schwarzhans, W., & Nielsen, J. G. (2016). A new classification of viviparous brotulas (Bythitidae) –with family status for Dinematichthyidae– based on molecular, morphological and fossil data. Molecular Phylogenetics and Evolution, 100, 391–408. https://doi.org/10.1016/j.ympev.2016.04.008
Mosier, D. (1984). Cave dwelling populations of Rhamdia (Pimelodidae). Association for Mexican Cave Studies Activities Newsletter, 14, 40–44.
Nicholas, G. (1962). Checklist of troglobitic organisms of Middle America. The American Midland Naturalist, 68, 165–188.
Nickum, J., Bart Jr, H. L., Bowser, P. R., Greer, I. E., Hubbs, C., Jenkins, J. A. et al. (2004). Guidelines for the use of fishes in research. Fisheries-Bethesda, 29, 26–26.
Ornelas-García, C. P., Domínguez-Domínguez, O., & Doadrio, I. (2008). Evolutionary history of the fish genus Astyanax Baird and Girard (1854) (Actinopterygii, Characidae) in Mesoamerica reveals multiple morphological homoplasies. BMC Evolutionary Biology, 8, 340. https://dx.doi.org/10.1186%2F1471-2148-8-340
Padial, J. M., Miralles, A., De la Riva, I., & Vences, M. (2010). The integrative future of taxonomy. Frontiers in Zoology, 7, 16. https://doi.org/10.1186/1742-9994-7-16
Paradis, E., Claude, J., & Strimmer, K. (2004). APE: Analyses of Phylogenetics and Evolution in R language. Bioinformatics, 20, 289–290. https://doi.org/10.1093/bioinformatics/btg412
Perdices, A., Bermingham, E., Montilla, A., & Doadrio, I. (2002). Evolutionary history of the genus Rhamdia (Teleostei: Pimelodidae) in Central America. Molecular Phylogenetics and Evolution, 25, 172–189. https://doi.org/10.1016/S1055-7903(02)00224-5
Pereira, L. H., Hanner, R., Foresti, F., & Oliveira, C. (2013). Can DNA barcoding accurately discriminate megadiverse Neotropical freshwater fish fauna? BMC Genetics, 14, 20. https://doi.org/10.1186/1471-2156-14-20
Pereira, L. H. G., Maia, G. M. G., Hanner, R., Foresti, F., & Oliveira, C. (2011). DNA barcodes discriminate freshwater fishes from the Paraíba do Sul River Basin, São Paulo, Brazil. Mitochondrial DNA, 22, 71–79. https://doi.org/10.3109/19401736.2010.532213
Popescu, A., Huber, K. T., & Paradis, E. (2012). ape 3.0: New tools for distance-based phylogenetics and evolutionary analysis in R. Bioinformatics, 28, 1536–1537. https://doi.org/10.1093/bioinformatics/bts184
Porter, M. L., & Crandall, K. A. (2003). Lost along the way: the significance of evolution in reverse. Trends in Ecology and Evolution, 18, 541–547. https://doi.org/10.1016/S0169-5347(03)00244-1
Poulson, T. L. (1963). Cave adaptation in amblyopsid fishes. American Midland Naturalist, 70, 257–290. https://doi.org/10.2307/2423056
Proudlove, G. S. (2006). Subterranean fishes of the world: an account of the subterranean (hypogean) fishes described up to 2003 with a bibliography 1541-2004. Moulis, France: International Society for Subterranean Biology.
QGIS Development Team. (2020). QGIS Geographic Information System. Open Source Geospatial Foundation Project.
Ratnasingham, S., & Hebert, P. D. N. (2007). bold: the Barcode of Life Data System (http://www.barcodinglife.org). Molecular Ecology Notes, 7, 355–364. https://dx.doi.org/10.1111%2Fj.1471-8286.2007.01678.x
Reddell, J. R. (1981). A review of the cavernicole fauna of Mexico, Guatemala, and Belize. Texas Memorial Museum, The University of Texas at Austin. No. 591 R4.
Rétaux, S., & Casane, D. (2013). Evolution of eye development in the darkness of caves: adaptation, drift, or both? EvoDevo, 4, 1–12. https://doi.org/10.1186/2041-9139-4-26
Ríos, N., Bouza, C., Gutiérrez, V., & García, G. (2017). Species complex delimitation and patterns of population structure at different geographic scales in Neotropical silver catfish (Rhamdia: Heptapteridae). Environmental Biology of Fishes, 100, 1047–1067. https://doi.org/10.1007/s10641-017-0622-1
Robertson, S. (1983). Zongolica Project 1983. Association for Mexican Cave Studies Activities Newsletter, 13, 36–41.
Rohlf, F. J. (2015). The tps series of software. Hystrix, the Italian Journal of Mammalogy, 26, 9–12. https://doi.org/10.4404/hystrix-26.1-11264
Romero, A., & Green, S. M. (2005). The end of regressive evolution: examining and interpreting the evidence from cave fishes. Journal of Fish Biology, 67, 3–32. https://doi.org/10.1111/j.0022-1112.2005.00776.x
Romero, A., & Paulson, K. M. (2001). It’s a wonderful hypogean life: a guide to the troglomorphic fishes of the world. In A. Romero (Ed.), The biology of hypogean fishes (pp. 13–41). Dordrecht: Springer.
Romero, A., Singh, A., McKie, A., Manna, M., Baker, R., Paulson, K. M. et al. (2002). Replacement of the troglomorphic population of Rhamdia quelen (Pisces: Pimelodidae) by an epigean population of the same species in the Cumaca Cave, Trinidad, West Indies. Copeia, 2002, 938–942. https://doi.org/10.1643/0045-8511(2002)002[0938:ROTTPO]2.0.CO;2
Sabaj, M. H. (2020). Codes for natural history collections in ichthyology and herpetology. Copeia, 108, 593–669. https://doi.org/10.1643/ASIHCODONS2020
Salinas, N. R., & Little, D. P. (2014). 2matrix: A utility for indel coding and phylogenetic matrix concatenation1. Applications in Plant Sciences, 2, 1300083. https://dx.doi.org/10.3732%2Fapps.1300083
Sbordoni, V., Argano, R., & Vomero, V. (1986). Relazione biologica sulle spedizioni “Malpaso” 1981–82 e 1984. Notiziario del Circolo Speleologico Romano, Nouva Serie, 1, 73–88.
Schliep, K. P. (2011). Phangorn: phylogenetic analysis in R. Bioinformatics, 27, 592–593. https://doi.org/10.1093/bioinformatics/btq706
Silfvergrip, A. M. C. (1996). A systematic revision of the Neotropical catfish genus Rhamdia (Teleostei, Pimelodidae) (PhD. Thesis). Stockholm University.
Soares, D., & Niemiller, M. L. (2020). Extreme adaptation in caves. The Anatomical Record, 303, 15–23. https://doi.org/10.1002/ar.24044
Sokal, R. R., & Michener, C. D. (1958). A statistical method for evaluating systematic relationship. University of Kansas Science Bulletin, 28, 1409–1438.
Stamatakis, A. (2006). RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics, 22, 2688–2690. https://doi.org/10.1093/bioinformatics/btl446
Stamatakis, A., Hoover, P., & Rougemont, J. (2008). A rapid bootstrap algorithm for the RAxML web servers. Systematic Biology, 57, 758–771. https://doi.org/10.1080/10635150802429642
Strecker, U., Hausdorf, B., & Wilkens, H. (2012). Parallel speciation in Astyanax cave fish (Teleostei) in Northern Mexico. Molecular Phylogenetics and Evolution, 62, 62–70. https://doi.org/10.1016/j.ympev.2011.09.005
Usso, M. C., Santos, A. R. D., Gouveia, J. G., Frantine-Silva, W., Araya-Jaime, C., Oliveira, M. L. M. D. et al. (2019). Genetic and chromosomal differentiation of Rhamdia quelen (Siluriformes, Heptapteridae) revealed by repetitive molecular markers and DNA barcoding. Zebrafish, 16, 87–97. https://doi.org/10.1089/zeb.2018.1576
Valdez-Moreno, M., Ivanova, N. V., Elías-Gutiérrez, M., Contreras-Balderas, S., & Hebert, P. D. N. (2009). Probing diversity in freshwater fishes from Mexico and Guatemala with DNA barcodes. Journal of Fish Biology, 74, 377–402. https://doi.org/10.1111/j.1095-8649.2008.02077.x
Venarsky, M. P., Huntsman, B. M., Huryn, A. D., Benstead, J. P., & Kuhajda, B. R. (2014). Quantitative food web analysis supports the energy-limitation hypothesis in cave stream ecosystems. Oecologia, 176, 859–869. https://doi.org/10.1007/s00442-014-3042-3
Walsh, S. J., & Chakrabarty, P. (2016). A new genus and species of blind sleeper (Teleostei: Eleotridae) from Oaxaca, Mexico: first obligate cave gobiiform in the Western Hemisphere. Copeia, 104, 506–517. https://doi.org/10.1643/CI-15-275
Ward, R. D. (2009). DNA barcode divergence among species and genera of birds and fishes. Molecular Ecology Resources, 9, 1077–1085. https://doi.org/10.1111/j.1755-0998.2009.02541.x
Weber, A., Allegrucci, G., & Sbordoni, V. (2003). Rhamdia laluchensis, a new species of troglobitic catfish (Siluriformes: Pimelodidae) from Chiapas, Mexico. Ichthyological Exploration of Freshwaters, 14, 273–280.
Weber, A., & Wilkens, H. (1998). Rhamdia macuspanensis: A Nnew Sspecies of troglobitic pimelodid catfish (Siluriformes; Pimelodidae) from a cave in Tabasco, Mexico. Copeia, 1998, 998–1004. https://doi.org/10.2307/1447347
Wiens, J. J., & Penkrot, T. A. (2002). Delimiting species using DNA and morphological variation and discordant species limits in spiny lizards (Sceloporus). Systematic Biology, 51, 69–91. https://doi.org/10.1080/106351502753475880
Wilkens, H. (1993). A new species of Rhamdia (Pisces: Pimelodidae) from a cave in the Sierra de Zongolica (Veracruz, México). Mitteilungen aus dem Hamburgischen Zoologischen Museum und Institut, 90, 375–378.
Wilkens, H. (2001). Convergent adaptations to cave life in the Rhamdia laticauda catfish group (Pimelodidae, Teleostei). Environmental Biology of Fishes, 62, 251–261. https://doi.org/10.1023/A:1011897805681
Wilkens, H. (2016). Genetics and hybridization in surface and cave Astyanax (Teleostei): a comparison of regressive and constructive traits. Biological Journal of the Linnean Society, 118, 911–928. https://doi.org/10.1111/bij.12773
Wilkens, H., & Strecker, U. (2017). Evolution in the dark: Darwin’s loss without selection. Berlin: Springer. https://doi.org/10.1007/978-3-662-54512-6