Discussion :

[pontarose]

Bonjour,

J'ai écris un script (Bio::SeqIO) qui me permet d' obtenir la taille de chacune de mes séquences. Mais voilà, j'ai à priori 104763 séquences mais mon script m'en compte 85673!!!!
Quand je teste mon script avec la ligne
$hash{$seq->primary_id} = $seq->seq; j'ai le bon nombre de séquences (104763) mais avec
$hash{$seq->seq} = $seq->primary_id.$seq->desc ; j'ai le bon 85673 séquences

je ne vois d'ou proviens le problème.

MON FICHIER : ftp://ftp.ensembl.org/pub/release-74....cds.all.fa.gz

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
 
cat Homo_sapiens.GRCh37.74.cds.all.fa |awk ' { if($0 ~ /^>/){print $1;}}' | wc -l
   104763

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
 
my $in  = Bio::SeqIO->new(-file => "Homo_sapiens.GRCh37.74.cds.all.fa", '-format' => 'Fasta');
 
my %hash;
my $count; 
 
while ( my $seq = $in->next_seq()){
 	$hash{$seq->seq} = $seq->primary_id.$seq->desc ;#pour récup la ligne entière 
}
 
foreach my $t (keys %hash){
	$count++;
        my $taille=length($t);
}
#print "nb de sequences: $count\n; # je trouve 85673

Quelqu'un aurait une suggestion?

[Gardyen]

tu as tout simplement des séquences identiques avec différents primary id

[pontarose]

Merci.
Tu veux dire des primary_id identique !!!

[Gardyen]

non, les primary id sont tous différents, vu que ce sont des identifiants. D'ailleurs tu obtiens bien toutes les séquences quand tu utilises $hash{$seq->primary_id}.
Ce sont les séquences ADN qui ne sont pas forcément différentes, donc quand tu fais $hash{$seq->seq}, tu écrases en fait l'entrée précédente.

Donc utilise plutôt les primary_id pour faire tes clés, c'est leur rôle après tout...