Cuando se acumulan diferentes versiones del mismo genoma, como pasa con la cebada, a menudo necesitaremos proyectar anotaciones de una versión a otra. Esta operación se llama lift-over en la literatura en inglés y tiene sus complicaciones, como se ven en la figura:
 |
| fuente: https://doi.org/10.12688/f1000research.14148.2 |
En una entrada anterior explicaba cómo hacerlo para genes, por ejemplo con el software LiftOff. Sin embargo, a veces lo que queremos mapear son SNPs, que se habían definido sobre una versión del genoma, sobre la siguiente.
Una manera, para genomas que tengan precalculados alineamientos en UCSC o Ensembl (chain files), es usar el software BCFtools/liftover, que se puede descargar como binario o compilar, y requiere bcftools 1.20 o superior. Puedes leer más sobre esta opción en https://doi.org/10.1093/bioinformatics/btae038 y https://github.com/freeseek/score. Una importante limitación es que solamente hay chain files pare ciertas especies. Por ejemplo, para plantas puedes consultar https://ftp.ebi.ac.uk/ensemblgenomes/pub/plants/current/assembly_chain
Para cualquier pareja de genomas podemos usar una estrategia que usábamos en Ensembl Plants, consiste en cortar la secuencia flanqueante de cada SNP en el genoma1 y mapearla sobre el genoma2 con BWA mem. Esta estrategia tiene como limitación que se pierde una fracción de las variantes originales, aquellas cuyas secuencias no mapeen bien en genoma2, o que estén en regiones repetidas, pero eso no es necesariamente malo. La ventaja que tiene es que no necesitas calcular alineamientos de dos genomas completos, lo cual es complejo y puede requerir grandes cantidades de RAM. Además en todo momento controlas lo que estás haciendo y si algo sale mal lo puedes ver y tratar de corregir. Esta estrategia se describe paso a paso en: https://github.com/eead-csic-compbio/eead-csic-compbio.github.io
Como resultado produce texto separado por tabuladores (TSV) cómo este (ver fichero completo):
1 51976 - LR890096.1 77101 C G 1 51988 - LR890096.1 77089 C G 1 51995 - LR890096.1 77082 G C 1 52015 - LR890096.1 77062 C G 1 263632 + LR890096.1 148230 G G 1 263634 + LR890096.1 148232 A A 1 263635 + LR890096.1 148233 A A 1 263637 + LR890096.1 148235 T T 1 263638 + LR890096.1 148236 G G 1 263646 + LR890096.1 148244 C C 1 263654 + LR890096.1 148252 C C 1 263699 + LR890096.1 148297 C C 1 263706 + LR890096.1 148304 A A 1 270084 + LR890096.1 154681 C C 1 270087 + LR890096.1 154684 G G
Un control de calidad posible es comprobar que la base de ambos genomas es la misma, aunque a veces estará un el reverso complementario, como se ve en el ejemplo para dos regiones de los cromomas 1 (genoma1) y LR890096.1 (genoma2).
Hasta pronto,
Bruno
