2009-10-09 13 views
15

tengo un dato en que siempre viene en bloque de cuatro en el siguiente formato (llamado FASTQ):conversión de FASTQ a FASTA con SED/AWK

@SRR018006.2016 GA2:6:1:20:650 length=36 
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNGN 
+SRR018006.2016 GA2:6:1:20:650 length=36 
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!+! 
@SRR018006.19405469 GA2:6:100:1793:611 length=36 
ACCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
+SRR018006.19405469 GA2:6:100:1793:611 length=36 
7);;).;);;/;*.2>/@@7;@77<..;)58)5/>/ 

¿Hay un simple awk/modo SED/bash para convertirlos en este formato (llamado FASTA):

>SRR018006.2016 GA2:6:1:20:650 length=36 
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNGN 
>SRR018006.19405469 GA2:6:100:1793:611 length=36 
ACCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 

en principio queremos extraer las dos primeras líneas de cada bloque de 4- y reemplazar @ con >.

+0

Bien, acabo de recibir un dolor de cabeza. – Homework

Respuesta

19

Esta es una vieja pregunta, y se han ofrecido muchas soluciones diferentes. Como la respuesta aceptada usa sed pero tiene un problema evidente (que es que reemplazará a @ con> cuando el signo @ aparece como la primera letra de la línea de calidad), me siento obligado a ofrecer una solución sencilla basada en sed que realmente funcione :

sed -n '1~4s/^@/>/p;2~4p' 

el único supuesto hecho es que cada lectura ocupa exactamente 4 líneas en el archivo FASTQ, pero que parece bastante seguro, en mi experiencia.

El script fastq_to_fasta en el kit de herramientas fastx también funciona. (Vale la pena mencionar que necesitas especificar la opción -Q33 para acomodar las codificaciones de Phred + 33 qual ahora comunes. Lo cual es gracioso, ya que está desperdiciando los datos de calidad de todos modos!)

+2

¡El convertidor más rápido que he visto hasta ahora! ¡10 millones de lecturas en 14s! –

+2

Esto es oro puro. – darxsys

+1

¡Gran solución! –

1
awk 'BEGIN{P=1}{if(P==1||P==2){gsub(/^[@]/,">");print}; if(P==4)P=0; P++}' data 

>SRR018006.2016 GA2:6:1:20:650 length=36 
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNGN 
>SRR018006.19405469 GA2:6:100:1793:611 length=36 
ACCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 

continuación

awk '{gsub(/^[@]/,">"); print}' data 

donde los datos es el archivo de datos. que he recibido:

>SRR018006.2016 GA2:6:1:20:650 length=36 
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNGN 
+SRR018006.2016 GA2:6:1:20:650 length=36 
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!+! 
>SRR018006.19405469 GA2:6:100:1793:611 length=36 
ACCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
+SRR018006.19405469 GA2:6:100:1793:611 length=36 
7);;).;);;/;*.2>/@@7;@77<..;)58)5/>/ 
+1

¡esto es correcto! – bua

1

Aquí está la solución a la parte de "saltar cada dos líneas" del problema que acabo de aprender de SO:

while read line 
do 
    # print two lines 
    echo "$line" 
    read line_to_print 
    echo "$line_to_print" 

    # and skip two lines 
    read line_to_skip 
    read line_to_skip 
done 

Si todo lo que hay que hacer es cambiar una @->, entonces supongo

while read line 
do 
    echo "$line" | sed 's/@/>/' 
    read line 
    echo "$line" 

    read line_to_skip 
    read line_to_skip 
done 

hará el trabajo.

+0

debe haber una redirección de entrada para el archivo de entrada. para reemplazar personajes en bash, $ {line/@ /}} debería ser suficiente. no es necesario sed – ghostdog74

1

Algo así como:

awk 'BEGIN{a=0}{if(a==1){print;a=0}}/^@/{print;a=1}' myFastqFile | sed 's/^@/>/' 

