FACULTAD DE CIENCIAS VETERINARIAS Y BIOLÓGICAS CARRERA PROFESIONAL DE BIOLOGÍA MARINA “Diversidad genética de la anchoveta Engraulis ringens a lo largo de la costa peruana basada en secuencias de ADN mitocondrial” Tesis para optar el título profesional de: BIÓLOGO MARINO Presentado por: Deivis Alexander Cueva Shapiama (0000-0002-1937-7452) Asesor: Giovanna Elizabeth Sotil Caycho (0000-0002-1583-8821) Lima – Perú 2021 FACULTAD DE CIENCIAS VETERINARIAS Y BIOLÓGICAS CARRERA PROFESIONAL DE BIOLOGÍA MARINA “Diversidad genética de la anchoveta Engraulis ringens a lo largo de la costa peruana basada en secuencias de ADN mitocondrial” Tesis para optar el título profesional de: BIÓLOGO MARINO Presentado por: Deivis Alexander Cueva Shapiama (0000-0002-1937-7452) Asesor: Giovanna Elizabeth Sotil Caycho (0000-0002-1583-8821) Lima – Perú 2021 ACTA DE SUSTENTACIÓN DE TESIS Lima, 07 de mayo de 2020. Los integrantes del Jurado de tesis: Presidente: Dr. PEPE ANTONIO ESPINOZA SILVERA Miembro: Dr. HECTOR ALONSO APONTE UBILLUS Miembro: Dr. ROGER IZIGA GOICOCHEA Se reúnen para evaluar la tesis titulada: “GENETIC DIVERSITY OF THE ANCHOVY ENGRAULIS RINGENS ALONG PERUVIAN COAST BASED ON MITOCHONDRIAL DNA SEQUENCES” Presentada por el estudiante/bachiller: ● DEIVIS ALEXANDER CUEVA SHAPIAMA Para optar el Título Profesional de Biólogo Marino Asesorado por: Mg. GIOVANNA ELIZABETH SOTIL CAYCHO Luego de haber evaluado el informe final de tesis y evaluado el desempeño de los estudiantes de la carrera de Biología Marina en la sustentación, concluyen de manera unánime ( X ) por mayoría simple ( ) calificar a: Tesista: DEIVIS ALEXANDER CUEVA SHAPIAMA Nota (en letras): Dieciocho Aprobado ( ) Aprobado - Muy buena ( X ) Aprobado - Sobresaliente ( ) Desaprobado ( ) Los miembros del jurado firman en señal de conformidad. Dr. PEPE ANTONIO ESPINOZA SILVERA Mg. GIOVANNA ELIZABETH SOTIL CAYCHO Presidente Asesora Dr. HECTOR ALONSO APONTE UBILLUS Miembro 1 Dr. ROGER IZIGA GOICOCHEA Miembro 2 Índice CARÁTULA………………………………………………………………………………………………………………………1 Acta de sustentación………………………………………………………………………………………………………4 RESUMEN……………………………………………………………………………………………………………………….6 ABSTRACT………………………………………………………………………………………………………………………7 CONTENIDO EN ESPAÑOL………………………………………………………………………………………………8 Introducción…………………………………………………………………………………………………………………10 Materiales y Métodos…………………………………………………………………………………………………..12 Resultados……………………………………………………………………………………………………………………16 Discusión………………………………………………………………………………………………………………………20 Figuras y tablas……………………………………………………………………………………………………………..27 Material suplementario………………………………………………………………………………………………..37 CONTENIDO EN INGLÉS………………………………………………………………………………………………..52 Introduction…………………………………………………………………………………………………………………54 Materials and Methods…………………………………………………………………………………………………56 Results…………………………………………………………………………………………………………..……………..60 Discussion…………………………………………………………………………………………………………………….63 Figures and tables…………………………………………………………………………………………………………69 Supplementary material.………………………………………………………………………………………………79 REFERENCIAS BIBLIOGRÁFICAS…………………………………………………………………………………….94 Resumen La anchoveta Engraulis ringens exhibe niveles abrumadores de biomasa y experimenta oscilaciones demográficas drásticas. Dada la inestabilidad poblacional observada, es crucial explorar la variabilidad genética, especialmente por la importante relación propuesta por estudios teóricos y observaciones empíricas entre los niveles de variabilidad intraespecífica y el destino de las poblaciones locales. A corto plazo, esta relación tiene implicancias en la sostenibilidad de las poblaciones explotadas; a largo plazo, en la perpetuación de las especies. Para la anchoveta, una de las especies más estudiadas en la parte Norte del Sistema de la Corriente de Humboldt, la variabilidad genética poblacional apenas se ha abordado. En el siguiente trabajo se describe la diversidad genética mitocondrial de E. ringens mediante evaluaciones temporales y la comparación con estimaciones reportadas para otros peces engraúlidos que habitan en entornos oceanográficos similares (sistemas de surgencia costera). Los resultados sugieren que E. ringens muestra los niveles más bajos de diversidad genética en relación con sus congéneres. Esta señal genética ya se había descrito en individuos del extremo sur de su rango geográfico, ahora ha sido confirmada en la parte norte del rango, que es el área donde ocurren los mayores niveles de biomasa de este recurso pesquero (costa norte y centro de Perú). Teniendo en cuenta que las inferencias basadas en ADN mitocondrial podrían reflejar tanto el efecto de la dinámica demográfica histórica o de presiones deterministas como la selección natural, se enfatiza la necesidad de generar más estudios que aporten con información basada en otros tipos de marcadores y bajo un muestreo exhaustivo para verificar la consistencia espacial y temporal de la señal genética observada. Abstract The anchoveta Engraulis ringens exhibits overwhelming biomass and experiences drastic demographic oscillations. Given this population instability, it is crucial to explore genetic variability, especially for the important relationship proposed by theoretical studies and empirical observations between levels of intraspecific variability and the fate of local populations. In the short term, this relationship has consequences for the sustainability of harvested populations; in the long term, for the perpetuation of the species. For anchoveta, one of the most studied species in the Northern Humboldt Current System, population genetic variability has been scarcely addressed. Here, anchoveta's mitochondrial genetic diversity is described through several temporal assessments and the comparison with diversity estimates reported for other engraulid fish inhabiting similar oceanographic settings (coastal upwelling systems). Results suggest that E. ringens displays the lowest levels of genetic diversity in relation to their congeners. This genetic signature had already been described for individuals from the southern end of the species range; now the same signal has been confirmed to the northern of the geographic range; that is, the area where the highest biomass of this fishery resource occurs (northern and central coasts of Peru). Considering that mtDNA-based inferences could reflect both historical demographic dynamics or deterministic pressures such as natural selection, the need to generate more studies is emphasized in order to provide information based on other types of markers and under a comprehensive sampling to verify the spatial and temporal consistency of the observed genetic signal. Diversidad genética de la anchoveta Engraulis ringens a lo largo de la costa peruana basada en secuencias de ADN mitocondrial Deivis Alexander Cueva Shapiama. Laboratorio de Genética Molecular, Instituto del Mar del Perú, Esquina Gamarra y General Valle S/N Chucuito Callao. Correo electrónico: dcueva@imarpe.gob.pe Giovanna Sotil. Laboratorio de Genética Molecular, Instituto del Mar del Perú, Esquina Gamarra y General Valle S/N Chucuito Callao. Correo electrónico: gsotil@imarpe.gob.pe Correspondencia. Correo electrónico: gsotil@imarpe.gob.pe Título corto: Diversidad genética de la anchoveta peruana Número de palabras para el texto principal: 4465 Resumen La anchoveta Engraulis ringens exhibe niveles abrumadores de biomasa y experimenta oscilaciones demográficas drásticas. Dada la inestabilidad poblacional observada, es crucial explorar la variabilidad genética, especialmente por la importante relación propuesta por estudios teóricos y observaciones empíricas entre los niveles de variabilidad intraespecífica y el destino de las poblaciones locales. A corto plazo, esta relación tiene implicancias en la sostenibilidad de las poblaciones explotadas; a largo plazo, en la perpetuación de las especies. Para la anchoveta, una de las especies más estudiadas en la parte Norte del Sistema de la Corriente de Humboldt, la variabilidad genética poblacional apenas se ha abordado. En el siguiente trabajo se describe la diversidad genética mitocondrial de E. ringens mediante evaluaciones temporales y la comparación con estimaciones reportadas para otros peces engraúlidos que habitan en entornos oceanográficos similares (sistemas de surgencia costera). Los resultados sugieren que E. ringens muestra los niveles más bajos de diversidad genética en relación con sus congéneres. Esta señal genética ya se había descrito en individuos del extremo sur de su rango geográfico, ahora ha sido confirmada en la parte norte del rango, que es el área donde ocurren los mayores niveles de biomasa de este recurso pesquero (costa norte y centro de Perú). Teniendo en cuenta que las inferencias basadas en ADN mitocondrial podrían reflejar tanto el efecto de la dinámica demográfica histórica o de presiones deterministas como la selección natural, se enfatiza la necesidad de generar más estudios que aporten con información basada en otros tipos de marcadores y bajo un muestreo exhaustivo para verificar la consistencia espacial y temporal de la señal genética observada. Palabras clave: anchoveta, genética de poblaciones, diversidad genética, ADN mitocondrial, Sistema de la Corriente de Humboldt. Introducción En la búsqueda de un manejo sostenible de las poblaciones explotadas, la variación genética es un atributo importante para explorar y comprender porque permite que las poblaciones muestren flexibilidad y se adapten frente a escenarios que desafían su supervivencia y perpetuación (Frankham et al. 2002; Reed y Frankham 2003). La Anchoveta Engraulis ringens es una de las especies marinas más estudiadas en Perú, con varios aspectos biológicos bien registrados para respaldar su manejo pesquero. Sin embargo, el conocimiento sobre la genética poblacional de la especie es notablemente escaso y aún se desconoce cómo se distribuye su variabilidad genética a lo largo de su rango geográfico en el Sistema de la Corriente de Humboldt (HCS). Al principio del estudio de Anchoveta, Tsukayama (1966) y Rojas 1971 (como se cita en Pauly et al., 1987, pág. 5) informaron diferencias merísticas y morfométricas entre individuos del centro norte y sur de Perú. Además, Jordan y Málaga (1972) mostraron desplazamientos desiguales (mediante experimentos de captura-recaptura) de individuos a lo largo de la costa peruana. A partir de estas observaciones se sugirió el manejo de dos unidades poblacionales, una correspondiente al centro-norte de Perú y la segunda al sur de Perú-norte de Chile (IMARPE 1973). Posteriormente, Mathisen (1989) propuso una organización alternativa compuesta por varias unidades locales basada en un patrón de distribución irregular que observó en mapas de densidad. El autor hipotetizó que tal agregación especial podría estar respaldada por diferencias genéticas. Posteriormente, Mendo (1991) evaluó estas hipótesis sobre stocks de anchoveta y respaldó la existencia de subpoblaciones centro-norte de Perú y sur de Perú-norte de Chile, una conclusión que extrajo de los datos de captura-recaptura y variabilidad de aloenzimas, mientras que los otros tipos de datos que utilizó no mostraron una diferenciación espacial. En el trabajo de Mendo (1991) se evaluaron ocho alozimas pero solo tres loci resultaron informativos, lo que contrasta con el número habitual de marcadores utilizados cuando se explora la variación genética nuclear por este medio (Gyllensten 1985; Waples 1987; Hedgecock et al. 1994; Tudela et al. al. 1999). En cuanto a la anchoveta del sur del rango geográfico, Ferrada et al. (2002) evaluaron la diferenciación genética entre muestras del norte y sur de Chile mediante el uso de marcadores de ADN nuclear y electroforesis, pero no detectaron una señal de diferenciación, posiblemente debido al bajo polimorfismo de los marcadores genéticos utilizados (Ferrada-Fuentes et al. 2018). Recientemente, se han realizado esfuerzos para desarrollar marcadores microsatélites para el estudio poblacional de E. ringens, pero pocos loci se han desarrollado con éxito a partir de muestras de Perú (Rojas 2011) y otros, más prometedores a partir de muestras de Chile (Ferrada-Fuentes et al. 2018), aún no han sido reportados bajo un diseño de estudio de poblaciones que incluya muestras de Perú. De similar importancia que el enfoque de diferenciación genética es el estudio de la historia demográfica de E. ringens en el HCS, un ecosistema marino conocido por exhibir una alta variabilidad oceanográfica (Chavez y Messié 2009). Se ha notado que la distribución y biomasa de la anchoveta es variable en diferentes escalas de tiempo y aparentemente en respuesta a factores ambientales (Gutiérrez et al 2012; Chávez et al. 2003). La inferencia sobre la historia demográfica a partir de información genética podría revelar patrones demográficos pasados que ayuden a comprender la dinámica poblacional contemporánea. Por ejemplo, Grant y Bowen (1998) estimaron la diversidad genética mitocondrial en varias especies de anchoveta, incluida E. ringens, y propusieron una dinámica de metapoblación para explicar la “arquitectura genética poblacional superficial” detectada en la anchoveta. El estudio de Grant y Bowen (1998) se basó en muestras de hace casi 30 años del extremo sur del rango de la especie. Por lo tanto, una actualización y verificación de esta señal genética en la parte norte del área de distribución (centro-norte de Perú) es pertinente, especialmente porque los registros de biomasa más grandes ocurren en esta área (Jordan 1971; Pauly y Tsukayama 1987; Ñiquen et al. 2000). Aunque el estudio de la variabilidad genética de la población de anchoveta comenzó hace casi 30 años, la información generada todavía es escasa para generar inferencias poblacionales robustas porque los estudios se han basado en un número reducido de marcadores o un muestreo restringido en tiempo y/o espacio. Estas limitaciones resaltan la necesidad de incrementar los esfuerzos para cubrir este vacío de información. En este sentido, en el siguiente estudio se ha explorado la diversidad genética mitocondrial de E. ringens a lo largo de la costa peruana durante tres años con el fin de describir una señal genética consistente basada en réplicas temporales y selección de marcadores moleculares que nos permitan comparar la diversidad genética de la anchoveta con las reportadas para otras especies del género Engraulis que habitan en sistemas de surgencia costera en otras partes del mundo. Materiales y Métodos Sitios de colecta y material biológico Se colectaron un total de 445 especímenes adultos de E. ringens frente a diferentes lugares de la costa norte, centro y sur del Perú (figura 1). El muestreo incluyó ejemplares de los dos stocks pesqueros (norte-centro y sur) gestionados en la zona económica exclusiva (ZEE) del Perú. Con la excepción de la muestra AS1016 (Tacna, octubre del 2016), todos los muestreos fueron realizados por la Dirección General de Investigación de Recursos Pelágicos (DGIRP) del Instituto del Mar del Perú (IMARPE) a través de cruceros de investigación para la evaluación hidroacústica de recursos pelágicos desde el 2014 hasta el 2016. Las campañas de muestreo se realizaron durante los dos principales períodos de mayor actividad reproductiva habitualmente observados en el Perú (febrero-marzo y septiembre-octubre). Las muestras de Tacna provinieron de un varamiento natural de peces ocurrido en octubre de 2016. En total se obtuvieron 13 muestras (tabla 1), cada una compuesta de 27 a 41 ejemplares de E. ringens. Los peces colectados se congelaron a -20 ° C y luego se almacenaron para su posterior análisis. Una vez en el laboratorio, se cortó un pequeño trozo (0,5 x 0,5 cm3) de tejido muscular dorsal de cada individuo y se conservó en etanol al 96%. Los tejidos se mantuvieron a 4 ° C hasta su posterior procesamiento. Extracción de ADN, amplificación por PCR y secuenciación de marcadores mitocondriales La extracción de ADN se realize aplicando el método CTAB2X/cloroformo-alcohol isoamílico. Se amplificaron dos regiones del ADN mitocondrial (mtDNA) mediante PCR: el dominio hipervariable 1 (HVR1) de 520 pares de bases (pb) de la región de control con los iniciadores que Lee et al. (1995) propusieron para peces teleósteos, y un fragmento de 913 pb del gen citocromo b (Cytb) con los iniciadores utilizados por Silva et al. (2014) para la anchoveta de Europa Engraulis encrasicolus. Se utilizó el kit HotStart Taq Polymerase Master Mix (Qiagen) para realizar reacciones de PCR de 10 µL, manteniendo las concentraciones de solución buffer, enzima ADN polimerasa y desoxinucleótidos trifosfato (dNTPs) de acuerdo con las recomendaciones del fabricante, mientras que las concentraciones finales de MgCl2, iniciadores y ADN molde se ajustaron para cada marcador. Los productos de PCR de cada muestra se secuenciaron bidireccionalmente en un analizador genético ABI 3500 (Applied Biosystems). Revisión y edición de electroferogramas Los electroferogramas fueron revisados (verificación de un buen proceso de llamada de bases) y editados (recorte de extremos) con FinchTv 1.4.0 (Equipo de investigación de Geospiza). Las secuencias de nucleótidos de HVR1 y Cytb se alinearon de forma independiente utilizando el algoritmo ClustalW implementado en Bioedit 7.0.0 (Larkin et al. 2007). Para ambos marcadores se generó una secuencia consenso por muestra basada en secuencias de las hebras sentido y antisentido. Polimorfismos de ADN y estimaciones de diversidad genética Se utilizó el software DnaSP (Librado y Rozas 2009) para explorar el polimorfismo de ADN (número y tipos de mutación y sitios variables) y generar datos de haplotipos. Los índices de diversidad de nucleótidos (π) y de haplotipos (Hd) se calcularon utilizando el software Arlequin 3.5.1.2 (Excoffier y Lischer 2010). Además, se comparó la diversidad genética de E. ringens con estimaciones reportadas para otros miembros del género Engraulis que habitan en regiones de surgencia costera. Las especies nominales consideradas fueron Engraulis encrasicolus, Engraulis capensis, Engraulis eurystole, Engraulis japonicus, Engraulis australis, Engraulis mordax, Engraulis anchoita y Engraulis ringens. En los estudios recopilados la diversidad genética se calculó a partir de secuencias parciales y casi completas del gen citocromo b (de 400 a 1000 pb) y la Región Control (de 400 a 600 pb), y tamaños de muestra entre 15 y 30 individuos (tablas complementarias S4 y S5). Diferenciación y estructura genética La diferenciación genética se evaluó mediante comparaciones por pares del estadístico Φst y pruebas exactas de diferenciación, utilizando Arlequin 3.5.1.2 (Excoffier y Lischer 2010). Los posibles escenarios de estructura genética (espacial y temporal) se evaluaron mediante el análisis jerárquico de la varianza molecular utilizando la prueba AMOVA implementada en Arlequin 3.5.1.2 (Excoffier y Lischer 2010). La organización de muestras por zonas frente a Perú (Norte, Centro y Sur), stocks pesqueros (norte- centro de Perú y sur de Perú-norte de Chile), campañas de muestreo (meses) y años fueron considerados como los niveles jerárquicos de comparación de la variación dentro y entre niveles. Además, se estimó la genealogía de ambos marcadores mitocondriales mediante la construcción de redes de haplotipos utilizando el método Median-Joining (Bandelt et al. 1999) y el software PopArt (Leigh y Bryant 2015). Una red de haplotipos es una forma gráfica de mostrar las relaciones genealógicas entre variantes del mtDNA y la distribución de sus frecuencias en relación con un diseño de muestreo en particular. Pruebas de neutralidad y distribución Mismatch Se utilizó las pruebas D de Tajima (Tajima 1989) y Fs de Fu (Fu 1997) para evaluar la hipótesis nula de evolución neutral. Además, realizamos analisis de mismatch para cada marcador. EL análisis mismatch permite inferir posibles eventos demográficos que la población pudo experimentar en el pasado. Una distribución unimodal refleja un rápido crecimiento de la población a partir de un tamaño poblacional efectivo pequeño; mientras que, una distribución multimodal refleja la estabilidad de la población a través de un largo periodo (Slatkin y Hudson 1991, Rogers y Harpending 1992). Resultados Polimorfismos de ADN y diversidad genética en HVR1 y Cytb Se encontraron 54 mutaciones en HVR1 (520 pb) de 445 individuos. El mayor número de sustituciones consistió en transiciones (46 transiciones frente a 8 transversiones). El número total de sitios variables fue 52 (50 con dos variantes de nucleótidos por sitio y dos con tres variantes por sitio). El número de singletons y sitios informativos de parsimonia fue casi proporcional (28 y 24 respectivamente). Todos estos cambios a nivel de nucleótidos generaron 84 haplotipos, de los cuales 29 (35%) estaban presentes en más de una muestra (haplotipos compartidos), mientras que 55 (65%) solo ocurrieron en sitios y fechas específicas (haplotipos privados). Los haplotipos privados presentaron baja frecuencia (tabla complementaria S1). En la región parcial (913 pb) del gen Cytb, se encontraron 47 mutaciones en 236 individuos. Al igual que en HVR1, la mayor parte de sustituciones fueron transiciones (44 transiciones frente a 3 transversiones). Todos los sitios variables fueron polimórficos por presentar solo dos variantes de nucleótidos por sitio, por lo que su número fue igual al número total de mutaciones. A diferencia de HVR1, el número de singletons (35) fue mayor que el número de sitios informativos de parsimonia (12). Todos los sitios variables definieron 44 haplotipos, siendo 9 (20%) compartidos y 35 (80%) privados. La mayoría de estos haplotipos también mostraron una frecuencia extremadamente baja (Tabla complementaria S1). En general, se evidenció un mayor polimorfismo en HVR1 que en Cytb (Tabla complementaria S1). Por ejemplo, HVR1 mostró una mayor proporción de transversiones (un tipo de sustitución menos común que las transiciones), presencia de sitios variables con más de dos variantes y una menor proporción de singletons. Nivel de diversidad genética mitocondrial en E. ringens La diversidad total de haplotipos basada en HVR1 indicó una alta diversidad (Hd = 0,88), mientras que la estimación con Cytb mostró un nivel regular (Hd = 0,50). La diversidad total de nucleótidos siguió un patrón similar, siendo mayor en HVR1 (π = 0,0046) en comparación con Cytb (π = 0,0008). En general, ambas regiones mitocondriales indican una cantidad de regular a alta de variantes de mtDNA, la mayoría con baja frecuencia y niveles de divergencia respaldada por pocos pasos mutacionales (Tabla 2). Entre los diferentes sitios y fechas de muestreo, la diversidad de haplotipos en HVR1 se mantuvo alta, variando entre 0,82 y 0,93; mientras que los valores de diversidad de nucleótidos oscilaron entre 0,0040 y 0,0050. Con respecto a Cytb, la diversidad de haplotipos osciló entre 0,37 y 0,68; mientras que la diversidad de nucleótidos varió entre 0,0005 y 0,0013 (Tabla 2). No se detectó diferenciación significativa basada en las pruebas de permutación (10000 interacciones) en los estimados de diversidad genética entre sitios o fechas de colecta (Tabla complementaria S2 y S3). En general, la diversidad genética en ambos marcadores presentó rangos de variación estrechos sin señales de diferencias significativas. Incluso, al agrupar los índices de diversidad por campaña de muestreo o año resultó en una reducción mucho mayor de las diferencias (Tabla 2). Los niveles de diversidad genética de E. ringens calculados en este studio difieren notablemente de los reportados para las anchovetas del viejo mundo (OWA) y, en menor medida, de los valores reportados para la anchoveta de California (Figura 2). La anchoveta peruana mostró un Hd entre 0.37 y 0.68 mientras que las OWA en general mostraron valores superiores a 0.90. Las estimaciones de Hd para E. mordax también fueron más altas (Hd = 0,74-0,95) que las de E. ringens, aunque ligeramente más bajas que aquellas de las OWA. Las estimaciones para E. anchoita (Hd = 0,45) fueron muy similares a las de E. ringens; sin embargo, se basan en un solo muestreo. En términos de diversidad de nucleótidos, las OWA exhiben un amplio rango de variación que va desde 0,006 a 0,017; mientras que, el rango de π en E. mordax fue más estrecho (0,004-0,009). En contraste, E. ringens exhibe valores extremadamente bajos de π (~ 0.001) con una variabilidad mínima entre diferentes sitios del norte, centro y sur de Perú, y similar a las estimaciones calculadas a partir de la muestra de Talcahuano (Chile) (Tabla complementaria S4). Diferenciación y estructura genética poblacional Las comparaciones por pares del estadístico Φst indicaron niveles bajos de diferenciación genética, con rangos de 0,000-0,036 usando Cytb y 0,000-0,066 usando HVR1 (Tabla complementaria S6 y S8). En ambos marcadores, algunas comparaciones por pares mostraron una diferenciación baja pero estadísticamente significativa. Para HVR1, las diferencias significativas no fueron respaldadas por las pruebas exactas de diferenciación, mientras que para Cytb algunas diferencias genéticas están respaldadas por las dos pruebas (Tabla complementaria S7 y S9, respectivamente). En todos los escenarios de estructura genética evaluados mediante pruebas AMOVA, el porcentaje de variación siempre fue mayor dentro de poblaciones hipotéticas que entre ellas, sugiriendo ausencia de una estructura espacial (zonas frente a Perú: Norte, Centro y Sur; Φst = 0.001, P- valor = 0,292 / Poblaciones pesqueras; Φst = -0,008, valor P = 0,988) o temporal (campañas de muestreo; Φst = -0,007, valor P = 0,903 / año, Φst = -0,006, valor P = 0,982). Debido al resultado general del análisis de la varianza molecular, solo se muestran resultados de niveles jerárquicos de comparación de las siguientes organizaciones: 1) Norte, Centro y Sur a lo largo de la costa peruana; y 2) stocks pesqueros (Tabla 3 y 4). En general, las redes de haplotipos mostraron una distribución espacial y temporal aleatoria de las variantes de mtDNA. La distribución de frecuencia de los haplotipos más comunes tanto en HVR1 como en Cytb no mostró cambios drásticos en todas los sitios y fechas de colecta (Figura. 4). Genealogía de haplotipos, pruebas de neutralidad y análisis mismatch La topología de la red Median-Joining basada en HVR1 fue reticulada y compleja, y consistió de un haplotipo dominante (30,8%) con varios haplotipos satélites (normalmente 0,2%) a su alrededor, separados por uno o dos pasos mutacionales; y un grupo de cuatro haplotipos de frecuencia intermedia (5,4-9,4%) relacionados con el haplotipo principal a través de muchos haplotipos de baja frecuencia y otros hipotéticos. Por otro lado, la red basada en Cytb exhibió claramente una topología de tipo estrella, con un solo haplotipo dominante y varios haplotipos de satélite de baja frecuencia (Figure 4). Ambos marcadores muestraron una genealogía superficial de los linajes del mtDNA (Figura 4). Las pruebas de neutralidad generaron valores significativamente negativos del estadístico D (-2.6418, valor P <0.05) y Fs (-29.0335, valor P <0.01), rechazando la hipótesis de neutralidad (Tabla 5). Además, el análisis mismatch mostró distribuciones unimodales en la frecuencia de las diferencias por pares tanto en HVR1 como en Cytb (Figura 3), que se ajusta al patrón esperado en un escenario de expansión de la población. Discusión Nivel de diversidad genética mitocondrial en E. ringens La diversidad genética total estimada en este estudio fue similar a aquellas reportadas por Grant y Bowen (1998), cuya muestra provino de Talcahuano, Chile, el área más austral de la distribución de E. ringens. Esto indica una estimación consistente espacial y temporalmente de la diversidad genética de E. ringens, basada en regiones codificante y no codificante del mtDNA. Las comparaciones de la diversidad de nucleótidos y haplotipos entre los peces engraúlidos sugirieron que E. ringens es una de las especies con los niveles más bajos de diversidad genética entre estos pequeños peces pelágicos que habitan en los sistemas de surgencia costera. Grant y Bowen (1998) ya habían descrito esta señal genética para algunas anchovetas de Talcahuano, Chile. Ellos notaron que una coalescencia superficial (escasa divergencia entre haplotipos) de los principales linajes del mtDNA en E. ringens contrastaba con la profunda divergencia genética entre la anchoveta y sus taxones hermanos (E. anchoita de Argentina y E. mordax de California). Este hallazgo les llevó a descartar una colonización reciente de las poblaciones regionales de Argentina o California. En cambio, sugirieron que esta señal genética sería el resultado de eventos recurrentes de extinción y colonización a nivel intrarregional, invocando así una estructura de metapoblación. Esta posible explicación supone un alto dinamismo en la abundancia y distribución de la anchoveta y una organización poblacional compleja ¿Qué evidencia existe de esto? Es conocido que varias poblaciones de peces clupeidos se caracterizan por exhibir una biomasa abrumadora, pero sujeta a fluctuaciones drásticas (Lluch-Belda et al 1989; Schwartzlose et al. 1999). La anchoveta no es un caso diferente, ya que muestra una alta variabilidad demográfica en diferentes escalas de tiempo, incluyendo cambios interanuales (Gutiérrez et al 2012) e interdecadales (Chávez et al. 2003). Esta variabilidad tiene sentido en la interpretación de cambios demográficos históricos inferidos a partir de secuencias de ADNmt. En la mayoría de los estudios filogeográficos, los cambios climáticos a gran escala, como las glaciaciones durante el Pleistoceno, generalmente se sugieren como factores que han dado forma a la arquitectura genética de las poblaciones de peces (Díaz- Viloria et al.2012; Silva et al.2014; Silva et al., 2017; Delrieu-Trottin et al.2017). En la HCS, se sabe que otros procesos climáticos a escala oceánica y global pueden explicar los cambios drásticos a escala interdecadal en la biomasa de anchoveta (Chavez et al. 2003). Además, a corto plazo, los eventos de El Niño pueden reducir drásticamente la abundancia de anchoveta (Barber y Chavez 1983; Gutiérrez et al. 1999; Ñiquen y Bouchón 2004) y hay evidencia que sugiere que estos eventos han ocurrido episódicamente en los últimos años 17000 años (Rein et al. 2005). Por lo tanto, los determinantes ambientales que actúan en diferentes escalas de tiempo pueden inducir cambios demográficos radicales en E. ringens. Volviendo a la propuesta de Grant y Bowen (1998), la evidencia sobre la dinámica de la población de E. ringens respalda cambios recurrentes en la abundancia en lugar de una estabilidad demográfica, y la teoría predice que la inestabilidad demográfica puede causar una reducción en la diversidad genética (Pannell y Charlesworth 1999). Sin embargo, la existencia de una organización poblacional compleja (metapoblación intrarregional) aún no es clara en E. ringens. Es posible que la dinámica de la población de anchoveta no se ajuste estrictamente a un modelo clásico de extinción-colonización porque no se han reportado casos de extinción local, incluso en condiciones oceanográficas que desafían su supervivencia. Por ejemplo, pocos meses después de uno de los eventos de El Niño más fuertes en las últimas décadas (El Niño 97/98), ocurrió una recuperación sorprendentemente rápida de la población centro-norte de Perú, lo que sugiere que las anchovetas adultas pudieron sobrevivir a una reducción del habitat (Bertrand et al.2004). En lugar del exterminio local, una explicación más plausible describiría un drástico colapso de la población y posterior recuperación con individuos de otras partes del área de distribución, pero también con sobrevivientes de la misma área local. Por tanto, esta posible explicación sugiere importantes niveles de conectividad en todo el rango de especies. Esta dinámica está más en línea con la definición de metapoblación que Kritzer y Sale (2004) propusieron para las poblaciones marinas, donde se enfatiza el grado de conectividad y su relevancia en el destino de la subpoblación más que los eventos de extinción/colonización. Diferenciación y estructura genética de la población En el HCS hay tres pesquerías importantes de anchoveta: un componente norte, centro y sur (Figura 1). En la ZEE peruana se gestiona el componente norte (stock norte-centro de Perú) y parte del componente central (stock sur Perú-norte Chile). El tercer componente, menor en términos de biomasa, ocurre en el extremo sur del rango de distribución, en el mar jurisdiccional de Chile (Stock de Chile Centro-Sur) (Cubillos et al. 2007). En este trabajo se probaron diferentes escenarios de estructura poblacional, considerando varios niveles jerárquicos de comparación (temporal y espacial), y no se observó señal de estructuración en congruencia con los componentes pesqueros. Este resultado contrasta con el estudio de Mendo (1991), donde se describió una restricción en la dispersión entre las anchovetas del centro-norte y sur a partir de los datos de captura-recaptura. La restricción fue más pronunciada durante el invierno- primavera austral, cuando ocurre el período de mayor actividad reproductiva, lo que sugiere que el apareamiento aleatorio podría no ser siempre el caso para la anchoveta y conlleva a la idea de posibles diferencias genéticas. De hecho, Mendo (1991) también apoyó los resultados de diferenciación con datos sobre la variación genética mediante el uso de marcadores de aloenzimas. Aunque solo tres loci polimórficos fueron informativos, mostraron diferencias significativas en la frecuencia de los alelos entre las regiones norte-centro y sur. Las diferencias entre los resultados de mtDNA y aloenzimas no solo se deben a las diferencias en el tipo de herencia de ambos caracteres moleculares, sino que también podría estar involucrado algún tipo de dinámica de población. Grant y Leslie (1993) hicieron simulaciones sobre el efecto de los eventos periódicos de extinción- colonización sobre la variabilidad mitocondrial y de aloenzimas y surgió un patrón general: después de un número dado de generaciones, los eventos de extinción- colonización producen una erosión drástica en la variabilidad genética mitocondrial con valores de Fst muy cerca a cero; seguido de una lenta reducción de la variabilidad genética nuclear con valores de Fst manteniéndose alejado de cero. Grant y Leslie (1993) explicaron que en el contexto de metapoblación, los linajes mitocondriales son más propensos a la extinción en comparación con los nucleares debido a su naturaleza haploide y herencia uniparental (lo que reduce el tamaño efectivo de la población aproximadamente una cuarta parte con relación al del AND nuclear). Los resultados contrastantes entre el estudio de Mendo (1991) y este trabajo se ajustan a esta descripción. Por lo tanto, si la población de anchoveta experimenta este tipo de dinámica, entonces el mtDNA no es un marcador apropiado para hacer inferencias sólidas sobre la estructura genética de la población. La evidencia hasta ahora proporciona algunas pistas sobre una posible estructura de la población de E. ringens, pero los estudios aún fallan en 1) cubrir completamente el rango de la especie o 2) emplear un número apropiado de marcadores genéticos. Este aspecto se encuentra aún más encriptado por la naturaleza pelágica de la anchoveta y su implicación en la dinámica de dispersión y agregación espacial. Por ejemplo, las perturbaciones del hábitat de la anchoveta (por ejemplo, eventos de El Niño o eventos cálidos similares a El Niño) hacen que los peces se muevan a diferentes escalas espaciales, desde movimientos verticales hasta desplazamientos longitudinales/latitudinales más grandes (Gutiérrez et al. 2012). En esta situación de mezcla, si los migrantes permanecen hasta el período reproductivo, la diferenciación genética sería un desafío debido al efecto amortiguador de combinar individuos de distintas subpoblaciones. En resumen, debido a que el grado de conectividad entre las poblaciones centro-norte de Perú y el sur de Perú-norte de Chile no se ha abordado de manera exhaustiva (aquí o en estudios anteriores), no es possible concluir si la evidencia genética respalda una metapoblación o una población irregular de E. ringens (es decir, estructurado espacialmente pero mostrando una dinámica de población interna igualmente en todo su rango; Kritzer y Sale 2004) .Sin embargo, la comparación interespecífica de la diversidad genética muestra claramente un patrón de niveles bajos en la anchoveta, y podrían invocarse explicaciones complejas para dar cuenta de eso, como por ejemplo una dinámica de metapoblación. Por lo tanto, será necesario incluir información de otros tipos de marcadores moleculares (por ejemplo, microsatélites o los generados por secuenciación masiva de ADN) así como diseños de muestreo adecuados para verificar si la baja diversidad genética en el mtDNA de la anchoveta refleja ciertamente una Agradecimientos Agradecemos a la Dirección General de Investigación de Recursos Pelágicos (DGIRP) del Instituto del Mar del Perú (IMARPE) por ayudarnos con la colecta del material biológico. Conflicto de intereses Los autores declaran no tener ningún conflicto de intereses FIGURAS Y TABLAS Diversidad genética de la anchoveta Engraulis ringens a lo largo de la costa peruana basada en secuencias de ADN mitochondrial D Cueva1, G Sotil1 1Laboratorio de Genética Molecular, Instituto del Mar del Peru, Esquina Gamarra y General Valle S/N Chucuito E-mail: gsotil@imarpe.gob.pe Figura 1. Sitios de muestreo en el mar jurisdiccional de Perú. En el panel interno (lado izquierdo inferior) se muestran el rango geográfico de E. ringens (gris oscuro) y los principals stocks de pesca (gris claro) en el Sistema de la Corriente de Humboldt. Table 1. Características de las 13 muestras colectadas en este estudio entre 2014 y 2016. Para simplificar, el código de muestra se compone de las dos primeras letras de la región o localidad de muestreo, seguidas de dos números que indican el mes y finalmente dos números del año. Zona Stock pesquero Localidad Fecha Código de muestra latitudinal Norte – Centro Perú Norte – Perú Paita Febrero 2014 AP2014 Paita Febrero 2015 AP0215 Salaverry Agosto 2015 AN0815 Paita Abril 2016 AN0416 Octubre Sechura AN1016 2016 Centro – Perú Callao Febrero 2014 AC2014 Callao Febrero 2015 AC0215 Callao Agosto 2015 AC0815 Pisco Agosto 2015 AS0815 Callao Abril 2016 AC0416 Octubre Callao AC1016 2016 Sur Perú – Norte Chile Sur – Perú Ilo Febrero 2014 AS2014 Octubre Tacna AS1016 2016 Tabla 2. Diversidad genética mitocondrial de E. ringens a lo largo de la costa peruana. Un guión (-) indica muestras en las que no se generaron secuencias de ADNmt. Para obtener una descripción general de la diversidad genética, las estimaciones también se calcularon mediante la combinación de secuencias de diferentes sitios de la misma campaña de muestreo (mes), año y como valor total. n Hd (± s.d.) π (± s.d.) Muestra HVR1 Cytb HVR1 Cytb HVR1 Cytb AP2014 33 - 0.85 ± 0.05 - 0.0046 ± 0.0005 - AP0215 39 30 0.87 ± 0.04 0.68 ± 0.09 0.0049 ± 0.0005 0.0013 ± 0.0003 AN0815 32 30 0.88 ± 0.04 0.51 ± 0.11 0.0046 ± 0.0004 0.0006 ± 0.0002 AN0416 39 30 0.84 ± 0.05 0.51 ± 0.11 0.0042 ± 0.0004 0.0009 ± 0.0003 AN1016 37 30 0.89 ± 0.04 0.51 ± 0.11 0.0044 ± 0.0003 0.0008 ± 0.0002 AC2014 27 - 0.92 ± 0.03 - 0.0045 ± 0.0005 - AC0215 32 28 0.86 ± 0.04 0.50 ± 0.12 0.0040 ± 0.0005 0.0007 ± 0.0002 AC0815 32 30 0.84 ± 0.07 0.37 ± 3.00 0.0047 ± 0.0005 0.0006 ± 0.0002 AC0416 41 28 0.90 ± 0.04 0.38 ± 0.11 0.0049 ± 0.0004 0.0005 ± 0.0002 AC1016 31 30 0.93 ± 0.03 0.52 ± 0.11 0.0047 ± 0.0004 0.0009 ± 0.0003 AS0815 40 - 0.92 ± 0.03 - 0.0050 ± 0.0003 - AS2014 30 - 0.90 ± 0.04 - 0.0043 ± 0.0004 - AS1016 32 - 0.82 ± 0.06 - 0.0045 ± 0.0003 - Febrero 2014 90 - 0.89 ± 0.02 - 0.0045 ± 0.0003 - Febrero 2015 71 58 0.87 ± 0.03 0.59 ± 0.08 0.0045 ± 0.0003 0.0010 ± 0.0002 Agosto 2015 104 60 0.89 ± 0.03 0.44 ± 0.08 0.0048 ± 0.0002 0.0006 ± 0.0001 Abril 2016 80 58 0.86 ± 0.03 0.45 ± 0.08 0.0045 ± 0.0003 0.0007 ± 0.0002 Octubre 2016 100 60 0.88 ± 0.02 0.51 ± 0.08 0.0045 ± 0.0002 0.0009 ± 0.0002 2015 175 118 0.88 ± 0.02 0.52 ± 0.06 0.0047 ± 0.0002 0.0008 ± 0.0001 2016 180 118 0.87 ± 0.02 0.48 ± 0.06 0.0045 ± 0.0002 0.0008 ± 0.0001 Total 445 236 0.88 ± 0.01 0.50 ± 0.04 0.0046 ± 0.0001 0.0008 ± 0.0001 0.020 E. ringens A 0.018 E. anchoita Anchovetas del nuevo mundo 0.016 E. mordax E. japonicus 0.014 E. australis 0.012 E. capensis Anchovetas del viejo mundo 0.010 E. encrasicolus E. eurystole 0.008 0.006 0.004 0.002 0.000 0.0 0.2 0.4 0.6 0.8 1.0 1.2 Diversidad de haplotipos (Hd) 0.040 B 0.035 E. ringens Anchoveta del nuevo mundo 0.030 E. japonicus Anchovetas del viejo mundo 0.025 E. encrasicolus 0.020 0.015 0.010 0.005 0.000 0.0 0.2 0.4 0.6 0.8 1.0 1.2 Diversidad de haplotipos (Hd) Figura 2. Diagrama de dispersión de la diversidad genética mitocondrial de especies marinas del género Engraulis que habitan en zonas de surgencia costera. Las estimaciones se calcularon a partir del gen Cytb (A) y las secuencias de la Región de Control (B). Consulte el material complementario para obtener detalles sobre la literatura utilizada. Diversidad de nucleótidos (π) Diversidad de nucleótidos (π) Tabla 3. Análisis de la varianza molecular (AMOVA) empleando secuencias de HVR1. Estructura genética evaluada: Norte, Centro y Sur del Perú. Pruebas de significancia basadas en 20022 permutaciones. Porcentaje Fuente de Grados de Suma de Componentes de variación Estadísticos Valor p variación libertad cuadrados de la varianza de los AMOVA componentes Entre las zonas Fct = 2 3.160 0.00200 0.17 0.21701 norte, 0.00168 centro y sur Entre las Fsc = muestras 10 13.116 0.00365 0.31 0.28673 0.00307 dentro de cada zona Dentro de Fst = cada 432 512.453 1.18623 99.53 0.21881 0.00475 muestra Total 444 528.730 1.19189 Tabla 4. Análisis de la varianza molecular (AMOVA) empleando secuencias de HVR1. Estructura genética evaluada: stock centro-norte Perú versus sur Perú-norte Chile. Pruebas de significancia basadas en 20022 permutaciones. Porcentaje Fuente de Grados de Suma de Componentes de variación Estadísticos Valor p variación libertad cuadrados de la varianza de los AMOVA componentes Fct = - Entre stocks 1 0.403 -0.00956 -0.81 0.98756 0.00807 pesqueros Entre muestras Fsc = 11 15.874 0.00746 0.63 0.15133 dentro de 0.00625 cada stock pesquero Dentro de Fst = - cada 432 512.453 1.18623 100.18 0.21881 .00177 muestra Total 444 528.730 1.18414 Tabla 5. Pruebas de neutralidad para las secuencias de HVR1 y Cytb. Los valores en negrita indican desviaciones significativas del modelo neutral (Nivel de significancia = 0.05). Tajima's D Fu's FS Muestra HVR1 Cytb HVR1 Cytb AP2014 -1.1744 - -9.4778* - AC2014 -0.3300 - -7.9154* - AS2014 -0.8683 - -9.1784* - AP0215 -0.9123 -2.1456* -10.7640* -7.6703* AC0215 -0.7352 -2.2798* -7.6521* -8.0920* AN0815 -0.8244 -1.9744* -6.5902* -6.5545* AC0815 -0.6032 -2.2355* -15.1188* -5.2786* AS0815 -1.0151 - -12.0377* - AN0416 -1.2314 -1.9535* -9.1137* -6.2841* AC0416 -1.2818 -1.7197 -10.3126* -2.6376* AN1016 -0.8442 -1.8828* -10.6453* -3.7452* AC1016 -0.6189 -2.2964* -9.7874* -8.3947* AS1016 0.1148 - -5.5473* - Total -1.94513* -2.6418* -26.37530* -29.0335* Figura 3. Análisis de la distribución de diferencias por pares (distribución mismatch) bajo escenarios de una historia poblacional constante y en crecimiento-declive. Figura 4. Redes de haplotipos de Cytb (A) y HVR1 (B) construídas con el algoritmo Median-Joining. MATERIAL SUPLEMENTARIO Diversidad genética de la anchoveta Engraulis ringens a lo largo de la costa peruana basada en secuencias de ADN mitochondrial D Cueva1, G Sotil1 1Laboratorio de Genética Molecular, Instituto del Mar del Peru, Esquina Gamarra y General Valle S/N Chucuito E-mail: gsotil@imarpe.gob.pe Tabla S1. Número y tipo de variaciones en las secuencias del mtDNA en Engraulis ringens. Marcador mitocondrial HVR1 Cytb Tamaño de muestra (n) 445 236 Número total de sitios de nucleótidos (tamaño de locus) 520 913 Número total de mutaciones (Eta) 54 47 Nro. de INDELs (inserciones-deleciones) 0 0 Nro. de sustituciones 54 47 Nro. de transiciones 46 44 Nro. de transversiones 8 3 Número total de sitios variables (frecuencia relativa) 52 (10%) 47 (5%) Nro. de singletons (frecuencia relativa) 28 (54%) 35 (74%) Nro. de singletons con dos variantes 27 35 Nro. de singletons con tres variantes 1 0 Nro. de sitios parsimonia-informativos (frecuencia relativa) 24 (46%) 12 (26%) Nro. de sitios parsimonia-informativos con dos variantes 23 12 Nro. de sitios parsimonia-informativos con tres variantes 1 0 Número total de haplotipos 84 44 Nro. de haplotipos privados (frecuencia relativa) 55 (65%) 35 (80%) Nro. de haplotipos compartidos (frecuencia relativa) 29 (35%) 9 (20%) Tabla S2a. Pruebas de permutación para diferencias significativas en la diversidad de nucleótidos (π) entre poblaciones basadas en datos de secuencia de HVR1. El procedimiento fue tomado de Alexander et al. (2016). El script R está disponible en https://github.com/laninsky/genetic_diversity_diffs 0.00 π 0.0046 0.0045 0.0043 0.0049 0.0040 0.0046 0.0047 0.0050 0.0042 0.0049 0.0044 0.0047 45 π diff / AN081 AN041 AN101 AS1 AP2014 AC2014 AS2014 AP0215 AC0215 AC0815 AS0815 AC0416 AC1016 p-values 5 6 6 016 0.84 AP2014 -- 0.883 0.609 0.646 0.394 0.937 0.902 0.560 0.538 0.584 0.787 0.888 5 0.0001 0.96 AC2014 -- 0.758 0.519 0.539 0.798 0.802 0.453 0.660 0.513 0.881 0.774 00176 8 0.0003 0.0002 0.76 AS2014 -- 0.307 0.698 0.580 0.543 0.277 0.910 0.286 0.871 0.549 06899 06723 8 0.0002 0.0003 0.0006 0.47 AP0215 -- 0.158 0.692 0.719 0.905 0.225 0.966 0.393 0.734 97185 97361 04084 6 0.0005 0.0004 0.0002 0.0008 0.53 AC0215 -- 0.353 0.306 0.141 0.795 0.155 0.562 0.364 45741 45565 38843 42926 0 0.0000 0.0001 0.0003 0.0002 0.0006 0.76 AN0815 -- 0.970 0.595 0.463 0.642 0.685 0.968 55220 55396 62119 41965 00962 6 0.0000 0.0001 0.0003 0.0002 0.0006 0.0000 0.75 AC0815 -- 0.618 0.456 0.674 0.665 0.977 74606 74782 81505 22579 20347 19386 9 0.0003 0.0004 0.0006 0.0000 0.0009 0.0003 0.0002 0.42 AS0815 -- 0.185 0.932 0.314 0.670 68814 68990 75712 71629 14555 13594 94208 2 0.0003 0.0002 0.0000 0.0006 0.0001 0.0004 0.0004 0.0007 0.66 AN0416 -- 0.202 0.726 0.441 82769 82593 75871 79954 62972 37990 57375 51583 8 0.0003 0.0004 0.0006 0.0000 0.0008 0.0002 0.0002 0.0000 0.0007 0.45 AC0416 -- 0.365 0.700 19624 19800 26523 22439 65366 64404 45018 49189 02394 2 0.0001 0.0000 0.0001 0.0004 0.0003 0.0002 0.0002 0.0005 0.0001 0.0005 0.91 AN1016 -- 0.650 89985 89809 16914 87170 55757 45205 64591 58798 92785 09609 4 0.0000 0.0001 0.0003 0.0002 0.0006 0.0000 0.0000 0.0002 0.0004 0.0002 2.8139 0.74 AC1016 -- 91407 91583 98306 05778 37148 36187 16801 77407 74176 28217 E-04 0 0.0001 0.0000 0.0001 0.0004 0.0004 0.0001 0.0001 0.0004 0.0002 0.0004 6.6855 0.0002 AS1016 -- 23130 22954 83769 20315 22612 78350 97736 91943 59640 42754 E-05 14537 } Tabla S2b. Pruebas de permutación para diferencias significativas en la diversidad de haplotipos (Hd) entre poblaciones basadas en datos de secuencia de HVR1. El procedimiento fue tomado de Alexander et al. (2016). El script R está disponible en https://github.com/laninsky/genetic_diversity_diffs 0.8 Hd 0.850 0.923 0.903 0.873 0.859 0.883 0.841 0.919 0.837 0.895 0.892 0.931 21 Hd diff AP201 AC201 AS201 AP021 AC021 AN081 AC081 AS081 AN041 AC041 AN10 AC101 AS1 /P- 4 4 4 5 5 5 5 5 6 6 16 6 016 values 0.6 AP2014 -- 0.257 0.409 0.717 0.877 0.584 0.871 0.265 0.811 0.46 0.484 0.21 21 0.0726 0.1 AC2014 -- 0.759 0.436 0.344 0.539 0.232 0.942 0.192 0.642 0.60 0.902 98135 42 0.0530 0.0196 0.1 AS2014 -- 0.614 0.474 0.733 0.307 0.781 0.286 0.897 0.829 0.675 69488 28647 91 0.0227 0.0499 0.0303 0.3 AP0215 -- 0.794 0.867 0.584 0.389 0.543 0.704 0.729 0.344 65612 32524 03876 77 0.0084 0.0642 0.0445 0.0142 0.5 AC0215 -- 0.699 0.75 0.303 0.711 0.542 0.571 0.268 9218 05955 77308 73432 36 0.0326 0.0400 0.0203 0.0099 0.0241 0.3 AN0815 -- 0.485 0.513 0.425 0.817 0.894 0.430 85728 12407 8376 20117 93548 21 0.0096 0.0823 0.0627 0.0324 0.0181 0.0423 0.7 AC0815 -- 0.198 0.942 0.365 0.379 0.160 52981 51117 22469 18593 45161 3871 15 0.0688 0.0038 0.0157 0.0460 0.0603 0.0361 0.0785 0.1 AS0815 -- 0.16 0.648 0.623 0.843 51981 46154 82493 8637 59801 66253 04963 17 0.0136 0.0863 0.0667 0.0364 0.0221 0.0463 0.0040 0.0825 0.7 AN0416 -- 0.292 0.329 0.139 71635 69771 41123 37247 63815 57364 18654 23617 76 0.0447 0.0279 0.0083 0.0219 0.0362 0.0120 0.0543 0.0241 0.0584 0.2 AC0416 -- 0.952 0.531 43163 54972 26325 77552 50983 57435 96145 08818 14799 08 0.0415 0.0311 0.0115 0.0187 0.0330 0.0088 0.0511 0.0273 0.0551 0.0032 0.2 AN1016 -- 0.485 13104 85031 56384 47492 20924 27376 66085 38877 84739 30059 36 0.0808 0.0081 0.0277 0.0580 0.0723 0.0481 0.0904 0.0119 0.0944 0.0360 0.039 0.0 AC1016 -- 04008 05873 3452 38396 11828 1828 56989 52026 75643 60844 2909 87 0.0298 0.1025 0.0828 0.0525 0.0383 0.0201 0.0986 0.0161 0.0745 0.071 0.110 AS1016 0.0625 -- 14272 12407 8376 79883 06452 6129 66253 42636 57435 32738 61828 Tabla S3a. Pruebas de permutación para detectar diferencias significativas en la diversidad de nucleótidos (π) entre poblaciones según los datos de secuencia de Ctyb. El procedimiento fue tomado de Alexander et al. (2016). El script R está disponible en https://github.com/laninsky/genetic_diversity_diffs π 0.0013 0.0007 0.0006 0.0006 0.0009 0.0005 0.0008 0.0009 π diff / p-values AP0215 AC0215 AN0815 AC0815 AN0416 AC0416 AN1016 AC1016 AP0215 -- 0.13 0.08 0.06 0.35 0.05 0.24 0.27 AC0215 0.00055 -- 0.86 0.73 0.54 0.62 0.70 0.65 AN0815 0.00062 6.205E-05 -- 0.89 0.45 0.77 0.59 0.50 AC0815 0.00067 1.200E-04 5.791E-05 -- 0.35 0.89 0.47 0.42 AN0416 0.00033 2.200E-04 2.820E-04 0.00034 -- 0.28 0.80 0.90 AC0416 0.00072 1.681E-04 1.060E-04 0.00005 0.00039 -- 0.40 0.34 AN1016 0.00042 1.369E-04 1.989E-04 0.00026 0.00008 0.00030 -- 0.92 AC1016 0.00039 1.671E-04 2.291E-04 0.00029 0.00005 0.00034 3.021E-05 -- Tabla 3a. Pruebas de permutación para detectar diferencias significativas en la diversidad de haplotipos (Hd) entre poblaciones según los datos de secuencia de Ctyb. El procedimiento fue tomado de Alexander et al. (2016). El script R está disponible en https://github.com/laninsky/genetic_diversity_diffs Hd 0.68 0.50 0.51 0.37 0.51 0.38 0.51 0.52 Hd/P-values AP0215 AC0215 AN0815 AC0815 AN0416 AC0416 AN1016 AC1016 AP0215 -- 0.28 0.30 0.06 0.32 0.08 0.31 0.33 AC0215 0.1808 -- 0.98 0.41 0.96 0.49 0.99 0.94 AN0815 0.1678 0.0130 -- 0.38 0.92 0.45 0.94 0.90 AC0815 0.3126 0.1318 0.1448 -- 0.36 0.95 0.40 0.35 AN0416 0.1632 0.0176 0.0046 0.1494 -- 0.42 0.90 0.95 AC0416 0.2946 0.1138 0.1267 0.0181 0.1313 -- 0.47 0.42 AN1016 0.1701 0.0107 0.0023 0.1425 0.0069 0.1244 -- 0.89 AC1016 0.1609 0.0199 0.0069 0.1517 0.0023 0.1336 0.0092 -- Tabla S4. Características sobre las estimaciones compiladas de diversidad genética de diferentes miembros marinos de Engraulis. Las estimaciones se calcularon a partir de secuencias de nucleótidos del gen mitocondrial Cytb. Referencia Especie Localidad de muestreo Fecha de muestreo pb n Hd π Grant and Bowen 1998 E. anchoita Not shown Not shown 521 19 0.4 0.001 Grant et al 2010 E. anchoita Argentina Not shown 687 19 0.5 0.001 Grant and Bowen 1998 E. australis Not shown Not shown 521 16 0.9 0.009 Silva et al 2017 E. australis Australia 2008 1044 34 1.0 0.007 Silva et al 2017 E. australis Nueva Zelanda 2005 1044 35 0.7 0.001 Grant and Bowen 1998 E. capensis Southern Africa* Not shown 521 18 0.2 0.004 Silva et al 2017 E. capensis Namibia 2007 1044 24 0.9 0.010 Silva et al 2017 E. capensis Sudafrica 2007 1044 13 0.9 0.019 Grant and Bowen 1998 E. encrasicolus Not shown Not shown 521 16 0.9 0.015 Grant et al 2010 E. encrasicolus Europe Not shown 687 17 0.9 0.014 Borrel et al 2012 E. encrasicolus Bay of Biscay - CA Not shown 524 39 0.7 0.012 Borrel et al 2012 E. encrasicolus Bay of Biscay - LI Not shown 524 20 0.9 0.016 Borrel et al 2012 E. encrasicolus Bay of Biscay - BB Not shown 524 17 1.0 0.016 Borrel et al 2012 E. encrasicolus Bay of Biscay - GE Not shown 524 12 0.8 0.013 Borrel et al 2012 E. encrasicolus Bay of Biscay - FR Not shown 524 31 0.9 0.016 Borrel et al 2012 E. encrasicolus Mediterranean - Adriatic Sea Not shown 524 20 1.0 0.014 Silva et al 2014 E. encrasicolus Málaga, ML Not shown 1044 31 1.0 0.007 Silva et al 2014 E. encrasicolus Tangier, TG Not shown 1044 38 0.4 0.003 Silva et al 2017 E. encrasicolus Noruega 2007 1044 24 1.0 0.009 Silva et al 2017 E. encrasicolus Polonia 2008 1044 9 0.9 0.014 Silva et al 2017 E. encrasicolus El Canal Inglés 2007 1044 27 1.0 0.010 Silva et al 2017 E. encrasicolus El Golfo de Vizcaya 2007 1044 23 1.0 0.015 Silva et al 2017 E. encrasicolus Norte de Portugal 1998 1044 25 1.0 0.013 Silva et al 2017 E. encrasicolus Sur de Portugal 2007 1044 29 1.0 0.011 Silva et al 2017 E. encrasicolus Islas Canarias 1999 1044 24 1.0 0.007 Silva et al 2017 E. encrasicolus Senegal 1999 1044 25 1.0 0.006 Silva et al 2017 E. encrasicolus Guinea-Bisáu 2006 1044 20 1.0 0.006 Silva et al 2017 E. encrasicolus Ghana 2008 1044 25 1.0 0.006 Silva et al 2017 E. eurystole US 2006 1044 12 0.9 0.004 Grant and Bowen 1998 E. japonicus Not shown Not shown 521 20 0.9 0.010 Yu et al 2005 E. japonicus Global (China Sea-Yellow Sea) Not shown 420 46 1.0 0.006 Chen et al 2010 E. japonicus Taiwan Este - Spring 2001 740 44 1.0 0.014 Chen et al 2010 E. japonicus Taiwan Oeste - Spring 2001 740 44 1.0 0.013 Chen et al 2010 E. japonicus Taiwan Este - Autum 2001 740 45 1.0 0.014 Grant et al 2010 E. japonicus Japan Not shown 687 22 0.9 0.010 Zheng et al 2015 E. japonicus Yellow Sea - A 2011 402 18 1.0 0.007 Zheng et al 2015 E. japonicus Yellow Sea - B1 2011 402 16 0.8 0.003 Zheng et al 2015 E. japonicus Yellow Sea - B2 2011 402 19 0.9 0.007 Zheng et al 2015 E. japonicus Yellow Sea - C 2011 402 7 0.9 0.007 Zheng et al 2015 E. japonicus Yellow Sea - D 2010 402 17 0.5 0.001 Zheng et al 2015 E. japonicus East chine Sea - E1 2009 402 24 1.0 0.004 Zheng et al 2015 E. japonicus East chine Sea - E2 2010 402 24 0.9 0.006 Zheng et al 2015 E. japonicus Yellow Sea - WH 2008 402 15 1.0 0.009 Zheng et al 2015 E. japonicus Yellow Sea - JN 2008 402 18 0.7 0.004 Zheng et al 2015 E. japonicus Yellow Sea - K 2008 402 14 0.9 0.006 Silva et al 2017 E. japonicus Japon 2006 1044 24 1.0 0.011 Grant and Bowen 1998 E. mordax Not shown Not shown 521 14 0.9 0.007 Lecomte et al 2004 E. mordax Washintong - Grays Harbor 1995 539 17 0.7 0.004 Lecomte et al 2004 E. mordax North California - Eureka 1994 539 20 0.8 0.004 Lecomte et al 2004 E. mordax South California - San Diego 94 1994 539 14 0.9 0.004 Lecomte et al 2004 E. mordax South California - San Diego 01 2001 539 30 0.9 0.005 Lecomte et al 2004 E. mordax Ensenada 1995 539 23 0.8 0.004 Lecomte et al 2004 E. mordax Upper baja California Sur - South of Cedros Is. 1995 539 35 0.8 0.004 Lecomte et al 2004 E. mordax Lower Baja California Sur - Bahía Magdalena 1995 539 37 0.9 0.006 Lecomte et al 2004 E. mordax Gulf of California - Guaymas 1995 539 20 0.9 0.005 Grant et al 2010 E. mordax California-Mexico Not shown 687 14 0.9 0.006 Grant and Bowen 1998 E. ringens Chile-Peru* Not shown 521 17 0.4 0.001 Grant et al 2010 E. ringens Talcahuano, Chile Not shown 687 17 0.4 0.001 This study E. ringens North Peru - Paita AP0215 feb-15 913 30 0.7 0.001 This study E. ringens Center Peru - Callao AC0215 feb-15 913 28 0.5 0.001 This study E. ringens North Peru - Salaverry AN0815 ago-15 913 30 0.5 0.001 This study E. ringens Center Peru - Callao AC0815 ago-15 913 30 0.4 0.001 This study E. ringens North Peru - Paita AN0416 abr-16 913 30 0.5 0.001 This study E. ringens Center Peru - Callao AC0416 abr-16 913 28 0.4 0.001 This study E. ringens North Peru - Sechura AN1016 oct-16 913 30 0.5 0.001 This study E. ringens Center Peru - Callao AC1016 oct-16 913 30 0.5 0.001 Tabla S5. Caracteristicas de las estimaciones compiladas de diversidad genética de diferentes miembros marinos de Engraulis. Las estimaciones se calcularon a partir de secuencias de nucleótidos de la región control del ADN mitocondrial. Referencia Especie Localidad de muestreo Fecha de muestreo pb n Hd π Viña et al. 2004 E. encrasicolus Islas canarias 1997 402 44 1.0 0.013 Viña et al. 2004 E. encrasicolus Cádiz 1997 402 47 1.0 0.017 Viña et al. 2004 E. encrasicolus Bay of Biscay - Donostia 1996 402 51 1.0 0.030 Viña et al. 2004 E. encrasicolus Bay of Biscay - Arcachon 1996 402 52 1.0 0.034 Viña et al. 2004 E. encrasicolus Bay of Biscay - Garona 1996 402 51 1.0 0.034 Viña et al. 2004 E. encrasicolus Northern Alborán Sea 1997 402 50 1.0 0.027 Viña et al. 2004 E. encrasicolus Northern Spain - Cadaqués 1996 402 52 1.0 0.031 Viña et al. 2004 E. encrasicolus Gulf of Lions 1993 402 23 1.0 0.032 Viña et al. 2004 E. encrasicolus Elba Island 1993 402 23 1.0 0.029 Viña et al. 2004 E. encrasicolus Sicilia 1996 402 37 1.0 0.031 Viña et al. 2004 E. encrasicolus North Adriatic sea 1997 402 42 1.0 0.028 Viña et al. 2004 E. encrasicolus South Adriatic sea 1997 402 38 1.0 0.031 Viña et al. 2004 E. encrasicolus Aegean Sea 1996 402 49 1.0 0.025 Cuttitta et al. 2015 E. encrasicolus Sicilia 1_2 2010 608 28 0.9 0.023 Cuttitta et al. 2015 E. encrasicolus Sicilia 3 2010 608 17 1.0 0.018 Cuttitta et al. 2015 E. encrasicolus Sicilia 5 2010 608 13 1.0 0.024 Cuttitta et al. 2015 E. encrasicolus Malta 2010 608 16 0.9 0.017 Liu et al. 2006 E. japonicus East Chine Sea - P 2001 522 16 1.0 0.015 Liu et al. 2006 E. japonicus East Chine Sea - E 2001 522 22 1.0 0.020 Liu et al. 2006 E. japonicus East Chine Sea - F 2001 522 20 1.0 0.018 Liu et al. 2006 E. japonicus East Chine Sea - L 2001 522 17 1.0 0.017 Liu et al. 2006 E. japonicus East Chine Sea - J 2001 522 20 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - R 2005 522 20 1.0 0.018 Liu et al. 2006 E. japonicus Sea japan - M 1997 522 18 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - I 1997 522 18 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - N 1996 522 18 1.0 0.019 Liu et al. 2006 E. japonicus Sea Japan - K 1997 522 18 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - T 1997 522 18 1.0 0.020 Liu et al. 2006 E. japonicus Sea Japan - G 1997 522 18 1.0 0.016 Liu et al. 2006 E. japonicus Sea Japan - U 1997 522 18 1.0 0.020 This study E. ringens North-Peru (Paita) 02_2014 520 33 0.9 0.005 This study E. ringens Center-Peru (Callao) 02_2014 520 27 0.9 0.004 This study E. ringens South-Peru (Ilo) 02_2014 520 30 0.9 0.004 This study E. ringens North-Peru (Paita) 02_2015 520 39 0.9 0.005 This study E. ringens Center-Peru (Callao) 02_2015 520 32 0.9 0.004 This study E. ringens North-Peru (Salaverry) 08_2015 520 32 0.9 0.005 This study E. ringens Center-Peru (Callao) 08_2015 520 32 0.8 0.005 This study E. ringens Center-Peru (Pisco) 08_2015 520 40 0.9 0.005 This study E. ringens North-Peru (Paita) 04_2016 520 39 0.8 0.004 This study E. ringens Center-Peru (Callao) 04_2016 520 41 0.9 0.005 This study E. ringens North-Peru (Sechura) 04_2016 520 37 0.9 0.004 This study E. ringens Center-Peru (Callao) 10_2016 520 31 0.9 0.005 This study E. ringens South-Peru (Tacna) 10_2016 520 32 0.8 0.004 Tabla S6. Comparaciones por pares de Фst basadas en secuencias de HVR1. Los valores en negrita indican diferencias genéticas significativas (Nivel de significancia = 0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 0.046 AS2014 -0.009 0.033 AP0215 -0.017 0.014 -0.013 AC0215 0.003 0.001 0.011 -0.011 AN0815 -0.021 0.028 -0.001 -0.015 -0.009 AC0815 -0.002 0.066 0.009 -0.002 0.017 0.009 AS0815 0.022 -0.018 0.014 0.002 -0.004 0.013 0.045 AN0416 -0.021 0.059 -0.006 -0.008 0.011 -0.013 -0.005 0.036 AC0416 -0.014 0.018 -0.016 -0.016 -0.002 -0.011 0.015 0.003 -0.005 AN1016 -0.003 0.002 -0.008 -0.007 -0.003 -0.007 0.025 -0.004 0.004 -0.012 AC1016 0.024 -0.006 0.035 0.001 -0.020 0.007 0.032 0.000 0.035 0.014 0.011 AS1016 0.000 0.011 0.021 -0.015 -0.017 -0.008 0.011 0.002 0.013 0.004 0.009 -0.017 Tabla S7. Valores p de las pruebas exactas de diferenciación de poblaciones basadas en secuencias de HVR1 (longitud de la cadena de Markov: 200000 pasos). Los valores p en negrita indican diferencias genéticas significativas (Nivel de significancia = 0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 0.215 AS2014 0.499 0.672 AP0215 0.998 0.459 0.879 AC0215 0.111 0.756 0.672 0.545 AN0815 0.922 0.279 0.164 0.841 0.238 AC0815 0.880 0.103 0.622 0.843 0.086 0.095 AS0815 0.736 0.536 0.636 0.748 0.297 0.224 0.473 AN0416 0.915 0.492 0.443 0.888 0.335 0.595 0.744 0.321 AC0416 0.917 0.550 0.858 0.960 0.505 0.456 0.401 0.545 0.608 AN1016 0.386 0.560 0.791 0.745 0.734 0.331 0.355 0.476 0.446 0.432 AC1016 0.189 0.775 0.459 0.646 0.822 0.533 0.146 0.455 0.331 0.509 0.680 AS1016 0.215 0.093 0.125 0.524 0.252 0.143 0.564 0.613 0.259 0.144 0.500 0.577 Tabla S8. Comparaciones por pares de Фst basadas en secuencias de Cytb. Un guión (-) indica muestreos en los que no se generaron secuencias de mtDNA. Los valores en negrita indican diferencias genéticas significativas (Nivel de significancia = 0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 - AS2014 - - AP0215 - - AC0215 - - - 0.014 AN0815 - - - 0.036 0.004 AC0815 - - - 0.014 -0.008 0.000 AS0815 - - - - - - - AN0416 - - - 0.010 -0.001 0.016 0.010 - AC0416 - - - 0.016 0.005 0.022 -0.004 - 0.017 AN1016 - - - -0.013 0.002 0.027 0.001 - 0.001 0.011 AC1016 - - - 0.003 -0.004 0.013 -0.007 - 0.006 -0.004 -0.005 AS1016 - - - - - - - - - - - - Tabla S9. Valores p de las pruebas exactas de diferenciación poblacional basadas en secuencias de Cytb (longitud de la cadena de Markov: 200000 pasos). Los valores en negrita indican diferencias genéticas significativas (Nivel de significancia = 0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 - AS2014 - - AP0215 - - - AC0215 - - - 0.291 AN0815 - - - 0.016 0.470 AC0815 - - - 0.017 0.402 0.407 AS0815 - - - - - - - AN0416 - - - 0.381 1.000 0.428 0.287 - AC0416 - - - 0.265 0.407 0.066 0.533 - 0.319 AN1016 - - - 0.874 0.390 0.033 0.143 - 0.433 0.379 AC1016 - - - 0.218 1.000 0.252 0.549 - 0.741 0.504 0.394 AS1016 - - - - - - - - - - - - Tabla S10. Análisis de la varianza molecular (AMOVA) basado en secuencias de Cytb. Estructura genética evaluada: norte versus centro de Perú. Pruebas de significancia (20022 permutaciones). Porcentaje de Fuente de la Grados de Suma de Componentes de Estadísticos variación de los Valores p variación libertad cuadrados la varianza AMOVA componentes Entre norte y 1 0.532 0.00091 0.25 Fct = 0.00249 0.26091 centro de Perú Entre muestreas 6 2.544 0.00202 0.55 Fsc = 0.00551 0.11028 dentro del norte y centro de Perú Dentro de cada 228 83.086 0.36441 99.20 Fst = 0.00798 0.05184 muestra Total 235 86.161 0.36734 Genetic diversity of the anchovy Engraulis ringens along Peruvian coast based on mitochondrial DNA sequences Deivis Alexander Cueva Shapiama. Laboratorio de Genética Molecular, Instituto del Mar del Perú, Esquina Gamarra y General Valle S/N Chucuito Callao. E-mail: dcueva@imarpe.gob.pe Giovanna Sotil. Laboratorio de Genética Molecular, Instituto del Mar del Peru, Esquina Gamarra y General Valle S/N Chucuito Callao. E-mail: gsotil@imarpe.gob.pe Correspondence. E-mail: gsotil@imarpe.gob.pe Running title: Genetic diversity of Peruvian anchovy Word count for main text: 4465 Abstract The anchoveta Engraulis ringens exhibits overwhelming biomass and experiences drastic demographic oscillations. Given this population instability, it is crucial to explore genetic variability, especially for the important relationship proposed by theoretical studies and empirical observations between levels of intraspecific variability and the fate of local populations. In the short term, this relationship has consequences for the sustainability of harvested populations; in the long term, for the perpetuation of the species. For anchoveta, one of the most studied species in the Northern Humboldt Current System, population genetic variability has been scarcely addressed. Here, anchoveta's mitochondrial genetic diversity is described through several temporal assessments and the comparison with diversity estimates reported for other engraulid fish inhabiting similar oceanographic settings (coastal upwelling systems). Results suggest that E. ringens displays the lowest levels of genetic diversity in relation to their congeners. This genetic signature had already been described for individuals from the southern end of the species range; now the same signal has been confirmed to the northern of the geographic range; that is, the area where the highest biomass of this fishery resource occurs (northern and central coasts of Peru). Considering that mtDNA-based inferences could reflect both historical demographic dynamics or deterministic pressures such as natural selection, the need to generate more studies is emphasized in order to provide information based on other types of markers and under a comprehensive sampling to verify the spatial and temporal consistency of the observed genetic signal. Keywords: Anchoveta, population genetics, genetic diversity, mitochondrial DNA, Humboldt Current System. Introduction In the search for sustainable management of harvested populations, genetic variation is a crucial attribute to explore and understand because it allows populations to show flexibility and adapt to scenarios that challenge their survival and perpetuation (Frankham et al. 2002; Reed and Frankham 2003). The anchoveta Engraulis ringens is one of the most studied marine species in Peru, with several well-registered biological aspects to support its fishery management. However, knowledge about population genetics of the species is remarkably scarce and it is still unknown how its genetic variability is repartitioned throughout its geographic range in the Humboldt Current System (HCS). Early in the study of Anchoveta, Tsukayama (1966) and Rojas 1971 (as cited in Pauly et al., 1987, pp. 5) reported meristic and morphometric differences between individuals from North-Center and South of Peru. Furthermore, Jordan and Malaga (1972) showed uneven displacements (by using Tag-recovery experiments) of individuals along Peruvian coast. From these observations, the management of two population units, one of North-Central Peru and the second from South Peru-Northern Chile was suggested (IMARPE 1973). Subsequently, Mathisen (1989) proposed an alternative organization composed of several local units based on a patchy distribution pattern that was regularly observed in density maps. The author also hypothesized that genetic differences could underlie such spatial aggregation. Later, Mendo (1991) tested these hypotheses and supported the existence of North- Central Peru and South Peru-North Chile subpopulations, a conclusion that was drawn from tag-recovery experiments and allozymes variability, while other types of data he used did not show a spatial differentiation. In Mendo's work (1991), eight Allozymes were evaluated but only three loci resulted informative, contrasting with the usual number of markers used when nuclear genetic variation is explored by this means (Gyllensten 1985; Waples 1987; Hedgecock et al. 1994; Tudela et al. 1999). As for the anchoveta from the southern of the geographic range, Ferrada et al. (2002) assessed the genetic differentiation between samples from North and South of Chile by using nuclear DNA markers and electrophoresis, but no signal of differentiation was detected, possibly due to the low polymorphism of the genetic markers used (Ferrada- Fuentes et al. 2018). Recently, efforts have been made to develop microsatellite markers for the population study of E. ringens but few loci have been successfully developed from samples of Peru (Rojas 2011) and others, more promising from samples from Chile (Ferrada-Fuentes et al. 2018), have not yet been reported under a population study design including samples from Peru. Equally important as the genetic differentiation approach is the study of demographic history of E. ringens in the HCS, a marine ecosystem known for exhibiting high oceanographic variability (Chavez y Messié 2009). It has been noticed that Anchoveta’s distribution and biomass is variable at different time scales and in response to environmental forcing apparently (Gutiérrez et al 2012; Chávez et al. 2003). Inference about demographic history emerging from genetic information could reveal past demographic patterns that help to understand contemporary population dynamics. For example, Grant and Bowen (1998) estimated mitochondrial genetic diversity in several anchovy species, including E. ringens and, based on demographic history, proposed metapopulation dynamics to explain a “shallow population genetic architecture” detected in Anchoveta. Grant and Bowen’s study (1998) was based on samples from almost 30 years ago and from the southern end of the species range. Therefore, an update and verification of this genetic signal in the northern part of the range (North-Central Peru) is important, because the largest biomass records have historically occurred in this area (Jordan 1971; Pauly and Tsukayama 1987; Ñiquen et al. 2000). Although the study of Anchoveta population genetic variability started almost 30 years ago, the information generated is still scarce for population inferences because studies were based on reduced number of markers or temporal/spatial restricted sampling designs. These limitations highlight the need to increase efforts for fulfil this gap. In this study we have explored the mitochondrial genetic diversity of E. ringens along the Peruvian coast for three years in order to describe a consistent genetic signal based on repeated assessments and selection of molecular markers that allow us to compare Anchoveta’s genetic diversity with those reported for other Engraulis species inhabiting coastal upwelling systems around the world. Materials and Methods Samples collection A total of 445 adult specimens of E. ringens were collected in front of different locations on the north, central and southern coast of Peru (Figure 1). The sampling, therefore, includes samples from the two fishery stocks (North-Center and South) managed in Peruvian Exclusive Economic Zone (EEZ). Except for AS1016 collection (Tacna, October 2016), samplings were all performed by the Dirección General de Investigación de Recursos Pelágicos (DGIRP) from the Instituto del Mar del Perú (IMARPE) through research cruises for the hydroacoustic evaluation of pelagic resources from 2014 to 2016. Samplings were carried out during the two main periods of greater reproductive activity usually observed in Perú (February-March and September-October). Tacna samples came from a natural fish stranding occurred in October 2016. In total 13 samples were obtained (Table 1) each one composed of 27 to 41 specimens of E. ringens. The collected fish were frozen at -20 °C and then stored to later analysis. Once in the laboratory a small piece (0.5x0.5 cm3) of dorsal muscle tissue from each individual was collected and preserved in 96% ethanol. The collected tissues were kept at 4 °C until subsequent processing. DNA extraction, PCR amplification and sequencing of mitochondrial markers The DNA was extracted with the CTAB2X-Chloroform-Isoamyl Alcohol method. Two regions of mitochondrial DNA (mtDNA) were amplified by PCR: the hypervariable region 1 (HVR1) of 520 base pair (bp) from the mitochondrial Control Region with primers that Lee et al. (1995) proposed for teleost fish; and a 913 bp fragment of the cytochrome b gene (Cytb) with primers used by Silva et al. (2014) for European anchovy Engraulis encrasicolus. The HotStart Taq Polymerase Master Mix kit (Qiagen) was used to perform 10µL PCR reactions, keeping buffer, enzyme and dNTP's concentrations according to manufacturer recommendation, while the final concentrations of MgCl2, primers and DNA template were adjusted for each marker. PCR products from each sample were sequenced for both DNA strands on ABI 3500 Genetic Analyzer (Applied Biosystems). Reviewing and editing electropherograms Electropherograms were checked (verification of good base calling process) and edited (trimming unreliable ends regions) with FinchTv 1.4.0 (Geospiza Research Team). The nucleotide sequences of HVR1 and Cytb were independently aligned using the ClustalW algorithm implemented in Bioedit 7.0.0 (Larkin et al. 2007). For both markers a consensus sequence was generated per sample based on forward and reverse strand sequences. Polymorphisms and genetic diversity estimates DnaSP software (Librado and Rozas 2009) was used to explore DNA polymorphism (number and types of mutation and variable sites) and to generate haplotype data. Nucleotide (π) and haplotype (Hd) diversity estimators were computed using Arlequin 3.5.1.2 software (Excoffier and Lischer 2010). We compared E. ringens genetic diversity (haplotype and nucleotide diversity) with estimates reported for other members of Engraulis genus inhabiting regions of coastal upwelling. The nominal species considered were Engraulis encrasicolus, Engraulis capensis, Engraulis eurystole, Engraulis japonicus, Engraulis australis, Engraulis mordax, Engraulis anchoita and Engraulis ringens. Genetic diversity was generated from partial and almost complete sequences of cytochrome b gene (approximately 400-1000 bp) and mitochondrial Control Region (approximately 400-600 bp); and sample sizes between 15 and 30 individuals (Supplementary Table S4 and S5). Genetic differentiation and population genetic structure Genetic differentiation was assessed with Φst pairwise comparisons and exact tests of differentiation, using Arlequin 3.5.1.2 (Excoffier and Lischer 2010). Possible scenarios of genetic structure (spatial and temporal) were evaluated by a hierarchical analysis of molecular variance using AMOVA test implemented in Arlequin 3.5.1.2 (Excoffier and Lischer 2010). The arrange of samples by zones in front of Peru (North, Center and South), fishing stocks (North-Center Perú and South Peru-North Chile), sampling campaigns (months) and years were considered as hierarchical levels of comparison of the estimated variation within and between levels. Additionally, we explore the genealogy of both mitochondrial markers by building haplotype networks using the Median-Joining method (Bandelt et al. 1999) and PopArt software (Population Analysis with Reticulate Trees; Leigh and Bryant 2015). A haplotype network is a graphic way to display genealogical relationships between mtDNA variants and the distribution of their frequencies in relation to a particular sampling design. Neutrality tests and Mismatch distribution We use the Tajima's D (Tajima 1989) and Fu's Fs (Fu 1997) tests to evaluate the null hypothesis of neutrality. Additionally, we performed a mismatch distribution for each marker. A mismatch analysis allows to infer possible demographic events that population could undergone in the past. A unimodal distribution reflects rapid population growth from a small population; while a multimodal distribution reflects long-term population stability (Slatkin and Hudson 1991, Rogers and Harpending 1992). Results DNA polymorphisms and genetic diversity in HVR1 and Cytb Fifty-four mutations were found in HVR1 (520 bp) from 445 individuals. The largest number of substitutions consisted of transitions (46 transitions vs. 8 transversions). The total number of variable sites was 52 (50 with two nucleotide variants per site, and two with three variants per site). The number of singleton and parsimony-informative sites was almost proportional (28 and 24 respectively). All these changes at nucleotide level generated 84 haplotypes, of which 29 (35%) were present in more than one sample (sample-shared haplotypes), while 55 (65%) only occurred in specific locations and dates (sample-private haplotypes). Sample-private haplotypes were mostly comprised of low frequency haplotypes (Supplementary Table S1). In partial region (913 bp) Cytb gene, 47 mutations were found from 236 individuals. Similar to HVR1, the majority of substitutions were transitions (44 transitions vs. 3 transversions). All the variable sites were polymorphic because they presented only two nucleotide variants per site, so their number was equal to the total number of mutations. Unlike HVR1, the number of singletons (35) was greater than the number of parsimony-informative sites (12). All variable sites defined 44 haplotypes, being 9 (20%) of them shared and 35 (80%) private. Most of these haplotypes also showed extremely low frequency (Supplementary Table S1). In general, It was evident a greater incorporation of mutations in HVR1 than in Cytb (Supplementary Table S1). For example, HVR1 showed a higher proportion of transversions (a type of substitution less common than transitions), presence of variable sites with more than two variants, and a lower proportion of singletons. Level of the mitochondrial genetic diversity in E. ringens Total haplotype diversity based on HVR1 indicated high diversity (Hd = 0.88), while estimate with Cytb showed a regular level (Hd = 0.50). Total nucleotide diversity followed a similar pattern, being higher in HVR1 (π = 0.0046) compared to Cytb (π = 0.0008). In general, both mitochondrial regions indicate regular to high amount of mtDNA variants, but most of them display low frequency and divergences supported by few mutational steps (Table 2). Among different sampling sites and dates the haplotype diversity in HVR1 remained high, varying between 0.82 and 0.93; while nucleotide diversity values ranged from 0.0040 to 0.0050. Regarding to Cytb, the haplotype diversity ranged between 0.37 and 0.68; while nucleotide diversity varied between 0.0005 and 0.0013 (Table 2). No significant differentiation was detected based on permutation tests (10000 interactions) for significance differences of genetic diversity (Supplementary Table S2 and S3). In general, genetic diversity in both markers presented narrow ranges of variation with no signals of significant differences. What is more, grouping diversity levels by sampling campaign or year resulted in a much more reduction of differences (Table 2). Levels of genetic diversity of E. ringens estimated here sharply differ from those reported for Old World Anchovies (OWA) and to a lesser extent from values reported for the California anchovy (Figure 2). Anchoveta showed Hd between 0.37 and 0.68 while the OWA generally showed values greater than 0.90. The estimates of Hd for E. mordax were also higher (Hd = 0.74-0.95) than those for E. ringens, although slightly lower than the OWA. Estimates for E. anchoita (Hd = 0.45) were very similar to E. ringens; however, they are based on a single sampling. In terms of nucleotide diversity, OWA exhibits a wide range of variation ranging from 0.006 to 0.017, while the range of π in E. mordax was narrower (0.004-0.009). In contrast, E. ringens exhibits extremely low values of π (~ 0.001) with minimal variability in different areas of northern, central and south of Peru, and quite similar to estimates compute from Talcahuano sample (Chile) (Supplementary Table S4). Genetic differentiation and population genetic structure Φst pairwise comparisons indicated low levels of genetic differentiation, with ranges of 0.000-0.036 using Cytb and 0.000-0.066 using HVR1 (Supplementary Table S6 and S8). In both markers, some pairwise comparisons showed low but statistically significant differentiation. For HVR1 the significant differences were not supported by the exact tests of differentiation, while for Cytb some genetic differences are supported by the two tests (Supplementary Table S7 and S9 respectively). In all scenarios of genetic structure, evaluated by AMOVA tests, the percentage of variation was always greater within hypothetical populations than among them, suggesting absence of a spatial (zones in front of Peru: North, Central and South; Φst = 0.001, P-value = 0.292 / Fishing stocks; Φst = -0.008, P-value = 0.988) or temporary structure (sampling campaigns; Φst = -0.007, P-value = 0.903 / years, Φst = -0.006, P- value = 0.982). Due to the general result of the analysis of the molecular variance, only results of hierarchical levels of comparison of the following organizations are shown: 1) North, Central and South along the Peruvian coast; and 2) Fishing stocks (Table 3 and 4). In general, haplotype networks showed a random spatial and temporal distribution of mtDNA variants. The frequency distribution of the most common haplotypes in both HVR1 and Cytb showed no drastic changes across all sampling locations and dates (Figure 4). Haplotypes genealogy, Neutrality tests and mismatch analysis The topology of Median-Joining network based on HVR1 was reticulated and complex, consisting of a dominant haplotype (30.8%) with several satellite haplotypes (usually 0.2%) around it, separated by one or two mutational steps; and a group of four intermediate frequency haplotypes (5.4-9.4%) related to the main haplotype through many low frequency haplotypes and other hypothetical ones. On the other hand, the Cytb-based network clearly exhibited a star-type topology, with a single dominant haplotype and several low-frequency satellite haplotypes (Fig. 4). Both markers show shallow genealogy of mtDNA lineages (Figure 4). The neutrality tests generated significantly negative values of the D statistic (-2.6418, P-value <0.05) and Fs (-29.0335, P-value <0.01), rejecting the neutrality hypothesis (Table 5). Furthermore, the mismatch distribution analysis showed unimodal distributions of pairwise differences frequency in both HVR1 and Cytb (Figure 3), which fits to the expected pattern under a population expansion scenario. Discussion Level of the mitochondrial genetic diversity in E. ringens The total genetic diversity of this study was similar to estimates reported by Grant and Bowen (1998), which sample came from Talcahuano, Chile, the southernmost area of Anchovetas distribution. This suggest a spatial and temporal consistent estimate of E. ringens genetic diversity, based on coding and no coding regions of mtDNA. Comparisons of nucleotide and haplotype diversity among engraulid fish suggested that E. ringens is one of the species with the lowest levels of genetic diversity among these small pelagic fish inhabiting coastal upwelling systems. Grant and Bowen (1998) had already been described this genetic signal for some Anchovetas from Talcahuano, Chile. They noticed that a shallow coalescence of main mtDNA lineages (i.e. a genealogy with scarce inter-haplotype divergence) in E. ringens contrasted sharply with deep genetic divergence between Anchoveta and their sister taxa (E. anchoita from Argentina and E. mordax from California). This finding led them to rule out a recent colonization from the regional populations of Argentina or California. Instead, they suggested this genetic signal would be the result of recurrent intraregional extinction and colonization events, invoking a metapopulation structure. This possible explanation supposes a high dynamism of the Anchoveta’s abundance and distribution and a complex population organization. What evidence exists of this? It is known that several populations of clupeids fish are characterized by exhibiting overwhelming biomass, but subject to drastic fluctuations (Lluch-Belda et al 1989; Schwartzlose et al. 1999). Anchoveta is not a different case, as it shows a high demographic variability at different time scales, including interannual (Gutiérrez et al 2012) and interdecadal changes (Chávez et al. 2003). This variability makes sense in the interpretation of historical demographic changes inferred from mtDNA sequences. In most phylogeographic studies, large-scale climatic changes, such as glaciations during the Pleistocene, are usually suggested as factors that have shaped the genetic architecture of fish populations (Díaz-Viloria et al. 2012; Silva et al. 2014; Silva et al., 2017; Delrieu-Trottin et al. 2017). In the HCS other climatic processes at ocean and global scale are known to be congruent with drastic changes at interdecadal scale in Anchoveta’s biomass (Chavez et al. 2003). In addition, at short-time scale, El Niño events can dramatically reduce the abundance of anchovy (Barber and Chavez 1983; Gutiérrez et al. 1999; Ñiquen and Bouchón 2004) and there is evidence to suggest these events have episodically occurred in the last 17000 years (Rein et al. 2005). Therefore, environmental determinants, acting at different time scales, can induce radical demographic changes in E. ringens. Back to the proposal by Grant and Bowen (1998), evidence on population dynamics supports recurrent changes in abundance rather than demographic stability, and theory predicts that demographic instability may cause a reduction in genetic diversity (Pannell and Charlesworth 1999). However, the existence of a complex population organization (intra-regional metapopulation) is not yet clear in E. ringens. The anchoveta population dynamics may not strictly fit to a classic extinction- colonization model because local extinction has not been reported even under oceanographic conditions that challenge its survival. For example, few months after one of the strongest El Niño events in recent decades (El Niño 97/98), a surprisingly rapid recovery of the North-Central Peru Stock occurred, suggesting adult anchovies were able to survive a range reduction of habitat (Bertrand et al. 2004). Instead of local extermination, a more plausible explanation would describe drastic population crash and subsequent recovery with individuals from other parts of the range but also with survivors from the same local area. This possible explanation thus suggests important levels of connectivity throughout the species range. This dynamic is more in line with the definition of metapopulation that Kritzer and Sale (2004) proposed for marine populations, where the degree of connectivity and its relevance in the fate of the subpopulation is emphasized rather than the extinction/colonization events. Genetic Differentiation and Population Structure In the HCS there are three important anchovy fisheries; a north, center and south component (Figure 1). In Peruvian EEZ, the north component (Stock North-Central Peru) and part of the center component (South Stock Peru-North Chile) are managed. The third component, smaller in terms of biomass, occurs at the southern end of the distribution range, in the jurisdictional sea of Chile (Central-South Chile Stock) (Cubillos et al. 2007). In this work, different scenarios of population structure were tested, considering several hierarchical levels of comparison (temporal and spatial), and no structuring signal was observed in congruence with the fishery components. This result contrasts with Mendo's study (1991), where a restriction in dispersion between North-Central and South anchovies was described from Tag-Recovery data. The restriction was more pronounced in the southern Winter-Spring, when the period of greatest reproductive activity occurs, suggesting that random mating may not always be the case for Anchoveta and leading to the idea of possible genetic differences. Indeed Mendo (1991) also supported the previous results with data on genetic variation by using allozyme markers. Although only three polymorphic loci were informative, they showed significant differences in allele frequency between North-Central and South regions. Differences between mtDNA and Allozymes results are not only due to the obvious differences in inheritance of both molecular traits, but also some kind of population dynamics could be involved. Grant and Leslie (1993) made simulations on the effect of periodic extinction-colonization events on mitochondrial and Allozyme variability and a general pattern emerged: after a given number of generations, extinction-colonization events produce drastic erosion in mitochondrial genetic variability with Fst very close to zero; followed by slow reduction in nuclear genetic variability with Fst keeping away from zero. Grant and Leslie (1993) explained that, under metapopulation context, mitochondrial lineages are more prone to extinction compared to nuclear ones because mtDNA is haploid and has uniparental inheritance (reducing the effective population size approximately one-quarter than nuclear one). The contrasting results between Mendo’s study (1991) and this work fit this description. Therefore, if the Anchoveta population experiences this type of dynamics, then mtDNA is not an appropriate marker to make robust inferences about population genetic structure. The evidence so far provides some clues about a possible population structure of E. ringens, but studies still failed to 1) cover entirely the species range or 2) employ appropriate number of genetic markers. This aspect is further encrypted by the pelagic nature of Anchoveta and its implication on dispersal and spatial aggregation dynamics. For example, disturbs of Anchoveta´s habitat (e.g. El Niño events or El Niño-like warm events) lead fish to move at different spatial scales, from vertical movements to larger longitudinal/latitudinal displacements (Gutiérrez et al. 2012). In this mixing situation, if migrants stay until the reproductive period, genetic differentiation would be challenge because the buffering effect of combining individual from distinct subpopulations. To summarize, because the degree of connectivity between the North-Central Peru and South Peru-North Chile stocks has not been comprehensively addressed (here or in previous studies), we could not conclude whether genetic evidence support a metapopulation or patchy population of E. ringens; that is, spatially structured but displaying an internal population dynamic equally throughout their range (Kritzer and Sale 2004).However, the interspecific comparison of genetic diversity sharply shows a pattern of low levels in Anchoveta, and several complex explanations could be invoked to account for it. One of them suggesting metapopulation dynamics. Therefore, it will be necessary to include information from other types of molecular markers (e. g. microsatellites or those generated by high throughput DNA sequencing) and appropriate sampling designs to verify if the low genetic diversity in Anchoveta’s mtDNA certainly reflect a genomic feature of this specie. Acknowledgements We would like to thank Dirección General de Investigación de Recursos Pelágicos (DGIRP) from Instituto del mar del Perú (IMARPE) for help with sample collections. Conflict of Interest The authors declare no conflict of interest. FIGURES AND TABLES Genetic diversity of the anchovy Engraulis ringens along Peruvian coast based on mitochondrial DNA sequences D Cueva1, G Sotil1 1Laboratorio de Genética Molecular, Instituto del Mar del Peru, Esquina Gamarra y General Valle S/N Chucuito E-mail: gsotil@imarpe.gob.pe Figure 1. Sampling locations in Peruvian EEZ. The geographic range of E. ringens and the main fishery components in the Humboldt Current System are also showed. Table 1. Details about the 13 sampling campaigns carried out in this study from 2014 to 2016. Fishery Stock and Peruvian coast zone which samples come from is indicated. To simplify, the sample code is made up of the first two letters of the sampling region or locality, followed by two numbers indicating the month and finally two numbers for the year. Peruvian coast Fishery stock Location Date Sample code zone February North-Center Peru North-Peru Paita AP2014 2014 February Paita AP0215 2015 Salaverry August 2015 AN0815 Paita April 2016 AN0416 October Sechura AN1016 2016 February North-Center Peru Center-Peru Callao AC2014 2014 February Callao AC0215 2015 Callao August 2015 AC0815 Pisco August 2015 AS0815 Callao April 2016 AC0416 October Callao AC1016 2016 South Peru – North February South-Peru Ilo AS2014 Chile 2014 October Tacna AS1016 2016 Table 2. Mitochondrial genetic diversity of E. ringens along Peruvian coast. A hyphen (-) indicates samples where no mtDNA sequences were generated. For an overview of genetic diversity estimates were also compute by pooling sequences from different locations from the same sampling campaign (month), year and as a total value. n Hd (± s.d.) π (± s.d.) Samples HVR1 Cytb HVR1 Cytb HVR1 Cytb AP2014 33 - 0.85 ± 0.05 - 0.0046 ± 0.0005 - AP0215 39 30 0.87 ± 0.04 0.68 ± 0.09 0.0049 ± 0.0005 0.0013 ± 0.0003 AN0815 32 30 0.88 ± 0.04 0.51 ± 0.11 0.0046 ± 0.0004 0.0006 ± 0.0002 AN0416 39 30 0.84 ± 0.05 0.51 ± 0.11 0.0042 ± 0.0004 0.0009 ± 0.0003 AN1016 37 30 0.89 ± 0.04 0.51 ± 0.11 0.0044 ± 0.0003 0.0008 ± 0.0002 AC2014 27 - 0.92 ± 0.03 - 0.0045 ± 0.0005 - AC0215 32 28 0.86 ± 0.04 0.50 ± 0.12 0.0040 ± 0.0005 0.0007 ± 0.0002 AC0815 32 30 0.84 ± 0.07 0.37 ± 3.00 0.0047 ± 0.0005 0.0006 ± 0.0002 AC0416 41 28 0.90 ± 0.04 0.38 ± 0.11 0.0049 ± 0.0004 0.0005 ± 0.0002 AC1016 31 30 0.93 ± 0.03 0.52 ± 0.11 0.0047 ± 0.0004 0.0009 ± 0.0003 AS0815 40 - 0.92 ± 0.03 - 0.0050 ± 0.0003 - AS2014 30 - 0.90 ± 0.04 - 0.0043 ± 0.0004 - AS1016 32 - 0.82 ± 0.06 - 0.0045 ± 0.0003 - February 2014 90 - 0.89 ± 0.02 - 0.0045 ± 0.0003 - Febreary 2015 71 58 0.87 ± 0.03 0.59 ± 0.08 0.0045 ± 0.0003 0.0010 ± 0.0002 August 2015 104 60 0.89 ± 0.03 0.44 ± 0.08 0.0048 ± 0.0002 0.0006 ± 0.0001 April 2016 80 58 0.86 ± 0.03 0.45 ± 0.08 0.0045 ± 0.0003 0.0007 ± 0.0002 October 2016 100 60 0.88 ± 0.02 0.51 ± 0.08 0.0045 ± 0.0002 0.0009 ± 0.0002 2015 175 118 0.88 ± 0.02 0.52 ± 0.06 0.0047 ± 0.0002 0.0008 ± 0.0001 2016 180 118 0.87 ± 0.02 0.48 ± 0.06 0.0045 ± 0.0002 0.0008 ± 0.0001 Total 445 236 0.88 ± 0.01 0.50 ± 0.04 0.0046 ± 0.0001 0.0008 ± 0.0001 0.020 E. ringens A 0.018 E. anchoita New Worl Anchovies 0.016 E. mordax E. japonicus 0.014 E. australis 0.012 E. capensis Old World Anchovies 0.010 E. encrasicolus E. eurystole 0.008 0.006 0.004 0.002 0.000 0.0 0.2 0.4 0.6 0.8 1.0 1.2 Haplotype diversity (Hd) 0.040 B 0.035 E. ringens New World Anchovy 0.030 E. japonicus Old World Anchovies 0.025 E. encrasicolus 0.020 0.015 0.010 0.005 0.000 0.0 0.2 0.4 0.6 0.8 1.0 1.2 Haplotype diversity (Hd) Figure 2. Scatter plot of mitochondrial genetic diversity of marine anchovy species inhabiting coastal upwelling zones. Estimates were computed from Cytb gene (A) and mitochondrial Control Region (B) sequences. See supplementary material for details about the literature used. nucleotide diversity (π) nucleotide diversity (π) Table 3. Analysis of Molecular Variance (AMOVA) based on HVR1 sequences. Genetic structure to test: North, Central and South of Peru. Significance tests based on 20022 permutations. Percentage of Source of Degree of Sum of Variance AMOVA variance P Value variation freedom squares components statistics components Among North, Fct = Central and 2 3.160 0.00200 0.17 0.21701 0.00168 South of Peru groups Among samples within Fsc = 10 13.116 0.00365 0.31 0.28673 North, 0.00307 center and South groups Within Fst = 432 512.453 1.18623 99.53 0.21881 samples 0.00475 Total 444 528.730 1.19189 Table 4. Analysis of Molecular Variance (AMOVA) based on HVR1 sequences. Genetic structure to test: Stock North-Central Peru versus South Peru-North Chile. Significance tests based on 20022 permutations. Percentage of Source of Degree of Sum of Variance AMOVA variance P Value variation freedom squares components statistics components Among Fct = - 1 0.403 -0.00956 -0.81 0.98756 fishery 0.00807 stock Among samples Fsc = 11 15.874 0.00746 0.63 0.15133 within 0.00625 fishery stocks Within Fst = - 432 512.453 1.18623 100.18 0.21881 samples .00177 Total 444 528.730 1.18414 Table 5. Neutrality tests for HVR1 and Cytb sequences. Bold values indicate significant deviations from neutral model (Significance Level=0.05). Tajima’s D Fu’s FS Samples HVR1 Cytb HVR1 Cytb AP2014 -1.1744 - -9.4778* - AC2014 -0.3300 - -7.9154* - AS2014 -0.8683 - -9.1784* - AP0215 -0.9123 -2.1456* -10.7640* -7.6703* AC0215 -0.7352 -2.2798* -7.6521* -8.0920* AN0815 -0.8244 -1.9744* -6.5902* -6.5545* AC0815 -0.6032 -2.2355* -15.1188* -5.2786* AS0815 -1.0151 - -12.0377* - AN0416 -1.2314 -1.9535* -9.1137* -6.2841* AC0416 -1.2818 -1.7197 -10.3126* -2.6376* AN1016 -0.8442 -1.8828* -10.6453* -3.7452* AC1016 -0.6189 -2.2964* -9.7874* -8.3947* AS1016 0.1148 - -5.5473* - Total -1.94513* -2.6418* -26.37530* -29.0335* Figure 3. Analysis of the pairwise differences distribution (mismatch distribution) under sceneries of constant and growth-decline population history. Figure 4. Haplotype networks based on Cytb (A) and HVR1 (B) sequences SUPPLEMENTARY MATERIAL Genetic diversity of the anchovy Engraulis ringens along Peruvian coast based on mitochondrial DNA sequences D Cueva1, G Sotil1 1Laboratorio de Genética Molecular, Instituto del Mar del Peru, Esquina Gamarra y General Valle S/N Chucuito E-mail: gsotil@imarpe.gob.pe Table S1. Number and types of mtDNA sequence variations in Engraulis ringens. Mitochondrial markers HVR1 Cytb Sample size (n) 445 236 Total number of nucleotide sites (locus size) 520 913 Total number of mutation (Eta) 54 47 No. of INDELs (insertions-deletions) 0 0 No. of nucleotide Substitutions 54 47 No. of transition mutations 46 44 No. of transversion mutations 8 3 Total number of variable sites (relative frequency) 52 (10%) 47 (5%) No. of singleton sites (relative frequency) 28 (54%) 35 (74%) No. of singleton with two variants 27 35 No. of singleton with three variants 1 0 No. of parsimony-informative sites (relative frequency) 24 (46%) 12 (26%) No. of parsimony-informative with two variants 23 12 No. of parsimony-informative sites with three variants 1 0 Total number of haplotypes 84 44 No. of sample-private haplotypes (relative frequency) 55 (65%) 35 (80%) No. of sample-shared haplotypes (relative frequency) 29 (35%) 9 (20%) Table S2a. Permutation tests for significant differences in nucleotide diversity (π) between populations based on HVR1 sequence data. The procedure was taken from Alexander et al. (2016). The R script is available in https://github.com/laninsky/genetic_diversity_diffs 0.00 π 0.0046 0.0045 0.0043 0.0049 0.0040 0.0046 0.0047 0.0050 0.0042 0.0049 0.0044 0.0047 45 π diff / p- AN101 AS1 AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AC1016 values 6 016 0.84 AP2014 -- 0.883 0.609 0.646 0.394 0.937 0.902 0.560 0.538 0.584 0.787 0.888 5 0.00010 0.96 AC2014 -- 0.758 0.519 0.539 0.798 0.802 0.453 0.660 0.513 0.881 0.774 0176 8 0.00030 0.00020 0.76 AS2014 -- 0.307 0.698 0.580 0.543 0.277 0.910 0.286 0.871 0.549 6899 6723 8 0.00029 0.00039 0.00060 0.47 AP0215 -- 0.158 0.692 0.719 0.905 0.225 0.966 0.393 0.734 7185 7361 4084 6 0.00054 0.00044 0.00023 0.00084 0.53 AC0215 -- 0.353 0.306 0.141 0.795 0.155 0.562 0.364 5741 5565 8843 2926 0 0.00005 0.00015 0.00036 0.00024 0.00060 0.76 AN0815 -- 0.970 0.595 0.463 0.642 0.685 0.968 5220 5396 2119 1965 0962 6 0.00007 0.00017 0.00038 0.00022 0.00062 0.00001 0.75 AC0815 -- 0.618 0.456 0.674 0.665 0.977 4606 4782 1505 2579 0347 9386 9 0.00036 0.00046 0.00067 0.00007 0.00091 0.00031 0.00029 0.42 AS0815 -- 0.185 0.932 0.314 0.670 8814 8990 5712 1629 4555 3594 4208 2 0.00038 0.00028 0.00007 0.00067 0.00016 0.00043 0.00045 0.00075 0.66 AN0416 -- 0.202 0.726 0.441 2769 2593 5871 9954 2972 7990 7375 1583 8 0.00031 0.00041 0.00062 0.00002 0.00086 0.00026 0.00024 0.00004 0.00070 0.45 AC0416 -- 0.365 0.700 9624 9800 6523 2439 5366 4404 5018 9189 2394 2 0.00018 0.00008 0.00011 0.00048 0.00035 0.00024 0.00026 0.00055 0.00019 0.00050 0.91 AN1016 -- 0.650 9985 9809 6914 7170 5757 5205 4591 8798 2785 9609 4 0.00009 0.00019 0.00039 0.00020 0.00063 0.00003 0.00001 0.00027 0.00047 0.00022 2.8139 0.74 AC1016 -- 1407 1583 8306 5778 7148 6187 6801 7407 4176 8217 E-04 0 0.00012 0.00002 0.00018 0.00042 0.00042 0.00017 0.00019 0.00049 0.00025 0.00044 6.6855 0.00021 AS1016 -- 3130 2954 3769 0315 2612 8350 7736 1943 9640 2754 E-05 4537 Table S2b. Permutation tests for significant differences in Haplotype diversity (Hd) between populations based on HVR1 sequence data. The procedure was taken from Alexander et al. (2016). The R script is available in https://github.com/laninsky/genetic_diversity_diffs 0.82 Hd 0.850 0.923 0.903 0.873 0.859 0.883 0.841 0.919 0.837 0.895 0.892 0.931 1 Hd diff /P- AN101 AS1 AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AC1016 values 6 016 0.62 AP2014 -- 0.257 0.409 0.717 0.877 0.584 0.871 0.265 0.811 0.46 0.484 0.21 1 0.07269 0.14 AC2014 -- 0.759 0.436 0.344 0.539 0.232 0.942 0.192 0.642 0.60 0.902 8135 2 0.05306 0.01962 0.19 AS2014 -- 0.614 0.474 0.733 0.307 0.781 0.286 0.897 0.829 0.675 9488 8647 1 0.02276 0.04993 0.03030 0.37 AP0215 -- 0.794 0.867 0.584 0.389 0.543 0.704 0.729 0.344 5612 2524 3876 7 0.00849 0.06420 0.04457 0.01427 0.53 AC0215 -- 0.699 0.75 0.303 0.711 0.542 0.571 0.268 218 5955 7308 3432 6 0.03268 0.04001 0.02038 0.00992 0.02419 0.32 AN0815 -- 0.485 0.513 0.425 0.817 0.894 0.430 5728 2407 376 0117 3548 1 0.00965 0.08235 0.06272 0.03241 0.01814 0.04233 0.71 AC0815 -- 0.198 0.942 0.365 0.379 0.160 2981 1117 2469 8593 5161 871 5 0.06885 0.00384 0.01578 0.04608 0.06035 0.03616 0.07850 0.11 AS0815 -- 0.16 0.648 0.623 0.843 1981 6154 2493 637 9801 6253 4963 7 0.01367 0.08636 0.06674 0.03643 0.02216 0.04635 0.00401 0.08252 0.77 AN0416 -- 0.292 0.329 0.139 1635 9771 1123 7247 3815 7364 8654 3617 6 0.04474 0.02795 0.00832 0.02197 0.03625 0.01205 0.05439 0.02410 0.05841 0.20 AC0416 -- 0.952 0.531 3163 4972 6325 7552 0983 7435 6145 8818 4799 8 0.04151 0.03118 0.01155 0.01874 0.03302 0.00882 0.05116 0.02733 0.05518 0.00323 0.23 AN1016 -- 0.485 3104 5031 6384 7492 0924 7376 6085 8877 4739 0059 6 0.08080 0.00810 0.02773 0.05803 0.07231 0.04811 0.09045 0.01195 0.09447 0.03606 0.0392 0.08 AC1016 -- 4008 5873 452 8396 1828 828 6989 2026 5643 0844 909 7 0.02981 0.10251 0.08288 0.05257 0.03830 0.02016 0.09866 0.01614 0.07455 0.0713 0.1106 AS1016 0.0625 -- 4272 2407 376 9883 6452 129 6253 2636 7435 2738 1828 Table 3a. Permutation tests for significant differences in nucleotide diversity (π) between populations based on Ctyb sequence data. The procedure was taken from Alexander et al. (2016). The R script is available in https://github.com/laninsky/genetic_diversity_diffs π 0.0013 0.0007 0.0006 0.0006 0.0009 0.0005 0.0008 0.0009 π diff / p-values AP0215 AC0215 AN0815 AC0815 AN0416 AC0416 AN1016 AC1016 AP0215 -- 0.13 0.08 0.06 0.35 0.05 0.24 0.27 AC0215 0.00055 -- 0.86 0.73 0.54 0.62 0.70 0.65 AN0815 0.00062 6.205E-05 -- 0.89 0.45 0.77 0.59 0.50 AC0815 0.00067 1.200E-04 5.791E-05 -- 0.35 0.89 0.47 0.42 AN0416 0.00033 2.200E-04 2.820E-04 0.00034 -- 0.28 0.80 0.90 AC0416 0.00072 1.681E-04 1.060E-04 0.00005 0.00039 -- 0.40 0.34 AN1016 0.00042 1.369E-04 1.989E-04 0.00026 0.00008 0.00030 -- 0.92 AC1016 0.00039 1.671E-04 2.291E-04 0.00029 0.00005 0.00034 3.021E-05 -- Table 3a. Permutation tests for significant differences in haplotype diversity (Hd) between populations based on Ctyb sequence data. The procedure was taken from Alexander et al. (2016). The R script is available in https://github.com/laninsky/genetic_diversity_diffs Hd 0.68 0.50 0.51 0.37 0.51 0.38 0.51 0.52 Hd/P-values AP0215 AC0215 AN0815 AC0815 AN0416 AC0416 AN1016 AC1016 AP0215 -- 0.28 0.30 0.06 0.32 0.08 0.31 0.33 AC0215 0.1808 -- 0.98 0.41 0.96 0.49 0.99 0.94 AN0815 0.1678 0.0130 -- 0.38 0.92 0.45 0.94 0.90 AC0815 0.3126 0.1318 0.1448 -- 0.36 0.95 0.40 0.35 AN0416 0.1632 0.0176 0.0046 0.1494 -- 0.42 0.90 0.95 AC0416 0.2946 0.1138 0.1267 0.0181 0.1313 -- 0.47 0.42 AN1016 0.1701 0.0107 0.0023 0.1425 0.0069 0.1244 -- 0.89 AC1016 0.1609 0.0199 0.0069 0.1517 0.0023 0.1336 0.0092 -- Table S4. Details about the compilated estimates of genetic diversity from different marine Engraulis members. Estimates were calculated from nucleotide sequences of mitochondrial gene Cytochrome b. Sampling Refference species Sampling location pb n Hd π date Grant and Bowen 1 E. anchoita Not shown Not shown 521 0.4 0.001 1998 9 1 Grant et al 2010 E. anchoita Argentina Not shown 687 0.5 0.001 9 Grant and Bowen 1 E. australis Not shown Not shown 521 0.9 0.009 1998 6 104 3 Silva et al 2017 E. australis Australia 2008 1.0 0.007 4 4 104 3 Silva et al 2017 E. australis Nueva Zelanda 2005 0.7 0.001 4 5 Grant and Bowen 1 E. capensis Southern Africa* Not shown 521 0.2 0.004 1998 8 104 2 Silva et al 2017 E. capensis Namibia 2007 0.9 0.010 4 4 104 1 Silva et al 2017 E. capensis Sudafrica 2007 0.9 0.019 4 3 Grant and Bowen E. 1 Not shown Not shown 521 0.9 0.015 1998 encrasicolus 6 E. 1 Grant et al 2010 Europe Not shown 687 0.9 0.014 encrasicolus 7 E. 3 Borrel et al 2012 Bay of Biscay - CA Not shown 524 0.7 0.012 encrasicolus 9 E. 2 Borrel et al 2012 Bay of Biscay - LI Not shown 524 0.9 0.016 encrasicolus 0 E. 1 Borrel et al 2012 Bay of Biscay - BB Not shown 524 1.0 0.016 encrasicolus 7 E. 1 Borrel et al 2012 Bay of Biscay - GE Not shown 524 0.8 0.013 encrasicolus 2 E. 3 Borrel et al 2012 Bay of Biscay - FR Not shown 524 0.9 0.016 encrasicolus 1 E. 2 Borrel et al 2012 Mediterranean - Adriatic Sea Not shown 524 1.0 0.014 encrasicolus 0 E. 104 3 Silva et al 2014 Málaga, ML Not shown 1.0 0.007 encrasicolus 4 1 E. 104 3 Silva et al 2014 Tangier, TG Not shown 0.4 0.003 encrasicolus 4 8 E. 104 2 Silva et al 2017 Noruega 2007 1.0 0.009 encrasicolus 4 4 E. 104 Silva et al 2017 Polonia 2008 9 0.9 0.014 encrasicolus 4 E. 104 2 Silva et al 2017 El Canal Inglés 2007 1.0 0.010 encrasicolus 4 7 E. 104 2 Silva et al 2017 El Golfo de Vizcaya 2007 1.0 0.015 encrasicolus 4 3 E. 104 2 Silva et al 2017 Norte de Portugal 1998 1.0 0.013 encrasicolus 4 5 E. 104 2 Silva et al 2017 Sur de Portugal 2007 1.0 0.011 encrasicolus 4 9 E. 104 2 Silva et al 2017 Islas Canarias 1999 1.0 0.007 encrasicolus 4 4 E. 104 2 Silva et al 2017 Senegal 1999 1.0 0.006 encrasicolus 4 5 E. 104 2 Silva et al 2017 Guinea-Bisáu 2006 1.0 0.006 encrasicolus 4 0 E. 104 2 Silva et al 2017 Ghana 2008 1.0 0.006 encrasicolus 4 5 104 1 Silva et al 2017 E. eurystole US 2006 0.9 0.004 4 2 Grant and Bowen 2 E. japonicus Not shown Not shown 521 0.9 0.010 1998 0 4 Yu et al 2005 E. japonicus Global (China Sea-Yellow Sea) Not shown 420 1.0 0.006 6 4 Chen et al 2010 E. japonicus Taiwan Este - Spring 2001 740 1.0 0.014 4 4 Chen et al 2010 E. japonicus Taiwan Oeste - Spring 2001 740 1.0 0.013 4 4 Chen et al 2010 E. japonicus Taiwan Este - Autum 2001 740 1.0 0.014 5 2 Grant et al 2010 E. japonicus Japan Not shown 687 0.9 0.010 2 1 Zheng et al 2015 E. japonicus Yellow Sea - A 2011 402 1.0 0.007 8 1 Zheng et al 2015 E. japonicus Yellow Sea - B1 2011 402 0.8 0.003 6 1 Zheng et al 2015 E. japonicus Yellow Sea - B2 2011 402 0.9 0.007 9 Zheng et al 2015 E. japonicus Yellow Sea - C 2011 402 7 0.9 0.007 1 Zheng et al 2015 E. japonicus Yellow Sea - D 2010 402 0.5 0.001 7 2 Zheng et al 2015 E. japonicus East chine Sea - E1 2009 402 1.0 0.004 4 2 Zheng et al 2015 E. japonicus East chine Sea - E2 2010 402 0.9 0.006 4 1 Zheng et al 2015 E. japonicus Yellow Sea - WH 2008 402 1.0 0.009 5 1 Zheng et al 2015 E. japonicus Yellow Sea - JN 2008 402 0.7 0.004 8 1 Zheng et al 2015 E. japonicus Yellow Sea - K 2008 402 0.9 0.006 4 104 2 Silva et al 2017 E. japonicus Japon 2006 1.0 0.011 4 4 Grant and Bowen 1 E. mordax Not shown Not shown 521 0.9 0.007 1998 4 1 Lecomte et al 2004 E. mordax Washintong - Grays Harbor 1995 539 0.7 0.004 7 2 Lecomte et al 2004 E. mordax North California - Eureka 1994 539 0.8 0.004 0 1 Lecomte et al 2004 E. mordax South California - San Diego 94 1994 539 0.9 0.004 4 3 Lecomte et al 2004 E. mordax South California - San Diego 01 2001 539 0.9 0.005 0 2 Lecomte et al 2004 E. mordax Ensenada 1995 539 0.8 0.004 3 Upper baja California Sur - South of Cedros 3 Lecomte et al 2004 E. mordax 1995 539 0.8 0.004 Is. 5 3 Lecomte et al 2004 E. mordax Lower Baja California Sur - Bahía Magdalena 1995 539 0.9 0.006 7 2 Lecomte et al 2004 E. mordax Gulf of California - Guaymas 1995 539 0.9 0.005 0 1 Grant et al 2010 E. mordax California-Mexico Not shown 687 0.9 0.006 4 Grant and Bowen 1 0. 0.00 E. ringens Chile-Peru* Not shown 521 1998 7 4 1 1 0. 0.00 Grant et al 2010 E. ringens Talcahuano, Chile Not shown 687 7 4 1 3 This study E. ringens North Peru - Paita AP0215 feb-15 913 0.7 0.001 0 2 This study E. ringens Center Peru - Callao AC0215 feb-15 913 0.5 0.001 8 3 This study E. ringens North Peru - Salaverry AN0815 ago-15 913 0.5 0.001 0 3 This study E. ringens Center Peru - Callao AC0815 ago-15 913 0.4 0.001 0 3 This study E. ringens North Peru - Paita AN0416 abr-16 913 0.5 0.001 0 2 This study E. ringens Center Peru - Callao AC0416 abr-16 913 0.4 0.001 8 3 This study E. ringens North Peru - Sechura AN1016 oct-16 913 0.5 0.001 0 3 This study E. ringens Center Peru - Callao AC1016 oct-16 913 0.5 0.001 0 Table S5. Details about the compilated estimates of genetic diversity from different marine Engraulis members. Estimates were calculated from nucleotide sequences of mitochondrial Control Region. Refference species Sampling location Sampling date pb n Hd π Viña et al. 2004 E. encrasicolus Islas canarias 1997 402 44 1.0 0.013 Viña et al. 2004 E. encrasicolus Cádiz 1997 402 47 1.0 0.017 Viña et al. 2004 E. encrasicolus Bay of Biscay - Donostia 1996 402 51 1.0 0.030 Viña et al. 2004 E. encrasicolus Bay of Biscay - Arcachon 1996 402 52 1.0 0.034 Viña et al. 2004 E. encrasicolus Bay of Biscay - Garona 1996 402 51 1.0 0.034 Viña et al. 2004 E. encrasicolus Northern Alborán Sea 1997 402 50 1.0 0.027 Viña et al. 2004 E. encrasicolus Northern Spain - Cadaqués 1996 402 52 1.0 0.031 Viña et al. 2004 E. encrasicolus Gulf of Lions 1993 402 23 1.0 0.032 Viña et al. 2004 E. encrasicolus Elba Island 1993 402 23 1.0 0.029 Viña et al. 2004 E. encrasicolus Sicilia 1996 402 37 1.0 0.031 Viña et al. 2004 E. encrasicolus North Adriatic sea 1997 402 42 1.0 0.028 Viña et al. 2004 E. encrasicolus South Adriatic sea 1997 402 38 1.0 0.031 Viña et al. 2004 E. encrasicolus Aegean Sea 1996 402 49 1.0 0.025 Cuttitta et al. 2015 E. encrasicolus Sicilia 1_2 2010 608 28 0.9 0.023 Cuttitta et al. 2015 E. encrasicolus Sicilia 3 2010 608 17 1.0 0.018 Cuttitta et al. 2015 E. encrasicolus Sicilia 5 2010 608 13 1.0 0.024 Cuttitta et al. 2015 E. encrasicolus Malta 2010 608 16 0.9 0.017 Liu et al. 2006 E. japonicus East Chine Sea - P 2001 522 16 1.0 0.015 Liu et al. 2006 E. japonicus East Chine Sea - E 2001 522 22 1.0 0.020 Liu et al. 2006 E. japonicus East Chine Sea - F 2001 522 20 1.0 0.018 Liu et al. 2006 E. japonicus East Chine Sea - L 2001 522 17 1.0 0.017 Liu et al. 2006 E. japonicus East Chine Sea - J 2001 522 20 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - R 2005 522 20 1.0 0.018 Liu et al. 2006 E. japonicus Sea japan - M 1997 522 18 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - I 1997 522 18 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - N 1996 522 18 1.0 0.019 Liu et al. 2006 E. japonicus Sea Japan - K 1997 522 18 1.0 0.018 Liu et al. 2006 E. japonicus Sea Japan - T 1997 522 18 1.0 0.020 Liu et al. 2006 E. japonicus Sea Japan - G 1997 522 18 1.0 0.016 Liu et al. 2006 E. japonicus Sea Japan - U 1997 522 18 1.0 0.020 This study E. ringens North-Peru (Paita) 02_2014 520 33 0.9 0.005 This study E. ringens Center-Peru (Callao) 02_2014 520 27 0.9 0.004 This study E. ringens South-Peru (Ilo) 02_2014 520 30 0.9 0.004 This study E. ringens North-Peru (Paita) 02_2015 520 39 0.9 0.005 This study E. ringens Center-Peru (Callao) 02_2015 520 32 0.9 0.004 This study E. ringens North-Peru (Salaverry) 08_2015 520 32 0.9 0.005 This study E. ringens Center-Peru (Callao) 08_2015 520 32 0.8 0.005 This study E. ringens Center-Peru (Pisco) 08_2015 520 40 0.9 0.005 This study E. ringens North-Peru (Paita) 04_2016 520 39 0.8 0.004 This study E. ringens Center-Peru (Callao) 04_2016 520 41 0.9 0.005 This study E. ringens North-Peru (Sechura) 04_2016 520 37 0.9 0.004 This study E. ringens Center-Peru (Callao) 10_2016 520 31 0.9 0.005 This study E. ringens South-Peru (Tacna) 10_2016 520 32 0.8 0.004 Table S6. Фst pairwise comparisons based on HVR1 sequences. Bold Фst values indicate significant genetic differences (Significance Level=0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 0.046* AS2014 -0.009 0.033 AP0215 -0.017 0.014 -0.013 AC0215 0.003 0.001 0.011 -0.011 AN0815 -0.021 0.028 -0.001 -0.015 -0.009 AC0815 -0.002 0.066* 0.009 -0.002 0.017 0.009 AS0815 0.022 -0.018 0.014 0.002 -0.004 0.013 0.045* AN0416 -0.021 0.059* -0.006 -0.008 0.011 -0.013 -0.005 0.036* AC0416 -0.014 0.018 -0.016 -0.016 -0.002 -0.011 0.015 0.003 -0.005 AN1016 -0.003 0.002 -0.008 -0.007 -0.003 -0.007 0.025 -0.004 0.004 -0.012 AC1016 0.024 -0.006 0.035 0.001 -0.020 0.007 0.032 0.000 0.035 0.014 0.011 AS1016 0.000 0.011 0.021 -0.015 -0.017 -0.008 0.011 0.002 0.013 0.004 0.009 -0.017 Table S7. P Values for exact tests of population differentiation based on HVR1 sequences (Markov chain length: 200000 steps). Bold P values indicate significant genetic differences (Significance Level=0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 0.215 AS2014 0.499 0.672 AP0215 0.998 0.459 0.879 AC0215 0.111 0.756 0.672 0.545 AN0815 0.922 0.279 0.164 0.841 0.238 AC0815 0.880 0.103 0.622 0.843 0.086 0.095 AS0815 0.736 0.536 0.636 0.748 0.297 0.224 0.473 AN0416 0.915 0.492 0.443 0.888 0.335 0.595 0.744 0.321 AC0416 0.917 0.550 0.858 0.960 0.505 0.456 0.401 0.545 0.608 AN1016 0.386 0.560 0.791 0.745 0.734 0.331 0.355 0.476 0.446 0.432 AC1016 0.189 0.775 0.459 0.646 0.822 0.533 0.146 0.455 0.331 0.509 0.680 AS1016 0.215 0.093 0.125 0.524 0.252 0.143 0.564 0.613 0.259 0.144 0.500 0.577 Table S8. Фst pairwise comparisons based on Cytb sequences. A hyphen (-) indicates samplings where no mtDNA sequences were generated. Bold Фst values indicate significant genetic differences (Significance Level=0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 - AS2014 - - AP0215 - - AC0215 - - - 0.014 AN0815 - - - 0.036* 0.004 AC0815 - - - 0.014 -0.008 0.000 AS0815 - - - - - - - AN0416 - - - 0.010 -0.001 0.016 0.010 - AC0416 - - - 0.016 0.005 0.022* -0.004 - 0.017 AN1016 - - - -0.013 0.002 0.027* 0.001 - 0.001 0.011 AC1016 - - - 0.003 -0.004 0.013 -0.007 - 0.006 -0.004 -0.005 AS1016 - - - - - - - - - - - - Table S9. P Values for exact tests of population differentiation based on Cytb sequences (Markov chain length: 200000 steps). Bold P values indicate significant genetic differences (Significance Level=0.05). AP2014 AC2014 AS2014 AP0215 AC0215 AN0815 AC0815 AS0815 AN0416 AC0416 AN1016 AC1016 AS1016 AP2014 AC2014 - AS2014 - - AP0215 - - - AC0215 - - - 0.291 AN0815 - - - 0.016* 0.470 AC0815 - - - 0.017* 0.402 0.407 AS0815 - - - - - - - AN0416 - - - 0.381 1.000 0.428 0.287 - AC0416 - - - 0.265 0.407 0.066 0.533 - 0.319 AN1016 - - - 0.874 0.390 0.033* 0.143 - 0.433 0.379 AC1016 - - - 0.218 1.000 0.252 0.549 - 0.741 0.504 0.394 AS1016 - - - - - - - - - - - - Table S10. Analysis of Molecular Variance (AMOVA) based on Cytb sequences. Genetic structure to test: North versus Central of Peru. Significance tests (20022 permutations). Percentage of Source of Degree of Sum of Variance AMOVA variance P Value variation freedom squares components statistics components Among North and 1 0.532 0.00091 0.25 Fct = 0.00249 0.26091 Central of Peru Among samples 6 2.544 0.00202 0.55 Fsc = 0.00551 0.11028 within North and center of Peru Within 228 83.086 0.36441 99.20 F samples st = 0.00798 0.05184 Total 235 86.161 0.36734 Referencias Bibliográficas Bandelt HJ, Forster P, Röhl A (1999) Median-joining networks for inferring intraspecific phylogenies. Molecular biology and evolution. Vol 16(1), pp 37-48. Barber RT, Chavez FP (1983) Biological consequences of el niño. Science. Vol 222(4629), pp 1203-1210. Bertrand A, Segura M, Gutiérrez M, Vásquez L (2004) From small‐scale habitat loopholes to decadal cycles: a habitat‐based hypothesis explaining fluctuation in pelagic fish populations off Peru. Fish and fisheries. Vol 5(4), pp 296-316. Chavez FP, Ryan J, Lluch-Cota SE, Ñiquen M (2003) From anchovies to sardines and back: multidecadal change in the Pacific Ocean. Science. Vol 299(5604), pp 217-221. Chavez FP, Messié M (2009). A comparison of eastern boundary upwelling ecosystems. Progress in Oceanography, 83(1-4), 80-96. Cubillos LA, Serra R, Fréon P (2007). Synchronous pattern of fluctuation in three anchovy fisheries in the Humboldt Current System. Aquatic Living Resources, 20(1), 69-75. Delrieu-Trottin E, Mona S, Maynard J, Neglia V, Veuille M, Planes S (2017). Population expansions dominate demographic histories of endemic and widespread Pacific reef fishes. Sci Rep 7: 40519. Díaz-Viloria N, Sánchez-Velasco L, Pérez-Enríquez R (2012). Recent population expansion in the evolutionary history of the Californian anchovy Engraulis mordax. Hidrobiológica, 22(3), 258-266. Excoffier L, Lischer HEL (2010) Arlequin suite ver 3.5: A new series of programs to perform population genetics analyses under Linux and Windows. Molecular Ecology Resources. Vol 10, pp 564-567. Ferrada S, Hernández K, Montoya R, Galleguillos R (2002) Estudio poblacional del recurso anchoveta (Engraulis ringens Jenyns 1842) (Clupeiformes, Engraulidae), mediante análisis de ADN. Gayana (Concepción). Vol 66(2), pp 243-248. Ferrada S, Galleguillos R, Canales CB, Herrera V (2018) Development and characterization of thirty-two microsatellite markers for the anchovy, Engraulis ringens Jenyns, 1842 (Clupeiformes, Engraulidae) via 454 pyrosequencing. Latin american journal of aquatic research. Vol 46(2), pp 452-456. Frankham R, Briscoe DA, Ballou JD (2002) Introduction to conservation genetics. Cambridge university press. Fu YX (1997) Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics. Vol 147(2), pp 915-925. Gyllensten U (1985). The genetic structure of fish: differences in the intraspecific distribution of biochemical genetic variation between marine, anadromous, and freshwater species. Journal of Fish Biology, 26(6), 691-699. Gutiérrez M, Castillo R, Segura M, Peraltilla S, Flores M. (2012) Trends in spatio-temporal distribution of Peruvian anchovy and other small pelagic fish biomass from 1966-2009. Latin American Journal of Aquatic Research. Vol 40(3 Esp.), pp 633-648. Gutiérrez M (1999) Distribución y abundancia de anchoveta y otras especies pelágicas entre los eventos El Nino 1982-83 y 1997-98. Bol. Inst. Mar Perú. Vol 18, pp 77-88. Grant, WS, & Leslie, RW (1993). Effect of metapopulation structure on nuclear and organellar DNA variability in senli-arid environments of southern Africa. South African Journal of Science, 89(6), 287-293. Grant WS, Bowen BW (1998) Shallow population histories in deep evolutionary lineages of marine fishes: insights from sardines and anchovies and lessons for conservation. Journal of heredity. Vol 89(5), pp 415-426. Hedgecock D, Hutchinson E, Sly FL, Nelson K (1994). The central stock of northern anchovy (Engraulis mordax) is not a randomly mating population. CalCoFi rep, 35, 121- 136. IMARPE (1973). Report of the Third Session of the Panel of Experts on the population dynamics of Peruvian anchoveta. Bol. Inst. Mar Peru-Callao 2(9):525-599. Jordán, R. 1971. Distribution of anchovy (Engraulis ringens) in relation to the environment. Inv. Pesq. Vol 35(1), pp 113-12 Jordán R, Málaga A (1972) Resultado de la primera marcación experimental de anchoveta (Engraulis ringens J.) En el mar. Abril 1970 - Setiembre 1971. Kritzer JP, Sale PF (2004) Metapopulation ecology in the sea: from Levins' model to marine ecology and fisheries science. Fish and Fisheries. Vol 5(2), pp 131-140. Larkin MA, Blackshields G, Brown NP, Chenna R, McGettigan PA, McWilliam H, et al. (2007) Clustal W and Clustal X version 2.0. Bioinformatics. Vol 23(21), pp 2947-2948. Lee WJ, Conroy J, Howell WH, Kocher TD (1995) Structure and evolution of teleost mitochondrial control regions. Journal of Molecular Evolution. Vol 41(1), pp 54-66. Leigh JW, Bryant D (2015) popart: full‐feature software for haplotype network construction. Methods in Ecology and Evolution. Vol 6(9), pp 1110-1116. Librado P, Rozas J (2009) DnaSP v5: A software for comprehensive analysis of DNA polymorphism data. Bioinformatics. Vol 25, pp 1451-1452. Lluch-Belda DRJM, Crawford RJ, Kawasaki T, MacCall AD, Parrish RH, Schwartzlose RA, et al. (1989) World-wide fluctuations of sardine and anchovy stocks: the regime problem. South African Journal of Marine Science. Vol 8(1), pp 195-205. Mathisen OA (1989) Adaptation of the anchoveta (Engraulis ringens) to the Peruvian upwelling system. In The Peruvian upwelling ecosystem: dynamics and interactions. ICLARM Conf Proc, Instituto del Mar del Peru (IMARPE). Vol 18, pp 220-234. Mendo J (1991) Stock identification of Peruvian anchoveta (Engraulis ringens): morphometric, tagging/recapture, electrophoretic and ecological studies. Ph.D. thesis, University of Bremen. Ñiquen M, Bouchon M, Cahuín S, Díaz E (2000) Pesquería de anchoveta en el mar peruano 1950-1999. Boletín IMARPE. Vol 19, pp 117-123. Ñiquen M, Bouchon M (2004) Impact of El Niño events on pelagic fisheries in Peruvian waters. Deep sea research part II: topical studies in oceanography, Vol 51(6-9), pp 563- 574. Pauly D (1987) The Peruvian anchoveta and its upwelling ecosystem: three decades of change WorldFish. Vol 15. Pannell JR, Charlesworth B (1999) Neutral genetic diversity in a metapopulation with recurrent local extinction and recolonization. Evolution. Vol 53(3), pp 664-676. Reed DH, Frankham R (2003) Correlation between fitness and genetic diversity. Conservation biology. Vol 17(1), pp 230-237. Rein B, Lückge A, Reinhardt L, Sirocko F, Wolf A, Dullo WC (2005) El Niño variability off Peru during the last 20,000 years. Paleoceanography. Vol 20(4). Rogers AR, Harpending H (1992) Population growth makes waves in the distribution of pairwise genetic differences. Molecular biology and evolution. Vol 9(3), pp 552-569. Rojas DE (2011) Evaluación de marcadores moleculares Ilps y Strs heterólogos en la anchoveta peruana Engraulis Ringens para estudios poblacionales. Tesis de pregrado. Universidad Nacional Mayor de San Marcos. Silva G, Horne JB, Castilho R (2014) Anchovies go north and west without losing diversity: Post‐glacial range expansions in a small pelagic fish. Journal of Biogeography. Vol 41(6), pp 1171-1182. Silva G, Cunha RL, Ramos A, Castilho R (2017). Wandering behaviour prevents inter and intra oceanic speciation in a coastal pelagic fish. Scientific reports, 7(1), 1-12. Slatkin M, Hudson RR (1991) Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations. Genetics. Vol 129(2), pp 555-562. Schwartzlose RA, Alheit J (1999) Worldwide large-scale fluctuations of sardine and anchovy populations. African Journal of Marine Science. Vol 21. Tajima F (1989) Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. Vol 123, pp 585-595. Tudela S, Garcıa-Marın JL, Pla C (1999). Genetic structure of the European anchovy, Engraulis encrasicolus L., in the north-west Mediterranean. Journal of Experimental Marine Biology and Ecology, 234(1), 95-109. Waples RS (1987). A multispecies approach to the analysis of gene flow in marine shore fishes. Evolution, 41(2), 385-400