[ensembl-dev] Bio::EnsEMBL::Feature::strand doesn't match output strand in plugin

[Will McLaren]

Hi Matt,

The strand in $line_hash->{Extra}->{STRAND} is the strand of the
overlapping genomic feature (transcript, regulatory feature etc), not the
variation feature itself.

It should match $tva->feature->strand in your example code.

If your input is VCF, you will find that the strand of the variation
feature is _always_ 1, as all variants encoded in VCF are on the forward
strand.

If the variation feature overlaps a transcript on the negative strand, then
before any sequence annotation is done on the variant allele, it is reverse
complemented; you can see the difference between the alleles in these
situations by calling:

# variant allele sequence relative to transcript
$tva->feature_seq()

# "original" variant allele sequence from variation feature
$tva->variation_feature_seq()

HTH

Will McLaren
Ensembl Variation

On 21 October 2014 21:34, Matt Wood <matt.wood at codifiedgenomics.com> wrote:

> I'm writing a plugin and am running into an issue with strand. The strand
> I get from Bio::EnsEMBL::Feature::strand doesn't match the strand I'm
> seeing in the output or the strand that is in the line hash provided to my
> run method. Specifically I'm having an issue with the negative strand.
>
> This is my test plugin:
>
> package StrandTest;
> use strict;
> use warnings;
> use Bio::EnsEMBL::Variation::Utils::BaseVepPlugin;
> use base qw(Bio::EnsEMBL::Variation::Utils::BaseVepPlugin);
>
> sub run {
>     my ($self, $tva, $line_hash) = @_;
>
>     my $strand = $tva->variation_feature->strand;
>     my $extra_strand = $line_hash->{Extra}->{STRAND};
>
>     print("Strand: $strand/$extra_strand\n");
>
>     return {};
> }
>
> 1;
>
> My output shows this: "Strand: 1/-1"
>
> The strand written to my output file is -1, and the strand in the line
> hash is -1, but the strand from the VariationFeature is 1.
>
> My VCF:
>
> ##fileformat=VCFv4.0
> #CHROM    POS    ID    REF    ALT    QUAL    FILTER    INFO    FORMAT
> PAT
> 6    133078809    .    C    A    60    PASS    .    GT:VR:RR:DP:GQ
> 0/1:12:15:27:.
>
> I'm running VEP with these options: --plugin StrandTest --offline --vcf
> --force_overwrite --pick -i test.vcf -o test.vep
>
> Any idea what's going on? For the time being I'm going to use
> "$line_hash->{Extra}->{STRAND}" for my strand info. Is that safe to do? Are
> there going to be cases where that hash is undef?
>
> Thanks,
> Matt Wood
>
>
>
> _______________________________________________
> Dev mailing list    Dev at ensembl.org
> Posting guidelines and subscribe/unsubscribe info:
> http://lists.ensembl.org/mailman/listinfo/dev
> Ensembl Blog: http://www.ensembl.info/
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.ensembl.org/pipermail/dev_ensembl.org/attachments/20141022/ca1d079a/attachment.html>