Reappraisal of the systematics of Microglanis cottoides (Siluriformes, Pseudopimelodidae), a catfish from southern Brazil

The southern region of Brazil is characterized by high species diversity and endemism of freshwater fishes distributed across geographically isolated river basins. Microglanis cottoides has a widespread range across these river basins and occurs in sympatry with other endemic species of the genus (e.g. M. cibelae, M. eurystoma, and M. malabarbai). Herein we tested the monophyly of M. cottoides and presented for the first time information about the molecular phylogeny of species in the genus. The results suggest that M. cottoides currently forms a non-monophyletic group which includes populations endemic to the Uruguay River basin that are more closely related to M. malabarbai, and excludes M. cibelae, found to be nested within M. cottoides. Based on an integrative approach using morphological and molecular data, we propose M. cibelae as a junior synonym of M. cottoides, and the populations of the Uruguay River basin previously assigned to M. cottoides in fact belong to M. malabarbai. Our molecular phylogeny shows that M. cottoides is sister to M. parahybae, which is also a coastal species, and M. malabarbai is sister of M. garavelloi, both endemic to inland river basins. The time-calibrated phylogeny indicates that the separation between inland and the coastal clades occurred in the Tertiary period, and that the species within the coastal basins diverged in the Pliocene, which overlaps with the diversification times estimated for the two inland species as well. This pattern of diversification corroborates some previous studies with other fishes from the same region.


Introduction
In southern South America, between 6 and 10 Kya, sea level oscillations in the Atlantic formed a series of coastal watersheds geographically isolated from each other, which drain directly a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 into the ocean [1], [2]. In addition, this region comprises the Uruguay River basin, an inland drainage isolated from the coastal ones by the steep mountains on the eastern shore of the Brazilian crystalline shield [1]. As results of geomorphological and sea level variations, studies have identified smaller ecoregions [3] with high levels of endemism across coastal [4], [5], [6], [7] and inland drainages [8], [9], [10], [11] in southern South America.
Microglanis Eigenmann, 1912 (Pseudopimelodidae) are a group of fishes with small body size (standard length of less than 11 cm [12]) comprising 28 species distributed across South America [13], [14], [15], [16]. The geographic distribution of these fishes is affected by physioecological constraints since Microglanis do not tolerate salt water, which prevent species from dispersing through the sea, sensu Myers  America have not yet been revised on the basis of genetics and morphological data using a statistical framework to test whether current taxonomic decisions would receive support.
Identifications and descriptions of species solely based on morphological characters may present limitations when the characters used are masked by the effect of phenotypic plasticity and/or genetic variability [25]. This scenario would bring instability to species recognition and produce biased reports of species distribution and diversity. For this reason, in recent decades, several molecular tools have been used collectively with morphological data to study biodiversity. For instance, molecular tools have contributed to understanding the ontogenetic pattern of melanin in the lateral region of cyprinid fishes and led the authors to synonymize 12 nominal species [26].
The modern integrative taxonomy was formally introduced in 2005 to increase accuracy and reliability to delimit and describe taxa by integrating information from different kind of data and methodologies, such as molecular markers, ecological and morphological characters [27], [28]. In the last years, the increase of articles using multiple lines of evidence corroborates the success of the modern integrative approach, indicating a clear renewal of the taxonomy [29].
DNA barcoding is among the molecular tools applied in biodiversity studies, identification of cryptic species, and taxonomy [ Vogler et al. [50] argue against the usefulness of DNA barcoding in the delimitation and identification of species/lineages based only on genetic distances, since it does not take into account the differences in the time of divergence between them. However, improved statistical methods were proposed to delimit species with DNA barcoding data, such as the General Mixed Yule Coalescent Method (GMYC) [51], which is based on an evolutionary model incorporating the topology of a tree, which allows researches to differentiate between the interspecific (''diversification'') and intraspecific (''coalescence'') processes of the branching processes of the lineages [51], [52].
Here, we tested whether 1) Microglanis cottoides may comprise cryptic species across coastal drainages given the high level of endemism of other fish species; 2) endemic species of Microglanis in the coastal drainages from São Paulo to Rio Grande do Sul are synonyms of M. cottoides; and 3) populations from the Uruguay River basin belong to M. cottoides.

