Hola de nuevo, estimado explorer y resto de la comunidad.
Después de haber corregido mi código en base a tus sugerencias y a nueva información que me han dado para acabar la tarea creo que casi-casi lo tengo.
Los cambios son de sintaxis (todo aquello que entendía de tu código, explorer, no entendía por qué a un
filehandler se le asignaba una variable con "$", por ejemplo), y el más radical e importante es que hay que hacer varios archivos y varias consultas con
blastp, ya que no se puede pasar un multiFASTA a
blastp. De esto último me he ido dando cuenta con las pruebas.
Pues bien, mi código casi hace lo que le pido, es decir: me almacena para cada secuencia su ID, y la secuencia aminoacídica, pero no lo hace bien, ya que es justo a continuación, sin retorno de línea como le indico hasta en dos ocasiones.
Aún así parece que el
blastp en algunas de ellas sí lo acepta aunque no he conseguido corregirlo.
Y por otro lado no consigo que pase por el
blastp los archivos que he generado, dentro del bucle while.
Me dice esto:
Illegal division by zero at ./capturador2.pl line 29, <QUERY> line 2.Lo cual es muy desconcertante. El código que me queda al final es algo así, en la ruta de la base de datos (es una ruta local, del
cluster que utilizamos en el curso).
Using perl Syntax Highlighting
#!/usr/bin/perl
use autodie;
die "USAGE: capturador.pl file.fa\n" unless ( $ARGV[0] and -s $ARGV[0] );
my $id;
my $name;
my @seq;
open QUERY, '<', $ARGV[O];
while (<QUERY>) {
next if /^\s*$/;
chomp;
if (/^>(\S+)/) {
$id = $1;
$name = $id;
$name =~ s/[|]/_/g;
open INT, ">workfile$name.fa";
$id =~ s/[\n]//g;
print INT "$id\n";
}
else {
push( @seq, $_ );
$juntos = join( "", @seq );
$juntos =~ s/[*]//g;
$juntos =~ s/[\n]//g;
print INT "\n$juntos";
undef @seq;
system(blastp -query workfile$name.fa -db /mnt/beegfs/bioinfo_db/Blast/swissprot_db -out result $name.fas);
}
# system (blastp -query workfile$name.fa -db /mnt/beegfs/bioinfo_db/Blast/swissprot_db -out result$name.fas);
}
close(QUERY);
close(INT);
Coloreado en 0.001 segundos, usando
GeSHi 1.0.8.4
Si os dais cuenta he intentado ejecutar el
blastp, dentro y fuera del else{}, pero el fallo es el mismo.
A modo de resumen, mis problemas son dos:
- No consigo que los archivos generados queden de la forma "ID \n secuencia", sino que quedan "IDsecuencia".
- El blastp, no consigo que lo ejecute.
Muchas gracias, explorer. He aprendido mucho de todo lo que me has sugerido/indicado.