Unrecognized species diversity and endemism in the cichlid genus Bujurquina (Teleostei: Cichlidae) together with a molecular phylogeny document large-scale transformation of the western Amazonian river network and reveal complex paleogeography of the Ecuadorian Amazon

The Amazon basin and especially the western Amazon with the adjacent Tropical Andes are the most biodiverse areas of our planet (Reis et al., 2003; Hubbell et al. 2008; Hoorn & Wesselingh, 2010a, b; Albert & Reis, 2011). The basin is also an important provider of global ecosystem services, providing 20% of the world’s oxygen supply and responsible for 10% of the net primary productivity of the whole terrestrial biosphere (http://earthobservatory.nasa.gov; Subramaniam, 2008). Furthermore, the basin is a significant source of nutrients for oceanic life that sequesters globally relevant amounts of carbon (Subramaniam, 2008).

The Amazon basin hosts the most diverse and extreme terrestrial and aquatic ecosystems on Earth. Paramount among its megadiverse groups are the freshwater fishes with the Amazon hosting the largest community of freshwater fishes in the world. The whole Neotropical region with its diversity centre in the Amazon hosts more than 5,600 valid species of freshwater fishes, with current estimates far exceeding 7,000 species, and these represent over 10% of all living vertebrate species (Reis et al., 2003; Albert & Reis, 2011; Van der Sleen & Albert, 2018).

Our knowledge of the Amazonian species megadiversity is still incomplete, in part due to insufficient sampling across the vast complex region, especially among the freshwater fishes (Hoorn & Wesselingh, 2010a, b; Van der Sleen & Albert, 2018). Even more limited is our understanding of the evolutionary mechanisms that have generated this diversity (Hoorn & Wesselingh, 2010a, b). Of paramount importance for Amazonian faunal evolution appears to have been the geological processes which have shaped the Andes and the river network of western Amazonia. Collectively, these appear to be major determinants of Amazonian diversity, both for the aquatic and terrestrial faunas (Hoorn & Wesselingh, 2010a, b; Albert & Reis, 2011; Albert et al., 2018). The originally south–north-oriented drainage of the paleo-Amazon changed with progressive Andean elevation. This led to an eastward propagation, and sedimentation, that resulted in the present west–east-oriented river network. This large-scale reorientation of the Amazon river network may have heavily reshaped the biodiversity patterns and promoted further diversification, both in the foothills as well as in the lowlands. The present mosaic of muddy, clear, and black-water rivers was established from a previous still little known river network, that was itself preceded in the western Amazon by a mega-wetland, the Pebas formation (Hoorn & Wesselingh, 2010a, b; Wesselingh & Hoorn, 2011).

The general outlines of the large-scale transformation of the western Amazon river system described above were derived predominantly from geological and paleontological data. The input from extant faunal data and phylogenies has been very limited, with no well-sampled studies performed, and none with groups endemic to the area of the western Amazon that are highly philopatric, with small-scale endemism and a large species diversity, i.e. groups that would have the largest explanatory power.

Freshwater fishes are among the best model species groups for use in reconstructing landscape evolution, given they can be highly dispersal limited (Wallace, 1876; Darlington, 1957). One group of Amazonian freshwater fishes exhibiting high degrees of endemism, philopatry and strict habitat preference is the cichlid fishes, Cichlidae, which is one of the most important groups of Neotropical fishes. In terms of species richness, cichlids are the dominant group of larger sized fishes in Middle America (Myers, 1966; Bussing, 1976, 1985; Říčan et al., 2013, 2016) and the third most species-rich family of fishes in South America (Van der Sleen & Albert, 2018).

Among the cichlid fishes, the best model group for the study of Andean and western Amazonian landscape evolution, and its direct influence on the formation of biodiversity, is Bujurquina Kullander, 1986. This genus has its centre of geographical distribution and highest diversity in the western Amazon, and its distribution and endemism follow the putative course of the paleo-Amazon (Kullander, 1986; Říčan, 2017). Bujurquina thus has a distribution unlike any other large genera of South American cichlids, which generally have highest diversity and endemism in the older eastern Amazon in the Brazilian and Guiana shield areas (www.cichlidae.com; Stawikowski & Werner, 1998, 2004; Kullander et al., 2018; Van der Sleen & Albert, 2018; Říčan et al., 2021).

Bujurquina follows an arc from the south in Argentina to Paraguay, Brazil and Bolivia. It then has its highest currently known diversity and endemism in Peru, and continues to the north through Ecuador, Colombia and Venezuela (Kullander, 1986). Bujurquina species are highly philopatric, show small-scale endemism and strong adaptation lo local environmental conditions (Říčan, 2017). Molecular-based dating of Bujurquina diversification (Musilová et al., 2015) supports the hypothesis that most of its evolution occurred during the period of major Amazonian river network reorganization described above.

Due to its time frame of evolution, geographical distribution, species diversity and narrow endemism, Bujurquina is the best candidate group among the cichlid fishes, and most likely among all fish groups, to provide information on general patterns and processes in the evolution of western Amazonia. Understanding diversity patterns and evolution of Bujurquina will not only strengthen our understanding of the processes of cichlid and fish diversification in this area, but will also provide a groundwork for the whole faunal and biogeographical evolution of the western Amazon.

The biogeography and endemism of Bujurquina play out in two biogeographic zones. There appears to be a separate band of endemism in Bujurquina in the Andean foothills, which has a distinct fauna from the supposedly more widely distributed species in the Amazonian lowlands (Kullander, 1986). The two zones appear to host different faunas of Bujurquina since the lowland species are usually not known from the adjacent foothills and vice versa (Kullander, 1986). Foothill areas of the same river basins thus appear to have a different Bujurquina fauna from lowland areas of the same basins, at least in Peru where Kullander (1986) proposed this hypothesis (Říčan, 2017).

Bujurquina currently includes 18 valid species, but this species diversity is clearly significantly underestimated as will be outlined and demonstrated in this study. If we generalize the biogeographic pattern to all described species of Bujurquina (i.e. not only those found in Peru), the foothill band of narrowly distributed endemic species includes B. oenolaemus Kullander, 1987 (Paraguay basin, Bolivia), B. eurhinus Kullander, 1986 (Madeira basin, Peru), B. robusta Kullander, 1986, B. labiosa Kullander, 1986, B. megalospilus Kullander, 1986, B. apoparuana Kullander, 1986, B. hophrys Kullander, 1986 (Ucayali basin, Peru), B. huallagae Kullander, 1986, B. ortegai Kullander, 1986, B. zamorensis (Regan, 1905), B. pardus Arbour et al., 2014 (Marañón basin, Peru and Ecuador) and the northernmost species B. mariae (Eigenmann, 1923) (Orinoco basin in Colombia and Venezuela). The lowland group of supposedly more widely distributed species includes B. vittata (Heckel, 1840) (Paraguay basin, Argentina, Paraguay, Brazil and Bolivia), B. cordemadi Kullander, 1986, B. tambopatae Kullander, 1986 (both with types in Madeira basin, Peru, ranging possibly into Bolivia and Brazil), B. moriorum Kullander, 1986 (type locality in lower Ucayali basin, Peru), B. peregrinabunda Kullander, 1986, B. syspilus (Cope, 1872) (both with type localities in the Amazon proper in Peru, ranging into Colombia and westernmost Brazil).

In this study, we use a newly obtained collection of Bujurquina specimens which, for the first time, provides a representative sampling throughout the distribution area of the genus in the western Amazon and which reveals numerous new species and a strong localized endemism. The majority of the collected material and putatively new species come from foothill areas, while the expansive western Amazonian lowland areas still remain sparsely sampled. Even from here, however, putative new species are known (Říčan, unpublished data; Vriesendorp et al., 2007; Gilmore et al., 2010; Pitman et al., 2013) and one is included in the present study. We use morphological and molecular characters, and molecular phylogenetic and biogeographic analyses based on these specimens to ascertain the species diversity of the genus and to study the phylogeography and historical biogeography of the genus. Our results reveal that the species diversity of the genus in general, and particularly in Ecuador, is severely underestimated. We further reconstruct the biogeographic history of the genus in the western Amazon and in Ecuador in particular. From the detailed sampling in Ecuador, we derive a detailed biogeographic history of the genus and the paleogeography of the rivers in the Ecuadorian Amazon which we compare with the geological composition and known geological history of the Ecuadorian Amazon.

Species determination

We combine morphological species determination with post hoc species delimitation using the molecular mtDNA cytochrome b (cytb) marker. Specimens were identified to species with the use of original descriptions, identification keys and comparative material of all presently recognized Bujurquina species. The majority of specimens of Bujurquina in Ecuador are based on our analyses identified as putatively new species. These are recognized for the first time in the present study predominantly based on a number of colour/pigmentation traits (Figs. 1, 3, S3) and are reported as “Bujurquina sp.” (possible new/unidentified species). The candidate new species are here provisionally named using their postulated areas of endemism, which are given by river basins (or their portions) within the context of the Ecuadorian territory. Representative specimens of most studied species-level taxa of Ecuadorian Bujurquina are shown in Fig. 1 and the remaining described species are shown in Fig. S1.

Representative specimens of the majority of studied species-level taxa of Bujurquina in Ecuador. The coloration pattern of the spinous portion of the dorsal fin together with the coloration of the head (chiefly the development of the suborbital stripe-blotch and the pattern of the opalescent markings) provides a unique diagnosis for each species (see Text). Note that each species is diagnosed by a unique combinations of coloration characters as shown by the pictographs (the upper portion shows coloration pattern in the spinous portion of the dorsal fin, the lower left portion of the pattern of the suborbital stripe in juveniles, the lower right in adults, and the letters C and D show the course of the midlateral stripe on the body, to the caudal peduncle, or to the posterior insertion of the dorsal fin, respectively). Scale bar in the inset phylogeny shows 4 My

Full size image

Specimen sampling

We have collected a representative sampling of specimens and tissue samples covering the whole distribution of Bujurquina with particular focus on the Ecuadorian Amazon, especially in the more diverse and accessible foothills (Fig. 2). The lowlands of the Ecuadorian Amazon as well as in other western Amazonian countries remain undersampled (Figs. 2, 3). Our molecular sampling in this study includes 13 of the 18 described species (Figs. 1, 2, 3, 4, 5, S2, Table S1), with species found in and around the territory of Ecuador sampled from their type localities (B. zamorensis, B. moriorum, B. peregrinabunda, B syspilus). In addition to the described species, our sampling includes twenty putative undescribed species of Bujurquina, i.e. 33 species-level taxa in total and thus more than double the currently recognized species diversity.

Distribution of identified Bujurquina species in Ecuador. Dots show localities, solid lines between dots show phylogeography within species based on Figs. 4, 5, S2 (phylogeography not available for all localities since only a representative sample of each putative species was sequenced). Colours and large same coloured numbers identify each species (species are in text identified by these numbers and also by names derived from distribution). Small coloured numbers below species numbers show basal divergence within species in Ma (see Figs. 4, 5). Same sized black numbers show dating in Ma of main geographic divides based on biogeographic analysis in Fig. 5. Rivers are identified in white letters. The three main basins of Ecuadorian western Amazon (Curaray-Napo-Putumayo plus extralimital Caquetá; Santiago-Morona-Pastaza which are tributaries of the Marañón; and the Tigre basin) are identified and their drainage divides are shown by dotted lines (these principal drainage divides are also identified in the biogeographic analysis in Fig. 5). White interrupted arrows identify the inferred paleodrainage flowing towards the north based on geological information (see Fig. 6) which finds support in the biogeographic reconstruction in Fig. 5. Yellow circle identifies the Sangay volcano, that is potentially important in the formation of the modern Pastaza river and species found in the Pastaza basin (see text)

Full size image

Distribution and the two main groups of described species of Bujurquina (plus the long known B. sp. Bolivia). In red are main geographic drainage basin divides of the western Amazon. All species were included in morphological analyses, species with an asterisk were also available for molecular analyses. Note that colour codes of species in this figure do not correspond to colour codes used in other maps and in other figures

Full size image

Phylogenetic relationships and species delimitation within Bujurquina based on the cytb marker from BEAST analysis. Outgroup taxa are omitted from the figure (see Fig. S2). Species identified by colour and numbers are those occurring in Ecuador. Numbers at nodes show age in Ma, grey bars at nodes show confidence intervals of age estimates (95% HPD). Bayesian posterior probabilities above 0.95 are shown as black points on nodes. Species delimitations using the bPTP and GMYC methods are shown as coloured boxes. Note that species delimitations using the molecular marker delimit 25 species in total (8 in Ecuador), morphology delimits 32 species (12 in Ecuador), and molecular phylogeny delimits 30 species (11 in Ecuador). Scale bar shows 5 My

Full size image

Phylogenetic relationships and biogeographic reconstruction within Bujurquina based on the cytb marker from BEAST analysis. Outgroup taxa are omitted from the figure (see Fig. S2). The two main phylogenetic groups of Bujurquina are identified. Species identified by colour and numbers are those occurring in Ecuador. Numbers at nodes show age in Ma, grey bars at nodes show confidence intervals of age estimates (95% HPD). Bayesian posterior probabilities above 0.95 are shown as black points on nodes

Full size image

The fish specimens and tissue samples were collected under permanent permits to Cecilia Elizabeth Rodriguez Haro (IKIAM, Tena, and UEA, Puyo, Ecuador), Javier Maldonado Ocampo (Universidad Javeriana, Bogota, Colombia), and Hernán Ortega (Universidad Nacional Mayor de San Marcos, Lima, Perú). Specimens were obtained by fishing using gill nets, cast nets, hook and line, and spear fishing in all available microhabitats at all localities (following the recommendations of Říčan, 2017). Localities were chosen specifically so that they capture the full regional ecological habitat variation of the studied area. At localities with sufficient visibility (all upland and rainforest streams and some lowland locations), direct observation of the fish fauna was carried out using mask and snorkel.

Collected specimens (Table S1) were individually and permanently tagged and a DNA tissue sample taken from their right pelvic fin. DNA tissue samples were preserved in 1.5-ml Eppendorf Safe lock tubes in absolute ethanol and deposited at the fish collection of the Department of Zoology, Faculty of Science, University of South Bohemia in České Budějovice, Czech Republic. Specimens were preserved in 10% formalin and later transferred to 70% ethanol. Specimens are currently deposited in the fish collection of the Department of Zoology, Faculty of Science, University of South Bohemia in České Budějovice, Czech Republic. Upon formal description of the species will be transferred to publically available museum collections for permanent storage and study. This study is based on 446 specimens.

Molecular marker

For molecular phylogenetic analyses and delimitation of putative species, we have used the mitochondrial (mtDNA) cytochrome b (cytb) gene. The cytb gene is commonly used as a molecular marker in cichlid phylogenetic studies, and is able to provide well-resolved phylogenies of the Neotropical cichlids that have subsequently been confirmed by later studies, including those using phylogenomic approaches (e.g. Říčan et al., 2016; Ilves et al., 2017; Burress et al., 2017, 2018). The cytb gene has recently been used successfully in similar studies focussing on diversity of other South American cichlid genera (e.g. Piálek et al., 2019a, b; Říčan et al., 2019, 2021).

For the cytb marker, we have sequenced 190 specimens from the 33 ingroup species and species-level taxa. We also include sequences from our previous studies (Musilová et al., 2008, 2009), as well as other sequences obtained from GenBank (Table S1). The total number of ingroup sequences is 201 plus 50 outgroup sequences including one heroine genus (Hypselecara) as the primary outgroup and five cichlasomatine genera (same tribe as Bujurquina), including the two successive sister genera of Bujurquina, namely Tahuantinsuyoa and Andinoacara (Table S1).

Laboratory methods

Genomic DNA was extracted from ethanol-preserved fin tissue using the JETQUICK Tissue DNA Spin Kit (Genomed, Germany) following standard protocol. The primers and reaction conditions of polymerase chain reaction (PCR) amplification are as in Říčan et al. (2008). The products were analysed in an ABI 3730XL automated sequencer (Applied Biosystems; Macrogen Inc., Korea). Contiguous sequences of the gene segments were created by assembling DNA strands (forward and reverse) using GENEIOUS v. 11.0.2 (http://geneious.com, Kearse et al., 2012). Nucleotide coding sequences were also translated into protein sequences to check for possible stop codons or other ambiguities. All newly generated sequences were submitted to GenBank under Accession numbers OP436067–OP436260 (Table S1). Sequences were aligned using MUSCLE v. 3.8 (Edgar, 2004), using the default settings.

Phylogenetic methods

Maximum parsimony (MP) analysis in PAUP* 4b.10 (Swofford, 2003) and Bayesian inference analyses (BI) in MrBayes v. 3.1.2 (Huelsenbeck & Ronquist, 2001; Ronquist & Huelsenbeck, 2003) and BEAST v.1.8.4 (Drummond & Rambaut, 2007) were used for phylogenetic inference. The MP phylogenetic analyses in PAUP* were run with 500 random sequence additions, 10 trees kept per addition, and a hs (heuristic) search on the saved trees to find all the shortest trees. Bootstrap analyses were done using the same approach, with 5 random sequence additions per one replication. Bootstrap analyses were run with 1,000 replications. The BI analyses in MrBayes (and BEAST, see below) were run with haplotype dataset and were partitioning into codon positions (1st + 2nd vs. 3rd). The haplotype dataset was created in Fabox (Villesen, 2007). An optimal model of evolution according to Akaike criterion was selected using MrModeltest 2.2 (Nylander, 2004) and PAUP* v. 4.0b10 (Swofford, 2003). The BI analysis using the Markov chain Monte Carlo (MCMC) simulation was run for 2 million generations with trees sampled and saved every 1,000 generations. Two independent analyses, each comprising two runs with eight chains, were performed to compare results of independent analyses. The analyses were run at the freely available Cipres server (https://www.phylo.org/). The first 10% of trees from each run were discarded as burnin. Convergence of the runs was estimated with the use of graphical visualization and diagnostics (especially the effective sample size; ESS) in Tracer v. 1.8.4 (Rambaut & Drummond, 2007). The remaining trees were used for reconstruction of the 50% majority-rule consensus tree with posterior probability (PP) values of the branches. Trees were rooted with a representative sampling of heroine and cichlasomatine cichlids including the most closely related genera Tahuantinsuyoa and Andinoacara.

Molecular-clock dating analyses in BEAST

For divergence time estimation, we used the Bayesian Evolutionary Analysis by Sampling Trees (BEAST) software package version v.1.8.4 (Drummond & Rambaut, 2007) with parameters as in MrBayes analyses (except for 30 million generations). We used the relaxed molecular-clock model with lognormal distribution of rates and for tree prior the coalescent model with constant size. The calibration of the molecular clock was performed using secondary calibration from the study of Musilová et al. (2015) which was focussed on the dating of Neotropical cichlids in general and the Cichlasomatine cichlids in particular, with focus of Andinoacara, the sister group of the Bujurquina-Tahuantinsuyoa clade. Musilová et al. (2015) employed a calibration using a set of Neotropical fossil cichlid species and is so far the only study that includes Bujurquina and its most closely related genera in a dated phylogeny. For the calibration, we have used three nodes from the study (Fig. S2), (1) the Andinoacara basal node, estimated by Musilová et al. (2015) at 19.33 Ma (95% HPD 15.03–24.25 Ma), (2) the Bujurquina-Tahuantinsuyoa-Andinoacara node, estimated at 31.51 Ma (95% HPD 26.22–37.17 Ma) and (3) the Bujurquina-Tahuantinsuyoa-Andinoacara vs. rest of cichlasomatine cichlids node, estimated at 46.47 Ma (95% HPD 40.96–52.62 MA).

The analyses were also run on the Cipres server (https://www.phylo.org/). Runs were checked for convergence with Tracer v.1.8.4 (Rambaut & Drummond, 2007). Four well converged runs were combined in LogCombiner v.1.8.4 with a burnin of 10% for each of the data partition schemes. The final tree for each data partition scheme was produced from these data with TreeAnnotator v.1.8.4.

Species delimitation analyses using GMYC and PTP

We employed the General Mixed Yule Coalescent (GMYC) and Poisson tree processes (PTP) analyses for molecular species delimitation using the cytb marker. Both methods were designed for delimiting species based primarily on single molecular markers (hence where multilocus coalescent-based methods are not applicable).

The General Mixed Yule Coalescent (GMYC) model (Pons et al., 2006; Fujisawa & Barraclough, 2013) is frequently used in empirical studies (Fontaneto et al., 2007; Monaghan et al., 2009; Carstens & Dewey, 2010; Vuataz et al., 2011; Powell, 2012) and the newer Poisson tree processes (PTP) model (Zhang et al., 2013) has been shown to even outperform the GYMC method where distances between species are small. Both methods outperform OTU-picking methods (relying on simple sequence similarity thresholds) and are more robust to cases where the barcoding gap is absent (Zhang et al., 2013). Both methods model speciation (among-species branching events) via a pure birth process and within-species branching events as neutral coalescent processes. The methods identify the transition points between inter- and intra-species branching rates by maximizing the likelihood score of the model. While the GMYC method uses time to identify branching rate transition points (hence only on a time-calibrated ultrametric tree), the PTP method directly uses the number of substitutions and does not require a time-calibrated ultrametric tree. Both methods assume that all lineages leading from the root to the transition points are different species.

PTP and GMYC analyses were run on the web interface (http://species.h-its.org/). Prior to analyses, a haplotype dataset was created in Fabox (Villesen, 2007). The bPTP analysis was run both on the MrBayes tree and on the ultrametric BEAST tree in order to compare results (see Zhang et al., 2013) and in the latter case to make the bPTP result directly comparable to the GMYC result.

Biogeographic analyses

Reconstruction of ancestral areas for all nodes in the phylogenetic tree was carried out using the event-based Bayesian statistical dispersal-vicariance analysis (S-DIVA; implemented in RASP 2.0; Yu et al., 2015). Biogeographic units used for analyses were present river basins, when a species distribution was not restricted to single river basin, the species range was coded as a sum of the river basins. The analyses were carried out using a number of different ‘maxareas’ options in S-DIVA up to the maximum number of areas in the analysis. When results are the same for all ‘maxareas’ analyses, then just this one result is reported.

Species diversity of Bujurquina in the Ecuadorian Amazon based on morphological delimitation

Our determination of species based on morphological data suggests fifteen species of Bujurquina in the Ecuadorian Amazon (Fig. 1), including the two species with type localities in Ecuador, B. zamorensis and B. pardus. All Ecuadorian Bujurquina species share two putative synapomorphies, a patterned dorsal fin and a blue lower lip. Both character states are however not exclusive to Ecuadorian Bujurquina. They are shared with all known Colombian species, including one that is described, B. mariae (Eigenmann, 1923) and several that are not described (OŘ, unpublished data; two are included in the present work). These states are also shared with a subset of Peruvian species. The Peruvian species are (1) those found in the Marañón basin including Huallaga, including two described species, B. huallagae and B. ortegai, plus several undescribed species (OŘ, unpublished data; three of which are included in the present study), and (2) in the Amazon basin per se, i.e. downstream from the confluence of the Marañón and Ucayali, including B. moriorum, the type species of the genus, B. peregrinabunda, B. syspilus and several undescribed species (OŘ, unpublished data, one of which is included in the present study). The remaining species found in the Ucayali basin (B. apoparuana, B. hophrys, B. labiosa, B. robusta, B. megalospilus) and all species of southern Bujurquina (B. eurhinus, B. cordemadi, B. tambopatae, B. oenolaemus, B. vittata, and at least one undescribed species that is included in the present study) have an immaculate (unpatterned) dorsal fin and do not have a blue lower lip.

A short characterization featuring the main and most readily observed identification characters of all included described and putative species in this study from the territory of Ecuador follows below and we also include a pictogram featuring the most important subset of these main distinguishing characters for each species alongside the overview molecular phylogeny in Figs. 1 and S3.

Throughout the text, we will refer to the species by their provisory names that are derived from their so far known distribution areas and/or by their sp. numbers (which are found on the distribution maps) in situations where it is more convenient to do so.

Bujurquina zamorensis (sp. no. 1) is diagnosed by a combination of (1) dorsal fin with circular spots, (2) fully ocellated cheek blotch, (3) violet snout (in other species less violet and more maroon), (4) fine and numerous opalescent spots and complete absence of lines on cheek and operculum. Bujurquina zamorensis was described from the Rio Zamora (without precise type locality), a tributary of the Upper Santiago and based on our analyses also includes specimens from the Upper Pastaza and portion of the Upper Napo.

Bujurquina pardus (sp. 11) as judged from a single pictured preserved specimen is diagnosed by a combination of (1) dorsal fin with circular spots, (2) white lappets on dorsal fin spines, (3) distinctly blotched flank scales in the zones of the vertical body stripes, (4) continuation of vertical bars onto the dorsal fin and (5) a fully developed suborbital stripe that is distinctly anteriorly bend and starts behind the eye. Development of opalescent cheek coloration is not known. Character 3 is however not unique to B. pardus contrary to Arbour et al. (2014) because it is present in most Ecuadorian putative species and we are thus not using it in the overview of species below where it is almost uninformative.

Bujurquina sp. Santiago (sp. 14) is most similar to B. zamorensis, being however distinct from it in (1) not having a fully ocellated suborbital blotch (mostly not ocellated at all), (2) in having an almost immaculate dorsal fin with only slight traces of the spotting and (3) in having in juveniles a bent posteriorly inclined suborbital stripe vs. a straight posteriorly inclined suborbital stripe. Bujurquina sp. Santiago is so far only known from one tributary of the lowermost Ecuadorian portion of the Santiago basin.

Bujurquina sp. Middle Pastaza-Morona (sp. 12) is most similar to B. pardus in sharing (1) a continuation of vertical bars onto the dorsal fin, and 2) having a fully developed suborbital stripe. It is however distinct from B. pardus in (1) lacking such distinct spots on flank scales, (2) in having red lappets on dorsal fin spines, (3) in having the suborbital stripe vertical vs. anteriorly bend, and starting below the centre of the eye vs. behind the eye and (4) has a dorsal fin with a horizontal stripe vs. dorsal fin with circular spots. Also, (5) the cheek and opercular opalescent coloration in B. pardus is not apparent in preserved specimens, while B. sp. Middle Pastaza-Morona has the most extensively developed cheek and opercular coloration among Ecuadorian Bujurquina, with thick and long lines (dorsally) and spots (ventrally) on the cheek and spots on the opercular series. Bujurquina sp. Middle Pastaza-Morona is so far known from the Middle Pastaza basin (Upper Bobonaza, Upper Gualino) and the Upper Morona basin (Middle Macuma, Lower Shaime). We additionally identify one juvenile specimen from the Upper Villano (Upper Curaray basin) as also conspecific with Bujurquina sp. Middle Pastaza-Morona.

Bujurquina sp. Lower Pastaza-Morona (sp. 13) has a suborbital stripe shape and pattern development most similar to B. zamorensis and B. sp. Upper Napo, differing from the latter and similar to the former in that in the cheek blotch is always accented and always present as the only blotched part of the suborbital stripe. Bujurquina sp. Lower Pastaza-Morona is distinct from all other Ecuadorian Bujurquina species in that it has the smallest area of opalescent pigmentation on the head. Interestingly, the opposite extreme is found in the parapatric species B. sp. Middle Pastaza-Morona found in the same river basin and apparently only separated along an elevational gradient. The opalescent pigmentation on the head is in B. sp. Lower Pastaza-Morona only found well posteriorly from the mouth angle, while in all other species in Ecuador, it is found at least from the mouth angle posteriorly (B. zamorensis and B. sp. Santiago have the next smallest extent) or is found even more anteriorly on the head. Another distinguishing character of B. sp. Lower Pastaza-Morona from all other Ecuadorian Bujurquina species is the pigmentation of the dorsal fin, which is nearly immaculate dorsally and has a longitudinal line in the ventral portion that branches into spotting, i.e. in the lower part similar to e.g. B. zamorensis or B. sp. Upper Napo. Bujurquina sp. Lower Pastaza-Morona is so far known from the Lower Pastaza basin (uppermost Lower Bobonaza) and the Lower Morona.

Bujurquina sp. Upper Napo (sp. 2) and B. sp. Upper Aguarico (sp. 6) are similar to B. peregrinabunda (from the Amazon basin below the Napo mouth in Peru) in (1) the position of the suborbital stripe that is shifted posteriorly to the preopercular margin and starting behind the eye, and (2) in having a dorsal fin with circular spots. Bujurquina sp. Upper Napo (sp. 2) is however distinct from B. peregrinabunda in (1) having the suborbital stripe anteriorly bend vs. straight vertical, and in juveniles posteriorly inclined vs. straight vertical, (2) by having the suborbital stripe usually fragmented into an irregular series of blotches or just a cheek blotch, (3) by having a richly ornamented cheek and opercular series by opalescent markings vs. few ornamentations, mostly absent from opercular series and by (4) having a distinctly differently (maroon) coloured snout vs. same coloration as rest of head and body. Bujurquina sp. Upper Napo is so far known from the Upper Napo tributaries. We have additionally to the above clearly diagnosable specimens also collected specimens from the Middle Napo and Middle Aguarico tributaries that have intermediate characters between B. sp. Upper Napo and B. sp. Aguarico-Napo and B. sp. Upper Aguarico. We interpret these specimens as hybrids between B. sp. Upper Napo and the latter two species.

Bujurquina sp. Upper Aguarico (sp. 6) is more similar to B. peregrinabunda in that its suborbital stripe is not anteriorly bent and that it is complete (not fragmented into an irregular series of blotches or just a cheek blotch) (distinguishing it from B. sp. Upper Napo) but it is distinct from B. peregrinabunda in that its suborbital stripe is posteriorly inclined (including in juveniles as in B. sp. Upper Napo) vs. straight vertical. It is additionally distinguished from B. peregrinabunda by the same remaining characters as B. sp. Upper Napo. Bujurquina sp. Upper Aguarico is so far only known from tributaries of the uppermost Aguarico.

Bujurquina sp. Aguarico-Napo (sp. 4) is similar to B. sp. Upper Napo (sp. 2) in having the same pattern and shape of the suborbital stripe (posteriorly shifted to start behind the eye, posteriorly inclined and anteriorly bent, usually fragmented into an irregular series of blotches or just a cheek blotch) but is distinguished from it by (1) even more developed opalescent markings on the cheek and opercular series that do not form only spots but mostly are in the form of long vermiculated lines interspersed with spots on the cheek (and spots on the opercular series), and (2) more importantly, by that its dorsal fin does not have circular spots but is patterned by series of distinct diagonal lines. Bujurquina sp. Aguarico-Napo is so far known from the Upper Indillana (Middle Napo), Upper and Lower Tiputini (tributary of Lower Napo), and the Upper Cuyabeno (Lower Aguarico tributary, a tributary of the Lower Napo) and there is an extensive number of specimens with intermediate coloration between B. sp. 4 (sp. Aguarico-Napo) and B. sp. 2 (sp. Upper Napo) which we interpret as hybrid specimens in a hybrid zone. There are also other intermediate specimens between other species which we also interpret as possible hybrids and these are indicated in the distribution maps and phylogenetic analyses.

Bujurquina sp. Putumayo (sp. 7) is the northernmost species known from Ecuador. Bujurquina sp. Putumayo is similar to the another northern species B. sp. Coca in that they belong among those species with a dorsal fin with circular spots and in particular in that they both share a cheek blotch instead of a suborbital stripe. This character combination they thus share with the three southernmost Ecuadorian species, B. zamorensis, B. sp. Santiago and B. sp. Lower Pastaza-Morona. Bujurquina sp. Putumayo and B. zamorensis are the only Ecuadorian species that in adult specimens never show any traces of the suborbital stripe and only have a distinctly developed cheek blotch. Bujurquina sp. Putumayo is distinguished from B. zamorensis and B. sp. Santiago in (1) having (retaining from the juvenile phase) opalescent lines on the cheek below the snout area vs. cheeks spotted without any lines, and (2) in having a distinctly continuous midlateral stripe (vs. discontinuous, fragmented into midlateral blotches) that continues distinctly onto the caudal peduncle vs. into the axil of the dorsal in and caudal peduncle. Bujurquina sp. Putumayo is known from two Ecuadorian localities in the tributaries of the Rio San Miguel, Rio Putumayo basin.

Bujurquina sp. Upper Coca (sp. 5) is similar to Bujurquina sp. Putumayo except having traces of a suborbital stripe, in having the suborbital stripe in juveniles straight vs. bend, and in body proportions. Bujurquina sp. Coca is known only from two localities in the Rio Dashiño basin, Upper Coca basin.

Bujurquina sp. Upper Curaray (sp. 8) is most similar to B. sp. Upper Napo and to B. sp. Villano. It lacks the extensions of the anterior body bars onto the dorsal fin of the latter and hence is similar to B. sp. Upper Napo, and is intermediate in the shape of the suborbital stripe. Bujurquina sp. Upper Curaray has the bluest lower lip and the most extensive and the best developed blue coloration of the lower head of all Ecuadorian Bujurquina species. This blue coloration also extends onto the cheek and opercle in the form of well-developed rather thick lines, in contrast to predominantly small spots or very short and narrow lines in the other two similar species. Bujurquina sp. Upper Curaray is known from the uppermost Curaray and Upper Nushiño (also Curaray basin).

Bujurquina sp. Villano (sp. 9) has a distinctive very well-developed suborbital stripe that is straight, posteriorly inclined, wide and distinctly wider dorsally than ventrally. No other species of Bujurquina has such a suborbital stripe; the most similar is found in B. sp. Aguarico, where it is however posteriorly shifted and is not distinctly wider dorsally. Additionally, B. sp. Villano has completely uniquely coloured dorsal fin in adult females, where two anterior body bars continue onto the dorsal fin. This condition is somewhat similar to B. sp. Middle Pastaza-Morona and B. pardus, but in B. sp. Villano, the coloration is limited to only two vertical bars, reaches much higher and is much more intense and is only found in females (vs. in both sexes). Bujurquina sp. Villano is so far only known from a single locality in the Villano river, Upper Curaray basin.

Bujurquina sp. Lower Curaray (sp. 10) is distinct from all other Ecuadorian Bujurquina species by having very dominant, wide and widely spaced cheek stripes that cover the whole of the cheek, with the absence of any spots or thin lines in between the stripes. The vertical body bars are also unique among the Ecuadorian Bujurquina species in showing the largest degree of fusions and most specimens only show four vertical bars as adults (vs. five to six). Its suborbital stripe is also distinctive from most other Ecuadorian Bujurquina species in being covered (i.e. interrupted) by the cheek stripes and not being sharply bordered, similarly as in B. sp. Middle Pastaza-Morona, but unlike most other species, where the suborbital stripe interrupts the cheek lines and is sharply bordered (or has other forms, e.g. is reduced to a suborbital blotch). Its dorsal fin coloration is also distinct from most species in being diagonally striped (hence similar only to Bujurquina sp. Aguarico-Napo; sp. 4, among Ecuadorian species) vs. with circular spots in most other species, but the striping is more coarse in juveniles and more fine and faint in adults than e.g. in B. sp. 4. Bujurquina sp. Lower Curaray is so far only known from several localites in upland areas close to Lorocachi on both sides of the Curaray river.

Bujurquina sp. Cononaco (sp. 15) has a well-developed suborbital stripe that is straight, posteriorly inclined and rather wide. Similar suborbital stripe is found in B. sp. Villano, B. sp. Aguarico and B. moriorum. In adults of B. sp. Cononaco, the suborbital stripe has however a tendency to become shortened so that only its ventral portion is well developed, a situation that is most similar to B. moriorum. The dorsal fin coloration of B. sp. Cononaco is different from all these three species and in fact from all known Bujurquina species in that its dorsal portion has short thick diagonal stripes (similar to B. sp. Oran, which however has two series of thin short stripes and unlike the whole dorsal fin and with much thinner stripes in B. sp. Aguarico-Napo, B. sp. Lower Curaray, B. ortegai or B. mariae) and the bottom portion of the dorsal fin has very irregular blotches that sometimes approach circles as seen in other species, where they are however found on the whole dorsal fin (e.g. B. sp. Upper Napo, B. sp. Upper Aguarico, B. zamoresis among the geographically closest species).

Bujurquina sp. Napo (sp. 3) is known only from a few juvenile specimens from a single locality and morphological characterization is thus far difficult. B. sp. 3 is however distinct form the surrounding B. sp. Upper Napo (sp. 4) by having even in small specimens a suborbital stripe reduced to a blotch (vs. present as a full stripe) that is located much more ventral and posterior than in adults of B. sp. 4. The suborbital stripe course in B. sp. 3 is also straight from the eye postero-ventrad (vs. bend from a vertical position starting behind the eye and not below it). The dorsal fin in B. sp. 3 has circular blotches as in B. sp. 4.

Species diversity of Bujurquina in the Ecuadorian Amazon based on molecular data

Morphological determination suggests fifteen species-level taxa of Bujurquina in the Ecuadorian Amazon (Fig. 1). Our molecular analyses included all but two (B. pardus and Bujurquina sp. Lower Curaray, i.e. species 11 and 10) and the analyses find eleven monophyletic clades, that correspond to all thirteen included putative species (species 1, 2, 3, 4, 5 + 6, 7, 8 + 9, 12, 13, 14, 15; Figs. 1, 4, 5, S2). Four putative morphological species are recovered by the molecular analyses in two clades, each grouping two species without reciprocal monophyly (species 5 + 6, and 8 + 9; Figs. 1, 4, 5, S2).

In total, both bPTP and GMYC molecular delimitation analyses delimit the same total number (25; Fig. 4) of Bujurquina species throughout our whole sampling (i.e. including also species outside of Ecuador). Molecular phylogenetic analyses find 30 monophyletic species-level groups, while morphological characters delimit 33 species in this sampling (Figs. 4, 5, S2).

Endemism of Bujurquina in Ecuador

The following candidate species identified in this study are endemic to the foothill areas of Ecuador: B. sp. Upper Napo (= B. sp. Ecuador 2), The little known B. sp. Napo (= B. sp. Ecuador 3), B. sp. Upper Coca (= B. sp. Ecuador 5), B. sp. Upper Aguarico (= B. sp. Ecuador 6), B. sp. Upper Curaray (= B. sp. Ecuador 8), B. sp. Villano (= B. sp. Ecuador 9), B. sp. Middle Pastaza-Morona (= B. sp. Ecuador 12) and B. sp. Santiago (= B. sp. Ecuador 14).

The following candidate species from this study are endemic to the lowland areas of Ecuador: B. sp. Bolivia, B. sp. Aguarico-Napo (= B. sp. Ecuador 4), B. sp. Putumayo (= B. sp. Ecuador 7), B. sp. Lower Curaray (= B. sp. Ecuador 10),= and B. sp. Lower Pastaza-Morona (= B. sp. Ecuador 13).

The other striking feature of the endemism of Bujurquina in Ecuador is that there are six species (species 1, 2, 4, 12, 13, 15) which are endemic to areas not corresponding to presently existing river basins, two of these are even endemic to parts of presently completely isolated river basins (sp. 1, sp. 13), and that most Ecuadorian species are distributed parapatrically in the individual basins, with likely no basin hosting just a single species (see below for the Putumayo and Tigre basins; Figs. 2, 4, 5, S2).

Bujurquina sp. 1 is endemic to the upper Santiago, upper Pastaza, and upper Napo, B. sp. 2 to in upper Napo, middle Coca and middle Aguarico, B. sp. 4 to the middle Napo and middle-lower Aguarico, B. sp. 12 to the lower parts of Pastaza and Morona, B. sp. 13 to the middle parts of Pastaza and Morona, and also upper Villano, part of the Curaray and hence Napo basin, and B. sp. 15 to the upper Cononaco, i.e. Curaray basin and hence lower Napo basin, and to the upper Napo.

Within Ecuador from south to north, the individual river basins host the following numbers of species: the Santiago basin two species (sp. 1 and 14), the Morona also two species (sp. 12 and 13), the Pastaza three species (sp. 1, 12 and 13), the Curaray hosts four morphospecies (sp. 8, 9, 10 and 15), the upper Napo per se hosts three species (sp. 1, 2, 3) plus hybrids with a fourth (sp. 4), including the Coca hosts five species (sp. 5), the Aguarico hosts two species (sp. 4 and 6) plus hybrids with a third (sp. 2). Only the Putumayo (sp. 7) and Tigre (B. pardus) are currently known to host a single species within Ecuadorian territory.

The Napo basin portion of the Ecuadorian Amazon hosts ten Bujurquina species (sp. 1–6, 8–10, 15), while the Marañón basin, covering the remaining portion of the Ecuadorian Amazon, hosts five species (sp. 1, 11–14).

Bujurquina zamorensis (sp. 1) is found both in the Zamora river (type locality) and in other areas of the upper Santiago basin, as well as in the adjacent but non-communicating upper Pastaza basin, and as well in the uppermost Napo basin (Fig. 2). The Santiago and Pastaza basins are now completely separated by the Sangay volcano foothill area. The only sample we have obtained and analysed from the Zamora river basin shares the same haplotype as some specimens from the upper Pastaza river. Bujurquina zamorensis is additionally found in the upper Pastaza plateau area also in upper reaches of rivers that flow into the upper Napo. Bujurquina zamorensis is (together with the parapatric sp. 12) the only species in Ecuador found in three of its largest non-communicating river basins (Santiago, Pastaza, Napo).

Bujurquina sp. 12 and sp. 13 have distributions shared between the Pastaza and Morona river basins and the species are distributed parapatrically in the basins (Fig. 2). Bujurquina sp. 12 is found in the foothill areas and B. sp. 13 in the lowlands, the first predominantly on stony bottoms with faster current and the second on soft bottoms with low currents. Both species continue in this biogeographic pattern all the way south in Peruvian territory, at least just south of the Pongo de Manseriche (or Rentema) on the Marañón river (OR, unpublished data). Bujurquina sp. 12 has a basal phylogeographic dichotomy caused by the Pastaza river dated at 0.75 Ma (Figs. 2, 4, 5), and species 13 has the same age (0.75 Ma) of divergence from its sister species (sp. 1, i.e. B. zamorensis). Bujurquina sp. 12 is found in three separate large river basins, apart from the foothills of the Pastaza and Morona also in the upper Villano river, which is the southernmost tributary of the upper Curaray. The Villano sequenced specimen shares the same haplotype with the adjacent Bobonaza specimens (middle Pastaza tributary).

Bujurquina sp. 15 is known only from the upper Cononaco basin, a tributary of the Curaray, the largest tributary of the Napo, which however joins the Napo faraway deep inside Peruvian territory. The sequenced specimen from the upper Cononaco basin shares the same haplotype with a specimen morphologically rather clearly determined as B. sp. Upper Napo.

Bujurquina sp. 4 has the least intriguing distribution among those species not found in a single river basin. It is known from the middle Napo and its tributaries in Ecuador, specifically in the Tiputini and Indillana basins south of the Napo and in the Cuyabeno basin north of the Napo (which is a tributary of the Aguarico, itself a tributary of the Napo) (Figs. 2, 4, 5, S2).

Bujurquina sp. 2 has morphologically clear-cut specimens found only in the upper Napo, but its phylogeography finds its haplotypes also in specimens in the middle Aguarico, in the Jivino, in the Cuyabeno and even in the Putumayo (Figs. 2, 4, 5, S2). The dating of divergence of these haplotypes falls between 0.14 and 0.33 Ma. Bujurquina sp. 2 (upper Napo) is separated in distribution from the similar sp. 6 (upper Aguarico) and from the dissimilar sp. 5 (upper Coca) by the outlying eastern volcanic massifs of the Sumaco (3,990 m) and Reventador (3,562 m) volcanoes.

Phylogeny, time frame of evolution and biogeography of Bujurquina based on mtDNA analyses

Phylogenetic analyses of the 1055 characters cytb data matrix performed with NJ and MP in PAUP* (performed with the total 251 sequences including 201 ingroup sequences) and BI analyses in MrBayes and BEAST (performed with the corresponding haplotypes data matrix) provided robust and very similar results (Figs. 1, 4, 5, S2). BI runs in independent analyses with the best-fitting model of evolution (GTR + I + G) converged well (ESS > 200 for all parameters) and led to identical topologies under the same settings. The NJ and MP analyses were done with the complete dataset (Fig. S2), the BI and BEAST analyses with the haplotype dataset (Figs. 1, 4, 5).

At the species level, only one difference was found between the NJ and MP vs. BI and BEAST analyses in the position of samples E17-25-74, and E17-29-92 and 93 which are in the first analyses (Fig. S2) placed as an unsupported sister clade of sp. 5 + 6, while in the BI and BEAST analyses (Figs. 1, 4, 5), they are placed within sp. 2.

No differences were found at the supra-specific level of relationships between the MP, BI and BEAST analyses.

The Ecuadorian Bujurquina species do not form a monophyletic group, since they are interspersed by species from Colombia and northern Peru. Most Ecuadorian Bujurquina species are, however, found in one clade (all but species 12). Species 12 is thus the clear phylogenetic outlier among the Ecuadorian species, but biogeographically this is hard to justify since species 12 is nested geographically between other Ecuadorian species (spp. 1, 9, 11, 13, 14) and has a parapatric distribution between species 1 and 13 in the same river basin (the Pastaza). The relationships of the remaining Ecuadorian Bujurquina species have a clear biogeographic foundation, with the phylogeny having a south to north topology. The northern terminal clade includes species in the Upper Napo (spp. 2, 4, 5, 6), with successive sister groups in the adjacent the Putumayo (sp. 7), in the main tributary of the Napo, the Curaray (spp. 8 + 9 and 15; and species 3 also in the Upper Napo) and in the more distant Caquetá in Colombia. This northern clade has a basal dating of 3.2 Ma (Figs. 4, 5), the Upper Napo clade (spp. 2, 4, 5, 6) of 1.1 Ma. The sister group of the northern clade is species 14 in the non-adjacent lower Santiago (plus two other species further south in Peru) with a divergence of 3.6 Ma. The geographical gap between the northern clade and species 14 is filled by the next clade, species 1 + 13, with a shallow divergence of 0.75 between its species and a divergence of 4.6 Ma from its sister group. The divergence of all these Ecuadorian species from the distantly related yet geographically embedded species 12 is the considerable 13.9 Ma.

Our phylogenetic reconstruction divides Bujurquina into two main clades which can be diagnosed by patterned vs. unpatterned dorsal fin, corresponding in geographical terms to Northern vs. Southern group (Figs. 3, 4, 5).

Biogeographical analyses reconstruct a predominantly vicariant history of the genus with increasing number of ancestral areas at nodes with increasing number of ‘maxareas’ reconstruction allowed at each node. The reported results are thus given from the analysis allowing the maximum ‘maxareas’ setting. Details of the reconstruction are shown in Fig. 5 and are briefly described below.

The main biogeographic dichotomy within the genus is between the Southern and Northern groups, separated within the present upper Ucayali basin (Figs. 3, 4, 5). The Southern group has a wide ancestral area that became fragmented by vicariance first into Paraguay-Paraná vs. Madre de Dios-upper Ucayali at 10.2 Ma, and second into the upper Ucayali vs. Madre de Dios at 8.4 Ma.

The Northern group also has a wide ancestral area, which includes all of the northern western Amazon plus Orinoco, and which also became fragmented by vicariance (Figs. 3, 5), but in a more complicated pattern than the Southern group, and the pattern at the basal nodes requires better areal sampling within the Northern group. The oldest clearly identifiable vicariant event is the separation of the Orinoco basin from the northern western Amazon at 7.9 Ma. The oldest node in the Northern group is much older at 13.9 Ma. The clade that diverged at 13.9 Ma includes one species from Ecuador (sp. 12; Fig. 3), while all other Ecuadorian species are found in a much younger clade, dated at 4.6 Ma, that diverges as the next node below the Orinoco-northwestern Amazon divergence. This clade had a wide ancestral area that included most of northern western Amazonia, namely (from S to N) the Marañón, Santiago, Pastaza, Curaray, Napo, Putumayo and Caquetá river basins, and their connecting paleo-stem of the main proto-Amazon. This clade became fragmented by vicariance between 4.6 and 3.6 Ma, with the separation of the southern river basins of Marañón, Santiago and Pastaza. The next to separate river basins from the rest were the southernmost Curaray and the northernmost Caquetá, which both separated at 1.5–1.6 Ma. Separation of the Napo and Putumayo occurred much later at 0.14 Ma.

Paleogeography of Ecuadorian Amazonian rivers based on Bujurquina phylogeography

History of connections between river basins in Ecuadorian Amazonia

The species boundaries of identified Ecuadorian species are often not corresponding to present river basins or limited to only their portions (Fig. 2) suggesting relatively recent reconfigurations of the drainage divides as hinted at above, while phylogenetic relationships of all species suggest more ancient configurations of the river basins (Fig. 5).

The biogeography of Bujurquina suggests the following paleogeographic reconstruction of the river network of Ecuador, presented following the topology of the phylogeny, i.e. from the youngest changes in the configuration in the northern areas, to the oldest changes in the southern areas (Figs. 2, 5, S2).

From ca 0.2 Ma (0.14 Ma is the youngest dated event), there is the present configuration of the river network.

Up until 0.33–0.14 Ma, there appeared to have been connection between the upper Napo and the upper Putumayo. Based on the biogeography and phylogeography of species 2, 5, 6 (diverging at 0.7 Ma) and species 4 (diverging at 1.1 Ma and with a N–S and not W–E phylogeography), this configuration of a NE flowing upper Napo existed already prior to 1.1 Ma.

At 1.5 Ma and 1.6 Ma, there is divergence of upper Curaray and Caquetá from upper Napo. Previously to these dates, the upper Curaray thus appears to have been connected with the upper Napo and the Caquetá with northern Ecuador through the Putumayo. Prior to 1.5 Ma, we thus postulate connection between all northern basins of Ecuador (whole Napo basin plus Putumayo plus the extralimital Caquetá).

At 3.6 Ma, we find the main paleogeographic event in Ecuador, the divergence of the northern Ecuador clade (all of Napo basin, Putumayo, plus Caquetá) from the south (i.e. the Marañón basin that includes its Ecuadorian basins of Santiago, Morona, and Pastaza). This areal division is occupied by the geological formation known as the Guayusa platform (see section “Discussion”).

An older divergence at 4.6 Ma occurred within the present Marañón basin, between the Santiago basin and the Morona-Pastaza basin.

Prior to 4.6 Ma, our biogeography of Bujurquina thus reconstructs a single river system in Ecuador, that is dividing from south to north. This S–N paleo-river system is thus markedly different from the presently dominant W–E configuration of the river basins. The reconfiguration to the present W–E systems thus occurred in present day Ecuadorian territory between 4.6 and 3.6 Ma in the south and 1.6, 1.5, 1.1, 0.7 Ma in the north, with persistence of connection (likely a secondary connection) between upper Napo and upper Putumayo till 0.2 Ma.

The modern Pastaza and Tigre river basins occupy the divide between the northern clade and southern Ecuador on the central portion of the Guayusa platform. The Tigre basin was not represented in the present study and thus remains enigmatic but the Pastaza basin has a strange configuration of Bujurquina fauna. It hosts three vertically distributed species (spp. 1, 12 and 13), with species 1 and 13 closely related and tied in their origin to the origin of the Pastaza river (0.75 Ma), and with the sp. 12 found in between them. Bujurquina sp. 12 occupies its distribution area for longer time than the existence of the modern Pastaza (since it is the most basal Ecuadorian species), but shows its major phylogeographic dichotomy (0.75 Ma) tied to the origin of the Pastaza river.

S–N to W–E transformation of Ecuadorian Amazonian river network

A S–N flowing river network is indicated in our biogeographic reconstruction as dominant throughout the history of evolution of Bujurquina in Ecuador (Figs. 2, 5) and outside Ecuador (Figs. 3, 5). Within Ecuador, this S–N flowing river network persisted the longest in the north, with the connection of Napo and Putumayo persisting until 0.2–0.1 Ma. The next southern tributary (Curaray) and the northernmost tributary (Caquetá) were separated from tributaries of the S–N river into W–E rivers between 3.2 and 1.5 Ma, and the southern basins (Santiago, Pastaza) separated from the N–S river still flowing to the N (Putumayo, Napo, Curaray) into W–E basins between 4.6 and 3.6 Ma.

The above-described biogeographic events successively transforming the S–N paleo-river into an W–E system within Ecuador are paralleled outside of Ecuador by the phylogenetic and biogeographic relationships of the lowland Peruvian species in the E (B. sp. Oran, B. moriorum, B. peregrinabunda, and B. syspilus) which have diverged at 3.85 Ma, 2.58 Ma, 2.44 Ma and 2.16 Ma, respectively (Fig. 5), from their Ecuadorian counterparts to the W.

The W–E transformation of the river network in the area of Ecuador and adjacent Peru thus occurred between 4.6 and 1.5 Ma, while the actual modern river beds (the modern rivers) are much younger around 0.7–0.1 Ma (Figs. 3, 5).

Unrecognized species diversity of Bujurquina in the Ecuadorian Amazon

The most striking feature of the species diversity of Bujurquina in Ecuador is its unrecognized richness and endemism. Our results clearly demonstrate that the species diversity of Bujurquina in Ecuador is significantly underestimated. Our morphological determination distinguishes fifteen putative species of Bujurquina in Ecuador (Fig. 1). Thirteen of these have been included in molecular analyses in this study and eleven of them correspond to clades in the molecular analyses (and the remaining two putative species are also contained within these clades; species 1, 2, 3, 4, 5 + 6, 7, 8 + 9, 12, 13, 14, 15; Figs. 1, 4, 5, S2). Among these eleven molecular clades two used molecular delimitation analyses delimited eight and nine species-level groups (Fig. 4). Our results thus show that Ecuador hosts at least eight and more likely eleven to fifteen species based on a still incomplete area survey of the country.

Based on the morphological characters, we consider it plausible that the unsupported species 8 (nested within sp. 9) are hybrids between the adjacent sp. 9 and sp. 2, while sp. 5 and sp. 6 we interpret as having only very recently diverged and not yet having gained reciprocal monophyly or possibly as a case of mtDNA sweep (Figs. 4, 5, S2). Our analyses find other instances of likely hybridization or introgression between species 2, 4, 5, 7 (Figs. 4, 5, S2).

The molecular delimitation analyses recover all the molecularly highly divergent species (species 3, 4, 7, 9, 12, 14, 15) but do not delimit the shallowly divergent species 2 vs. 5 + 6 and 1 vs. 13 (Fig. 4). The bPTP analysis also does not delimit species 4 from species 2 and 5 + 6 (Fig. 4). While species 2 vs. 5 are rather similar in morphological characters and could represent a single species, this is unlikely in the case of species 4 vs. 2 and 5 and species 1 vs. 13, which are highly distinct in morphological characters (Fig. 1).

In our total sample analysed in this study that includes Bujurquina material from throughout its distribution, our morphological determination distinguishes 33 species, while molecular phylogeny and molecular delimitations delimit 28 and 25 species, respectively (Fig. 4). Presently, there are 18 valid species in Bujurquina and our results thus clearly indicate that the species diversity of the genus is severely underestimated not only within Ecuador, but within its whole distribution area. Our results (delimited vs. included species) suggest a total diversity of at least 34–38 species. The total diversity will however most likely be significantly higher, since endemism in Bujurquina is based on our results unparalleled among cichlasomatine cichlids and since large areas of the Amazon remain completely unexplored for Bujurquina diversity.

Present classification of the genus on the other hand only treats four species (Kullander, 1986) as occurring in Ecuador, of which two are not found in Ecuador based on our results (B. moriorum and B. peregrinabunda). These two species have mostly been considered as extending in their distribution from their type localities in Peru into Ecuador (Kullander, 1986). Barriga (2012) considered two other species described by Kullander (1986) as clear Peruvian endemics as also occurring In Ecuador (B. huallagae and B. syspilus). Morphological determination as well as our molecular phylogeny based on topotypes of these four species clearly demonstrate that neither is found in the so far studied areas of Ecuador (Figs. 1, 2, 3, S1). Only two currently recognized species have type localities in Ecuador (B. zamorensis and B. pardus) and both are here considered endemic to Ecuador.

The only recently described species of Bujurquina from Ecuador, B. pardus, described based on material collected in 1989, was not included in the present study because we have not been able to collect new specimens. We have collected specimens 8 and 14 air km from the GPS coordinates given in original description in tributaries of the Rio Shionoyacu and Rio Maratiyacu (tributaries of the Rio Tigre basin in this area), but these specimens do not agree in diagnostic characters with the diagnosis of B. pardus nor with the one preserved photographed specimen. In preserved state (where many coloration pattern characters are lost), B. pardus is among the here studied species most similar to B. sp. Middle Pastaza-Morona (B. sp. 12) in suborbital stripe, dorsal fin coloration, continuation of vertical bars onto dorsal fin; supported also in phylogenetic analysis based on coloration pattern data as presented in Říčan and Říčanová, under review, differing in white lappets of dorsal fin spines (vs. red), and in better developed blotches on flank scales. Among all material of Bujurquina known to us, B. pardus is most similar to Bujurquina sp. Panguana as judged from the photograph in Vriesendorp et al. (2007). Arbour et al. (2014) stated that B. pardus is most similar to B. huallagae, but these are rather dissimilar species as based on the diagnosis and description and the single available photo of B. pardus (Fig. 2 in Arbour et al., 2014 and cf. Fig. S2, a specimen of B. huallagae collected from close to the type locality in the uppermost Huallaga river) and on the phylogenetic analysis of the coloration pattern data (Říčan & Říčanová, under review), and also based on the distribution of B. huallagae, which is endemic to the Upper Huallaga river basin above the first pongo (Kullander, 1986; Říčan, unpublished data).

Despite our wide ranging sampling and potential for understanding of Bujurquina diversity, our present study has one potentially limiting factor and that is the use of a single uniparentally inherited molecular marker (mtDNA cytb). The potential limitations for species delimitation using mtDNA have been discussed in literature (e.g. Rannala & Zang, 2003; Yang & Rannala, 2010; Dupuis et al., 2012; Fujita et al., 2012). The use of additional nuclear markers is generally recommended as the use of these unlinked markers has the potential to improve the accuracy of phylogenetic reconstructions and species delimitation. In general, studies comparing species-level diversity based on both mtDNA and nDNA molecular phylogenies are still rare among Neotropical riverine cichlids. This is especially so within groups that have a good specimen sampling, and have comparable species diversity and ecology to Bujurquina. Such studies are limited to just two Neotropical cichlid groups; (1) the largest group of cichlids in Middle America (the heroine cichlids sensu Kullander, 1998, a junior synonym of the therapsine cichlids sensu Allgayer, 1989), and (2) the largest South American cichlid genus, Crenicichla.

In both the Middle American cichlids (Říčan et al., 2008, 2016) and in Crenicichla (Piálek et al., 2012, 2019a, b; Burress et al., 2018; Říčan et al., 2021), conflicts between the mtDNA and nDNA phylogenies have been found above the species level, but at the species level (i.e. in species delimitation), both mtDNA and nDNA are very congruent delimiting the same species diversity with very few exceptions. Hence, there is no clear indication from other similar cichlid groups that the species diversity of Bujurquina delimited solely by mtDNA, as in this study, should be questionable. This is especially so since all the species are diagnosable, and have been a priori delimited using coloration patterns, morphological characters and geography and only a posteriori using the mtDNA marker.

Underestimated Neotropical fish biodiversity

The minimum percentage increase of unrecognised species observed in this study in Bujurquina (83%; 33 species identified, 18 species valid) significantly surpasses the conclusions of Reis et al. (2016), who estimated that 34–42% of Neotropical freshwater fishes remain undescribed, and approaches the values in the recently studied largest Neotropical cichlid genus Crenicichla with a diversity centred also on the Amazon basin, where minimum percentage increase of unrecognised species was 65–121% (Říčan et al., 2021).

Reis et al. (2016) gave as the main explanations for the unrecognised diversity of Neotropical fishes: (1) the lack of systematic sampling throughout the distribution of the given group, in particular areas with difficult access; (2) widespread taxa or heterogeneous taxa with insufficient or overwhelming amounts of museum material or improperly preserved material and information; (3) cryptic or pseudocryptic (morphological differences apparent but overlooked) diversity in widespread species. In the cichlid fishes, an important factor contributing to last point is the limitation of available diagnostic characters retained in preserved specimens. Cichlids, being colorful and visually oriented fishes, lose almost all diagnostic coloration pattern characters when preserved. The study of live fishes and their photographs is therefore especially important, but detailed images of sufficient numbers of live specimens in the field are still very rarely taken in general field surveys. Existing collections typically lack properly documented live coloration of the preserved specimens, which greatly limits their usefulness.

The most importing factors impeding the detection of the actual species diversity in Bujurquina are encompassed in points 1 and 3. For the foreseeable future, the lack of sampling in areas of difficult access is likely to remain an impediment to our understanding. Meanwhile, modification and destruction of the cichlid habitats is likely to continue, as has occurred throughout the Neotropical region.

Morphological diagnosability of Bujurquina species

The best diagnostic characters of Bujurquina species are coloration pattern characters, in particular the suborbital stripe pattern, the pattern of the cheek stripes, dorsal fin coloration and the orientation of the midlateral stripe (Figs. 1, S1). All these characters have already been used since the first descriptions of Bujurquina species (Kullander, 1986), but only with limited success since all Bujurquina species descriptions have been based only on preserved material, while only live specimens show the full development and diversity of the coloration patterns.

Live specimens are thus critical for the evaluation of these characters, especially the pattern of the cheek stripes (which are completely lost upon preservation), and the dorsal fin coloration, which is difficult to see in preserved specimens. The patterns of the suborbital stripe and the midlateral stripe also change markedly between live and preserved specimens, and must be treated in descriptions entirely separately. Because of these differences between live and preserved coloration, Kullander (1986) and all subsequent publications that were only based on preserved specimens could see only a portion of the actual species diversity in Bujurquina. Notably, Kullander (1986) described B. peregrinabunda as a widespread species also occurring in Ecuador, while based on our data, B. peregrinabunda is not found in the studied area of Ecuador, and is not closely related to the here delimited species which could have been confused with it (B. sp. 2, 4, 6). The same confusion has been around B. moriorum, the type species of the genus, because the descriptions based solely on preserved material are simply not sufficient to diagnose species of Bujurquina, and species of cichlids in general (e.g. Piálek et al., 2015; Říčan et al., 2021).

Here, we presented the main differences between Bujurquina species found in Ecuador. As demonstrated in Fig. 1, each so far studied species (valid and newly proposed in this study) of Bujurquina in Ecuador and in the Northern group of Bujurquina can be diagnosed by a unique combination of just three main coloration pattern character complexes (Figs. 1, S3). These three character complexes can be studied in preserved specimens, but their actual diversity is only discernible in live specimens. These three character complexes are (1) the pattern of the suborbital stripe development from juveniles to adults and its adult configuration, (2) the pattern of the coloration of the dorsal fin and (3) the course of the midlateral stripe. For the true understanding of the suborbital stripe pattern, it is necessary to study both juvenile and adult specimens since only a developmental series shows the development of the particular species-specific pattern and enables to understand the sources and degrees of variability of the character in the specimens (Říčan et al., 2005). The diversity of the coloration patterns of the dorsal fin is much higher than the dichotomy given by Kullander (1986) and subsequently used by all studies (i.e. immaculate vs. spotted), as various patterns of spotting, lines, blotches and combinations of these can be found in the various species. The true diversity of these patterns can only be seen in live specimens where the dorsal fin is still pliable. Preservation can leave the dorsal fin stiff, impossible to open, and with a smeared patterning. The course of the midlateral stripe is the easiest character to study in preserved specimens, and is actually the only character that is better seen in preserved specimens than in live specimens due to the loss of most other coloration patterns on the body. There are dozens of other coloration patterns that have additional detailed diagnostic value in Bujurquina species. These need to be included into any well-performed species diagnoses, but for quick orientation and determination, our results show that careful study of just three main characters is sufficient to identify all so far studied species.

Biogeography and dating of Bujurquina diversification

The endemism of Bujurquina plays out in two main biogeographic zones since there is a separate band of endemism in the Andean foothills, which has a distinct fauna from the supposedly more widely distributed species in the Amazonian lowlands. Foothill areas of the same river basins thus have a different Bujurquina fauna from lowland areas of the same basins not only in Peru where originally proposed (Kullander, 1986), but also in Ecuador as found in this study.

Our sampling of Bujurquina in this study includes the majority of described species plus numerous putatively new species. Despite this major increase in species diversity in the genus, the biogeographic history of the genus can still be ascertained only to a limited level, both in the context of the whole western Amazon and of the country of Ecuador studied here. A significant collecting effort is still needed before this very species-rich cichlid genus comprising narrow endemic species will be able to provide clear insight into the evolution of the western Amazonian biogeography and landscape.

There are several large-scale events in the Amazonian evolution that can already be seen in our Bujurquina biogeographic analyses. These include first and foremost the major structural arches of western Amazonia that determine its present river basin configuration (Figs. 3, 5, 6). These are the Vaupes arch, Fitzcarrald arch and Michicola arch. There is also a major role of the Guyusa platform on Ecuadorian biogeography (see next section).

Main geographic barriers and drainage basin divides in the western Amazon with large-scale transformation of the western Amazonian river network based on the Bujurquina reconstruction. a Main geographic barriers and drainage basin divides in the western Amazon. b Putative paleo-course of the Amazon. c Large-scale transformation of the western Amazonian river network between 10 and 8 Ma. d Transformation of the western Amazonian river network in the area of Ecuador between 4.6 Ma and recent. See text for additional details

Full size image

The Vaupes arch is a complex geomorphological feature that separates the Orinoco and western Amazonian basins and the separation of the basins is dated in the literature at 8 Ma (Mora et al., 2010; Winemiller & Willis, 2011). The Orinoco basin has only one described species of Bujurquina (B. mariae), which was included in our study. In our analyses calibrated independently of this event, the dating of the separation of the two basins agrees with the generally accepted date very precisely (7.9 Ma; Fig. 5).

The Fitzcarrald arch is also a complex bulge (400,000 km², rising to 600 m), which divides the western Amazonian Basin into three subbasins: northern Amazonian foreland basin, southern Amazonian foreland basin and the eastern (or Acre) Amazonian foreland basin (Espurt et al., 2007, 2010; Mora et al, 2010; Roddaz et al., 2010; Hoorn et al., 2010). Hence, it separates the Ucayali headwaters (the longest axis of the Amazon) from the Madre de Dios headwaters (the longest and largest tributary of the Amazon). It also separates the Purus and Juruá basins from the Andes and each other. The Fitzcarrald arch is generally assumed to have started forming 4 Ma (e.g. Espurt et al., 2007; Hubert et al., 2007; Regard et al., 2009; Hoorn et al., 2010). But the Nazca Ridge, which causes the Fitzcarrald arch by a flat slab subduction below South America, began subducting much earlier, approximately 11.2 Ma at 11° S, with the current subduction location at 15° S (Hampel, 2002). In our analysis of Bujurquina, we have included species from the Ucayali and Madre de Dios basins and their divergence is dated at 8.4 Ma (Fig. 5). Our dating thus does not agree with the generally proposed 4 Ma dating from the literature, but it does agree with the dating of the Vaupes arch and probably also of the intervening Guayusa platform (see next section). The difference between the generally assumed 4 Ma and our dating of 8 Ma can have several different explanations. One is that divergence across a complex geomorphological structure can cause vicariance of different organisms with different ecologies at different times. A second explanation is that incomplete taxon sampling will cause different divergence dates. No study across the Fitzcarrald arch so far can claim to have had a complete species and area sampling and our study is no different, since we have only sampled known species from two basins, and have not sampled in the other two basins (Purus, Juruá). A third explanation is that the geological dating of the Fitzcarrald arch is not precise, being stated as not older than Pliocene (2.6–5.3 Ma) based on a combination of dating of the cessation of the volcanic activity in the adjacent Andes (4 Ma) and of sediments and geomorphic markers which are Pliocene to Pleistocene (Espurt et al., 2007, 2011). A fourth possibility is that our molecular-clock calibration suggests a date that is too early for the event, but we note correspondence of our results to the Vaupes and Michicola arch datings in the literature.

The Michicola arch separates the southern Amazonian foreland basin from the Paraguay river basin to the south and is in our analyses of Bujurquina dated at 10.2 Ma (Fig. 5). The Michicola arch has been interpreted in the literature as much more porous to biogeographic patterns than the other here mentioned arches, but our dating is compatible with the suggested dates of this structure (Lovejoy et al., 2010; Carvalho, 2011).

The 8–4.5 Ma time frame is considered a pivotal point in the history of evolution of the Amazon basin since it marks the turning point between a long existing S–N flowing western paleo-Amazon and a W–E flowing modern Amazon (Hoorn & Wesselingh, 2010a, b; Albert & Reis, 2011; Albert et al., 2018). In our biogeographic reconstructions, this pivotal point is visible in the correspondence of biogeography and dating of Bujurquina to the emergence of the arches that disrupted the S–N paleo-Amazon (Figs. 3, 5, 6). Our results for the territory centred on Ecuador and also including adjacent Peruvian and Colombian samples suggest that a very different paleogeography of river basins compared to the present one existed during the period between 8 Ma and until very recently (0.1 Ma) when the modern river basin pattern finally emerged. Hence, based on our results, the S–N to W–E transformation of the Amazon was not an instantaneous event, but an event that spanned most of the last 8 Ma.

Paleogeography of Ecuadorian rivers based on Bujurquina biogeography

The densest sampling of Bujurquina so far is for Ecuador in this study and hence possible reconstruction of paleogeography is more plausible for Ecuador at the moment than for any other area of the western Amazon. Our results demonstrate that Ecuadorian species diversity of Bujurquina is severely underestimated and that Bujurquina has highly complex biogeographic history in Ecuador. Half of the putative species discovered in our study in Ecuador do have distribution areas that do not correspond to a single major river basin. This observation coupled with the phylogenetic relationships of most of the species and the overall topology of the phylogeny of the genus in Ecuador suggest a markedly different paleogeography of Ecuadorian rivers compared to the present river network (Figs. 2, 5, S2). Both the degree of Bujurquina endemism in Ecuador and the detail of resolution of spatial biogeography attained by this genus are in our knowledge without parallel in the published ichthyological literature of Ecuador, and also not paralleled by any other cichlid genus.

While the present river network of Ecuadorian Amazonia has a clear W–E pattern of flow, our reconstruction of the paleogeography based on Bujurquina phylogeography is dominated by an S–N pattern. Based on molecular-clock dating, the restructuring of the originally S–N trending river into a W–E river network occurred after 8 Ma outside of Ecuador, and after 4.6 Ma within Ecuador. Final establishment of the modern rivers and basins, including the largest ones in Ecuadorian territory (i.e. the Pastaza, Napo, Putumayo), occurred only very recently between 0.75 and 0.1 Ma (Figs. 2, 5, S2).

The main cause of the transformation of the Ecuadorian Amazonian river network is clearly the propagation and rise of the Andes in general (Hoorn & Wesselingh, 2010a, b; Albert and Reis, 2011; Albert et al., 2018) and of the Guayusa platform in Ecuador in particular (Longo & Baldock, 1982; Figs. 2, 3, 5). The Guayusa platform is the result of subduction of the Carnegie Ridge (running E from the Galapagos rise) beneath western South America and beneath the Andes in an area that basically defines the territory of Ecuador (Figs. 2, 3). The Guayusa platform shares the same mechanism of origin and is hence a similar feature as the better known Fitzcarrald Arch in Peru (Espurt et al., 2007; Regard et al., 2009; Hoorn et al., 2010). The Carnegie Ridge is an aseismic ridge on the Nazca Plate that is thought to be a result of the passage of the Nazca Plate over the Galapagos hotspot (Gutscher et al., 1999; Michaud et al., 2009). The dating of the subduction of the Carnegie Ridge beneath South America is not precise, dated from about mid-Miocene (15 Ma; Spikings et al., 2001) to Pleistocene (2 Ma; Gutscher et al., 1999), and we have not been able to find any precise dating for the origin of the Guayusa platform itself in the literature, but the sediments that are found on top of the Guayusa platform are Miocene/Pliocene (5 Ma) to Pleistocene, as opposed to the surrounding Miocene sediments of the adjacent lowlands. The sedimentological situation is thus the same as in the Fitzcarrald arch (Espurt et al., 2007).

Similar to the better known Fitzcarrald or Vaupes arches, the Guayusa platform also separates the headwaters of major tributaries and source areas of the Amazon (Longo & Baldock, 1982; Figs. 2, 3). It separates the Napo basin from the Pastaza basin, and it also separates the Morona basin from the Pastaza basin, and likely is also the cause of the separation of the Morona basin from the Santiago basin by causing the rise of the narrow mountain spur between them. The platform also separates the largest tributary basin of the Napo, the Curaray, from its parental basin and is also the place where another large source of the Amazon finds its headwaters, the Tigre basin. The influence of the rise of the platform possibly continues so far north that it may be responsible for the separation of the Napo basin from the Putumayo basin (the Miocene/Pliocene sediments of the Guayusa platform continue north to the Aguarico basin), and even of the Putumayo basin from the Caquetá basin. If that is the case than the Guayusa platform, or the subduction of the Carnegie ridge in more general terms, is the most important biogeographic driver in the whole western Amazon since it separates the largest number of major Amazonian tributaries from each other.

Our dating of vicariant events on the Guayusa platform (Figs. 2, 5) indicates that the oldest events in the southern area of the platform are dated between 4.6 and 3.6 Ma; separating the basins of Santiago, Morona and Pastaza in the south from each other and from the Napo (including Curaray) in the north. The separation of the basins in the north took place progressively northwards, as is expected from the subduction of the Carnegie Ridge towards the NE (Gutscher et al., 1999). Towards the NE, the elevation of the extension of the Guayusa platform is progressively lower, and the subduction younger, and hence the influence on the area must be relatively recent. The separation of the Curaray from the Napo occurred, based on our dating, between 3.2 and 1.5 Ma. The restructuring of the Napo headwaters further north is younger still, dated between 1.1 and 0.1 Ma. The latter date marks the final separation and origin of present configurations of the Napo and Putumayo basins in the northernmost potential reaches of influence of the Guayusa platform.

More recent events in the south-central Guayusa platform are also identified by our data, and include the formation of the present Pastaza river (at 0.75 Ma), probably assisted in its upper reaches by the birth of the Sangay volcano, dated by geological evidence to the same time period (Monzier et al., 1999). The rise of the eastern slope of this massive volcanic structure and its periodical massive landslides may have diverted the upper Pastaza into its present course and forced it to breach through the elevated rim of the Pastaza plateau into the eastern lowlands. Based on this interpretation the upper Pastaza, previously to this change, probably emptied through the Upano basin into the Santiago basin, as suggested by the distribution area and young phylogeographic structure of B. zamorensis. Bujurquina zamorensis from the Zamora river basin shares the same haplotype as some specimens from the upper Pastaza river and this is hard to explain paleogeographically except in case the drainage divide between the two basins is extremely young. Bujurquina zamorensis is additionally found in the upper Pastaza plateau area also in upper reaches of rivers that flow into the upper Napo suggesting capture of the upper Pastaza plateau by the upper Napo tributaries.

A very recent separation of the upper Villano and upper Bobonaza basins is suggested by Bujurquina sp. 12, where the Villano sequenced specimen shares the same haplotype with the adjacent Bobonaza specimens (middle Pastaza tributary). Based on present configuration of the landscape, it appears quite likely that the upper Villano was until very recently a tributary of the upper Bobonaza, only recently captured by the rest of the Villano (and hence Curaray) basin.

A secondary young contact between the Upper Napo and the upper Cononaco (a tributary of the Curaray river) is suggested by haplotype sharing between B. sp. 15 (sp. Cononaco) and B. sp. 2 (sp. Upper Napo). This situation makes the phylogeography of this species in conflict with the apparently much older separation of the river basins (see the reconstructed paleogeography of the river basins in Ecuador). This phylogeographic pattern thus suggests secondary contact between these river basins, and a human mediation cannot be ruled out. The phylogeography of Bujurquina sp. 2 (sp. Upper Napo) also suggests past to recent connections between the upper Napo and the Putumayo dated between 0.14 and 0.33 Ma.

Reconstructed paleogeography of Ecuadorian rivers compared to geology of Ecuadorian Amazonia

The Ecuadorian Amazonian foreland basin, also referred to as the Oriente basin, is filled with sediments ranging from the upper Cretaceous to the present, which are divided by their type of origin into several main components (Longo & Baldock, 1982; Roddaz et al., 2010, 2012; Fig. S4). They include late Oligocene to late Pliocene river sediments (Chalcana, Arajuno, Chambira formations), early to late Miocene tidal sediments (Curaray formation), and Pliocene to recent lahars and alluvial sediments (Mesa and Mera formations).

The reconstructed paleogeography of Ecuadorian Amazonian rivers derived from Bujurquina biogeography (see above) is very precisely reflected in the sedimentary composition of the Oriente basin (Figs. 2, S4). The reconstructed main axis paleo-river of the basin running S–N corresponds to the late Oligocene to late Pliocene river sediments of the Chalcana, Arajuno and Chambira formations (Figs. 2, S4). These riverine sediments show deposition in both S–N and W–E directions, suggesting a main S–N trunk river and W–E sources of the river from the Andes. These riverine sediments are oldest in the W and progressively younger to the E, with increasing sedimentation rates and coarsening of sediments to younger ages, and with a westward/upward (younger) transition from meandering to braided, which all indicate an eastward progradation likely related to an increase in slope during the early to middle Miocene (Roddaz et al., 2012).

The Pliocene to recent lahars (Mesa and Mera formations) and alluvial sediments reflect the final transformation of the S–N river system into the present W–E river system. The Mesa and Mera formations constitute a Pliocene–Pleistocene large-scale fluvial fan (megafan; Bes de Berc et al., 2005), which is deposited in a W–E direction and which overlies, except where eroded, the youngest of the S–N riverine formations (the Chambira formation). The Pleistocene Mera megafan is overlain by the modern Holocene Pastaza megafan (Bernal et al., 2011) which demonstrates the recent origin of the modern Pastaza river. The Pliocene Mesa formation is the reason for the disruption of the originally S–N flowing paleo-river (the Chalcana, Arajuno and Chambira formations) into a SE flowing paleo-Pastaza-Morona and NE flowing paleo-Napo-Curaray. The Mesa formation is today isolated from its Andean sedimentary sources by an erosional belt that excavated the older underlying Arajuno formation (hence the name of the formation). This erosional excavation occurred as a direct result of the disruption of the S–N paleo-river along the present Pastaza-Curaray-Napo triple boundary into the SE flowing paleo-Pastaza-Morona and NE flowing paleo-Napo-Curaray. This oldest disruption of the S–N river in Ecuador and hence the start of the erosional excavation is dated at 3.6 Ma based on our Bujurquina biogeography.

The Mesa formation megafan caused not only the major S/N divide in the Oriente basin, but also defines the boundaries between the Napo and Curaray basins (dated in the Bujurquina biogeography between 1.6 and 1.5 Ma), and between the Cononaco and Curaray proper (dated in the Bujurquina biogeography at 3.2 Ma), and between the Curaray, Tigre and Pastaza basins (where we have so far not been able to sample sufficient numbers of Bujurquina localities and specimens).

The present W–E river system is clearly indicated in the geological landscape by the Holocene alluvial sediments (Fig. S4), and the westernmost locations of these Holocene sediments exactly correspond to the final reconstructed river flow direction changes from the original S–N system to the present W–E system in the Bujurquina phylogeography. These Holocene alluvial deposits transect the paleo-riverine sediments of the Arajuno and Chambira formations (mostly only the younger Chambira formation) and the cotemporary tidal sediments of the Curaray formation deposited to the east of the S–N paleo-river and to the east of the influence of the Guayusa platform elevation.

The distribution areas of Bujurquina species when compared with the geological composition of the Oriente basin (Figs. 2, S4) are highly congruent and reveal some interesting observations. The oldest species of Bujurquina in Ecuador (sp. 12, with a divergence of 13.9 Ma from the other species) has, as the only species in Ecuador, a distribution area that corresponds to the oldest S–N paleo-river sediments in the south of the Oriente. These include the Chalcana and Arajuno formations deposited from the late Oligocene to the middle of the late Miocene (23.3–7.3 Ma). This suggests that the species is a long-term endemic of Ecuador, and it explains its peculiar distribution area that runs perpendicular to the present river network, but follows the course of the ancient river network. The northern clade in Ecuador, comprised spp. 2–9 with a divergence starting at 3.6 Ma, corresponds to the northern branch of the S–N paleo-river (the Chalcana to Chambira sediments). The apparent sympatry of the long separated (3.22 Ma) sp. 2 and sp. 3 within this clade in the upper Napo is probably caused by secondary contact. The older sp. 3 plausibly occupied the paleo-Napo basin that dates back to the Chambira formation, as its sister group relationship with sp. 8 + 9 in the upper Curaray at 1.5 Ma testifies. Meanwhile, sp. 2 is a young species (0.7 Ma) reflecting a recent phase of the building of the upper Napo basin. One of the youngest species of Bujurquina in Ecuador is sp. 13, with a divergence of 0.75 Ma from its sister species. This species is found parapatrically to sp. 12 to the SE, and is only found in the area corresponding to the Pleistocene (Mera formation) to Holocene megafan sediments of the lower Pastaza megafan. This geological setting of sp. 13 provides direct evidence supporting our hypothesis derived from biogeography and dating of the Bujurquina phylogeny that sp. 13 (found in the lower part of the Pastaza megafan) is recently derived from sp. 1, or from a common ancestor with sp. 1 (B. zamorensis, found in the upper Pastaza part of the megafan), having become isolated from it by the origin of the modern Pastaza river at 0.75 Ma.

Based on our interpretation of our data, it seems plausible that Bujurquina species distributions and ages of divergence are correlated and can to a certain degree be predicted from the type, age, boundaries and topography of the sedimentary formations, and this probably is applicable not just to Ecuadorian Amazonia, but to the whole western Amazon.

The Bujurquina fauna of Ecuador presently includes only two supposedly endemic species (B. zamorensis, B. pardus) plus based on some studies two to three Peruvian species, which we suggest is erroneous. Our results on the other hand identify in Ecuador at least 12 species, 11 of which are also identified in the molecular phylogeny. The total number of species identified in Bujurquina in this study is 32 based on morphology and 30 based on molecular phylogeny. Bujurquina, based on both molecular and morphological (coloration pattern) data, is divided into two main groups that have a geographical foundation, termed here the Northern and Southern groups, which based on molecular-clock dating diverged around 16.6 Ma. The large-scale phylogeny and biogeography of Bujurquina in the western Amazon are found to be in very good agreement with geological evolution of the western Amazon, and on the smaller scale of the Ecuadorian Amazon reveal a complex paleogeography reflecting river network reconfiguration. Foremost among the large-scale events in the geological evolution of the western Amazon are the major structural arches that determine its past and present river basin configuration. All of these found in the western Amazon (the Vaupes arch, Fitzcarrald arch, Michicola arch) are identified as major events in the evolution of Bujurquina, and their geological dates are found to agree with the reconstructed and independently calibrated biogeographic events in the Bujurquina phylogeny. Based on the biogeographic relationships and endemism of the newly identified species in Ecuador, we introduce another major geological feature of the western Amazon, the Guyusa platform, which is interlinked with the Vaupes arch in Colombia. While the present river network of western Amazonia has a clear W–E pattern of flow, our biogeographic reconstruction based on the molecular phylogeny of Bujurquina is dominated by an S–N pattern. Based on molecular-clock dating, the restructuring of the originally S–N trending river (the Paleo-Amazon) into a W–E river network (the modern Amazon) occurred after 8 Ma outside of Ecuador, and between 4.6 and 3.6 Ma within Ecuador. Establishment of the modern rivers and basins, including the largest ones in Ecuadorian territory (i.e. the Pastaza, Napo, Putumayo), occurred only very recently between 0.75 and 0.1 Ma. The 8–4.5 Ma date is considered a pivotal point in the history of evolution of the Amazon basin, since it marks the turning point between a long existing S–N flowing western paleo-Amazon to a W–E flowing modern Amazon. In our biogeographic reconstructions, this pivotal point is visible in the correspondence of biogeography and dating of Bujurquina to the emergence of the arches, and the Guayusa platform, that disrupted the S–N paleo-Amazon and transformed it into a W–E flowing modern Amazon.