Foro - Perl en Español

por **Diego3D** » 2012-04-14 01:14 @093

Hola a todos

¡Espero que estén bien!

Resulta que estoy programando un script para anotar un genoma. Para eso busco los ORF en cada marco, los traduzco y los guardo en un archivo FASTA.

El problema es que el marco 1 y 2 me dan las mismas proteínas, y no entiendo por qué, ya que el marco 0 me funciona sin problemas. Lo mismo me pasa con los reversos complementarios.

Aquí está mi código:

Sintáxis: [ Descargar ] [ Ocultar ]

Using perl Syntax Highlighting

#!/usr/bin/perl
use strict;
use warnings;
use Bio::Seq;
use Bio::SeqIO;
 
 
my $obj_fasta = Bio::SeqIO->new(-file => "NC_008512.ffa",
                                -format => 'fasta');
my $ob = $obj_fasta->next_seq;
 
my $seq_obj= Bio::Seq->new(-display_id  =>      $ob->id,
                            -seq        =>      $ob->seq,
                            -alphabet   =>      'dna'
                      );
 
my $j=0;
my $frame;
 
for(my $i=0; $i<6;$i++){ # Busca proteínas en los 6 marcos de lectura
 
        if($i==3){ 
                $j=0;
        }
 
        if($i<=2){ # obtiene los tres primeros marcos directos
                $frame = substr($seq_obj->seq, $j,$seq_obj->length);
                open(MIFICH,">frame_$j.txt"); # se crea un archivo para cada frame
        }
 
        else{ # obtiene los tres marcos reversos complementarios
                $frame = substr($seq_obj->revcom->seq, $j,$seq_obj->length);
                open(MIFICH,">frame_-$j.txt"); # se crea un archivo para cada frame
        }
 
        $j++;
        my @orfs = $frame =~ /(ATG(?:.{3})*?(?:TAA|TAG|TGA|.{1,3}$))/g; # hace match con todos los ORF que encuentre
        
        foreach my $orf (@orfs) { 
                if((length($orf))>=240){ # proteínas con largo >= a 80 aminoácidos
 
                        my $prot_ob= Bio::Seq->new( # se crea un objeto secuencia con el ORF
                                -seq            =>      $orf,
                                -alphabet       =>      'dna'
                        );
 
                        if($i<=2){ # posición dentro del genoma del ORF directo
                                my $start = index($seq_obj->seq, $prot_ob->seq);
                                my $end = $start + length($prot_ob->seq);
                                print MIFICH "$start   $end \n";
                        }
 
                        else{ # posición dentro del genoma del ORF reverso
                                my $start = index($seq_obj->seq, $prot_ob->revcom->seq);
                                my $end = $start + length($prot_ob->seq);
                                print MIFICH "$end   $start \n";
                                
                        }
                        print MIFICH $prot_ob->seq->translate(-condotable=>11)->seq."\n\n"; # imprime la proteína en el archivo fasta generado para su respectivo frame
                }
        }
}
Coloreado en 0.003 segundos,  usando GeSHi 1.0.8.4

Estas secuencias, si los dividimos en codones, darían lugar a secuencias proteicas diferentes, ya que los codones serían distintos (cada secuencia es diferente, por estar desplazada una posición, luego se generan codones distintos).

Bueno... a cada uno de estos marcos le aplicas la expresión regular: /(ATG(?:.{3})*?(?:TAA|TAG|TGA|.{1,3}$))/gi .

Y aquí está el problema: el patrón que tienes NO tiene en cuenta la posición desplazada de cada secuencia, sino que busca, de forma lineal, dentro de la secuencia, por lo que, es muy, muy posible, que estés encontrando los mismos marcos de lectura abierta.

Ejemplo: tenemos esta secuencia, con un ORF dentro de ella:
aactgcatgacgtaacgtca

Ese ORF es encontrado sin problemas por la expresión regular:

Sintáxis: [ Descargar ] [ Ocultar ]

Using bash Syntax Highlighting

> perl -E '$seq = "aactgcatgacgtaacgtca"; say $1 while $seq =~ /(ATG(?:.{3})*?(?:TAA|TAG|TGA|.{1,3}$))/gi;'
atgacgtaa
Coloreado en 0.003 segundos,  usando GeSHi 1.0.8.4

Ahora, desplazamos la secuencia una posición más a la derecha, y vemos qué ocurre:
_aactgcatgacgtaacgtca

Sintáxis: [ Descargar ] [ Ocultar ]

