[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