Domanda:
Rimappare le coordinate genomiche per tenere conto degli indel
Greg Slodkowicz
2017-05-30 19:35:08 UTC
view on stackexchange narkive permalink

Sono interessato a ottenere sequenze codificanti del mio gene preferito in tutti gli individui dai 1000Genomes (e progetti simili). Uso GATK per ottenere il giusto sottoinsieme di varianti, vcf-consensus per mappare queste varianti sul genoma di riferimento e infine samtools per estrarre i singoli esoni. Funziona bene se le varianti sono SNP ma se ci sono indel, questo cambia le coordinate degli esoni e finisco per ottenere la regione sbagliata. Esiste un modo generico per rimappare le coordinate genomiche per tenere conto dei cambiamenti creati da indel?

Due risposte:
#1
+8
Karel Brinda
2017-05-30 20:23:34 UTC
view on stackexchange narkive permalink

Penso che tu abbia bisogno di un file LiftOver Chain per trasformare le tue coordinate. Puoi ottenere un file di questo tipo usando bcftools consensus con il parametro -c :

  -c, --chain <file> scrivi un file a catena for liftover  

Quindi puoi usarlo per trasformare le coordinate in vari formati genomici usando CrossMap.

Bello, +1. Sarebbe fantastico se potessi espandere questo per includere un comando di esempio che l'OP potrebbe eseguire, i passaggi necessari ecc. I siti Stack Exchange hanno linee guida molto rigide su cosa sia una risposta (ecco perché ti ho infastidito così) e idealmente , vogliamo che una risposta fornisca tutte le informazioni necessarie per risolvere la domanda posta. In altre parole, per favore non pubblicare suggerimenti (nemmeno quelli molto utili come questo) come risposte. Pubblica semplicemente un commento o trasforma il puntatore in una risposta completa.
Forse sono solo io, ma la sua risposta sembra a posto, anche per lo standard SE, ma sono solo io.
@nuin è d'accordo, motivo per cui ho votato. Penso solo che sarebbe ancora meglio con più dettagli, tutto qui.
Senza rancore.
#2
+2
finswimmer
2019-01-16 01:17:46 UTC
view on stackexchange narkive permalink

Se l'obiettivo è estrarre la sequenza di consenso per determinate regioni, è così che funziona nel 2019.

bgzip e indicizza il tuo vcf file.

   $ bgzip -c input.vcf > input.vcf.gz $  tabix input.vcf.gz  

Crea un region.txt che contiene una regione per riga nel formato chr: from-to . Se hai già un file bed , puoi utilizzare questo piccolo script awk per crearlo:

   $ awk '{print $ 1": "$ 2 + 1" - "$  3}' input.bed > region.txt  

Ora esegui questa combinazione di samtools e bcftools:

  $ samtools faidx -r region.txt genome.fa | bcftools consensus input.vcf.gz -o consensus. fa  


Questa domanda e risposta è stata tradotta automaticamente dalla lingua inglese. Il contenuto originale è disponibile su stackexchange, che ringraziamo per la licenza cc by-sa 3.0 con cui è distribuito.
Loading...