• Publicidad

Archivo multifasta

Perl aplicado a la bioinformática

Archivo multifasta

Notapor eliavalls » 2017-10-13 06:49 @325

Buenos días.

Partimos de un archivo de texto con cadenas de DNA en formato FASTA:

Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
>Rosalind_6404
CCTGCGGAAGATCGGCACTAGAATAGCCAGAACCGTTTCTCTGAGGCTTCCGGCCTTCCC
TCCCACTAATAATTCTGAGG
>Rosalind_5959
CCATCGGTAGCGCATCCTTAGTCCAATTAAGTCCCTATCCAGGCGCTCCGCCGAAGGTCT
ATATCCATTTGTCAGCAGACACGC
>Rosalind_0808
CCACCCTCGTGGTATGGCTAGGCATTCAGGAACCGGAGAACGCTTCAGACCAGCCCGGAC
TGGGAACCTGCGGGCAGTAGGTGGAAT
Coloreado en 0.001 segundos, usando GeSHi 1.0.8.4


Y queremos que nos devuelva el ID de la cadena con el porcentaje de CG más elevado, éste incluido.

La salida sería:
Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
 Rosalind_0808
 60.919540
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


¿Alguien sabría cómo hacerlo?

Gracias.
eliavalls
Perlero nuevo
Perlero nuevo
 
Mensajes: 5
Registrado: 2017-10-13 06:26 @310

Publicidad

Re: Archivo multifasta

Notapor explorer » 2017-10-13 07:37 @359

Bienvenida a los foros de Perl en Español, eliavalls.

Yo creo que se puede resolver muy fácilmente.

En estos foros encontrarás muchos ejemplos de cómo leer un archivo así. Busca por 'fasta' o 'multifasta'.

Una solución rápida sería:

Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
Abrir archivo
Por cada línea del archivo
        Si la línea comienza por un '>', (estamos en una cabecera)
                Hacemos la suma de las letras C y G que tengamos almacenadas hasta ahora,
                        y lo dividimos por el total de letras para obtener el porcentaje
                Si el porcentaje es superior al récord
                        Recordamos el ID
                        Recordamos el récord
                Ponemos a cero el almacén de letras
                Extraemos el ID de la nueva cabecera
        Si no,
                Agregamos la línea en el almacén de letras
Repetir
Cerrar archivo
Procesar el último almacén, igual que antes
Imprimir el valor de récord.
Coloreado en 0.001 segundos, usando GeSHi 1.0.8.4


Esta solución también aparece en varios hilos de este subforo de Bioinformática. Naturalmente, hay otras formas.

Prueba a ver si lo sacas. Si te atascas, vuelve por aquí y te daremos más pistas.
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Archivo multifasta

Notapor eliavalls » 2017-10-14 12:01 @542

¡Muchísimas gracias por la ayuda! Tenía la idea pero no sabía cómo hacerlo.
eliavalls
Perlero nuevo
Perlero nuevo
 
Mensajes: 5
Registrado: 2017-10-13 06:26 @310

Re: Archivo multifasta

Notapor eliavalls » 2017-10-17 14:20 @638

Y si, por ejemplo, también quisiera que me calculara el porcentaje de GC por ventanas (con un tamaño de ventana de 10 y un tamaño de paso de 5) de cada secuencia. ¿Cómo se haría? Uno de los resultados sería:

70
50
50
60
60
60
50
50
60
70
80
70
60
80
70
50
42.8571428571429
eliavalls
Perlero nuevo
Perlero nuevo
 
Mensajes: 5
Registrado: 2017-10-13 06:26 @310

Re: Archivo multifasta

Notapor explorer » 2017-10-17 15:33 @689

Humm... no me entero de cómo quieres hacerlo. ¿Puedes explicarlo de forma más detallada?
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Archivo multifasta

Notapor eliavalls » 2017-10-17 16:02 @710

En el archivo de texto hay 3 secuencias. Si, por ejemplo, cogemos esta:

Sintáxis: [ Descargar ] [ Ocultar ]
Using text Syntax Highlighting
>Rosalind_0808
CCACCCTCGTGGTATGGCTAGGCATTCAGGAACCGGAGAACGCTTCAGACCAGCCCGGAC
TGGGAACCTGCGGGCAGTAGGTGGAAT
Coloreado en 0.000 segundos, usando GeSHi 1.0.8.4


El porcentaje de GC de toda la secuencia es de:

60.919540

Pero si queremos calcular el porcentaje de GC por ventanas deslizantes (sliding windows) con un tamaño de ventana igual a 10 y un paso de ventana igual a 5 (coger 10 bases y calcula el porcentaje GC de estas 10; después avanza 5 bases y vuelve a coger 10; así sucesivamente). Nos devuelve los porcentajes de cada ventana por separado:

70
50
50
60
60
60
50
50
60
70
80
70
60
80
70
50
42.8571428571429

Yo he conseguido calcular el porcentaje de cada secuencia y que me devuelva el más grande, pero no sé cómo hacer que me calcule el porcentaje de una secuencia por partes.
eliavalls
Perlero nuevo
Perlero nuevo
 
Mensajes: 5
Registrado: 2017-10-13 06:26 @310

Re: Archivo multifasta

Notapor explorer » 2017-10-18 08:52 @411

Si la secuencia a analizar la tienes en una variable escalar puedes usar un bucle for para ir recorriendo las ventanas:

for (my $i = 0; $i < length $seq; $i += $winstep) {

y dentro del bucle, extraes las ventanas con substr():

my $window = substr $seq, $i, $winsize;

Y le aplicas el cálculo de porcentajes a $window.

Debes tener precaución con las ventanas que hay al final de la secuencia. Quiero decir que la última ventana es muy posible que no tenga un tamaño de 10 bases (salvo que la longitud de la secuencia principal sea múltiplo de 10). En el caso de que substr() te devuelva una ventana más pequeña, el porcentaje debería calcularse sobre el tamaño de esa ventana.

De forma programática, el caso normal y este caso especial se pueden resolver con las mismas líneas. Sería así:
Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
my $window_len = length $window;                # longitud de la ventana

my $gc = $window =~ tr/GC/GC/;                  # contar bases G y C

my $porcentaje = $gc / $windows_len * 100;      # cálculo del porcentaje
Coloreado en 0.010 segundos, usando GeSHi 1.0.8.4

y ya tienes calculado el porcentaje para esa ventana.
JF^D Perl Programming Language
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14102
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Re: Archivo multifasta

Notapor eliavalls » 2017-10-22 06:12 @300

¡¡¡Muchísimas gracias!!! Sí, ya tenía en cuenta que podría ser que el último substr() me devolviera una ventana más pequeña. Ahora ya está todo solucionado. Gracias de nuevo.
eliavalls
Perlero nuevo
Perlero nuevo
 
Mensajes: 5
Registrado: 2017-10-13 06:26 @310


Volver a Bioinformática

¿Quién está conectado?

Usuarios navegando por este Foro: No hay usuarios registrados visitando el Foro y 0 invitados