[Bioperl-l] What's wrong with the script about module Bio::SeqFeature::Gene::GeneStructure
Frank Schwach
fs5 at sanger.ac.uk
Fri Mar 18 15:05:29 UTC 2011
maybe I'm missing the point here but I think an "mRNA" feature should
not be a Bio::SeqFeature::Gene::GeneStructure and it should not have
exons - it would not make sense biologically.
Frank
On Fri, 2011-03-18 at 21:57 +0800, Tao Zhu wrote:
> Yes. But in the document
> http://doc.bioperl.org/releases/bioperl-current/bioperl-live/Bio/SeqFeature/Gene/GeneStructure.html
> It says that module Bio::SeqFeature::Gene::GeneStructure can inherit
> Bio::SeqFeature::Generic. I wonder how to do such inheritance.
>
> 在 2011-03-18五的 08:04 -0500,Chris Fields写道:
> > The SeqFeatures derived from a GenBank file are simple Bio::SeqFeature::Generic instances. I don't think there is a straightforward way to have the parser return Bio::SeqFeature::Gene::GeneStructure instances.
> >
> > chris
> >
> > On Mar 18, 2011, at 6:32 AM, Tao Zhu wrote:
> >
> > > Of course that's OK. But what I really want to know is how to use module
> > > Bio::SeqFeature::Gene::GeneStructure.
> > >
> > > 在 2011-03-18五的 10:17 +0100,Hans-Rudolf Hotz写道:
> > >> Hi Tao
> > >>
> > >> I don't fully understand your script, but I do see a major problem:
> > >>
> > >> why do you select for "primary_tag eq 'mRNA'" first?
> > >>
> > >>
> > >> this simple loop will probably do what you want:
> > >>
> > >> while( my $seq_obj = $catch_seq -> next_seq) {
> > >>
> > >> my @all_exon_features = grep {$_->primary_tag eq 'exon'}
> > >> $seq_obj ->get_SeqFeatures;
> > >> $exon_number = scalar(@all_exon_features);
> > >> print "$exon_number\n";
> > >> }
> > >>
> > >>
> > >> Alternatively, I recommend to read the HOWTO page:
> > >>
> > >> http://www.bioperl.org/wiki/HOWTO:Feature-Annotation
> > >>
> > >> which has a nice example to print out all 'primary_tag' (ie the 'feature
> > >> keys' in a GenBank formated file)
> > >>
> > >> I also recommend to read up on the DDBJ/EMBL/GenBank Feature Table
> > >> definition:
> > >> http://www.ncbi.nlm.nih.gov/projects/collab/FT/index.html
> > >>
> > >>
> > >> Regards, Hans
> > >>
> > >>
> > >>
> > >> On 03/18/2011 05:33 AM, Tao Zhu wrote:
> > >>> I wrote my script like this,
> > >>>
> > >>> #!/usr/bin/perl -w
> > >>> use Bio::SeqIO;
> > >>> my $catch_seq = Bio::SeqIO -> new(-file => 'test.gbk',-format=>
> > >>> 'genbank');
> > >>> while( my $seq_obj = $catch_seq -> next_seq)
> > >>> {
> > >>> my @all_mRNA_features = grep {$_->primary_tag eq 'mRNA'} $seq_obj ->
> > >>> get_SeqFeatures;
> > >>> for my $mRNA_feature (@all_mRNA_features)
> > >>> {
> > >>> if($mRNA_feature->isa('Bio::SeqFeature::Gene::GeneStructure'))
> > >>> {
> > >>> @exons=$mRNA_feature->exons;
> > >>> $exon_number = scalar(@exons);
> > >>> print "$exon_number\n";
> > >>> }
> > >>> }
> > >>> }
> > >>>
> > >>> I hope to count exon number in every mRNA. But it print nothing(You can
> > >>> arbitrarily get a genbank file to test it). What's wrong?
> > >>>
> > >>>
> > >
> > > --
> > > Tao Zhu, College of Life Sciences, Beijing Normal University, Beijing
> > > 100875, China
> > > Email: tzhu at mail.bnu.edu.cn
> > > Website: http://bnuzt.org (mainly written in Chinese)
> > >
> > >
> > > _______________________________________________
> > > Bioperl-l mailing list
> > > Bioperl-l at lists.open-bio.org
> > > http://lists.open-bio.org/mailman/listinfo/bioperl-l
> >
>
>
>
> _______________________________________________
> Bioperl-l mailing list
> Bioperl-l at lists.open-bio.org
> http://lists.open-bio.org/mailman/listinfo/bioperl-l
--
The Wellcome Trust Sanger Institute is operated by Genome Research
Limited, a charity registered in England with number 1021457 and a
company registered in England with number 2742969, whose registered
office is 215 Euston Road, London, NW1 2BE.
More information about the Bioperl-l
mailing list