Using bash Syntax Highlighting

> perl -E '$seq = " aactgcatgacgtaacgtca"; say $1 while $seq =~ /(ATG(?:.{3})*?(?:TAA|TAG|TGA|.{1,3}$))/gi;'
atgacgtaa
Coloreado en 0.001 segundos,  usando GeSHi 1.0.8.4

Sale lo mismo... porque a la expresión regular le ha dado lo mismo encontrar la secuencia en la posición seis que en la siete.

Ese es el problema: la expresión regular sí sabe de codones y sus posiciones múltiples -el subpatrón (?:.{3})*?-, pero no tiene en cuenta la posición inicial del codón de arranque 'ATG'.

Entonces... si el problema es encontrar los ORF en secuencias desplazadas (y sus inversas), lo que hay que hacer es, después del desplazamiento (y reversión) de las secuencias, pasar la secuencia a codones, y a partir de ellos, localizar los ORF.

Aquí tienes un ejemplo sencillo de cómo se podría resolver (sin BioPerl, que ahora mismo no lo tengo instalado):

Sintáxis: [ Descargar ] [ Ocultar ]

Using perl Syntax Highlighting

#!/usr/bin/perl
use v5.14;
 
my $seq = "aactgcatgacgtaacgtccatgacgtagaactcgatcgatcgatcgactgactagcttatagcattacacg";
 
my $regex_orf_mala  = qr/            (ATG(?:.{3})*?(?:TAA|TAG|TGA|.{1,3}$))/xi;
my $regex_orf_buena = qr/\G(?:.{3})*?(ATG(?:.{3})*?(?:TAA|TAG|TGA|.{1,3}$))/i;
 
(my $seqrev = $seq) =~ tr/atcg/tagc/;
$seqrev = reverse $seqrev;
 
say "Directa (mala):";
find_orf($seq, $regex_orf_mala);
 
say "Inversa (mala):";
find_orf($seqrev, $regex_orf_mala);
 
say;
 
say "Directa (buena):";
find_orf($seq, $regex_orf_buena);
 
say "Inversa (buena):";
find_orf($seqrev, $regex_orf_buena);
 
sub find_orf {
    my($seq, $regex) = @_;
    
    #say "Analizando [$seq] con [$regex]";
 
    for my $d (0..2) {
        my $s = ' ' x $d . $seq;
        say "$d : [$1]" while $s =~ /$regex/g;
    }
}
Coloreado en 0.001 segundos,  usando GeSHi 1.0.8.4

La salida es:

Sintáxis: [ Descargar ] [ Ocultar ]

Using text Syntax Highlighting

Directa (mala):

: [atgacgtaa]

: [atgacgtag]

: [atgacgtaa]

: [atgacgtag]

: [atgacgtaa]

: [atgacgtag]

Inversa (mala):

: [atgctataa]

: [atggacgttacgtcatgcagtt]

: [atgctataa]

: [atggacgttacgtcatgcagtt]

: [atgctataa]

: [atggacgttacgtcatgcagtt]

Directa (buena):

: [atgacgtaa]

: [atgacgtag]

Inversa (buena):

: [atgctataa]

: [atggacgttacgtcatgcagtt]

: [atgcagtt]Coloreado en 0.000 segundos,  usando GeSHi 1.0.8.4

Observa cómo la primera regex localiza la misma secuencia, tres veces.

Un detalle importante, también: esta solución no tiene en cuenta el caso de que quieras buscar los ORF que están incluidos unos dentro de otros (el caso de encontrarse varios codones de arranque que terminan en el mismo codón de terminación). Si ese fuera también el caso... habría que hacer un ajuste más fino para la búsqueda de codones.

por **Diego3D** » 2012-04-14 12:06 @546

Muchas gracias, explorer, por la gran respuesta que me diste

. Sospechaba que la expresión regular no discriminaba la distancia desde el inicio, cómo tu dices.

Lo que no entendí del código eso si era, ¿por qué hay una expresión regular buena y otra mala?

por **explorer** » 2012-04-14 12:23 @557

Es una forma de llamarlas, para distinguir entre el caso de un patrón que no tiene en cuenta la posición múltiplo de comienzo de los codones, de otra que sí.

Foro - Perl en Español

Problema con los frames y los ORF

Problema con los frames y los ORF

Publicidad

Re: Problema con los frames y los ORF

Re: Problema con los frames y los ORF

Re: Problema con los frames y los ORF

¿Quién está conectado?