Foro - Perl en Español

por **j_f_r_85** » 2010-12-06 06:37 @317

¡¡¡Buen día compañeros!!!

Estoy iniciando en la bioinformática, y pues Perl es el primer peldaño.

Tengo un pequeño problema con un script que comencé por curiosidad. Deseo generar secuencias aleatorias de ADN, con una longitud especificada con el usuario, así como el numero de secuencias. Igualmente el usuario puede definir el % de GC, es decir, qué porcentaje de C y G hay en la cadena.

He creado mi script, pero no me da la salida que deseo; yo deseo algo como:

Sintáxis: [ Descargar ] [ Ocultar ]

Using text Syntax Highlighting

CATCTAAAAG

TGCGTAGATC

GTGTACCAAC

AGGAGGCTTGColoreado en 0.000 segundos,  usando GeSHi 1.0.8.4

Pero me da algo así:

Sintáxis: [ Descargar ] [ Ocultar ]

Using text Syntax Highlighting

CATCTAAAAG

CATCTAAAAGTGCGTAGATC

CATCTAAAAGTGCGTAGATCGTGTACCAAC

CATCTAAAAGTGCGTAGATCGTGTACCAACAGGAGGCTTGColoreado en 0.000 segundos,  usando GeSHi 1.0.8.4

Sí, me crea las secuencias aleatorias, con el porcentaje de GC deseado, pero no en el patrón que deseo. ¿Me pueden ayudar a encontrar qué hice incorrectamente? ¡Gracias de antemano!

El script:

Sintáxis: [ Descargar ] [ Ocultar ]

Using perl Syntax Highlighting

#!/usr/bin/perl
#random_nt.pl
#Script generates "n" random DNA sequences given a %GC and length of sequence
 
use strict; use warnings;
 
die "usage: random_nt.pl <gc %> <length> <# sequences>" unless @ARGV == 3;
 
my ($gc, $length, $nseq) = @ARGV;
 
my $seq = '';
 
foreach (my $i = 0; $i < $nseq;   $i++){
        for (my $j = 0; $j < $length; $j++){
                my $newnt = rndmbp ($gc);
                $seq .= $newnt;
        }
        print "$seq\n";
}
 
sub rndmbp {
my ($gc)      = @_;
my $GC        = ($gc/2);
my $AT        = (0.5-$GC);
my ($base)    = '';
my $rndm      = rand;
if ($rndm     > ((2*$AT)+$GC)) 
        {$base   .= "G";}
elsif ($rndm  > (2*$AT))
        {$base   .= "C";}
elsif ($rndm  > $AT)
        {$base   .= "A";}
elsif ($rndm  > 0)
        {$base   .= "T";}
return $base;
}
Coloreado en 0.003 segundos,  usando GeSHi 1.0.8.4

¡Muchas gracias, wanako! Sabía que la respuesta seguro sería algo simple, pero la verdad no la encontraba, ¡y además ganamos 8 bytes! En serio, ¡muchas gracias!

No sé si se pueda optimizar el script, pero lo importante es que ahora funciona como quiero.

Ahora intentaré hacer algo con esa salida con un nuevo script. Gracias.

por **explorer** » 2010-12-06 13:02 @585

Perdona, j_f_r_85, pero creo que el programa no funciona muy bien.

Si ejecuto ./code_23733_2.pl .50 10 10, es decir, diez secuencias de diez bases, con un 50% de aparición de G y C, los resultados que me salen son:

Sintáxis: [ Descargar ] [ Ocultar ]

Using text Syntax Highlighting

CACAATGCTA 4/10

ATGTCTGATT 3/10

CCTACATCGC 6/10

GTTTGTCGGG 6/10

ATGAATTGGA 3/10

AGAAGACTTA 3/10

TCTAGCGACG 6/10

GAACTAAAAG 3/10

CAGGAAACTT 4/10

CCTTTCGCCT 6/10Coloreado en 0.000 segundos,  usando GeSHi 1.0.8.4

que, como ves, ninguno llega al 50% de aparición.

por **explorer** » 2010-12-06 15:53 @703

Esta es una solución que he preparado, que sí genera presencias de Cy G acorde a lo que nos pide el usuario.

Sintáxis: [ Descargar ] [ Ocultar ]

Using perl Syntax Highlighting

#!/usr/bin/perl
#
# Genera secuencias de ADN de una determinada longitud,
# en la que se indica el porcentaje de aparición de C y G
#
# JF 20101206.2153
 
use common::sense;              # siempre hay que tener sentido común
 
## Argumentos al programa
@ARGV == 3  or  die "Uso: $0 GC% <longitud de las secuencias> <número de secuencias>\n";
 
my ($GC_porcentaje, $longitud, $número_secuencias) = @ARGV;
 
