Muestreo
Se obtuvo el consentimiento informado por escrito del individuo saami cuyo genoma se analizó en este estudio que fue aprobado por el Comité Ético del Distrito Hospitalario de Helsinki y Uusimaa (Decisión 329/13/03/00/2013) y por el Comité Ético de la Universidad de Leipzig (Número de referencia de aprobación 398-13-16122013).
La toma de muestras y la extracción de ADN antiguo requiere un procedimiento estricto para evitar la contaminación introducida por material genético contemporáneo. En el caso de los 13 individuos de la Edad de Hierro de Finlandia de los que disponíamos, la toma de muestras tuvo lugar en unas instalaciones de sala blanca dedicadas al trabajo con ADN antiguo, en el Instituto de Ciencias Arqueológicas de Tübingen. El flujo de trabajo preliminar incluyó la documentación, la fotografía y el almacenamiento de las muestras en tubos de plástico y bolsas de plástico individuales con código de identificación. Como resultado de un primer estudio piloto, las muestras de dientes que utilizamos estaban fragmentadas y se eliminó parte de la dentina. El resto de la dentina se recogió separándola cuidadosamente del esmalte con un taladro de dentista y cabezas de broca de diamante enfriadas, giradas a una velocidad inferior a 15 rpm, para evitar posibles daños causados por el calor al ADN antiguo.
Para las muestras de los sitios de Bolshoy y Chalmny Varre, utilizamos polvo de diente sobrante que fue procesado originalmente en el Instituto de Antropología de la Universidad de Mainz con fines de réplica como se describe en la ref. 24. En resumen, los pasos de preparación de la muestra incluyeron la irradiación UV durante 30-45 minutos, seguida de un suave barrido de la superficie con blanqueadores comerciales diluidos. Los dientes fueron entonces chorreados con arena usando un abrasivo de óxido de aluminio (Harnisch & Rieth) y molidos hasta obtener un polvo fino usando un molino mezclador (Retsch).
Calibración de la fecha de radiocarbono
Calibramos la fecha de radiocarbono de Bolshoy, reportada en las refs 24,55 como 3473 ± 42 años BP, usando Intcal1356 como curva de calibración, usando OxCal 4.357.
Extracción de ADN y preparación de la biblioteca
El ADN de las seis muestras de Bolshoy y de las dos de Chalmny Varre fue extraído en las instalaciones de ADN antiguo del Instituto Max Planck para la Ciencia de la Historia Humana (MPI-SHH) en Jena, Alemania. La extracción de las muestras de Levänluhta se llevó a cabo de forma similar en las instalaciones de sala blanca del Instituto de Ciencias Arqueológicas de Tubinga. Para cada espécimen, se utilizaron ~ 50 mg de polvo de dentina para un procedimiento de extracción específicamente diseñado para la recuperación de ADN antiguo58. Se añadió al polvo de hueso un tampón de extracción que contenía 0,45 M de EDTA, pH 8,0 (Life Technologies) y 0,25 mg/ml de proteinasa K (Sigma-Aldrich) y se incubó a 37 °C con rotación durante toda la noche. El sobrenadante se separó del pellet de polvo de hueso por centrifugación (14.000 rpm). Se añadió al sobrenadante un tampón de unión consistente en GuHCL 5 M (Sigma Aldrich) y 40% de Isopropanol (Merck), junto con 400 μl de acetato de sodio 1 M (pH 5,5), y se purificó la solución haciéndola girar a través de una columna de purificación unida a un embudo de montaje de extensor de alta pureza (8 min. en 1500 rpm, con aceleración lenta). A continuación, se hizo girar la columna en un tubo de recogida (1 min. a 14.000 rpm) 1-2 veces para maximizar el rendimiento. A esto le siguieron dos pasos de lavado posteriores de 450 μl de tampón de lavado (High Pure Viral Nucleic Acid Large Volume Kit) y dos pasos de centrifugación en seco de 1 min a 14.000 rpm. El volumen total final de 100 μl de eluido se alcanzó mediante dos rondas de elución separadas de 50 μl de TET (10 mM Tris-HCL, 1 mM EDTA pH 8,0, 0,1% Tween20), cada una de ellas centrifugada durante 1 min a 14.000 rpm en un tubo Eppendorf nuevo de 1,5 ml. Los controles negativos (tampón en lugar de muestra) se procesaron en paralelo en una proporción de 1 control por cada 7 muestras.
De los 100 μl de extracto, se utilizaron 20 μl para inmortalizar el ADN de la muestra como una biblioteca de doble cadena. El procedimiento incluyó un paso de reparación del extremo romo, de ligación del adaptador y de relleno del adaptador, tal y como describen Meyer y Kircher59. Durante el paso de reparación del extremo romo, se utilizó una mezcla de 0,4 U/μl de T4 PNK (polinucleótido quinasa) y 0,024 U/μl de T4 ADN polimerasa, 1× NEB buffer 2 (NEB), 100 μM de mezcla de dNTP (Thermo Scientific), 1 mM de ATP (NEB) y 0.Se añadió 8 mg/ml de BSA (NEB) al ADN molde, seguido de la incubación en un termociclador (15 min 15 °C, 15 min 25 °C) y la purificación con un kit MinElute (QIAGEN). El producto se eluyó en 18 μl de tampón TET. El paso de ligación del adaptador incluyó una mezcla de 1× Quick Ligase Buffer (NEB), 250 nM de adaptadores Illumina (Sigma-Aldrich) y 0,125 U/μl de Quick Ligase (NEB), añadidos al eluido de 18 μl, seguido de una incubación de 20 min, y un segundo paso de purificación con columnas MinElute, esta vez en 20 μl de eluido. Para el paso de relleno, se añadió una mezcla de 0,4 U/μl de Bst-polimerasa y 125 μM de mezcla de dNTP, y a continuación se incubó la mezcla en un termociclador (30 min 37 °C, 10 min 80 °C). Se produjeron bibliotecas sin tratamiento con Uracil-DNA-glicosilasa (UDG) para los 13 extractos de Levänluhta. Además, se produjeron bibliotecas con tratamiento de UDG para siete de los 13 extractos originales de Levänluhta, y para todos los extractos de Bolshoy y Chalmny Varre. Para introducir el tratamiento UDG-half, se incluyó una etapa inicial en la preparación de la biblioteca, en la que se añadieron 250 U de enzima USER (NEB) en los 20 μl de extracto, seguida de una incubación a 37 °C durante 30 min, y luego a 12 °C durante 1 min. A esto le siguió de nuevo la adición de 200 UGI (inhibidor de la Uracil Glicosilasa, de NEB) y otra incubación idéntica para detener la escisión enzimática de los sitios desaminados, como se describe en60. Para cada biblioteca, se incorporó un par único de índices de ocho pb de longitud utilizando una Pfu Turbo Cx Hotstart DNA Polymerase y un programa de termociclaje con el siguiente perfil de temperatura: desnaturalización inicial (98 °C durante 30 s), ciclo de desnaturalización/reducción/extensión (98 °C durante 10 s/ 60 °C durante 20 s/ 72 °C durante 20 s) y extensión final a 72 °C durante 10 min61. El polvo de hueso de un oso de las cavernas se procesó en paralelo y sirvió como control positivo. Los controles negativos para las etapas de extracción y preparación de bibliotecas se mantuvieron junto a las muestras durante todo el flujo de trabajo.
La eficiencia del experimento se garantizó mediante la cuantificación de la concentración de las bibliotecas en qPCR (Roche) utilizando alícuotas de las bibliotecas antes y después de la indexación. El número de copias moleculares en las bibliotecas preindexadas osciló entre ~ 10E8 y ~ 10E9 copias/μl, lo que indica una preparación satisfactoria de la biblioteca, mientras que las bibliotecas indexadas oscilaron entre ~ 10E10 y ~ 10E12 copias/μl, lo que indica una eficiencia de indexación admisible. Los controles negativos mostraron una concentración de 4-5 órdenes de magnitud inferior a la de las muestras, lo que indica bajos niveles de contaminación procedentes de las etapas de procesamiento en el laboratorio.
Las bibliotecas se amplificaron con PCR, para la cantidad de ciclos correspondientes a las concentraciones de las bibliotecas indexadas, utilizando la polimerasa AccuPrime Pfx (5 μl de plantilla de biblioteca, 2 U de AccuPrime Pfx DNA polimerasa de Invitrogen, 1 U de mastermix 10× de PCR preparada y 0.3 μM de cebadores IS5 e IS6, por cada 100 μl de reacción) con un perfil térmico de 2 min de desnaturalización a 95 °C, 3-9 ciclos consistentes en 15 seg de desnaturalización a 95 °C, 30 seg de recocido a 60 °C, 2 min de elongación a 68 °C y 5 min de elongación a 68 °C. Las bibliotecas amplificadas se purificaron utilizando columnas de centrifugado MinElute con el protocolo estándar proporcionado por el fabricante (Qiagen), y se cuantificaron para la secuenciación utilizando un chip Agilent 2100 Bioanalyzer DNA 1000.
Para el individuo moderno saami, el ADN total fue extraído con fenol-cloroformo y esquilado físicamente utilizando la fragmentación COVARIS. Se realizó una preparación de biblioteca Illumina modificada utilizando la reparación del extremo romo seguida de la cola A del extremo 3′ y la ligación de adaptadores bifurcados. La PCR de indexación fue seguida por la escisión de fragmentos de entre 500 y 600 pb de un gel de agarosa al 2%.
Captura & secuenciación
Utilizamos el procedimiento de captura en solución de la ref. 62 para enriquecer nuestras bibliotecas en busca de fragmentos de ADN que se solaparan con 1.237.207 posiciones variables del genoma humano4. Las secuencias utilizadas como cebo, unidas a perlas magnéticas, se mezclaron con la plantilla de la muestra de ADN en solución, y se dejaron hibridar con el ADN objetivo durante una incubación de 24 horas a 60 °C en un horno de rotación. Se agruparon 4-6 muestras en proporciones de masa iguales con un total de 2000 ng de ADN. Las bibliotecas capturadas fueron secuenciadas (75 bp single-end, más paired-end adicional para las tres bibliotecas no-UDG de los individuos de Levänluhta) en una plataforma Illumina HiSeq 4000 en el Instituto Max Planck para la Ciencia de la Historia Humana en Jena. De las 13 muestras de la Edad de Hierro de Finlandia procesadas originalmente, siete resultaron ser de una calidad adecuada para ser utilizadas en los análisis posteriores. El genoma moderno de los saami fue secuenciado en un Genome Analyser II (8 carriles, 125 bp paired-end) en el Instituto Max Planck de Antropología Evolutiva en Leipzig.
Procesamiento de las lecturas secuenciadas
Utilizamos EAGER63 (versión 1.92.50) para procesar las lecturas secuenciadas, utilizando los parámetros por defecto (ver más abajo) para los datos de secuenciación de origen humano, tratados con UDG-half, de un solo extremo, al procesar las bibliotecas UDG-half de todos los individuos. Específicamente, se utilizó AdapterRemoval para recortar los adaptadores de secuenciación de nuestras lecturas, con un solapamiento mínimo de 1 pb, y utilizando una calidad de base mínima de 20 y una longitud de secuencia mínima de 30 pb. Se utilizó BWA aln (versión 0.7.12-r1039, https://sourceforge.net/projects/bio-bwa/files)64 para mapear las lecturas a la secuencia de referencia humana hs37d5, con una longitud de semilla (-l) de 32, un número máximo de diferencias (-n) de 0,01 sin hacer ningún filtro de calidad. La eliminación de duplicados se realizó con DeDup v0.12.1. El cálculo del daño por desaminación de las bases terminales se realizó con mapDamage, especificando una longitud (-l) de 100 pb (Tabla Suplementaria 1).
Para los análisis posteriores, utilizamos bamutils (versión 1.0.13, https://github.com/statgen/bamUtil.git) TrimBam para recortar dos bases al principio y al final de todas las lecturas. Este procedimiento elimina las posiciones que se ven afectadas por la desaminación, eliminando así los errores de genotipado que podrían surgir debido a un daño en el ADN antiguo.
Para tres individuos de Levänluhta que superaron el umbral de cobertura del 1% en el cribado preliminar, utilizamos las bibliotecas no tratadas con UDG para confirmar la autenticidad de los datos antiguos. Para estas bibliotecas no tratadas, se llevaron a cabo dos rondas de secuenciación, que se procesaron utilizando EAGER con los parámetros anteriores, pero especificando un modo de tratamiento no-UDG y estableciendo el tipo de secuenciación correcto entre las bibliotecas. Las lecturas fusionadas se extrajeron de los archivos bam resultantes, y se fusionaron con el archivo bam que contenía las lecturas del recorrido de la secuencia de extremo único utilizando samtools merge (versión 1.3)65.
El genoma moderno de Saami se generó utilizando Ibis para la llamada de bases y un script propio de recorte de adaptadores. Las lecturas resultantes se alinearon con el genoma humano de referencia hs37d5 utilizando bwa 0.5.9-r16 (parámetros -e 20 -o 2 -n 0.01).
Genotipado
Utilizamos un programa personalizado (pileupCaller) para genotipar los 15 individuos antiguos. Se generó un archivo pileup utilizando samtools mpileup con los parámetros -q 30 -Q 30 -B que contenía sólo los sitios que se solapaban con nuestro panel de captura. A partir de este archivo, para cada individuo y cada SNP del panel 1240K, se extrajo al azar una lectura que cubría el SNP, y se hizo una llamada pseudohaploide, es decir, se asumió que el individuo antiguo era homocigoto para el alelo en la lectura extraída al azar para el SNP en cuestión. PileupCaller está disponible en https://github.com/stschiff/sequenceTools.git.
Para las tres bibliotecas de Levänluhta que no fueron sometidas a tratamiento con UDG, sólo genotipificamos las transversiones, eliminando así los artefactos de daño post-mortem C- > T de los análisis posteriores.
El genoma shotgun del individuo saami moderno fue genotipado utilizando GATK (versión 1.3-25-g32cdef9) Unified Genotyper después de la realineación de indeles. Las llamadas de variantes se filtraron en busca de variantes con una puntuación de calidad superior a 30, y se utilizó un script personalizado para convertir las variantes en formato EigenStrat.
Los datos se fusionaron con un gran conjunto de datos que consistía en 3871 individuos antiguos y modernos genotipados en las matrices Human Origins y/o 1240K SNP, utilizando mergeit.
Determinación del sexo
Para determinar el sexo genético de cada individuo antiguo calculamos la cobertura en los autosomas así como en cada cromosoma sexual. Utilizamos un script personalizado (https://github.com/TCLamnidis/Sex.DetERRmine) para el cálculo de cada cobertura relativa así como sus barras de error asociadas (Figura Suplementaria 1, Nota Suplementaria 3 para más información sobre el cálculo de errores). Se espera que las hembras tengan una tasa x de 1 y una tasa y de 0, mientras que se espera que los machos tengan tanto una tasa x como una tasa y de 0,5 (ref. 49).
Autenticación
En primer lugar confirmamos que el patrón de desaminación en las bases terminales de las lecturas de ADN era característico del ADN antiguo (1.04-4,5% para las bibliotecas no-UDG, y 4,7-9,5% para las bibliotecas no-UDG, ver Tabla Suplementaria 1), utilizando mapDamage (versión 2.0.6)66. Realizamos una serie de pruebas diferentes para asegurar la autenticidad de nuestros datos antiguos. Para los individuos masculinos, investigamos los polimorfismos en el cromosoma X27 utilizando el paquete de software ANGSD (versión 0.910)67 . Esto reveló estimaciones robustas de contaminación para 2 individuos masculinos de Bolshoy, y 1 individuo masculino de Chalmny-Varre. Todos ellos estaban por debajo del 1,6% de contaminación (Tabla 1). Para los individuos femeninos de estos dos sitios, observamos que se proyectan cerca de los machos en el espacio PCA (Fig. 2a, Figura Suplementaria 3), sugiriendo efectos limitados de contaminación potencial. Además, generamos un conjunto de datos filtrados por PMD para todos los individuos utilizando pmdtools (versión 0.60)30. El filtrado PMD se realizó utilizando un genoma de referencia enmascarado para todas las posiciones del panel de captura 1240K, para evitar sesgos alélicos sistemáticos en las posiciones de SNP analizadas. Establecimos un umbral de pmd de 3, que, según la publicación original30, elimina eficazmente los posibles contaminantes modernos basándose en la ausencia de modificaciones de base consistentes con la desaminación.
Para proporcionar una estimación más cuantitativa de la posible contaminación en las mujeres, utilizamos el programa ContamMix (versión 1.0-10)29 para estimar la contaminación mitocondrial. Extrajimos las lecturas que mapeaban la referencia mitocondrial para cada uno de los individuos antiguos utilizando samtools (versión 1.3)65 . A continuación, generamos una secuencia consenso mitocondrial para cada uno de los individuos antiguos utilizando Geneious (versión 10.0.9, http://www.geneious.com,68), y llamando N a todos los sitios con una cobertura inferior a 5. Finalmente, todas las lecturas mitocondriales se alinearon con su respectiva secuencia consenso, utilizando bwa aln (versión 0.7.12-r1039)64 con un número máximo de diferencias en la semilla (-k) fijado en 5 y el número máximo de diferencias (-n) en 10, y bwa samse (versión 0.7.12-r1039)64. Se generó un alineamiento múltiple de la secuencia consenso y un conjunto de referencia de 311 genomas mitocondriales69 , utilizando mafft (versión v7.305b)70,71,72 con el parámetro –auto. La alineación de las lecturas, así como la alineación múltiple del consenso y de los 311 genomas mitocondriales de referencia se proporcionaron a ContamMix. Informamos aquí del modo de contaminación a posteriori, junto con los límites superior e inferior del intervalo posterior del 95 % (Tabla 1).
Para una autentificación adicional, ejecutamos ADMIXTURE28 supervisado (versión 1.3.0) para todas las muestras utilizando las seis poblaciones actuales (Atayal, French, Kalash, Karitiana, Mbuti y Papuan) como grupos genéticos definidos, para localizar cualquier diferencia grande en la agrupación genética entre los individuos del mismo sitio (Figura suplementaria 2). Probamos el poder de este método para detectar la contaminación y descubrimos que puede detectar la contaminación que está lejanamente relacionada con las ascendencias presentes dentro de los individuos de la prueba ya en tasas del 5-8%, pero carece del poder para identificar la contaminación estrechamente relacionada con los individuos de la prueba (véase la Nota Suplementaria 2). No observamos diferencias significativas (dentro de nuestra resolución) en los patrones de ascendencia entre los individuos antiguos del mismo sitio, con la excepción de Levänluhta, donde la muestra individual JK2065 parece derivar de una ascendencia diferente. Por lo tanto, le asignamos una etiqueta de población separada, Levänluhta_B en este estudio.
Por último, utilizando smartpca, proyectamos los conjuntos de datos filtrados por PMD y no filtrados en el mismo conjunto de componentes principales construidos en las poblaciones europeas modernas, para asegurar que los individuos antiguos permanezcan proyectados en las posiciones aproximadamente equivalentes independientemente del filtrado por PMD. Esto fue posible para todas las muestras con el tratamiento UDG-mitad, excepto para los individuos de Levänluhta, que representaban un daño demasiado pequeño para aplicar el filtrado PMD. En lo que respecta a este lugar, nos basamos por tanto en las bibliotecas sin UDG (utilizando únicamente las transversiones) que se generaron para los tres individuos utilizados en el análisis principal. Encontramos que, dentro del ruido esperado debido al bajo número de SNPs, todas las muestras muestran consistencia entre los conjuntos de datos filtrados y no filtrados, sugiriendo una baja cantidad de contaminación en todas las muestras (Figura Suplementaria 3a, b). Cuatro individuos adicionales de Levänluhta fueron excluidos del análisis principal y de esta prueba de autentificación debido a la baja cobertura (< 15.000 SNPs cubiertos) y a la falta de bibliotecas no-UDG.
Estadística F
Todos los programas utilizados para calcular la estadística F en este estudio se pueden encontrar como parte del paquete Admixtools (https://github.com/DReichLab/AdmixTools)2,42.
Utilizamos qp3Pop (versión 412) para todos los cálculos de F3.
Los estadísticos de F4 se calcularon utilizando qpDstat (versión 711), y qpAdm (versión 632)2 se utilizó para estimar las proporciones de la mezcla utilizando lo siguiente: Fuentes (poblaciones de la izquierda): Nganasan; WHG; EHG; Yamnaya_Samara; LBK_EN. Grupos externos (poblaciones de la derecha): OG1: Mbuti; CHG; Israel_Natufian; Onge; Villabruna; Ami; Mixe. OG2: Mbuti; CHG; Onge; Villabruna; Ami; Mixe. OG3: Mixe; CHG; Israel_Natufian; Villabruna; Onge; Ami. OG4: Mbuti; Israel_Natufian; Onge; Villabruna; Ami; Mixe. OG5: Mbuti; Samara_HG; CHG; Israel_Natufian; Villabruna; Ami.
Para garantizar que los conjuntos de grupos externos tuvieran suficiente poder para distinguir las ascendencias presentes en las fuentes, ejecutamos qpWave (versión 410) utilizando sólo las fuentes como poblaciones de la izquierda y cada conjunto de grupos externos como derechos. Todas estas ejecuciones de qpWave fueron consistentes sólo con el rango máximo, lo que significa que todos los conjuntos de grupos externos tenían suficiente poder para distinguir entre las cinco fuentes diferentes. Todos los modelos qpWave y qpAdm se ejecutaron utilizando la opción allsnps: SÍ. Cuando los modelos de mezcla de cinco vías proporcionados por qpAdm tenían valores p superiores a 0,05, pero incluían proporciones de mezcla no factibles y a una de las fuentes se le asignaba una proporción de mezcla negativa, ejecutamos el modelo de nuevo con esa fuente excluida. Para cada población de prueba, si el conjunto de grupos externos OG1 no producía un modelo completo que funcionara (p < 0,05), probamos conjuntos alternativos de grupos externos con una población derecha eliminada. Esto dio como resultado los conjuntos de grupos externos OG2-4. En el caso de Levänluhta, múltiples conjuntos de grupos externos produjeron modelos que funcionaron, los cuales están listados en los Datos Suplementarios 4. Siempre se muestra el modelo de mezcla que necesita el mínimo número de fuentes y que proporciona proporciones de mezcla factibles. En el caso de PWC de Suecia, donde ninguno de los conjuntos de grupos externos OG1-4 produjo un modelo de trabajo, se utilizó un conjunto revisado de poblaciones correctas (OG5) que incluye Samara_HG para proporcionar más poder para distinguir las ascendencias de cazadores-recolectores. Preferimos los modelos con OG1-4 sobre OG5 en general, porque OG5 contiene más genomas antiguos con potenciales sesgos en las poblaciones correctas, lo que causa más a menudo el fracaso de los modelos para las poblaciones modernas de Test. Las fuentes excluidas en los modelos mínimos se especificaron como N/A (Datos Suplementarios 4). Si se pudiera excluir cualquiera de los dos Yamnaya o EHG (como es el caso de Levänluhta), mostramos el modelo que es más consistente con publicaciones anteriores3,7,8,45 en la Fig. 4, pero mostramos ambos modelos en los Datos Suplementarios 4.
Análisis de componentes principales
Utilizamos smartpca (versión #16000)73 (https://github.com/DReichLab/EIG) para llevar a cabo el Análisis de Componentes Principales (PCA), utilizando el lsqproject: YES y shrinkmode: YES.
Para el PCA euroasiático (Fig. 2a), se utilizaron las siguientes poblaciones para construir componentes principales: Abkhasian, Adygei, Albanian, Altaian, Ami, Armenian, Atayal, Avar.SG, Azeri_WGA, Balkar, Balochi, Basque, BedouinA, BedouinB, Belarusian, Borneo, Brahui, Bulgarian, Buryat.SG, camboyano, canario, checheno, chuvash, croata, chipriota, checo, dai, daur, dolgan, druso, inglés, estonio, even, finlandés, francés, georgiano, griego, han, hazara, hezhen, húngaro, islandés, iraní, italiano_norte, italiano_sur, japonés, judío_shkenazi, judío_georgiano, judío_iraní, judío_iraquí, judío_libio, judío_marroquí, judío_tunecino, judío_turco, judío_yemenita, jordano, kalash, kalmyk, kinh, coreano, kumyk, kurd_WGA, kirguís, lahu, libanés, lezgin, lituano, makrani, mala, maltés, mansi, miao, mongola, mordoviano, naxi, nganasan, nogai, noruego.DG, noruego, orcadiano, oroqen, palestino, patán, ruso, saami.DG, saami_WGA, sardo, saudí, escocés, selkup, semende, she, sherpa.DG, siciliano, español, español_norte, sirio, tayiko, tailandés, tibetano.DG, Tu, Tubalar, Tujia, Turco, Turcomano, Tuviniano, Ucraniano, Ulchi, Uygur, Uzbeko, Xibo, Yakut, Yi, Yukagir.
Para el PCA de Eurasia Occidental (Figura Suplementaria 3a, b), se utilizaron las siguientes poblaciones para construir componentes principales: Abkhasian, Adygei, Albanian, Armenian, Balkar, Basque, BedouinA, BedouinB, Belarusian, Bulgarian, Canary_Islander, Chechen, Chuvash, Croatian, Cypriot, Czech, Druze, English, Estonian, Finnish, French, Georgian, Greek, Hungarian, Icelandic, Iranian, Italian_North, Italian_South, Jew_Ashkenazi, Jew_Georgian, judío_iraní, judío_iraquí, judío_libio, judío_marroquí, judío_tunecino, judío_turco, judío_yemenita, jordano, kumyk, libanés, lezgin, lituano, maltés, mordoviano, noruego, orcadiano, palestino, polaco, ruso, sardo, saudí, escocés, español, español_norte, sirio, turco, ucraniano.
Análisis ADMIXTURE
ADMIXTURE28 se ejecutó con la versión 1.3.0, tras la exclusión de las variantes con una frecuencia alélica menor de 0,01 y después de la poda LD utilizando plink (versión 1.90b3.29)74 con un tamaño de ventana de 200, un tamaño de paso de 5 y un umbral R2 de 0,5 (https://www.genetics.ucla.edu/software/admixture/download.html). Se ejecutaron cinco réplicas para cada valor de K, con valores de K que oscilaban entre 2 y 15. Las poblaciones utilizadas fueron: Ami, Ami.DG, armenio, Atayal, Atayal.DG, balochi, vasco, beduinoB, bielorruso, brahmán_tiwari, brahui, chuvash, croata, chipriota, checo, inglés, estonio, even, finlandés, finés.DG, francés, griego, gujaratiB, hadza, han, húngaro, islandés, kalash, karitiana, lituano, makrani, mala, mansi, mansi.DG, mari.SG, mbuti, mbuti.DG, mixe, mordoviano, nganasan, noruego, onge, orcadiano, papú, pima, ruso, saami.DG, ModernSaami, Sardinian, Scottish, Selkup, Spanish, Ukrainian, Ulchi, Yoruba, ALPC_Hungary_MN, Baalberge_MN, Baltic_BA, Baltic_CCC, Baltic_CWC, Baltic_LN, BolshoyOleniOstrov, Bu_kk_Culture_Hungary_MN, ChalmnyVarre, CHG, EHG, Esperstedt_MN, Ganj_Dareh_Iran_Neolítico, Hungary_MN, Hungary_Neolithic, Iran_Chalcolithic, JK2065, Koros_Hungary_EN, Kunda, Latvia_HG3, Latvia_MN1, LBK_EN, LBK_Hungary_EN, Levanluhta, Narva, PWC_Sweden_NHG.SG, Scandinavia_LNBA, SHG, Sweden_HG.SG, TRB, Ukraine_HG1, Ukraine_N1, WHG, Yamnaya_Samara.
Encontramos que K = 11 da como resultado el menor error de Validación Cruzada, como se muestra en la Figura Suplementaria 4b.
Haplotipos cromosómicos Y
Asignamos a los varones antiguos a los haplogrupos Y utilizando el programa yHaplo (https://github.com/23andMe/yhaplo)75. En resumen, este programa proporciona una búsqueda automatizada a través del árbol de haplogrupos Y (tal y como se proporciona dentro de yHaplo, tal y como se accedió desde ISOGG el 04 de enero de 2016) desde la raíz hasta la rama descendente basada en la presencia de alelos derivados y asigna el haplogrupo más descendente con alelos derivados. Para unos 15.000 SNP cromosómicos Y presentes tanto en nuestro panel de captura como en dos conjuntos de datos publicados76,77, muestreamos aleatoriamente una sola base y la utilizamos como genotipo haploide. Utilizamos un script personalizado para convertir los genotipos EigenStrat al formato yHaplo. Informamos del haplogrupo asignado más abajo proporcionado por el programa (Tabla 1). También comprobamos manualmente el estado derivado y la ausencia de mutaciones que definen el haplogrupo designado, ya que la falta de información podría conducir a una parada prematura en su búsqueda automatizada.
Haplotipado mitocondrial
Importamos las lecturas mitocondriales recortadas para cada individuo con calidad de mapeo >30 en Geneious (versión 10.0.9, https://www.geneious.com)68 y reensamblamos estas lecturas con el genoma de referencia RSRS78, utilizando el mapeador de Geneious, con sensibilidad media y 5 iteraciones. Utilizamos el llamador de variantes automatizado incorporado en Geneious para encontrar polimorfismos mitocondriales con una cobertura mínima de 3 y una frecuencia mínima de variantes de 0,67. Las variantes resultantes se exportaron a Excel y se compararon manualmente con los SNPs reportados en la filogenia de ADNmt en línea (árbol de ADNmt Build 17, 18 Feb 2016, http://www.phylotree.org/). Las posiciones de nucleótidos 309.1 C(C), 315.1C, indels AC en 515-522, 16182C, 16183C, 16193.1C(C) y 16519 fueron enmascaradas y no se incluyeron en nuestras llamadas de haplotipos.
SNPs fenotípicos
Utilizamos samtools mpileup (versión 1.3)65, filtrando la calidad del mapa (-Q) y de la base (-q) de 30, desactivando la calidad de la alineación por base (-B), en los archivos bam recortados, para generar un pileup de lecturas que se mapean a un conjunto de 43 SNPs fenotípicos4,40,41,79 que forman parte de nuestro panel de captura del genoma. Se utilizó un script personalizado de python para analizar la pileup en una tabla que contiene el número de lecturas que apoyan cada alelo (Datos Suplementarios 2).
Disponibilidad del código
Todo el software descrito por primera vez en este estudio está disponible de forma gratuita en los repositorios en línea. Sex.DetERRmine: https://github.com/TCLamnidis/Sex.DetERRmine
ContaminateGenotypes: https://github.com/TCLamnidis/ContaminateGenotypes