[Bio] / FortyEight / SeedExport.pm Repository:
ViewVC logotype

Diff of /FortyEight/SeedExport.pm

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.7, Thu Feb 15 22:41:39 2007 UTC revision 1.12, Fri Sep 25 13:40:30 2009 UTC
# Line 10  Line 10 
10  use Bio::FeatureIO;  use Bio::FeatureIO;
11  use Bio::SeqIO;  use Bio::SeqIO;
12  use Bio::Seq;  use Bio::Seq;
13    use Bio::Location::Split;
14    use Bio::Location::Simple;
15  use Bio::SeqFeature::Generic;  use Bio::SeqFeature::Generic;
16  use Bio::SeqFeature::Annotated;  use Bio::SeqFeature::Annotated;
17  use File::Basename;  use File::Basename;
# Line 22  Line 24 
24    my $genome          = $parameters->{'genome'};    my $genome          = $parameters->{'genome'};
25    my $directory       = $parameters->{'directory'};    my $directory       = $parameters->{'directory'};
26    my $format          = $parameters->{'export_format'};    my $format          = $parameters->{'export_format'};
27        my $strip_ec        = $parameters->{'strip_ec'};
28        my $filename        = $parameters->{'filename'};
29    
30    # set some variables    # set some variables
31    my $user = "master:master";    my $user = "master:master";
# Line 86  Line 90 
90      my $cloc = $contig.'_1_'.$fig->contig_ln($genome, $contig);      my $cloc = $contig.'_1_'.$fig->contig_ln($genome, $contig);
91      my $seq = $fig->dna_seq($genome, $cloc);      my $seq = $fig->dna_seq($genome, $cloc);
92      $seq =~ s/>//;      $seq =~ s/>//;
93      $bio->{$contig} = Bio::Seq->new(-id => "$contig", seq => $seq);          $bio->{$contig} = Bio::Seq->new( -id => $contig,
94                                            -seq => $seq,
95                                           );
96      $bio->{$contig}->desc("Contig $contig from $gs");      $bio->{$contig}->desc("Contig $contig from $gs");
97    
98      my $feature = Bio::SeqFeature::Generic->new(      my $feature = Bio::SeqFeature::Generic->new(
99                                                  -start  => 1,                                                  -start  => 1,
100                                                  -end    => $fig->contig_ln($genome, $contig),                                                  -end    => $fig->contig_ln($genome, $contig),
101                                                  -tag        => { dbxref     => "taxon: $taxid",                                                      -tag        => { db_xref     => "taxon: $taxid",
102                                                                   organism   => $gs,                                                                   organism   => $gs,
103                                                                   mol_type   => "genomic DNA",                                                                   mol_type   => "genomic DNA",
104                                                                   genome_id  => $genome,                                                                   genome_id  => $genome,
# Line 133  Line 140 
140      my %gos = map { $_ => 1 } @gos if (scalar(@gos));      my %gos = map { $_ => 1 } @gos if (scalar(@gos));
141    
142      # add GOs      # add GOs
143      push @{$note->{"Dbxref"}}, @gos;          push @{$note->{"db_xref"}}, @gos;
144    
145      # add ECs      # add ECs
146      push @{$note->{"EC_number"}}, keys(%ecs);      push @{$note->{"EC_number"}}, keys(%ecs);
# Line 143  Line 150 
150      my @rw_aliases = map { $fig->rewrite_db_xrefs($_->[0]) } $fig->mapped_prot_ids($pid);      my @rw_aliases = map { $fig->rewrite_db_xrefs($_->[0]) } $fig->mapped_prot_ids($pid);
151      my @aliases;      my @aliases;
152      foreach my $a (@rw_aliases) {      foreach my $a (@rw_aliases) {
153        push @{$note->{"Dbxref"}}, $a if ( $a );              push @{$note->{"db_xref"}}, $a if ( $a );
154      }      }
155    
156      # get the links      # get the links
# Line 162  Line 169 
169          print STDERR "Ignored link: $ln\n";          print STDERR "Ignored link: $ln\n";
170          next;          next;
171        }        }
172        push @{$note->{"Dbxref"}}, "$db:$id";              push @{$note->{"db_xref"}}, "$db:$id";
173      }      }
174    
175      # add FIG id as a note      # add FIG id as a note
176      push @{$note->{"Dbxref"}}, "FIG_ID:$peg";          # push @{$note->{"db_xref"}}, "FIG_ID:$peg";
177    
178      # get the features      # get the features
179    
180            my $loc_obj;
181      my @location = $fig->feature_location($peg);      my @location = $fig->feature_location($peg);
182            my @loc_info;
183            my $contig;
184      foreach my $loc (@location) {      foreach my $loc (@location) {
185                my($start, $stop);
186        $loc =~ /^(.*)\_(\d+)\_(\d+)$/;        $loc =~ /^(.*)\_(\d+)\_(\d+)$/;
187        my ($contig, $start, $stop) = ($1, $2, $3);              ($contig, $start, $stop) = ($1, $2, $3);
188        my $original_contig = $contig;        my $original_contig = $contig;
189        my $strand = '1';              my $strand = '+';
190        my $frame = $start % 3;        my $frame = $start % 3;
191        if ($start > $stop) {        if ($start > $stop) {
192          $frame = $stop % 3;          $frame = $stop % 3;
193          ($start, $stop, $strand) = ($stop, $start, '-1');                  ($start, $stop, $strand) = ($stop, $start, '-');
194        } elsif ($start == $stop) {        } elsif ($start == $stop) {
195          $strand = ".";          $strand = ".";
196          $frame = ".";          $frame = ".";
197        }        }
198    
199                push(@loc_info, [$contig, $start, $stop, $strand, $frame]);
200    
201                my $sloc = new Bio::Location::Simple(-start => $start,
202                                                     -end => $stop,
203                                                     -strand => $strand);
204                if (@location == 1)
205                {
206                    $loc_obj = $sloc;
207                }
208                else
209                {
210                    $loc_obj = new Bio::Location::Split() if !$loc_obj;
211                    $loc_obj->add_sub_Location($sloc);
212                }
213            }
214    
215        my $source = "FIG";        my $source = "FIG";
216        my $type = $fig->ftype($peg);        my $type = $fig->ftype($peg);
217        my $feature;        my $feature;
218    
219            # strip EC from function
220            my $func_ok = $func;
221            if ($strip_ec) {
222                $func_ok =~ s/\s+\(EC \d+\.(\d+|-)\.(\d+|-)\.(\d+|-)\)//g;
223                $func_ok =~ s/\s+\(TC \d+\.(\d+|-)\.(\d+|-)\.(\d+|-)\)//g;
224            }
225    
226        if ($type eq "peg") {        if ($type eq "peg") {
227          $feature = Bio::SeqFeature::Generic->new(              $feature = Bio::SeqFeature::Generic->new(-location => $loc_obj,
                                                  -start    => $start,  
                                                  -end      => $stop,  
                                                  -strand   => $strand,  
228                                                   -primary  => 'CDS',                                                   -primary  => 'CDS',
229                                                   -tag      => {                                                   -tag      => {
230                                                                 product     => $func,                                                           product     => $func_ok,
231                                                                 translation => $fig->get_translation($peg),                                                                 translation => $fig->get_translation($peg),
232                                                                },                                                                },
233                                                  );                                                  );
# Line 202  Line 236 
236            $feature->add_tag_value($tagtype, @{$note->{$tagtype}});            $feature->add_tag_value($tagtype, @{$note->{$tagtype}});
237          }          }
238    
         my $feature2 = Bio::SeqFeature::Annotated->new( -start    => $start,  
                                                         -end      => $stop,  
                                                         -strand   => $strand,  
                                                         -phase    => 0,  
                                                         -frame    => $frame,  
                                                         -source   => $source,  
                                                         -type     => "CDS",  
                                                         -seq_id   => $contig );  
   
         push(@$bio2, $feature2);  
   
239          # work around to get annotations into gff          # work around to get annotations into gff
240          push @$gff_export, "$contig\t$source\tCDS\t$start\t$stop\t.\t$strand\t$frame\t$func\n"              # this is probably still wrong for split locations.
241                $func_ok =~ s/ #.+//;
242                $func_ok =~ s/;/%3B/g;
243                $func_ok =~ s/,/2C/g;
244                $func_ok =~ s/=//g;
245                for my $l (@loc_info)
246                {
247                  my $ec = "";
248                  my @ecs = ($func =~ /[\(\[]*EC[\s:]?(\d+\.[\d-]+\.[\d-]+\.[\d-]+)[\)\]]*/ig);
249                  if (scalar(@ecs)) {
250                    $ec = ";Ontology_term=".join(',', map { "KEGG_ENZYME:" . $_ } @ecs);
251                  }
252                  my($contig, $start, $stop, $strand, $frame) = @$l;
253                  push @$gff_export, "$contig\t$source\tCDS\t$start\t$stop\t.\t$strand\t$frame\tID=".$peg.";Name=".$func_ok.$ec."\n";
254                }
255    
256    
257        } elsif ($type eq "rna") {        } elsif ($type eq "rna") {
258          $feature = Bio::SeqFeature::Generic->new(              my $primary;
259                                                   -start    => $start,              if ( $func =~ /tRNA/ ) {
260                                                   -end      => $stop,                  $primary = 'tRNA';
261                                                   -strand   => $strand,              } elsif ( $func =~ /(Ribosomal RNA|5S RNA)/ ) {
262                                                   -primary  => 'RNA',                  $primary = 'rRNA';
263                } else {
264                    $primary = 'RNA';
265                }
266    
267                $feature = Bio::SeqFeature::Generic->new(-location => $loc_obj,
268                                                         -primary  => $primary,
269                                                         -tag      => {
270                                                             product => $func,
271                                                         },
272    
273                                                  );                                                  );
274                $func_ok =~ s/ #.+//;
275                $func_ok =~ s/;/%3B/g;
276                $func_ok =~ s/,/2C/g;
277                $func_ok =~ s/=//g;
278          foreach my $tagtype (keys %$note) {          foreach my $tagtype (keys %$note) {
279            $feature->add_tag_value($tagtype, @{$note->{$tagtype}});            $feature->add_tag_value($tagtype, @{$note->{$tagtype}});
280    
281                    # work around to get annotations into gff
282                    for my $l (@loc_info)
283                    {
284                        my($contig, $start, $stop, $strand, $frame) = @$l;
285                        push @$gff_export, "$contig\t$source\t$primary\t$start\t$stop\t.\t$strand\t.\tID=$peg;Name=$func_ok\n";
286                    }
287          }          }
288    
289        } else {        } else {
290          print STDERR "unhandled feature type: $type\n";          print STDERR "unhandled feature type: $type\n";
291        }        }
   
292        $bio->{$contig}->add_SeqFeature($feature);        $bio->{$contig}->add_SeqFeature($feature);
293      }      }
294    }  
295    
296    # generate filename    # generate filename
297    my $filename = $directory . $genome . "." . $format_ending;      if (!$filename)
298        {
299            $filename = $directory . $genome . "." . $format_ending;
300        }
301    
302    # check for FeatureIO or SeqIO    # check for FeatureIO or SeqIO
303    if ($format eq "GTF") {    if ($format eq "GTF") {

Legend:
Removed from v.1.7  
changed lines
  Added in v.1.12

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3