Development of microsatellite loci and population genetics of the catfish Pimelodus yuma (Siluriformes: Pimelodidae)

Pimelodus yuma (formerly Pimelodus blochii) is a freshwater fish, endemic to the Colombian Magdalena-Cauca and Caribbean basins that experiences habitat disturbances resulting from anthropogenic activities. Due to the lack of information about the population genetics of this species, this study developed 14 species-specific microsatellite loci to assess the genetic diversity and population structure of samples from the lower section of the Cauca River. The studied species showed genetic diversity levels higher than the average values reported for Neotropical Siluriformes and significant inbreeding levels as was described for some congeners. Furthermore, P. yuma comprises two coexisting genetic groups that exhibit gene flow along the lower section of the Cauca River. This information constitutes a baseline for future monitoring of the genetic diversity and population structure in an anthropic influenced sector of the MagdalenaCauca basin.


INTRODUCTION
Seasonally migratory patterns may influence the gene flow or genetic structure in the natural populations of fishes of the family Pimelodidae. Since the gene flow is only possible if the fishes successfully reproduce once they have arrived at their new site (Freeland, 2020), both the migrations and reproductive cycles explain the gene flow in Pimelodus maculatus Lacepède, 1803 in the upper section of the Uruguay River, Tibagi River basin, or Tietê River (Almeida et al., 2001;Almeida et al., 2003;Ribolli et al., 2012) and in Pseudoplatystoma corruscans (Spix & Agassiz, 1829) in the São Francisco River (Dantas et al., 2013) and the Paraguay basins in Brazil (Prado et al., 2017).
Pimelodus yuma Villa-Navarro & Acero, 2017 is one of the 36 valid species of the genus Pimelodus Lacepède, 1803(Fricke et al., 2020 distributed in Cauca, Magdalena and Sinú River drainages in Colombia (Villa-Navarro et al., 2017). This species and their congener Pimelodus crypticus Villa-Navarro & Cala, 2017 were both considered as Pimelodus blochii Valenciennes, 1840 before 2017 (Villa-Navarro et al., 2017). Currently, there are no reports of any risk category for P. yuma in the Colombian red list of threatened freshwater fishes (Mojica et al., 2012) or in the red list of threatened species of the International Union for 3/15 ni.bio.br | scielo.br/ni Conservation of Nature, IUCN. Nevertheless, there is concern about the real status of the species considering the available information of fisheries in the Magdalena-Cauca basin.
Based on the bioeconomic model analysis, P. yuma is considered overexploited (Gutiérrez et al., 2011). Additionally, fisheries and the fish production modeling data show a decline in number of captures from 17,969 tons in 1975 to 8,370 tons in 2016, indicating the high fishing pressure exerted over this catfish (Barreto Reyes, 2017). Additionally, other factors that compromise this species are the pollution of water due to disposal materials, anthropogenic disturbance of hydrodynamics of the river (hydropower stations), deforestation, and introduction of alien fish species (Jiménez-Segura et al., 2016;Tognelli et al., 2016). Nonetheless, studies of population genetics of P. yuma remain absent, although they could provide information to estimate the degree of genetic vulnerability of this species.
The use of microsatellite loci permits a high-sensitivity evaluation of the population genetic diversity due to their polymorphism levels and wide distribution in the genome (Triantafyllidis et al., 2002). Additional to these advantages, considering the high mutation rates of microsatellite loci, they are also useful to examine recent events (Chistiakov et al., 2006). Despite of the available microsatellite loci for three species within the genus (Paiva, Kalapothakis, 2008;Restrepo-Escobar, Márquez, 2020;Savada et al., 2020) and some others within Pimelodidae (Batista et al., 2010;Carvalho, Beheregaray, 2011;Saulo-Machado et al., 2011;Souza et al., 2012;Prado et al., 2014), their use for studying the population genetics of P. yuma may be problematic. The limitations of cross-amplification of microsatellite loci, related to unsuccessful amplification in phylogenetically distant species (Barbará et al., 2007), include low levels of polymorphism, presence of null alleles (Rutkowski et al., 2009), allele size homoplasy (Estoup et al., 2002), and inability to evaluate properly orthologous loci (Yue et al., 2010).
Since Pimelodus yuma was catalogued as a medium-distance migrant (100-500 km) (Usma et al., 2013) and the sampling sites in this study are separated by less than 260 km with a geography lacking slopes, rapids, cascades, the existence of gene flow in P. yuma was hypothesized for three sites of the lower sections of the Cauca River. Furthermore, given the anthropogenic intervention in the Magdalena-Cauca basin and the decrease of captures, it is expected a loss in genetic diversity. To test these hypotheses, we developed a set of species-specific microsatellite loci that allow the study of the population genetics for the endemic Colombian catfish.

MATERIAL AND METHODS
This study analyzed 138 muscle tissues of Pimelodus yuma from the main channel and floodplain lakes in the lower sections of the Cauca River. Although a sampling effort using gillnet, and cast nets was made along the eight sections of the Cauca River described by Landínez-García, Márquez (2016), this species was found only in the area that corresponds to three sampling sections (S4/5, S6, and S7/8) (Fig. 1, S1). The S4/5 section consists of sites distributed along the main channel of the river and floodplain lakes, whereas S6 and S7/8 encompass floodplain lakes in the lower section of the Cauca River. The analyzed sections were selected according to the availability of samples and geographical distance. All the samples were collected by Integral S.A. from 2011 4/15 ni.bio.br | scielo.br/ni to 2014, before the construction of the Hydroelectric Ituango Project. The samples, preserved in ethanol 96%, were provided to the laboratory by Integral S.A. through two scientific cooperation agreements (September 19th, 2013;Grant CT-2013).
We followed the methodology described by Landínez-García, Marquez (2018) to develop the species-specific microsatellite loci for P. yuma. The DNA isolation was performed in one individual of P. yuma from S8, a site of the lower section of the Cauca River. The 454 GS FLX (+) (ROCHE) technology was used for pyrosequencing the previously obtained genomic library. Sequences reads were analyzed with the PRINSEQ-LITE v0.20.4 to assess the quality of the reads and to eliminate sequences with less than 100 bp in length. Then, PAL_FINDER v0.02.03 (Castoe et al., 2010) was used to identify and extract microsatellite loci with perfect di-, tri-, tetra-, penta-and hexa-nucleotide motifs. To design the primers, the flanking sequences of the microsatellite loci were analyzed using PRIMER3 v2.0 (Rozen, Skaletsky, 2000). An electronic PCR (primer-BLAST; available in https://ncbiinsights.ncbi.nlm.nih.gov/2017/06/28/e-pcr-is-retiring-useprimer-blast/) was carried out to evaluate the correct alignment of the selected primers.
Forty-one microsatellite loci were selected to evaluate their amplification capacity. To assess PCR conditions (Landínez-García, Márquez, 2016), DNA was extracted with the GeneJET Genomic DNA Purification Kit (ThermoScientific) following the manufacturer's instructions. The amplification capacity of the loci was tested using two random samples 5/15 ni.bio.br | scielo.br/ni and the products of the amplifications were separated through electrophoresis in 10% polyacrylamide gel setting 100 volts during 45 min in a Mini-PROTEAN Tetra Cell (Bio-Rad) and visualized with silver stain. Then, 28 loci that showed amplification capacity were tested using 12 random DNA samples from the three sections analyzed. Later, 22 loci were selected due to their level of polymorphism and band resolution within 100 and 350 bp in size. A set of 14 microsatellite loci that showed well-defined peaks and absence of stutter bands were considered for population genetic analysis.
The optimal conditions and concentrations for amplification reactions were as proposed by Landínez-García, Marquez (2018) with some modifications described below. In this three-primer strategy, 0.5 pmol/μL was used for each of the forward primer labeled in the 5' end with each one the adapters (Blacket et al., 2012), 1.0 pmol/μL of each reverse primer and 10-12 ng/μL of template DNA. Furthermore, the PCR amplifications were carried out with an initial denaturalization at 94°C for 3 min, followed by 45 cycles consisting of a denaturalization step of 90°C for 22 s and an annealing step of 56°C for 18 s without an extension step or final elongation. The products were separated by electrophoresis in an automatic sequencer ABI 3500 HD (Applied Biosystems), using the GeneScan™ 600 LIZ™ dye size standard (Applied Biosystems), and the allelic sizes were edited using GENEMAPPER ID-X v1.5. Potential genotyping errors were detected using Micro-Checker v.2.2.3 software (Van Oosterhout et al., 2004).
The above calculations were complemented with a Bayesian analysis in STRUCTURE v2.3.4 (Pritchard et al., 2000;Falush et al., 2002Falush et al., , 2007Hubisz et al., 2009). This analysis was performed with 200,000 Monte Carlo Markov Chain (MCMC) steps, 20,000 iterations as burn-in and setting the parameters admixture and non-admixture models, as well as correlated frequencies. The analysis was repeated 20 times per each simulated genetic group (K), which range between 1 and 6 groups. Finally, the webbased software STRUCTURESELECTOR (Li, Liu, 2018) was used to calculate the best estimated K, the estimators MEDMEANK, MAXMEANK, MEDMEDK, MAXMEDK (Puechmaille, 2016), and heuristic scores (Raj et al., 2014). The summarized results of the runs were plotted in a co-ancestry histogram of all the individuals using STRUCTURESELECTOR and the integrated software CLUMPAK (Kopelman et al., 2015). Finally, based on coancestry coefficients obtained with STRUCTURE and CLUMPP v1.1.1 (Jakobsson, Rosenberg, 2007), the individuals were rearranged by genetic stock and then were accordingly analyzed.

RESULTS
A total of 14 microsatellite loci were characterized in this study (Tab. 1); such loci did not show evidence of large allele dropout or stutter bands. Although two loci seemed to exhibit null alleles after Bonferroni correction, the departures from Hardy-Weinberg and heterozygote deficit, also observed in other loci in greater sample size, indicate a possible Wahlund effect (Tab. 2) as discussed below. The number of alleles per locus ranged between 8 and 23 with an average value of 13.071 alleles/locus. The allele sizes were in the expected range (100-350 bp) with a minimum size of 114 bp and a maximum of 347 bp. The polymorphic information content (PIC) showed values ranging from 0.508 to 0.937. Additionally, the observed and expected heterozygosities presented average values of 0.670 and 0.849, respectively. Five of the 14 developed microsatellite loci were excluded from the further genetic diversity analysis pending exploration in a greater sample size (Pyu6, Pyu7, Pyu11, Pyu22, and Pyu32).  The diversity estimators (Tab. 2) showed the lowest average expected heterozygosity in S4/5 (0.832) and the highest in S7/8 (0.873). The observed heterozygosity exhibited the lowest value in S6 (0.703) and the highest in S4/5 (0.744). Additionally, S4/5 (12.556 alleles/locus) showed the lowest average number of alleles per locus while S6 showed the highest average number (15.889 alleles/locus). Only two of the nine loci were departed from HWE in the three evaluated sites and the other were in the HWE in at least one site. Additionally, the across loci values of the inbreeding coefficients for all analyzed sites showed deficit of heterozygosity (0.109 < F IS < 0.186). These levels of genetic diversity and inbreeding were similar to those obtained for the genetic stocks (Stock1 and Stock2) supported by the genetic structure analysis.
Using the pairwise (S4/5-S6, S4/5-S7/8, S6-S7/8) standardized estimators F' ST (-0.014, 0.019, -0.004) and D EST (0.007, 0.011, 0.005) there was not statistical significance among the individuals of the different analyzed sites. Similarly, the AMOVA (F' ST = -0.003; P= 0.605) and the DAPC (Fig. 2) also showed that there were no differences among the three sections of the river. Hence, the individuals from the three sampling sites were genetically similar. The analysis in STRUCTURE showed two genetic stocks that coexist in the studied area (Fig. 3, S2; DeltaK = 2; MEDMEANK, MAXMEANK, MEDMEDK = 2; Mean LnP(K) = -5944.415 and -5930.270, for admixture and non-admixture models, respectively), without evidence of geographical structure. The analysis of diversity of the two stocks found in the study area evidenced that even without a geographical arrangement of the individuals, the heterozygosity deficit (0.000; 0.000) and F IS (0.157; 0.159) were statistically significant for both genetic stocks.

DISCUSSION
This study developed species-specific microsatellite loci to tests the hypothesis that Pimelodus yuma exhibits gene flow in three sites without geographical barriers, separated by less than 260 km in the middle and lower sections of the Cauca River. Among 41 microsatellite loci examined, 14 are considered highly informative showing values of PIC above the range proposed by Botstein et al. (1980). These loci exhibited values of average number alleles per locus (13.071) and expected heterozygosities (0.849) lower than those reported in the congener P. grosskopfii (Restrepo-Escobar, Márquez, 2020), but higher than those reported in P. microstoma Steindachner, 1877 (Savada et al., 2020) and were similar to those reported in P. maculatus (Paiva, Kalapothakis, 2008). Furthermore, these values were higher than those reported in  (Souza et al., 2012).
The genetic analysis of 138 samples of P. yuma showed high genetic diversity comparing values of alleles per locus (17.444) and expected heterozygosity (0.861) to those described in Neotropical Siluriformes (Na: 7.450, H E : 0.609; Hilsdorf, Hallerman, 2017). Additionally, P. yuma showed similar values to those reported for the congeners P. grosskopfii (Restrepo-Escobar, Márquez, 2020) and P. maculatus (Ribolli et al., 2012). In contrast, P. yuma showed higher diversity than the reported for some members of the family Pimelodidae such as Brachyplatystoma rousseauxii (Carvajal-Vallejos et al., 2014), P. corruscans (Vaini et al., 2016), and P. reticulatum (Prado et al., 2017). A high genetic diversity is a desirable characteristic for species, since it helps them to adapt to eventual adverse events (Li et al., 2017).
Despite the high genetic diversity, P. yuma exhibited significant deficit of observed heterozygosity that could be explained by several non-excluding alternatives. The presence of null alleles is invoked as a technical cause of the deficit of observed heterozygosity although the development of species-specific loci minimizes this possibility. As discussed below, P. yuma represents two genetic stocks that cohabitate the studied area, which may explain the deficit of observed heterozygosity (Wahlund effect). However, since the deficit of observed heterozygosity remained significant in the analysis by genetic stocks (not only by sampling site), other biological causes such as inbreeding may also explain this result considering the reduction in fishery production (Barreto Reyes, 2017), and the overexploitation exerted over the species since 2010 (Gutiérrez et al., 2011). Besides the fishing pressure, degradation of the Magdalena-Cauca basin caused by environmental, economic and demographic factors, also contribute to the observed decline in the number of catches of this fish species (Jiménez-Segura et al., 2016). Some Neotropical catfishes such as P. grosskopfii (Restrepo-Escobar et al., 2020), P. maculatus (Ribolli et al., 2012), and Pseudoplatystoma reticulatum (Abreu et al., 2009), have also shown evidences of inbreeding. Following Franklin (1980) and Soulé (1980), P. yuma requires a rigorous management of the populations since values above 10% indicate possible adverse effects over these populations.
According to the a priori expectation, this study showed evidences that the populations of P. yuma are not genetically structured in the middle and lower sections of the Cauca River, which is in line with the medium-range migration capacity of the species (100-500 km) and the absence of geographical barriers that facilitate the gene flow. This outcome was also found in some congeners like P. grosskopfii from the middle and lower sections of the Cauca River (Restrepo-Escobar et al., in press) and P. maculatus from upper section of the Uruguay River (Ribolli et al., 2012); although this species also showed genetic structure due to the multiple cascades in the rivers Uruguay and Paranapanema (Almeida et al., 2003;Ramella et al., 2006). In contrast to the study of P. grosskopfii (Restrepo-Escobar et al., in press), this work did not find individuals of P. yuma upstream from one of the slopes of the river, suggesting that this fish species is not capable to overcome this steeped geography.
Although P. yuma did not show genetic differences between the evaluated sections, it comprises two genetic stocks that cohabit. A similar genetic structure was found in its congener P. grosskopffi in the same area (Restrepo-Escobar, Márquez, 2020). The existence of those two genetic stocks might be related to spatial or seasonal reproductive isolation. Nevertheless, it remains to explore if the reproductive isolation by season can explain these differences since it was demonstrated that P. yuma (as P. blochii) from a floodplain lake in the medium section of the Magdalena River basin only reproduces during rainfall peaks, despite maintaining a constant condition factor at different environmental conditions (Lopez-Casas, Jimenez-Segura, 2007). Although our study was not designed to test the hypothesis of temporal genetic structure, the comparisons among available samples in the rain and drought periods per year (S3) did not show any tendency to seasonal or annual genetic structure.
This study led to stablish that the populations of P. yuma from the lower sections of the Cauca River do not exhibit spatial genetic structure, although there are two biological populations along the study area. Likewise, even though this species is affected by different anthropogenic activities, it was found that P. yuma shows a high genetic diversity, but with a deficit in the observed heterozygosity. The first 14 developed polymorphic microsatellite loci for this species are recommended for further studies of genetic diversity on population of P. yuma. Finally, these results are the baseline for future studies aiming to monitor the genetic diversity and structure of the populations of this endemic Colombian fish species in order to generate appropriate management plans.