$longitud          > 0  or  die "ERROR: la longitud debe ser mayor que 0\n";
$número_secuencias > 0  or  die "ERROR: el número de secuencias debe ser mayor que 0\n";
$GC_porcentaje     > 0  or  die "ERROR: el GC% debe estar entre 0 y 1\n";
$GC_porcentaje     < 1  or  die "ERROR: el GC% debe estar entre 0 y 1\n";
 
 
## Calculamos el número de caracteres de cada parte
my $GC_longitud_caracteres = int ($longitud * $GC_porcentaje);
my $AT_longitud_caracteres = $longitud - $GC_longitud_caracteres;
 
 
## Bucle por todas las secuencias
for (1 .. $número_secuencias) {
 
    ## Generamos la secuencia
    my @ATCG
        = ( (map { ('C','G')[rand 2] } 1 .. $GC_longitud_caracteres)
          , (map { ('A','T')[rand 2] } 1 .. $AT_longitud_caracteres)
          )
        ;  
 
    my $ATCG;
    while (@ATCG) {
        $ATCG .= splice @ATCG, rand @ATCG, 1;   # extraemos un elemento al azar
    }
 
    ## comprobación
    my $cuenta_CG = $ATCG =~ tr/CG/CG/;         # contamos las C y G en $ATCG
 
    say "[$ATCG] $cuenta_CG/$longitud";
}
Coloreado en 0.002 segundos,  usando GeSHi 1.0.8.4

que sale, para, por ejemplo ./code_23733_2.pl .32 100 10:

Sintáxis: [ Descargar ] [ Ocultar ]

Using text Syntax Highlighting

[ATTTGTAGCTTTCTAAAACTATATACAAGAGAGGTATTCATTTCCCATTCAAGCATAGGTGTAGCACTGAAACTTAACCGAGTTACCAAATTAATATTTG] 32/100

[CACGTAATCTGTAGTTATAGATCGTTCGTTGAGTTTACAGTGCAAAGGTTTACTAAGAATTCTAACGGGTTTCTCATTTAACTATCATTTTATAAGTTGT] 32/100

[TTAAAGATAAGAAGCTTAAGAATTTAAACTGGAGTCTTCGTAACATGTACGAAAAAACGCCGAGCTTATGTTTTTGTCATTTATGATCAAGGCTATGATT] 32/100

[GAAGTCTTCGGGCATGCCCGTATCTGAACGTGTTATAGAATTTATTAACTATAGACATCTATATTTATAACTAAAACTCGTCTATAATCTAGTCATTCAA] 32/100

[AACGATTTGCTGATAAAAGTGACCGGTTTCACATATTGCATACATTTTATTTTTAGACTATTCCCGTTGAATGACTTCGTGGATTGTTACTTTACTAATA] 32/100

[TAACATAGTAATCAAGGTCATGTTCAGGTAAGAATCGGTTTTCAAATTATTGAGGTTGTTATGGTCCTGTCCCTATACATCTTTTTATGAAATATCAGTT] 32/100

[AAAAGTTGATATTCAAAAAATACGAGCGGGGATACAATTTGTCAAATGAAAATAATCTACCATTGCTTTATCGTTAGACATACTTTAAGCGCTCCAGTTT] 32/100

[AAAACATAACTAAATCAATTCCATCTGTAGTAGGTATGTATTCACAGACACGACACAGCGATTTGAACGTCTATTTGTTATAAATCATTAAGGGTCAAAT] 32/100

[CAAGTAATCATGTACCTATAAGAAATAAAATTATTGACTTGTATTCAGTTTTACTATTCCCTAGCAGTCAAAAAATGGTGCGGCAAGTAATAACGTGAGA] 32/100

[TTTTTGCAAGTTATTAAGAATTAAGTAAGTGGTATAAAAGTGAACACGTTTTTATAACTTTAGAGTGAGCCTCGACTGCGTGATTAAATAGAAGGGAACA] 32/100Coloreado en 0.000 segundos,  usando GeSHi 1.0.8.4

La generación de la secuencia se hace con un array, lo cual no es lo recomendable si queremos generar secuencias de varios cientos de millones de bases. Pero bueno, para cosas pequeñas, puede valer

De todas maneras, creo que esto lo hace de forma directa el módulo Bio::Matrix::PSM::SiteMatrix.

Y me sonaba que había otro módulo de matemáticas que podía generar secuencias según la frecuencia que le indiquemos. Lo he intentando encontrar, pero nada...

por **j_f_r_85** » 2010-12-06 16:01 @709

Mi estimado explorer:

Inmediatamente se nota tu sapiencia al respecto del tema tu script luce muy profesional, y como tu dices ¡con sentido común!

¡Muchas gracias!

Tengo que seguir practicando.

por **explorer** » 2010-12-06 16:28 @728

common::sense es un módulo muy curioso.

Ya sé que los puristas del lenguaje me dirán que no hay que usarlo, y menos en unos foros, y menos en el foro básico, pero el caso es que tiene justo lo que quiero para todos mis programas: utf8, las nuevas características de Perl v5.10 y siguientes, strict, warnings, etc, etc.

Foro - Perl en Español

Secuencias aleatorias de ADN

Secuencias aleatorias de ADN

Publicidad

Re: Secuencias aleatorias de ADN

Re: Secuencias aleatorias de ADN

Re: Secuencias aleatorias de ADN

Re: Secuencias aleatorias de ADN

Re: Secuencias aleatorias de ADN

Re: Secuencias aleatorias de ADN

¿Quién está conectado?