debería funcionar.

+0

ya que está usando awk ya, no hay necesidad de perder un proceso extra llamando a sed. hacer la sustitución dentro de awk. – ghostdog74

+0

Tienes razón. He votado a favor su respuesta. – mouviciel

1

creo, con grep de GNU esto se podría hacer con esto:

grep -A 1 "^@" t.txt | grep -v "^--" | sed -e "s/^@/\>/" 
+0

si el archivo es muy grande, no tiene sentido encadenar greps y juntarse así. – ghostdog74

7

simplemente awk, sin necesidad de otras herramientas

# awk '/^@SR/{gsub(/^@/,">",$1);print;getline;print}' file 
>SRR018006.2016 GA2:6:1:20:650 length=36 
NNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNGN 
>SRR018006.19405469 GA2:6:100:1793:611 length=36 
ACCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
9

sed no está muerto. Si estamos golf:

sed '/^@/!d;s//>/;N' 

O, emulando http://www.ringtail.tsl.ac.uk/david-studholme/scripts/fastq2fasta.pl publicado por Pierre, que sólo imprime la primera palabra (el identificador) de la primera línea y lo hace (a) el tratamiento de errores:

#!/usr/bin/sed -f 
# Read a total of four lines 
$b error 
N;$b error 
N;$b error 
N 
# Parse the lines 
/^@\(\([^ ]*\).*\)\(\n[ACGTN]*\)\n+\1\n.*$/{ 
    # Output id and sequence for FASTA format. 
    s//>\2\3/ 
    b 
} 
:error 
i\ 
Error parsing input: 
q 

Parece que hay muchas herramientas existentes para convertir estos formatos; probablemente deberías usar estos en lugar de todo lo publicado aquí (incluido el anterior).

+1

sed está muy vivo, pero la solución sed ofrecida aquí probablemente hundirá su flujo de trabajo. No puede confiar en el carácter @ para indicar líneas de encabezado de manera exclusiva; las líneas de calidad también pueden comenzar con @. Por favor, mira mi solución a continuación. – Owen

9

Como se detalla en Cock, et al (2009) NAR, muchas de estas soluciones son incorrectas ya que "el carácter de marcador '@' (ASCII 64) puede aparecer en cualquier parte de la cadena de calidad. una línea que comience con '@' como indicando el inicio del próximo registro, sin verificar adicionalmente la longitud de la cadena de calidad hasta ahora coincide con la longitud de la secuencia.. "

Ver http://ukpmc.ac.uk/articlerender.cgi?accid=PMC2847217 para más detalles

+0

No es cierto para ninguna solución que especifique que el carácter @ está al comienzo de la línea con '^ @', que parece representar la mayoría de las respuestas. Cheers – Morlock

+1

En realidad, esta es una afirmación verdadera ya que un valor de calidad de @ puede aparecer en principio en cualquier parte de la cadena de calidad, incluido el primer carácter, no se garantiza que '^ @' capture las líneas de nombre solamente. –

+0

De hecho. Perdón por no haber tomado unos segundos más para pensar el problema correctamente. ¡Salud! – Morlock

1

Sé que estoy lejos en el futuro, pero en beneficio de los empleados de Google:

Es posible que desee utilizar fastq_to_fasta from the fastx toolkit que se mantendrá el signo @, sin embargo. . además, eliminará las líneas con Ns menos que le diga que no.

2

me gustaría escribir

awk ' 
    NR%4 == 1 {print ">" substr($0, 2)} 
    NR%4 == 2 {print} 
' fastq > fasta 
2

Este es el más rápido que tengo, una d me lo metió en mi archivo .bashrc:

alias fq2fa="awk '{print \">\" substr(\$0,2);getline;print;getline;getline}'" 

No deja en las líneas poco frecuentes pero no imposibles de calidad que comienzan con @ ... pero que fallan en FASTQ envuelto, si eso es incluso legal (que existe, sin embargo).