[Bioperl-l] What's wrong with the script about module Bio::SeqFeature::Gene::GeneStructure

[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