Unveiling biogeographical patterns of the ichthyofauna in the Tuichi basin, a biodiversity hotspot in the Bolivian Amazon, using environmental DNA

bolivia-soil-classification; https://catalog.data.gov/dataset/south-america-geologic-map- geo6ag) and from Schneider et al., 2017 [55] (S3 Table). Non-metric multidimensional scaling (NMDS) using Jaccard’s binary distance was per- formed and correlations between environmental variables and NMDS ordination were measured. Considering that some variables are not independent, such as longitude and altitude, the objective of this analysis is therefore to identify environmental variables which are not neces- sarily explanatory but that can be correlated with the structuring of the ichthyofauna. Correla- tions between environmental variables and NMDS ordination were computed with the function “envit” in the R package Vegan. A similarity analysis (ANOSIM, package Vegan 2.5.6) was performed to test whether the biogeographic entities identified by the HAC were significantly different in the NMDS projection. Data analyses were completed in R version 3.6.2 [56]. All command lines and scripts used for bioinformatic treatment, statistical analyses and carrying out graphs are available in S1 Text. Results Sequencing, assignation and controls results NGS sequencing and taxonomic assignment results. Overall for the 200 libraries, which included 154 samples, two controls, one alien sample, 42 negatives (6 PCR-blank, 22 negatives field controls, 14 negatives DNA extraction), sequencing produced 551.8 million raw reads (82 Gb). The mean number of reads per site was 14.7 million (SE = 12.6), see S4 Table for details.

Five libraries with abundance below 50 sequences assigned to fish have been excluded from the analyses: they included all libraries (3) from one site on the Tuichi (not reported), and one library each for the Chalala´n lake and the Tuichi site 14. Rarefaction curves done at the library and site level (S1 Fig) show that sequencing effort was enough for the remaining libraries. After cleaning reads, 157.6 million could be joined and were clustered in 210,035 OTUs. Of these, 5,246 OTUs made up from 12.59 million sequences were assigned to a fish COI refer- ence (median value per site = 170,742), 93.6% at the species level. OTUs (n = 8,502) with a match with a COI reference in the database, but with scores or percentages of identity values below the thresholds were discarded. Finally, 196,287 OTUs had no hit. After taxonomic assignment, 120599, 1418, 54 and 7 reads were aligned to Marcusenius monteiri (Bold: AMNHI782-12), Epinephelus marginatus (GenBank: JX124778.1), Zoarces gillii (GenBank: HM180941) and Sebastiscus marmoratus (GenBank: HM180872), respectively. As they correspond to marine or non-Amazonian taxa, we verified these 4 sequences by blasting them under BOLD / NCBI. The results clearly show that the names of these four accession numbers are erroneous as their sequences show high blast score and % identity with bacterial or human origin. They were, therefore, discarded. One hundred-two taxa were shared by the two primer sets (MK1 and MK2) and 90 taxa were unique to a single primer set (73 for MK1 and 17 for MK2). Although barcode MK1 is

smaller in size than MK2 (185 and 285 bp, respectively) their taxonomic resolutions were simi- lar, with more than 90.6% and 92.7% of the OTUs assigned at the species level. The percentage of identity between the OTU and the reference sequences of the database to which they are assigned was on average 98.71% (SE = 1.12). Taxa abundance per library are reported in S5 Table. Controls by negative, positive and alien samples. Of the 42 negative control libraries only 213 reads were assigned to fish taxa indicating a low level of contamination. Of 40 species in the mock control samples, 39 have been recovered by our primers, which validates their PLOS ONE Fish diversity and distribution in the Tuichi basin PLOS ONE | https://doi.org/10.1371/journal.pone.0262357 January 4, 2022 8 / 25