[Bioperl-l] Trying to get a conservation index out of a MSA

Téletchéa Stéphane stephane.teletchea at inserm.fr
Wed Feb 5 17:50:28 UTC 2014


Le 04/02/2014 20:48, Brian Osborne a écrit :
> Stephane,
>
> Can you send an example input file and complete script? You only showed a part of the script before, I think.
>
> Brian O.

Dear Brian,

Thank you for your patience, it turns out the problem was with my blast 
parsing:

instead of using the $hsp->hit_string to get the protein sequence of my 
target,
I created an alignment ($hsp->aln), and from this alignment, I took the 
sequence 1 using
get_sequence_by_pos(1), which I assumed was the second sequence (the hit 
one).

After thinking two seconds, I remembered that bioperl is intellingent, 
and human aware,
so it starts counting at one and not zero, so get_sequence_by_pos(1) is 
really the first sequence.

I have changed the parser to either use get_sequence_by_pos(2) or 
$hsp->hit_string and since
both are equal, I'm now using the later :-)

Bad brain, sorry for the noise...

Stéphane

PS: if you take one template sequence + 500 identical "hit" sequences,
this is then logical to get a 100% id for each position :-)

-- 
Equipe DSIMB - Dynamique des Structures et
des Interactions des Macromolécules Biologiques
INTS, INSERM-Paris-Diderot UMR_S 1134
6 rue Alexandre Cabanel - 75739 Paris cedex 15 - France
Tél : +33 144 493 057
Fax : +33 143 065 019
http://www.dsimb.inserm.fr / http://www.steletch.org




More information about the Bioperl-l mailing list