Ethical statement
All specimens used were collected in accordance with Brazilian laws, and the sampling was approved by the Sistema de Autorização e Informação em Biodiversidade (SISBIO number 12120-1) of the Instituto Chico Mendes de Conservação da Biodiversidade (ICMBio). After collection, the animals were anesthetized and sacrificed using 1% eugenol as approved by the Universidade Estadual de Londrina/UEL Ethics Committee on the Use of Animals (CEUA; protocol 37917-11) and accepted by the National Council for the Control of Animal Experimentation and Federal Board of Veterinary Medicine. The animals were preserved in 92.
In addition, we included samples of Microglanis species occurring in sympatry with Microglanis cottoides or distributed in the Southern Atlantic and Uruguay basins: M. cibelae (n = 18), Microglanis sp. (n = 4), and M. malabarbai (n = 6). The species identifications were made based on the external morphology, according to Shibatta & Benine [53], Malabarba & Mahler-Jr. [9], and Bertaco & Cardoso [54]. When only the tissues were available, we followed the identifications of collectors. For M. garavelloi, COI sequences were obtained from GenBank (JN989032, GU701443, GU701629, GU701627, and GU701626), and included only species from the Upper Paraná River basin. New sequences generated in this study were submitted to MF045829 to MF045833; M. parahybae: KP0603106 to KP063110). Tissue samples of Microglanis sp. from Uruguay River basin were obtained by donation and the species identification were maintained as appointed by the collector, because they were from other localities than the type-locality of M. malabarbai. These specimens were included in the analysis to verify if they could belong to the M. malabarbai species. Although M. eurystoma is also distributed in the Uruguay River basin, this species was not included in the molecular analyses as it was not possible to obtain samples from the type locality. Microglanis parahybae occurs north of the distribution of M. cottoides in coastal drainages of the Paraíba do Sul basin. Thus, we included five samples of this species to test the phylogenetic relationship between species sharing a similar environment but occurring in allopatry. Our analyses included a total of 119 sequences. All species used, as well as information on them, are summarized in S1 Table. DNA extraction from tissues followed a phenol/chloroform protocol [55]. Partial sequences of the mitochondrial gene COI were amplified using polymerase chain reaction (PCR) with the primers FishF1-5'TCA ACC AAC CAC AAA GAC ATT GGC AC-3' and FishR1-5'AGA CTT CTG GGT GGC CAA AGA ATC A-3' [56]. We used 15 μl as a total volume containing 7.5 μl GoTaq Green Master Mix (Promega), 0.15 μl of each primer at 20 μM, 2.0 μl genomic DNA (5 ng/μl), and 5.2 μl of double-distilled water. The PCR consisted of an initial denaturation (5', at 95˚C) followed by 35 cycles of chain denaturation (30s at 94˚C), primer hybridization (54˚C), and nucleotide extension (30s at 72˚C). PCR products were checked by electrophoresis in agarose gel, purified using EXOSAP (Exonuclease I and Shrimp Alkaline Phosphatase GE Healthcare1, Piscataway, USA) and sequenced in both directions using Big Dye Terminator v 3.1 (Applied Biosystems), with subsequent reading on an automatic sequencer ABI Prism 3500 XL (Applied Biosystems).
Alignment, genetic distance, phylogenetic estimation, and species delimitation analyses. Electropherogram Quality Analysis software [57] was used to produce consensus sequences for each individual, from the sequences of both primers. Subsequently, these sequences were edited and aligned using the ClustalW algorithm in MEGA 6.0 [58]. This same software was used to calculate the genetic distance within and between populations using the Kimura-2-Parameters (K2P) model [59]. Traditionally, a threshold of 2% genetic distance has   [62] was used to estimate a Bayesian phylogenetic tree. Fossil data show that the divergence time between Lophiosilurus alexandri and Cephalosilurus apurensis, both in the Pseudopimelodidae family, is estimated between 15.9 and 11.5 million years (Ma) [63]. We time-calibrated the phylogeny using a node-age approach based on divergence times of these two genera, since Microglanis species have not been included previously in timecalibrated phylogenies. We applied a strict clock model with a uniform prior, which is a generally well-justified prior within a species or among a few closely related species [64]. We used as tree prior the Birth-Death model, which is an extension of the Yule model and assumes that at any point in time each lineage may undergo speciation or extinction [65]. The data were analyzed as a single partition and the evolutionary model used was GTR+Γ+I, as specified by the program MrModeltest 2.3 [66]. In BEAST, we ran the analyses for ten million generations, sampling every 1,000th step. The convergence was assessed in Tracer 1.5 with 25% burn-in [67], and we summarized the MCMC samples using the maximum clade credibility topology using TreeAnnotator v1.5.3 [68]. The distribution of haplotypes and mutational step numbers was generated using the software package Network 4.1.0.8 (www.fluxus-engineering.com), with the median-joining method (MJN) [69].
From the results of the phylogenetic analysis, the species delimitation test was performed using the Generalized Mixed Yule Coalescent (GMYC) model [51], which is more suitable for data with one gene. This analysis was conducted using the "Species Limits by Threshold Statistics" package [70], implemented in R Core Team v3.0.1 [71]. By means of this package it is possible to calculate the number of clusters by classifying the bifurcation rates of a phylogeny as a result of interspecific or intra-specific branching processes [51].

Morphological analyses
A total of 21 linear measurements was obtained using a digital caliper according to Malabarba and Mahler-Jr [9], and Bertaco and Cardoso [54], with addition of the following variables: snout length (measured from the tip of the snout to the base of the anterior eye margin), pelvic fin length, posterior cleithral process length (measured from the origin of the elevated pectoral spine to its osseous tip), dorsal-fin to adipose-fin distance (measured between the posterior base of the dorsal fin and the anterior base of the adipose fin), anus to anterior anal fin base distance, posterior nostrils distance, anterior to posterior nostrils (measured between the posterior base of anterior nostril to anterior base of posterior nostril).
Specimens identified as Microglanis cottoides, M. cibelae and M. malabarbai, and belonging to the type-locality basins were included in the analysis in order to test the morphological similarity and its congruence with molecular analysis. Variables of body were presented as proportions of standard length (SL) and variables of head were presented as proportions of head length (HL). Differences in average body proportions were tested by One-way ANOVA. Principal components analysis (PCA) on covariance matrix of log transformed data was used to test for morphometrics differences between species, and obtain discriminant characters. Counts of lateral line pores were made on the left side of body, whenever possible. The box plot of counts was built with PAST [72], as well as the ANOVA and PCA were performed with this program.

Molecular analyses
Barcode and genetic distance. A total of 620 base pairs (bp) of the COI gene were analyzed. We did not find insertions, deletions, or stop codons in these sequences, indicating that all amplified regions correspond to a functional portion of the COI gene.
Genetic distance within the populations considered as Microglanis cottoides (RIB, GUA, PAR, ITA, MAD, ARA, PAT, and URU) ranges from 0% to 7%. Genetic distances within coastal drainage populations (RIB, GUA, PAR, ITA, MAD, ARA, PAT) ranges from 0 to 0.4% (Table 1). In contrast, the population of the Uruguay River basin presented the greatest genetic distance relative to the coastal drainages, varying from 5.7% to 7.3%. These distances were even greater than those found between M. cottoides (coastal drainage) and M. parahybae (from 3.7% to 4.2%). The latter species occurs north of the distribution of M. cottoides, in the Paraíba do Sul River basin. Populations of M. cottoides of the Uruguay River basin occur in sympatry with M. malabarbai and the genetic distance between them was 0.7%. Low genetic divergences were also observed between M. cibelae and M. cottoides (0.5% to 1.8%).
Haplotype network. We found 18 haplotypes among samples of Microglanis in Southern Brazil ( Fig 3A). Microglanis cottoides from the Uruguay River basin grouped with Microglanis sp., and M. malabarbai, whereas M. cottoides from the coastal drainages grouped with M.  cibelae. The network shows the presence of two distinct groups within M. cottoides, separated by 30 mutational steps (i.e. Uruguay River basin and coastal basins). The number of mutational steps separating individuals from these two regions was greater than that found between M. parahybae and M. cottoides (21 mutational steps). However, the number of mutational steps found between M. malabarbai and populations of M. cottoides from the Uruguay River basin (one to two mutational steps) was similar to the number of mutations between M. cottoides populations of the coastal basins ( Fig 3A).

Morphological analyses
The samples examined showed similar SL mean (Tab. 3; ANOVA: F = 0.3736, p = 0.6904). Regarding body proportions, the minimum and maximum values of all variables overlapped, preventing the prompt identification of species (Table 2) The Principal components analysis of the combined samples of Microglanis cottoides, M. cibelae and M. malabarbai showed better separation of the groups with components 2 and 3 (Fig 4). The first component, which retained most of the variation (89.4%), was the representative of size. The second and third components retained 3.8% and 1.0%, respectively, and best represented the shape. In PC2, Microglanis malabarbai differs from M. cibelae, for presenting longer maxillary barbel length, anterior internarial distance, prepelvic length, and mouth width (positive variable loadings, Table 3), and for smaller body width, dorsal-fin spine length, preanal length, posterior cleithral process length, and caudal peduncle length (negative variables loadings, Table 3). In this axis, M. cottoides is in intermediary position, not separating from M. cibelae, and evidencing greater similarity with M. cibelae than with M. malabarbai.

Molecular approach
Microglanis cottoides was originally described from the Camaquã River, a tributary of the Laguna dos Patos basin, but its distribution was extended to the Uruguay River basin [9], and the coastal rivers of São Paulo to Rio Grande do Sul States [23]. However, the molecular evidence of our study strongly suggests that M. cottoides, as currently defined, does not form a monophyletic group, since the population occurring in the Uruguay River basin is genetically distinct, and M. cibelae is nested within M. cottoides. The genetic distance between populations of Microglanis cottoides from the coastal basins and the Uruguay River basin (5.6% to 7.0%) corroborate the hypothesis that the Uruguay River population of M. cottoides has been misidentified. For instance, divergences of this magnitude are often found among congeners of other Neotropical fish [37], [39], [41]. In addition, our results also demonstrate no overlap between the intra and interspecific distances  [78], [79].
The haplotype network corroborates the results observed both in phylogenetic analysis and genetic distances. According to Hudson et al. [80], the haplotypes represent the nodes in a group of closely related taxa, linked to each other according to the similarity of the haplotype Reappraisal of the systematics of Microglanis cottoides sequences. Sequence divergence in Microglanis were lower within species than between closely related species, except for M. cibelae. Also, in the network is possible to visualize two very distinct haplogroups within M. cottoides (Uruguay River, corresponding to M. malabarbai, vs. coastal basins), separated by a high number of mutational steps, equivalent to those found among the other species of Microglanis added in the analysis. In fact, it is larger than that found between the haplogroups of M. cottoides from the coastal basin and M. parahybae from the Paraíba do Sul River basin.
The GMYC allows the calculation of the number of clusters resulting from interspecific or intraspecific branching processes [51] and this analysis indicated two clusters for Microglanis cottoides, containing the same clades recovered in the phylogenetic analysis, which reinforces that the species known as M. cibelae is a junior synonym of M. cottoides and that individuals from the Uruguay River basin, identified as M. cottoides, are actually M. malabarbai. In addition, this tool has a tendency to overestimate the number of lineages [51] in other words, to identify intraspecific lineages as distinct lineages, which did not occur with M. cibelae, reinforcing the hypothesis of synonymy.
The phylogenetic analysis showed that the coastal species Microglanis cottoides and M. parahybae are closely related and have diverged about 4.2 Ma; Microglanis garavelloi from the Upper Paraná River basin and M. malabarbai from the Uruguay River diverged also in the Pliocene and have been separated from the coastal species since the late Miocene (± 7.5 Ma). These results agree with Weitzman & Weitzman [81] and Lundberg [82], who suggest that diversification of fish up to the level of modern species predates the Pleistocene, as Albert & Reis [83] state that most Neotropical ichthyofauna diverged between 3.0 and 10 Ma.
The origin of many southern and southeastern Brazilian drainages that now run directly into the Atlantic Ocean was at the time of tectonic activity at the beginning of Tertiary (65-1.8 Ma), generating a complex system of 'failures' in the crystalline shield that culminated in erosions and subsequent events of headwater captures, as occurred between the Tietê River, of the Upper Paraná River basin, and the Paraíba do Sul River [1]. These events were also responsible for the distribution of some species into neighboring drainages, such as Ribeira de Iguape, as well as smaller coastal drainages [84], [85], [86], [87], [88] The differentiation of species by geographic isolation is the reflection of the geological past and the environmental changes occurring in the region [1], [89], and the degree of genetic divergence is strongly correlated with the age of physical isolation [90]. Thus, the tectonic activity that began in the Tertiary, pobably allowed the species of Microglanis of Southeast-South coast first diverge of the species of the Upper Paraná River system (± 7.5 Ma). Later, in a more recent episode (± 4.2 Ma) that likely involved the Serra do Mar, the diversification of coastal species occurred, separating M. cottoides (from São Paulo to Rio Grande do Sul States), from M. Parahybae (from Rio de Janeiro).
Species of fishes with similar distribution of Microglanis cottoides have already been observed south of Rio de Janeiro. For example, Cyphocharax santacatarinae (Fernández-Yépez 1948) occurs from Santa Catarina to São Paulo, in the Ribeira de Iguape River [91] and Oligosarcus hepsetus (Cuvier 1829) occurs in coastal rivers from Santa Catarina to Rio de Janeiro [92]. The comparison of inter and intraspecific patterns from co-distributed species allows the evaluation of how independent lineages have responded to the same historical processes in a given region [93], [94]. In this way, it is possible to suggest that this region served as an important ichthyofauna divisor of the rivers located to the north and to the south.
At the time of the separation between the coastal species and those of the adjacent plateau, the Uruguay River maintained a connection with what is now called the Upper Parana River and, according to Beurlen (apud [95]), only came to separate in the Miocene (24-5.3 Ma) on occasion of erosion in its basaltic cover. The close phylogenetic relationship between Microglanis garavelloi of the Upper Paraná River and M. malabarbai of the Uruguay River corroborates this fact.
The shared distribution of species among isolated basins may reflect recent vicariant events, such as changes in the drainage course or capture of rivers from one basin by the other, generated by geomorphological modifications [1]. Alternatively, in the case of Microglanis cottoides, dispersion between coastal rivers, due to the sea level fluctuations in the late Pleistocene, mainly by marine regression in the glacial periods, that allowed communication between drainages along the continental shelf [96] [97]. In the Brazilian continental shelf, recent studies have shown that sea level fluctuations during the Quaternary have left evidence, such as paleodrainages, that confirm successive stages of exposure and submersion [2], [98], [99], [100], [101]. Thus, geomorphological features associated with fluvial channels on continental shelves may indicate that coastal drainage, which now flows directly into the oceans, had communication with other nearby drainage basins [101], allowing ichthyofauna flow. This is one of the plausible explanations for the wide distribution of many species, including M. cottoides, in coastal drainages that are now isolated.

Morphological approach
Malabarba & Mahler-Jr. [9] diagnosed Microglanis cibelae for having a smaller head than M. cottoides (25.6-31.1% SL vs. 29.3-33.8%), and lower body width (25.4-29.8% SL vs. 28.5-33.9%). As can be observed there were overlaps in the proportions, which does not allow the species to be safely identified. In our analyses we also did not find intervals of these characteristics that could differentiate the samples as distinct species and, therefore, we used the PCA, because it was considered more suitable to measure the variation of the morphometric variables [102]. Nevertheless, there was a wide overlap between M. cottoides and M. cibelae on the second principal component axis. Malabarba & Mahler-Jr. ( [9], p.253) still presented a linear regression graph pointing to ontogenetic divergence in head size, which is relatively higher in M. cottoides than in M. cibelae in specimens larger than 30 mm SL. The same was observed in our study, corroborating those authors [9]. However, these variations in morphological characteristics may represent different populations, not different species. Malabarba & Mahler-Jr. [9] also include samples from the Uruguay River basin identified as M. cottoides. By our analysis, all the samples identified as M. cottoides of this basin correspond to the species M. malabarbai, whose differences are very large in relation to M. cibelae, and which may have influenced the decision of those authors.
Initially, the possibility of a cryptic species in sympatry with Microglanis malabarbai was considered, but was not confirmed, as the genetic divergence is too low to separate them into distinct species. However, a more complex situation seems to be involved, considering the two reviews of Microglanis from southern Brazil [9], [54]. Both studies, analyzing only the morphology, were not able to distinguish some populations of Microglanis from the Uruguay and Laguna dos Patos basins, even with high genetic divergence among them observed herein, raising the hypothesis that the morphological variation in populations of M. cottoides and M. malabarbai make them difficult to identify correctly. In addition, both papers describe differences in head length, body width, pectoral spine length, maxillary barbel length, and internareal distance to diagnose Based on our results, we propose that Microglanis cibelae is a junior synonym of M. cottoides and that individuals from the Uruguay River basin, previously identified as M. cottoides, are actually M. malabarbai, refuting the initial hypothesis of this work that there are cryptic species of Microglanis in the studied coastal drainages, and corroborating that only M. cottoides occurs in the coastal drainages from São Paulo to Rio Grande do Sul states and does not occur in the Uruguay River basin. Regarding the phylogenetic hypotheses proposed in this contribution, we are planning to test them with new analyses that include more genes.