BioPerl
view release on metacpan or search on metacpan
Bio/SeqIO/agave.pm view on Meta::CPAN
# line: $line\n\n");
}
} # close while loop
return;
}
# ==================================================================================
=head2 _process_sciobj
Title : _process_sciobj
Usage : $self->_process_sciobj
Function : Parses the data between the <sciobj></sciobj> tags.
Args : The string that holds the attributes for <sciobj>.
Returns : Data structure holding the values parsed between
the <sciobj></sciobj> tags.
Note : Method(s) that call(s) this method : _process
Method(s) that this method calls :
_helper_store_attribute_list , _process_contig
=cut
sub _process_sciobj {
my ($self, $attribute_line) = @_;
my $sciobj;
$self->_helper_store_attribute_list($attribute_line, \$sciobj);
my $line = $self->_readline;
# Zero or more <contig>
while ($line =~ /<contig\s?(.*?)\s?>/) {
my $contig = $self->_process_contig(\$line, $1);
push @{$sciobj->{'contig'}}, $contig;
# print "line in _process_sciobj: $line\n";
# $line changes value within the subs called in this sub (_process_contig).
}
return $sciobj;
}
# ==================================================================================
=head2 _process_contig
Title : _process_contig
Usage : $self->_process_contig
Function : Parses the data between the <contig></contig> tags.
Args : 2 scalars:
- reference to a scalar holding the line to be parsed.
- scalar holding the attributes for the <contig> tag
to be parsed.
Returns : Data structure holding the values parsed between
the <contig></contig> tags.
Note : Method(s) that call(s) this method : _process_sciobj
Method(s) that this method calls :
_helper_store_attribute_list, _one_tag , _process_fragment_order
=cut
sub _process_contig {
my ($self, $line, $attribute_line) = @_;
my $contig;
$self->_helper_store_attribute_list($attribute_line, \$contig);
$$line = $self->_readline;
# One <db_id>:
$self->_one_tag($line, \$contig, 'db_id');
# Zero or more <fragment_order>
$self->_process_fragment_order($line, \$contig);
return $contig;
}
# ==================================================================================
=head2 _process_fragment_order
Title : _process_fragment_order
Usage : $self->_process_fragment_order
Function : Parses the data between the <fragment_order></fragment_order> tags.
Args : 2 scalars:
- reference to a scalar holding the value of the line to be parsed.
- reference to a data structure to store the <fragment_order> data.
Returns : Nothing.
Note : Method(s) that call(s) this method : _process_contig
Method(s) that this method calls :
_helper_store_attribute_list , _process_fragment_orientation
=cut
sub _process_fragment_order {
my ($self, $line, $data_structure) = @_;
# Because I'm passing a reference to a data structure, I don't need to return it
# after values have been added.
while ($$line =~ /<fragment_order\s?(.*?)\s?>/) {
my $fragment_order;
$self->_helper_store_attribute_list($1, \$fragment_order);
# Store the attribute(s) for <fragment_order> into the
# $fragment_order data structure.
$$line = $self->_readline;
# One or more <fragment_orientation>
$self->_process_fragment_orientation($line, \$fragment_order);
# Don't forget: $line is a reference to a scalar.
push @{$$data_structure->{'fragment_order'}}, $fragment_order;
# Store the data between <fragment_order></fragment_order>
# in $$data_structure.
}
return;
}
# ==================================================================================
=head2 _process_fragment_orientation
Title : _process_fragment_orientation
Usage : $self->_process_fragment_orientation
Function : Parses the data between the <fragment_orientation> and
</fragment_orientation> tags.
Args : 2 scalars:
- reference to a scalar holding the value of the line to be parsed.
- reference to a data structure to store the <fragment_orientation> data.
Returns : Nothing.
Note : Method(s) that call(s) this method : _process_fragment_order
Method(s) that this method calls : _helper_store_attribute_list ,
_process_bio_sequence
=cut
sub _process_fragment_orientation {
my ($self, $line, $data_structure) = @_;
# counter to determine the number of iterations within this while loop.
my $count = 0;
# One or more <fragment_orientation>
while ($$line =~ /<fragment_orientation\s?(.*?)\s?>/) {
my $fragment_orientation;
$self->_helper_store_attribute_list($1, \$fragment_orientation);
$$line = $self->_readline;
# One <bio_sequence>
$$line =~ /<bio_sequence\s?(.*?)\s?>/;
# Process the data between <bio_sequence></bio_sequence>
my $bio_sequence = $self->_process_bio_sequence($line, $1);
$fragment_orientation->{'bio_sequence'} = $bio_sequence;
push @{$$data_structure->{'fragment_orientation'}}, $fragment_orientation;
++$count;
}
$self->throw("Error: Missing <fragment_orientation> tag. Got this: $$line\n\n")
if $count == 0;
return;
}
# ==================================================================================
=head2 _process_bio_sequence
Title : _process_bio_sequence
Usage : $self->_process_bio_sequence
Function : Parses the data between the <bio_sequence></bio_sequence> tags.
Args : 2 scalars:
- reference to a scalar holding the value of the line to be parsed.
- scalar holding the value of the attributes for <bio_sequence>
Returns : data structure holding the values between <bio_sequence></bio_sequence>
Note : Method(s) that call(s) this method : _process_fragment_orientation
Method(s) that this method calls : _helper_store_attribute_list ,
_one_tag , _question_mark_tag , _star_tag , _process_alt_ids ,
_process_xrefs , _process_sequence_map
=cut
sub _process_bio_sequence {
my ($self, $line, $attribute_line) = @_;
my $bio_sequence;
$self->_helper_store_attribute_list($attribute_line, \$bio_sequence);
$$line = $self->_readline;
# One <db_id>.
$self->_one_tag($line, \$bio_sequence, 'db_id');
# Zero or one <note>.
$self->_question_mark_tag($line, \$bio_sequence, 'note');
# Zero or more <description>
$self->_question_mark_tag($line, \$bio_sequence, 'description');
# Zero or more <keyword>
$self->_star_tag($line, \$bio_sequence, 'keyword');
# Zero or one <sequence>
$self->_question_mark_tag($line, \$bio_sequence, 'sequence');
# Zero or one <alt_ids>
# NOT IMPLEMENTED!!!!
#if ($line =~ /<alt_ids>/){ # NOT DONE YET!
# my $alt_ids;
# $bio_sequence->{'alt_ids'} = $self->_process_alt_ids(\$alt_ids);
#}
# Zero or one <xrefs>
if ($$line =~ /<xrefs\s?(.*?)\s?>/) {
my $xrefs = $self->_process_xrefs($line, \$bio_sequence);
$bio_sequence->{'xrefs'} = $xrefs || 'null';
}
# Zero or more <sequence_map>
if ($$line =~ /<sequence_map\s?(.*?)\s?>/) {
my $sequence_map = $self->_process_sequence_map($line);
push @{$bio_sequence->{'sequence_map'}}, $sequence_map;
}
# print Data::Dumper->Dump([$bio_sequence]); exit;
Bio/SeqIO/agave.pm view on Meta::CPAN
Usage : $self->_helper_store_attribute_list
Function : A helper method used to store the attributes from
the tags into the data structure.
Args : 2 scalars:
- scalar holding the attribute values to be parsed.
- reference to a data structure to store the data between the 2 tags.
Returns : Nothing.
Note : Method(s) that call(s) this method : Many.
Method(s) that this method call(s) : None.
=cut
sub _helper_store_attribute_list {
my ($self, $attribute_line, $data_structure) = @_;
my %attribs = ($attribute_line =~ /(\w+)\s*=\s*"([^"]*)"/g);
my $attribute_list;
for my $key (keys %attribs) {
# print "\tkey: $key , value: $attribs{$key}\n";
###$$data_structure->{$key} = $attribs{$key}; # <- The ORIGINAL.
push @{$$data_structure->{$key}}, $attribs{$key};
# Now, store them in an array because there may be > 1 tag, thus
# > 1 attribute of the same name.
# Doing this has made it necessary to change the _store_seqs method.
# ie: Change $bio_sequence->{'molecule_type'};
# to
# $bio_sequence->{'molecule_type'}->[0];
}
return;
}
# ==================================================================================
=head2 _store_seqs
Title : _store_seqs
Usage : $self->_store_seqs
Function : This method is called once in the life time of the script.
It stores the data parsed from the agave xml file into
the Bio::Seq object.
Args : None.
Returns : Nothing.
Note : Method(s) that call(s) this method : next_seq
Method(s) that this method calls : None.
=cut
sub _store_seqs {
my ($self) = @_;
for my $sciobj (@{$self->{'sciobj'}}) {
### $sciobj = $self->{'sciobj'}; # The root node.
for my $contig (@{$sciobj->{'contig'}}) { # Each contig has a fragment order.
for my $fragment_order (@{$contig->{'fragment_order'}}) { # Each fragment order has a fragment_orientation.
for my $fragment_orientation (@{$fragment_order->{'fragment_orientation'}}) {
# Each fragment_orientation contain 1 bio sequence.
my $bio_sequence = $fragment_orientation->{'bio_sequence'}; # <bio_sequence> contains all the
# interesting stuff:
my $sequence = $bio_sequence->{'sequence'};
my $accession_number = $bio_sequence->{'sequence_id'}->[0]; # also use for primary_id
my $organism = $bio_sequence->{'organism'};
my $description = $bio_sequence->{'description'};
my $molecule_type = $bio_sequence->{'molecule_type'}->[0];
my $primary_seq = Bio::PrimarySeq->new(
-id => $accession_number,
-alphabet => $molecule_type,
-seq => $sequence,
-desc => $description,
);
my $seq = Bio::Seq->new (
-display_id => $accession_number,
-accession_number => $accession_number,
-primary_seq => $primary_seq,
-seq => $sequence,
-description => $description,
);
my $organism_name = $bio_sequence->{organism_name}->[0];
if (defined $organism_name) {
my @classification = split(' ', $organism_name);
my $species = Bio::Species->new();
$species->classification(@classification);
$seq->species($species);
}
# Pull out the keywords: $keywords is an array ref.
my $keywords = $bio_sequence->{keyword};
my %key_to_value;
for my $keywords (@$keywords) {
# print "keywords: $keywords\n";
my @words = split(':', $keywords);
for (my $i = 0; $i < scalar @words - 1; $i++) {
if ($i % 2 == 0) {
my $j = $i; $j++;
# print "$words[$i] , $words[$j]\n";
$key_to_value{$words[$i]} = $words[$j];
}
}
# print Data::Dumper->Dump([%key_to_value]);
my $reference = Bio::Annotation::Reference->
new(-authors => $key_to_value{authors},
-title => $key_to_value{title},
-database => $key_to_value{database},
-pubmed => $key_to_value{pubmed},
);
$seq->annotation->add_Annotation('reference', $reference);
} # close for my $keywords
# print Data::Dumper->Dump([$bio_sequence]); print "here\n"; exit;
if (defined $bio_sequence->{'sequence_map'}) {
Bio/SeqIO/agave.pm view on Meta::CPAN
# Get the sequence features (ie genes, exons, etc) from this $sequence_map
for my $seq_feature (@{$sequence_map->{'annotations'}->{'seq_feature'}}) {
# print Data::Dumper->Dump([$seq_feature]); exit;
my $seq_location = $seq_feature->{'seq_location'};
my $start_coord = $seq_feature->{'least_start'}->[0];
my $feature_type = $seq_feature->{'feature_type'}->[0];
my $end_coord = $seq_feature->{'greatest_end'}->[0];
my $is_on_complement = $seq_feature->{'is_on_complement'}->[0];
# Specify the coordinates and the tag for this seq feature.
# print "Primary Tag for this SeqFeature: $feature_type\n";
my $feat = Bio::SeqFeature::Generic->
new(
-start => $start_coord,
-end => $end_coord,
-primary_tag => $feature_type,
);
if (defined $seq_feature->{'qualifier'} &&
ref($seq_feature->{'qualifier'}) eq 'ARRAY') {
for my $feature (@{$seq_feature->{'qualifier'}}) {
my $value = $feature->{'qualifier'};
my $feature_type = $feature->{'qualifier_type'};
for (my $i = 0;
$i < scalar @{$value};
$i++) {
$feat->add_tag_value(
$feature_type->[$i] => $value->[$i]
);
} # close the for loop
}
} # close if (defined $seq_feature->...
$seq->add_SeqFeature($feat);
} # close for my $seq_feature (@{$sequence_map->...
} # close if (defined $sequence_map->{annotations} &&
} # close for my $sequence_map (@{$bio_sequence->{'sequence_map'}}){
} # close if (defined $bio_sequence->{'sequence_map'}){
# This is where the Bio::Seq objects are stored:
push @{$self->{'sequence_objects'}}, $seq;
} # close for my $fragment_orientation
} # close for my $fragment_order
} # close for my $contig
} # close for my $sciobj
# Flag is set so that we know that the sequence objects are now stored in $self.
$self->{'seqs_stored'} = 1;
return;
}
# ==================================================================================
=head2 next_seq
Title : next_seq
Usage : $seq = $stream->next_seq()
Function : Returns the next sequence in the stream.
Args : None.
Returns : Bio::Seq object
Method is called from the script. Method(s) that this method calls:
_store_seqs (only once throughout the life time of script execution).
=cut
sub next_seq {
my ($self) = @_;
# convert agave to genbank/fasta/embl whatever.
$self->_store_seqs if $self->{'seqs_stored'} == 0;
$self->throw("Error: No Bio::Seq objects stored yet!\n\n")
if !defined $self->{'sequence_objects'}; # This should never occur...
if (scalar @{$self->{'sequence_objects'}} > 0) {
return shift @{$self->{'sequence_objects'}};
} else {
# All done. Nothing more to parse.
# print "returning nothing!\n";
return;
}
}
# ==================================================================================
=head2 next_primary_seq
Title : next_primary_seq
Usage : $seq = $stream->next_primary_seq()
Function: returns the next primary sequence (ie no seq_features) in the stream
Returns : Bio::PrimarySeq object
Args : NONE
Bio/SeqIO/agave.pm view on Meta::CPAN
sub write_seq {
# Convert the Bio::Seq object(s) to AGAVE xml file.
my ($self,@seqs) = @_;
foreach my $seq ( @seqs ) {
$self->_write_each_record( $seq ); # where most of the work actually takes place.
}
return;
}
# ==================================================================================
=head2 _write_each_record
Title : _write_each_record
Usage : $agave->_write_each_record( $seqI )
Function: change data into agave format
Returns : NONE
Args : Bio::SeqI object
=cut
sub _write_each_record {
my ($self,$seq) = @_;
# $self->{'file'} =~ s/>//g;
my $output = IO::File->new(">" . $self->{'file'});
my $writer = XML::Writer->new(OUTPUT => $output,
NAMESPACES => 0,
DATA_MODE => 1,
DATA_INDENT => 2 ) ;
$writer->xmlDecl("UTF-8");
$writer->doctype("sciobj", '', "sciobj.dtd");
$writer ->startTag('sciobj',
'version', '2',
'release', '2');
$writer->startTag('contig', 'length', $seq->length);
my $annotation = $seq ->annotation;
# print "annotation: $annotation\n"; exit; Bio::Annotation::Collection=HASH(0x8112e6c)
if ( $annotation->get_Annotations('dblink') ) {
# used to be $annotation->each_DBLink, but Bio::Annotation::Collection::each_DBLink
# is now replaced with get_Annotations('dblink')
my $dblink = $annotation->get_Annotations('dblink')->[0] ;
$writer ->startTag('db_id',
'id', $dblink->primary_id ,
'db_code', $dblink->database );
} else {
$writer ->startTag('db_id',
'id', $seq->display_id ,
'db_code', 'default' );
}
$writer ->endTag('db_id') ;
$writer->startTag('fragment_order');
$writer->startTag('fragment_orientation');
##start bio_sequence
####my $organism = $seq->species->genus . " " . $seq->species->species;
$writer ->startTag('bio_sequence',
'sequence_id', $seq->display_id,
'seq_length', $seq->length,
# 'molecule_type', $seq->moltype, # deprecated
'molecule_type', $self->alphabet,
#'organism_name', $organism
);
# my $desc = $seq->{primary_seq}->{desc};
# print "desc: $desc\n"; exit;
# print Data::Dumper->Dump([$seq]); exit;
##start db_id under bio_sequence
$annotation = $seq ->annotation;
# print "annotation: $annotation\n"; exit; Bio::Annotation::Collection=HASH(0x8112e6c)
if ( $annotation->get_Annotations('dblink') ) {
# used to be $annotation->each_DBLink, but Bio::Annotation::Collection::each_DBLink
# is now replaced with get_Annotations('dblink')
my $dblink = $annotation->get_Annotations('dblink')->[0] ;
$writer ->startTag('db_id',
'id', $dblink->primary_id ,
'db_code', $dblink->database );
} else {
$writer ->startTag('db_id',
'id', $seq->display_id ,
'db_code', 'default' );
}
$writer ->endTag('db_id') ;
##start note
my $note = "" ;
foreach my $comment ( $annotation->get_Annotations('comment') ) {
# used to be $annotations->each_Comment(), but that's now been replaced
# with get_Annotations()
# $comment is a Bio::Annotation::Comment object
$note .= $comment->text() . "\n";
}
$writer ->startTag('note');
$writer ->characters( $note ) ;
$writer ->endTag('note');
##start description
$writer ->startTag('description');
# $writer ->characters( $annotation->get_Annotations('description') ) ;
# used to be $annotations->each_description(), but that's now been
# replaced with get_Annotations.
# Simon added this: this is the primary_seq's desc (the DEFINITION tag in a genbank file)
$writer->characters($seq->{primary_seq}->{desc});
$writer ->endTag('description');
##start keywords
foreach my $genename ( $annotation->get_Annotations('gene_name') ) {
# used to be $annotations->each_gene_name, but that's now been
# replaced with get_Annotations()
$writer ->startTag('keyword');
Bio/SeqIO/agave.pm view on Meta::CPAN
my $medline = $ref->medline || 'null';
my $pubmed = $ref->pubmed || 'null';
my $database = $ref->database || 'null';
my $authors = $ref->authors || 'null';
my $title = $ref->title || 'null';
$writer ->characters( 'medline:' . "$medline" . ':' . 'pubmed:' .
"$pubmed" . ':' . 'database:' . "$database" .
':' .'authors:' . "$authors" . ':' . 'title:' . "$title" ) ;
$writer ->endTag('keyword');
}
## start sequence
$writer ->startTag('sequence');
$writer ->characters( $seq->seq ) ;
$writer ->endTag('sequence');
## start xrefs
$writer ->startTag('xrefs');
foreach my $link ( $annotation->get_Annotations('dblink') ) {
# link is a Bio::Annotation::DBLink object
$writer ->startTag('db_id',
'db_code', $link->database,
'id', $link->primary_id);
$writer ->characters( $link->comment ) ;
$writer ->endTag('db_id');
}
$writer ->endTag('xrefs') ;
##start sequence map
##we can not use : my @feats = $seq->all_SeqFeatures;
##rather, we use top_SeqFeatures() to keep the tree structure
my @feats = $seq->top_SeqFeatures ;
my $features;
##now we need cluster top level seqfeature by algorithm
my $maps;
foreach my $feature (@feats) {
my $map_type = $feature ->source_tag;
push (@{$maps->{ $map_type }}, $feature);
}
##now we enter each sequence_map
foreach my $map_type (keys %$maps ) {
$writer->startTag('sequence_map',
'label', $map_type );
$writer->startTag('annotations');
# the original author accidently entered 'annotation' instead of 'annotations'
foreach my $feature ( @{$maps->{ $map_type }} ) {
$self->_write_seqfeature( $feature, $writer ) ;
}
$writer->endTag('annotations');
$writer->endTag('sequence_map');
}
$writer->endTag('bio_sequence');
$writer->endTag('fragment_orientation');
$writer->endTag('fragment_order');
$writer->endTag('contig');
$writer->endTag('sciobj');
}
# ==================================================================================
=head2 _write_seqfeature
Usage : $agave->_write_each_record( $seqfeature, $write )
Function: change seeqfeature data into agave format
Returns : NONE
Args : Bio::SeqFeature object and XML::writer object
=cut
sub _write_seqfeature{
my ($self,$seqf, $writer) = @_;
##now enter seq feature
$writer ->startTag('seq_feature',
'feature_type', $seqf->primary_tag() );
my $strand = $seqf->strand();
$strand = 0 if !defined $strand;
# $strand == 1 ? 'false' : 'true';
my $is_on_complement;
if ($strand == 1) {
$is_on_complement = 'true';
} else {
$is_on_complement = 'false';
}
# die Data::Dumper->Dump([$seqf]) if !defined $strand;
$writer ->startTag('seq_location',
'lease_start', $seqf->start(),
'greatest_end', $seqf->end(),
# 'is_on_complement', $seqf->strand() == 1 ? 'false' : 'true') ;
'is_on_complement' , $is_on_complement);
# is_on_complement: is the feature found on the complementary
# strand (true) or not (false)?
$writer ->endTag('seq_location');
##enter qualifier
foreach my $tag ( $seqf->all_tags() ) {
$writer ->startTag('qualifier',
'qualifier_type', $tag);
$writer ->characters( $seqf->each_tag_value($tag) ) ;
$writer ->endTag('qualifier');
}
##now recursively travel the seqFeature
foreach my $subfeat ( $seqf->sub_SeqFeature ) {
$self->_write_seqfeature( $subfeat, $writer ) ;
}
$writer->endTag('seq_feature');
return;
( run in 1.040 second using v1.01-cache-2.11-cpan-364913b4093 )