[Bio] / FigKernelPackages / Observation.pm Repository:
ViewVC logotype

Diff of /FigKernelPackages/Observation.pm

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

revision 1.3, Tue Jun 12 20:17:44 2007 UTC revision 1.24, Tue Jul 10 20:11:38 2007 UTC
# Line 1  Line 1 
1  package Observation;  package Observation;
2    
3    use lib '/vol/ontologies';
4    use DBMaster;
5    
6  require Exporter;  require Exporter;
7  @EXPORT_OK = qw(get_objects);  @EXPORT_OK = qw(get_objects);
8    
9    use FIG_Config;
10  use strict;  use strict;
11  use warnings;  #use warnings;
12    use HTML;
13    
14  1;  1;
15    
# Line 21  Line 26 
26    
27  The data can be used to display information for a given sequence feature (protein or other, but primarily information is computed for proteins).  The data can be used to display information for a given sequence feature (protein or other, but primarily information is computed for proteins).
28    
 Example:  
   
 use FIG;  
 use Observation;  
   
 my $fig = new FIG;  
 my $fid = "fig|83333.1.peg.3";  
   
 my $observations = Observation::get_objects($fid);  
 foreach my $observation (@$observations) {  
     print "ID: " . $fid . "\n";  
     print "Start: " . $observation->start() . "\n";  
     ...  
 }  
   
 B<return an array of objects>  
   
   
 print "$Observation->acc\n" prints the Accession number if present for the Observation  
   
29  =cut  =cut
30    
31  =head1 BACKGROUND  =head1 BACKGROUND
# Line 64  Line 49 
49    
50  The public methods this package provides are listed below:  The public methods this package provides are listed below:
51    
 =head3 acc()  
52    
53  A valid accession or remote ID (in the style of a db_xref) or a valid local ID (FID) in case this is supported.  =head3 context()
54    
55    Returns close or diverse for purposes of displaying genomic context
56    
57  =cut  =cut
58    
59  sub acc {  sub context {
60    my ($self) = @_;    my ($self) = @_;
61    
62    return $self->{acc};    return $self->{context};
63  }  }
64    
65  =head3 description()  =head3 rows()
66    
67  The description of the hit. Taken from the data or from the our Ontology database for some cases e.g. IPR or PFAM.  each row in a displayed table
   
 B<Please note:>  
 Either remoteid or description is required.  
68    
69  =cut  =cut
70    
71  sub description {  sub rows {
72    my ($self) = @_;    my ($self) = @_;
73    
74      return $self->{rows};
75    }
76    
77    =head3 acc()
78    
79    A valid accession or remote ID (in the style of a db_xref) or a valid local ID (FID) in case this is supported.
80    
81    =cut
82    
83    sub acc {
84      my ($self) = @_;
85    return $self->{acc};    return $self->{acc};
86  }  }
87    
# Line 100  Line 94 
94    
95  =over 9  =over 9
96    
97    =item IDENTICAL (seq)
98    
99  =item SIM (seq)  =item SIM (seq)
100    
101  =item BBH (seq)  =item BBH (seq)
# Line 114  Line 110 
110    
111  =item PFAM (dom)  =item PFAM (dom)
112    
113  =item SIGNALP (dom)  =item SIGNALP_CELLO_TMPRED (loc)
114    
115  =item  CELLO(loc)  =item PDB (seq)
116    
117  =item TMHMM (loc)  =item TMHMM (loc)
118    
# Line 182  Line 178 
178    return $self->{stop};    return $self->{stop};
179  }  }
180    
181  =head3 evalue()  =head3 start()
182    
183  E-value or P-Value if present.  Start of hit in query sequence.
184    
185  =cut  =cut
186    
187  sub evalue {  sub qstart {
188    my ($self) = @_;    my ($self) = @_;
189    
190    return $self->{evalue};      return $self->{qstart};
191  }  }
192    
193  =head3 score()  =head3 qstop()
   
 Score if present.  
194    
195  B<Please note: >  End of the hit in query sequence.
 Either score or eval are required.  
196    
197  =cut  =cut
198    
199  sub score {  sub qstop {
200    my ($self) = @_;    my ($self) = @_;
201    
202    return $self->{score};      return $self->{qstop};
203  }  }
204    
205    =head3 hstart()
206    
207  =head3 display_method()  Start of hit in hit sequence.
208    
209  If available use the function specified here to display the "raw" observation.  =cut
 In the case of a BLAST alignment of fid1 and fid2 a cgi script  
 will be called to display the results of running the command "bl2seq fid1 fid2".  
210    
211  B<Please note> that URL linked to in display_method() is an external component and needs to added to the code for every class of evidence.  sub hstart {
212        my ($self) = @_;
213    
214        return $self->{hstart};
215    }
216    
217    =head3 end()
218    
219    End of the hit in hit sequence.
220    
221  =cut  =cut
222    
223  sub display_method {  sub hstop {
224    my ($self) = @_;    my ($self) = @_;
225    
226    # add code here      return $self->{hstop};
   
   return $self->{display_method};  
227  }  }
228    
229  =head3 rank()  =head3 qlength()
230    
231  Returns an integer from 1 - 10 indicating the importance of this observations.  length of the query sequence in similarities
   
 Currently always returns 1.  
232    
233  =cut  =cut
234    
235  sub rank {  sub qlength {
236    my ($self) = @_;    my ($self) = @_;
237    
238  #  return $self->{rank};      return $self->{qlength};
   
   return 1;  
239  }  }
240    
241  =head3 supports_annotation()  =head3 hlength()
   
 Does a this observation support the annotation of its feature?  
242    
243  Returns  length of the hit sequence in similarities
244    
245  =over 3  =cut
246    
247  =item 10, if feature annotation is identical to $self->description  sub hlength {
248        my ($self) = @_;
249    
250  =item 1, Feature annotation is similar to $self->annotation; this is computed using FIG::SameFunc()      return $self->{hlength};
251    }
252    
253  =item undef  =head3 evalue()
254    
255  =back  E-value or P-Value if present.
256    
257  =cut  =cut
258    
259  sub supports_annotation {  sub evalue {
260    my ($self) = @_;    my ($self) = @_;
261    
262    # no code here so far    return $self->{evalue};
   
   return $self->{supports_annotation};  
263  }  }
264    
265  =head3 url()  =head3 score()
266    
267  URL describing the subject. In case of a BLAST hit against a sequence, this URL will lead to a page displaying the sequence record for the sequence. In case of an HMM hit, the URL will be to the URL description.  Score if present.
268    
269  =cut  =cut
270    
271  sub url {  sub score {
272    my ($self) = @_;    my ($self) = @_;
273      return $self->{score};
274    }
275    
276    =head3 display()
277    
278    my $url = get_url($self->type, $self->acc);  will be different for each type
279    
280    =cut
281    
282    sub display {
283    
284      die "Abstract Method Called\n";
285    
   return $url;  
286  }  }
287    
288  =head3 get_objects()  =head3 display_table()
289    
290  This is the B<REAL WORKHORSE> method of this Package.  will be different for each type
291    
292    =cut
293    
294    sub display_table {
295    
296  It will probably have to:    die "Abstract Table Method Called\n";
297    
 - get all sims for the feature  
 - get all bbhs for the feature  
 - copy information from sim to bbh (bbh have no match location etc)  
 - get pchs (difficult)  
 - get attributes (there is code for this that in get_attribute_based_observations  
 - get_attributes_based_observations returns an array of arrays of hashes like this"  
   
   my $datasets =  
      [  
        [ { name => 'acc', value => '1234' },  
         { name => 'from', value => '4' },  
         { name => 'to', value => '400' },  
         ....  
        ],  
        [ { name => 'acc', value => '456' },  
         { name => 'from', value => '1' },  
         { name => 'to', value => '100' },  
         ....  
        ],  
        ...  
      ];  
    return $datasets;  
298   }   }
299    
300  It will invoke the required calls to the SEED API to retrieve the information required.  =head3 get_objects()
301    
302    This is the B<REAL WORKHORSE> method of this Package.
303    
304  =cut  =cut
305    
306  sub get_objects {  sub get_objects {
307      my ($self,$fid) = @_;      my ($self,$fid,$scope) = @_;
308    
309    my $objects = [];    my $objects = [];
310    my @matched_datasets=();    my @matched_datasets=();
311    
312    # call function that fetches attribute based observations    # call function that fetches attribute based observations
313    # returns an array of arrays of hashes    # returns an array of arrays of hashes
   #  
   get_attribute_based_observations($fid,\@matched_datasets);  
314    
315    # read sims      if($scope){
316            get_cluster_observations($fid,\@matched_datasets,$scope);
317        }
318        else{
319            my %domain_classes;
320            $domain_classes{'CDD'} = 1;
321            get_identical_proteins($fid,\@matched_datasets);
322            get_attribute_based_domain_observations($fid,\%domain_classes,\@matched_datasets);
323    get_sims_observations($fid,\@matched_datasets);    get_sims_observations($fid,\@matched_datasets);
324            get_functional_coupling($fid,\@matched_datasets);
325    # read sims + bbh (enrich BBHs with sims coordindates etc)          get_attribute_based_location_observations($fid,\@matched_datasets);
326    # read pchs          get_pdb_observations($fid,\@matched_datasets);
327    # read figfam match data from 48hr directory (BobO knows how do do this!)      }
   # what sources of evidence did I miss?  
328    
329    foreach my $dataset (@matched_datasets) {    foreach my $dataset (@matched_datasets) {
330      my $object = $self->new();          my $object;
331      foreach my $attribute (@$dataset) {          if($dataset->{'type'} eq "dom"){
332        $object->{$attribute->{'name'}} = $attribute->{'value'};              $object = Observation::Domain->new($dataset);
333      }      }
334  #    $object->{$attribute->{'feature_id'}} = $attribute->{$fid};          if($dataset->{'class'} eq "PCH"){
335      push (@$objects, $object);              $object = Observation::FC->new($dataset);
336            }
337            if ($dataset->{'class'} eq "IDENTICAL"){
338                $object = Observation::Identical->new($dataset);
339            }
340            if ($dataset->{'class'} eq "SIGNALP_CELLO_TMPRED"){
341                $object = Observation::Location->new($dataset);
342            }
343            if ($dataset->{'class'} eq "SIM"){
344                $object = Observation::Sims->new($dataset);
345            }
346            if ($dataset->{'class'} eq "CLUSTER"){
347                $object = Observation::Cluster->new($dataset);
348            }
349            if ($dataset->{'class'} eq "PDB"){
350                $object = Observation::PDB->new($dataset);
351      }      }
352    
353            push (@$objects, $object);
354        }
355    
356    return $objects;    return $objects;
357    
358  }  }
359    
360  =head1 Internal Methods  =head1 Internal Methods
# Line 358  Line 365 
365    
366  =cut  =cut
367    
368    sub get_attribute_based_domain_observations{
369    
370  =head3 get_url (internal)      # we read a FIG ID and a reference to an array (of arrays of hashes, see above)
371        my ($fid,$domain_classes,$datasets_ref) = (@_);
372    
373  get_url() return a valid URL or undef for any observation.      my $fig = new FIG;
374    
375  URLs are constructed by looking at the Accession acc()  and  name()      foreach my $attr_ref ($fig->get_attributes($fid)) {
376            my $key = @$attr_ref[1];
377            my @parts = split("::",$key);
378            my $class = $parts[0];
379    
380            if($domain_classes->{$parts[0]}){
381                my $val = @$attr_ref[2];
382                if($val =~/^(\d+\.\d+|0\.0);(\d+)-(\d+)/){
383                    my $raw_evalue = $1;
384                    my $from = $2;
385                    my $to = $3;
386                    my $evalue;
387                    if($raw_evalue =~/(\d+)\.(\d+)/){
388                        my $part2 = 1000 - $1;
389                        my $part1 = $2/100;
390                        $evalue = $part1."e-".$part2;
391                    }
392                    else{
393                        $evalue = "0.0";
394                    }
395    
396                    my $dataset = {'class' => $class,
397                                   'acc' => $key,
398                                   'type' => "dom" ,
399                                   'evalue' => $evalue,
400                                   'start' => $from,
401                                   'stop' => $to,
402                                   'fig_id' => $fid,
403                                   'score' => $raw_evalue
404                                   };
405    
406  Info from both attributes is combined with a table of base URLs stored in this function.                  push (@{$datasets_ref} ,$dataset);
407                }
408            }
409        }
410    }
411    
412    sub get_attribute_based_location_observations{
413    
414        my ($fid,$datasets_ref) = (@_);
415        my $fig = new FIG;
416    
417        my $location_attributes = ['SignalP','CELLO','TMPRED'];
418    
419        my $dataset = {'type' => "loc", 'class' => 'SIGNALP_CELLO_TMPRED','fig_id' => $fid};
420        foreach my $attr_ref ($fig->get_attributes($fid,$location_attributes)) {
421            my $key = @$attr_ref[1];
422            my @parts = split("::",$key);
423            my $sub_class = $parts[0];
424            my $sub_key = $parts[1];
425            my $value = @$attr_ref[2];
426            if($sub_class eq "SignalP"){
427                if($sub_key eq "cleavage_site"){
428                    my @value_parts = split(";",$value);
429                    $dataset->{'cleavage_prob'} = $value_parts[0];
430                    $dataset->{'cleavage_loc'} = $value_parts[1];
431                }
432                elsif($sub_key eq "signal_peptide"){
433                    $dataset->{'signal_peptide_score'} = $value;
434                }
435            }
436            elsif($sub_class eq "CELLO"){
437                $dataset->{'cello_location'} = $sub_key;
438                $dataset->{'cello_score'} = $value;
439            }
440            elsif($sub_class eq "TMPRED"){
441                my @value_parts = split(";",$value);
442                $dataset->{'tmpred_score'} = $value_parts[0];
443                $dataset->{'tmpred_locations'} = $value_parts[1];
444            }
445        }
446    
447        push (@{$datasets_ref} ,$dataset);
448    
449    }
450    
451    =head3 get_pdb_observations() (internal)
452    
453    This methods sets the type and class for pdb observations
454    
455  =cut  =cut
456    
457  sub get_url {  sub get_pdb_observations{
458        my ($fid,$datasets_ref) = (@_);
459    
460   my ($self) = @_;      my $fig = new FIG;
  my $url='';  
461    
462  # a hash with a URL for each observation; identified by name()      foreach my $attr_ref ($fig->get_attributes($fid,'PDB')) {
 #my $URL             => { 'PFAM' => "http://www.sanger.ac.uk/cgi-bin/Pfam/getacc?" ,\  
 #                       'IPR'    => "http://www.ebi.ac.uk/interpro/DisplayIproEntry?ac=" ,\  
 #                          'CDD' => "http://www.ncbi.nlm.nih.gov/Structure/cdd/cddsrv.cgi?uid=",\  
 #                       'PIR'    => "http://www.ncbi.nlm.nih.gov/Structure/cdd/cddsrv.cgi?uid=",\  
 #                       'FIGFAM' => '',\  
 #                          'sim'=> "http://www.theseed.org/linkin.cgi?id=",\  
 #                          'bbh'=> "http://www.theseed.org/linkin.cgi?id="  
 #};  
463    
464  # if (defined $URL{$self->name}) {          my $key = @$attr_ref[1];
465  #     $url = $URL{$self->name}.$self->acc;          my($key1,$key2) =split("::",$key);
466  #     return $url;          my $value = @$attr_ref[2];
467  # }          my ($evalue,$location) = split(";",$value);
468  # else  
469       return undef;          if($evalue =~/(\d+)\.(\d+)/){
470                my $part2 = 1000 - $1;
471                my $part1 = $2/100;
472                $evalue = $part1."e-".$part2;
473            }
474    
475            my($start,$stop) =split("-",$location);
476    
477            my $url = @$attr_ref[3];
478            my $dataset = {'class' => 'PDB',
479                           'type' => 'seq' ,
480                           'acc' => $key2,
481                           'evalue' => $evalue,
482                           'start' => $start,
483                           'stop' => $stop,
484                           'fig_id' => $fid
485                           };
486    
487            push (@{$datasets_ref} ,$dataset);
488  }  }
489    }
490    
491    =head3 get_cluster_observations() (internal)
492    
493  =head3 get_display_method (internal)  This methods sets the type and class for cluster observations
494    
495    =cut
496    
497  get_display_method() return a valid URL or undef for any observation.  sub get_cluster_observations{
498        my ($fid,$datasets_ref,$scope) = (@_);
499    
500        my $dataset = {'class' => 'CLUSTER',
501                       'type' => 'fc',
502                       'context' => $scope,
503                       'fig_id' => $fid
504                       };
505        push (@{$datasets_ref} ,$dataset);
506    }
507    
508    
509    =head3 get_sims_observations() (internal)
510    
511  URLs are constructed by looking at the Accession acc()  and  name()  This methods retrieves sims fills the internal data structures.
 and Info from both attributes is combined with a table of base URLs stored in this function.  
512    
513  =cut  =cut
514    
515  sub get_display_method {  sub get_sims_observations{
516    
517   my ($self) = @_;      my ($fid,$datasets_ref) = (@_);
518        my $fig = new FIG;
519        my @sims= $fig->nsims($fid,100,1e-20,"all");
520        my ($dataset);
521        foreach my $sim (@sims){
522            my $hit = $sim->[1];
523            my $percent = $sim->[2];
524            my $evalue = $sim->[10];
525            my $qfrom = $sim->[6];
526            my $qto = $sim->[7];
527            my $hfrom = $sim->[8];
528            my $hto = $sim->[9];
529            my $qlength = $sim->[12];
530            my $hlength = $sim->[13];
531            my $db = get_database($hit);
532            my $func = $fig->function_of($hit);
533            my $organism = $fig->org_of($hit);
534    
535            $dataset = {'class' => 'SIM',
536                        'acc' => $hit,
537                        'identity' => $percent,
538                        'type' => 'seq',
539                        'evalue' => $evalue,
540                        'qstart' => $qfrom,
541                        'qstop' => $qto,
542                        'hstart' => $hfrom,
543                        'hstop' => $hto,
544                        'database' => $db,
545                        'organism' => $organism,
546                        'function' => $func,
547                        'qlength' => $qlength,
548                        'hlength' => $hlength,
549                        'fig_id' => $fid
550                        };
551    
552  # a hash with a URL for each observation; identified by name()          push (@{$datasets_ref} ,$dataset);
553  #my $URL             => { 'sim'=> "http://www.theseed.org/featalign.cgi?id1=",\      }
554  #                        'bbh'=> "http://www.theseed.org/featalign.cgi?id1="  }
555  # };  
556    =head3 get_database (internal)
557    This method gets the database association from the sequence id
558    
559    =cut
560    
561    sub get_database{
562        my ($id) = (@_);
563    
564        my ($db);
565        if ($id =~ /^fig\|/)              { $db = "FIG" }
566        elsif ($id =~ /^gi\|/)            { $db = "NCBI" }
567        elsif ($id =~ /^^[NXYZA]P_/)      { $db = "RefSeq" }
568        elsif ($id =~ /^sp\|/)            { $db = "SwissProt" }
569        elsif ($id =~ /^uni\|/)           { $db = "UniProt" }
570        elsif ($id =~ /^tigr\|/)          { $db = "TIGR" }
571        elsif ($id =~ /^pir\|/)           { $db = "PIR" }
572        elsif ($id =~ /^kegg\|/)          { $db = "KEGG" }
573        elsif ($id =~ /^tr\|/)            { $db = "TrEMBL" }
574        elsif ($id =~ /^eric\|/)          { $db = "ASAP" }
575        elsif ($id =~ /^img\|/)           { $db = "JGI" }
576    
577        return ($db);
578    
 #if (defined $URL{$self->name}) {  
 #     $url = $URL{$self->name}.$self->feature_id."&id2=".$self->acc;  
 #     return $url;  
 # }  
 # else  
      return undef;  
579  }  }
580    
 =head3 get_attribute_based_evidence (internal)  
581    
582  This method retrieves evidence from the attribute server  =head3 get_identical_proteins() (internal)
583    
584    This methods retrieves sims fills the internal data structures.
585    
586  =cut  =cut
587    
588  sub get_attribute_based_observations{  sub get_identical_proteins{
589    
     # we read a FIG ID and a reference to an array (of arrays of hashes, see above)  
590      my ($fid,$datasets_ref) = (@_);      my ($fid,$datasets_ref) = (@_);
591        my $fig = new FIG;
592        my $funcs_ref;
593    
594      my $_myfig = new FIG;      my @maps_to = grep { $_ ne $fid and $_ !~ /^xxx/ } map { $_->[0] } $fig->mapped_prot_ids($fid);
595    
596      foreach my $attr_ref ($_myfig->get_attributes($fid)) {      foreach my $id (@maps_to) {
597            my ($tmp, $who);
598            if (($id ne $fid) && ($tmp = $fig->function_of($id))) {
599                $who = &get_database($id);
600                push(@$funcs_ref, [$id,$who,$tmp]);
601            }
602        }
603    
604          # convert the ref into a string for easier handling      my ($dataset);
605          my ($string) = "@$attr_ref";      my $dataset = {'class' => 'IDENTICAL',
606                       'type' => 'seq',
607                       'fig_id' => $fid,
608                       'rows' => $funcs_ref
609                       };
610    
611  #       print "S:$string\n";      push (@{$datasets_ref} ,$dataset);
         my ($key,$val) = ( $string =~ /\S+\s(\S+)\s(\S+)/);  
612    
         # THIS SHOULD BE DONE ANOTHER WAY FM->TD  
         # we need to do the right thing for each type, ie no evalue for CELLO and no coordinates, but a score, etc  
         # as fas as possible this should be configured so that the type of observation and the regexp are  
         # stored somewhere for easy expansion  
         #  
613    
614          if (($key =~ /PFAM::/) || ( $key =~ /IPR::/) || ( $key =~ /CDD::/) ) {  }
615    
616    =head3 get_functional_coupling() (internal)
617    
618              # some keys are composite CDD::1233244 or PFAM:PF1233  This methods retrieves the functional coupling of a protein given a peg ID
619    
620    =cut
621    
622    sub get_functional_coupling{
623    
624        my ($fid,$datasets_ref) = (@_);
625        my $fig = new FIG;
626        my @funcs = ();
627    
628        # initialize some variables
629        my($sc,$neigh);
630    
631        # set default parameters for coupling and evidence
632        my ($bound,$sim_cutoff,$coupling_cutoff) = (5000, 1.0e-10, 4);
633    
634        # get the fc data
635        my @fc_data = $fig->coupling_and_evidence($fid,$bound,$sim_cutoff,$coupling_cutoff,1);
636    
637        # retrieve data
638        my @rows = map { ($sc,$neigh) = @$_;
639                         [$sc,$neigh,scalar $fig->function_of($neigh)]
640                      } @fc_data;
641    
642        my ($dataset);
643        my $dataset = {'class' => 'PCH',
644                       'type' => 'fc',
645                       'fig_id' => $fid,
646                       'rows' => \@rows
647                       };
648    
649        push (@{$datasets_ref} ,$dataset);
650    
             if ( $key =~ /::/ ) {  
                 my ($firstkey,$restkey) = ( $key =~ /([a-zA-Z0-9]+)::(.*)/);  
                 $val=$restkey.";".$val;  
                 $key=$firstkey;  
651              }              }
652    
653              my ($acc,$raw_evalue, $from,$to) = ($val =~ /(\S+);(\S+);(\d+)-(\d+)/ );  =head3 new (internal)
654    
655              my $evalue= 255;  Instantiate a new object.
             if (defined $raw_evalue) { # some of the tool do not give us an evalue  
656    
657                  my ($k,$expo) = ( $raw_evalue =~ /(\d+).(\d+)/);  =cut
                 my ($new_k, $new_exp);  
658    
659                  #  sub new {
660                  #  THIS DOES NOT WORK PROPERLY    my ($class,$dataset) = @_;
661                  #  
662                  if($raw_evalue =~/(\d+).(\d+)/){    my $self = { class => $dataset->{'class'},
663                   type => $dataset->{'type'},
664                   fig_id => $dataset->{'fig_id'},
665                   score => $dataset->{'score'},
666               };
667    
668  #                   $new_exp = (1000+$expo);    bless($self,$class);
         #           $new_k = $k / 100;  
669    
670      return $self;
671                  }                  }
672                  $evalue = "0.01"#new_k."e-".$new_exp;  
673    =head3 identity (internal)
674    
675    Returns the % identity of the similar sequence
676    
677    =cut
678    
679    sub identity {
680        my ($self) = @_;
681    
682        return $self->{identity};
683              }              }
684    
685              # unroll it all into an array of hashes  =head3 fig_id (internal)
             # this needs to be done differently for different types of observations  
             my $dataset = [ { name => 'class', value => $key },  
                             { name => 'acc' , value => $acc},  
                             { name => 'type', value => "dom"} , # this clearly needs to be done properly FM->TD  
                             { name => 'evalue', value => $evalue },  
                             { name => 'start', value => $from},  
                             { name => 'stop' , value => $to}  
                             ];  
686    
687              push (@{$datasets_ref} ,$dataset);  =cut
688    
689    sub fig_id {
690      my ($self) = @_;
691      return $self->{fig_id};
692          }          }
693    
694    =head3 feature_id (internal)
695    
696    
697    =cut
698    
699    sub feature_id {
700      my ($self) = @_;
701    
702      return $self->{feature_id};
703      }      }
704    
705    =head3 id (internal)
706    
707    Returns the ID  of the identical sequence
708    
709    =cut
710    
711    sub id {
712        my ($self) = @_;
713    
714        return $self->{id};
715  }  }
716    
717  =head3 get_sims_observations() (internal)  =head3 organism (internal)
718    
719  This methods retrieves sims fills the internal data structures.  Returns the organism  of the identical sequence
720    
721  =cut  =cut
722    
723  sub get_sims_observations{  sub organism {
724        my ($self) = @_;
725    
726      my ($fid,$datasets_ref) = (@_);      return $self->{organism};
     my $fig =  
     my @sims= $fig->nsims($fid,100,1e-5);  
     foreach my $sim (@sims){  
         $hit = $sim->[1];  
         $evalue = $sim->[10];  
         $from = $sim->[8];  
         $to = $sim->[9];  
         my $dataset = [ { name => 'class', value => "SIM" },  
                         { name => 'acc' , value => $hit},  
                         { name => 'type', value => "seq"} ,  
                         { name => 'evalue', value => $evalue },  
                         { name => 'start', value => $from},  
                         { name => 'stop' , value => $to}  
                         ];  
727      }      }
     push (@{$datasets_ref} ,$dataset);  
728    
729    =head3 function (internal)
730    
731    Returns the function of the identical sequence
732    
733    =cut
734    
735    sub function {
736        my ($self) = @_;
737    
738        return $self->{function};
739  }  }
740    
741  =head3 get_sims_and_bbhs() (internal)  =head3 database (internal)
742    
743  This methods retrieves sims and also BBHs and fills the internal data structures.  Returns the database of the identical sequence
744    
745  =cut  =cut
746    
747  #     sub get_sims_and_bbhs{  sub database {
748        my ($self) = @_;
749    
750  #       # blast m8 output format      return $self->{database};
751  #       # id1, id2, %ident, align len, mismatches, gaps, q.start, q.stop, s. start, s.stop, eval, bit  }
   
 #       my $Sims=();  
 #       @sims_src = $fig->sims($fid,80,500,"fig",0);  
 #       print "found $#sims_src SIMs\n";  
 #       foreach $sims (@sims_src) {  
 #           my ($sims_string) = "@$sims";  
 # #       print "$sims_string\n";  
 #           my ($rfid,$start,$stop,$eval) = ( $sims_string =~ /\S+\s+(\S+)\s+\S+\s\S+\s+(\S+)\s+(\S+)\s+  
 #                                             \S+\s+\S+\s+\S+\s+\S+\s+(\S+)+.*/);  
 # #       print "ID: $rfid, E:$eval, Start:$start stop:$stop\n";  
 #           $Sims{$rfid}{'eval'}=$eval;  
 #           $Sims{$rfid}{'start'}=$start;  
 #           $Sims{$rfid}{'stop'}=$stop;  
 #           print "$rfid $Sims{$rfid}{'eval'}\n";  
 #       }  
   
 #       # BBHs  
 #       my $BBHs=();  
   
 #       @bbhs_src = $fig->bbhs($fid,1.0e-10);  
 #       print "found $#bbhs_src BBHs\n";  
 #       foreach $bbh (@bbhs_src) {  
 #           #print "@$bbh\n";  
 #           my ($bbh_string) = "@$bbh";  
 #           my ($rfid,$eval,$score) = ( $bbh_string =~ /(\S+)\s(\S+)\s(\S+)/);  
 #           #print "ID: $rfid, E:$eval, S:$score\n";  
 #           $BBHs{$rfid}{'eval'}=$eval;  
 #           $BBHs{$rfid}{'score'}=$score;  
 # #print "$rfid $BBHs{$rfid}{'eval'}\n";  
 #       }  
752    
753  #     }  sub score {
754      my ($self) = @_;
755    
756      return $self->{score};
757    }
758    
759    ############################################################
760    ############################################################
761    package Observation::PDB;
762    
763  =head3 new (internal)  use base qw(Observation);
764    
765  Instantiate a new object.  sub new {
766    
767        my ($class,$dataset) = @_;
768        my $self = $class->SUPER::new($dataset);
769        $self->{acc} = $dataset->{'acc'};
770        $self->{evalue} = $dataset->{'evalue'};
771        $self->{start} = $dataset->{'start'};
772        $self->{stop} = $dataset->{'stop'};
773        bless($self,$class);
774        return $self;
775    }
776    
777    =head3 display()
778    
779    displays data stored in best_PDB attribute and in Ontology server for given PDB id
780    
781  =cut  =cut
782    
783    sub display{
784        my ($self,$gd) = @_;
785    
786        my $fid = $self->fig_id;
787        my $dbmaster = DBMaster->new(-database =>'Ontology');
788    
789        my $acc = $self->acc;
790    
791        print STDERR "acc:$acc\n";
792        my ($pdb_description,$pdb_source,$pdb_ligand);
793        my $pdb_objs = $dbmaster->pdb->get_objects( { 'id' => $acc } );
794        if(!scalar(@$pdb_objs)){
795            $pdb_description = "not available";
796            $pdb_source = "not available";
797            $pdb_ligand = "not available";
798        }
799        else{
800            my $pdb_obj = $pdb_objs->[0];
801            $pdb_description = $pdb_obj->description;
802            $pdb_source = $pdb_obj->source;
803            $pdb_ligand = $pdb_obj->ligand;
804        }
805    
806        my $lines = [];
807        my $line_data = [];
808        my $line_config = { 'title' => "PDB hit for $fid",
809                            'short_title' => "best PDB",
810                            'basepair_offset' => '1' };
811    
812        my $fig = new FIG;
813        my $seq = $fig->get_translation($fid);
814        my $fid_stop = length($seq);
815    
816        my $fid_element_hash = {
817            "title" => $fid,
818            "start" => '1',
819            "end" =>  $fid_stop,
820            "color"=> '1',
821            "zlayer" => '1'
822            };
823    
824        push(@$line_data,$fid_element_hash);
825    
826        my $links_list = [];
827        my $descriptions = [];
828    
829        my $name;
830        $name = {"title" => 'id',
831                 "value" => $acc};
832        push(@$descriptions,$name);
833    
834        my $description;
835        $description = {"title" => 'pdb description',
836                        "value" => $pdb_description};
837        push(@$descriptions,$description);
838    
839        my $score;
840        $score = {"title" => "score",
841                  "value" => $self->evalue};
842        push(@$descriptions,$score);
843    
844        my $start_stop;
845        my $start_stop_value = $self->start."_".$self->stop;
846        $start_stop = {"title" => "start-stop",
847                       "value" => $start_stop_value};
848        push(@$descriptions,$start_stop);
849    
850        my $source;
851        $source = {"title" => "source",
852                  "value" => $pdb_source};
853        push(@$descriptions,$source);
854    
855        my $ligand;
856        $ligand = {"title" => "pdb ligand",
857                   "value" => $pdb_ligand};
858        push(@$descriptions,$ligand);
859    
860        my $link;
861        my $link_url ="http://www.rcsb.org/pdb/explore/explore.do?structureId=".$acc;
862    
863        $link = {"link_title" => $acc,
864                 "link" => $link_url};
865        push(@$links_list,$link);
866    
867        my $pdb_element_hash = {
868            "title" => "PDB homology",
869            "start" => $self->start,
870            "end" =>  $self->stop,
871            "color"=> '6',
872            "zlayer" => '3',
873            "links_list" => $links_list,
874            "description" => $descriptions};
875    
876        push(@$line_data,$pdb_element_hash);
877        $gd->add_line($line_data, $line_config);
878    
879        return $gd;
880    }
881    
882    1;
883    
884    ############################################################
885    ############################################################
886    package Observation::Identical;
887    
888    use base qw(Observation);
889    
890  sub new {  sub new {
891    
892        my ($class,$dataset) = @_;
893        my $self = $class->SUPER::new($dataset);
894        $self->{rows} = $dataset->{'rows'};
895    
896        bless($self,$class);
897        return $self;
898    }
899    
900    =head3 display_table()
901    
902    If available use the function specified here to display the "raw" observation.
903    This code will display a table for the identical protein
904    
905    
906    B<Please note> that URL linked to in display_method() is an external component and needs to added to the code for every class of evi
907    dence.
908    
909    =cut
910    
911    
912    sub display_table{
913    my ($self) = @_;    my ($self) = @_;
914    
915    $self = { acc => '',      my $fig = new FIG;
916              description => '',      my $fid = $self->fig_id;
917              class => '',      my $rows = $self->rows;
918              type => '',      my $cgi = new CGI;
919              start => '',      my $all_domains = [];
920              stop => '',      my $count_identical = 0;
921              evalue => '',      my $content;
922              score => '',      foreach my $row (@$rows) {
923              display_method => '',          my $id = $row->[0];
924              feature_id => '',          my $who = $row->[1];
925              rank => '',          my $assignment = $row->[2];
926              supports_annotation => ''          my $organism = $fig->org_of($fid);
927            };          my $single_domain = [];
928            push(@$single_domain,$who);
929            push(@$single_domain,&HTML::set_prot_links($cgi,$id));
930            push(@$single_domain,$organism);
931            push(@$single_domain,$assignment);
932            push(@$all_domains,$single_domain);
933            $count_identical++;
934        }
935    
936    bless($self, 'Observation');      if ($count_identical >0){
937            $content = $all_domains;
938        }
939        else{
940            $content = "<p>This PEG does not have any essentially identical proteins</p>";
941        }
942        return ($content);
943    }
944    
945    1;
946    
947    #########################################
948    #########################################
949    package Observation::FC;
950    1;
951    
952    use base qw(Observation);
953    
954    sub new {
955    
956        my ($class,$dataset) = @_;
957        my $self = $class->SUPER::new($dataset);
958        $self->{rows} = $dataset->{'rows'};
959    
960        bless($self,$class);
961    return $self;    return $self;
962  }  }
963    
964  =head3 feature_id (internal)  =head3 display_table()
965    
966    If available use the function specified here to display the "raw" observation.
967    This code will display a table for the identical protein
968    
969  Returns the ID  of the feature these Observations belong to.  
970    B<Please note> that URL linked to in display_method() is an external component and needs to added to the code for every class of evi
971    dence.
972    
973  =cut  =cut
974    
975  sub feature_id {  sub display_table {
976    
977        my ($self,$dataset) = @_;
978        my $fid = $self->fig_id;
979        my $rows = $self->rows;
980        my $cgi = new CGI;
981        my $functional_data = [];
982        my $count = 0;
983        my $content;
984    
985        foreach my $row (@$rows) {
986            my $single_domain = [];
987            $count++;
988    
989            # construct the score link
990            my $score = $row->[0];
991            my $toid = $row->[1];
992            my $link = $cgi->url(-relative => 1) . "?user=master&request=show_coupling_evidence&prot=$fid&to=$toid&SPROUT=";
993            my $sc_link = "<a href=$link>$score</a>";
994    
995            push(@$single_domain,$sc_link);
996            push(@$single_domain,$row->[1]);
997            push(@$single_domain,$row->[2]);
998            push(@$functional_data,$single_domain);
999        }
1000    
1001        if ($count >0){
1002            $content = $functional_data;
1003        }
1004        else
1005        {
1006            $content = "<p>This PEG does not have any functional coupling</p>";
1007        }
1008        return ($content);
1009    }
1010    
1011    
1012    #########################################
1013    #########################################
1014    package Observation::Domain;
1015    
1016    use base qw(Observation);
1017    
1018    sub new {
1019    
1020        my ($class,$dataset) = @_;
1021        my $self = $class->SUPER::new($dataset);
1022        $self->{evalue} = $dataset->{'evalue'};
1023        $self->{acc} = $dataset->{'acc'};
1024        $self->{start} = $dataset->{'start'};
1025        $self->{stop} = $dataset->{'stop'};
1026    
1027        bless($self,$class);
1028        return $self;
1029    }
1030    
1031    sub display {
1032        my ($thing,$gd) = @_;
1033        my $lines = [];
1034        my $line_config = { 'title' => $thing->acc,
1035                            'short_title' => $thing->type,
1036                            'basepair_offset' => '1' };
1037        my $color = "4";
1038    
1039        my $line_data = [];
1040        my $links_list = [];
1041        my $descriptions = [];
1042    
1043        my $db_and_id = $thing->acc;
1044        my ($db,$id) = split("::",$db_and_id);
1045    
1046        my $dbmaster = DBMaster->new(-database =>'Ontology');
1047    
1048        my ($name_title,$name_value,$description_title,$description_value);
1049        if($db eq "CDD"){
1050            my $cdd_objs = $dbmaster->cdd->get_objects( { 'id' => $id } );
1051            if(!scalar(@$cdd_objs)){
1052                $name_title = "name";
1053                $name_value = "not available";
1054                $description_title = "description";
1055                $description_value = "not available";
1056            }
1057            else{
1058                my $cdd_obj = $cdd_objs->[0];
1059                $name_title = "name";
1060                $name_value = $cdd_obj->term;
1061                $description_title = "description";
1062                $description_value = $cdd_obj->description;
1063            }
1064        }
1065    
1066        my $name;
1067        $name = {"title" => $name_title,
1068                 "value" => $name_value};
1069        push(@$descriptions,$name);
1070    
1071        my $description;
1072        $description = {"title" => $description_title,
1073                                 "value" => $description_value};
1074        push(@$descriptions,$description);
1075    
1076        my $score;
1077        $score = {"title" => "score",
1078                  "value" => $thing->evalue};
1079        push(@$descriptions,$score);
1080    
1081        my $link_id;
1082        if ($thing->acc =~/\w+::(\d+)/){
1083            $link_id = $1;
1084        }
1085    
1086        my $link;
1087        my $link_url;
1088        if ($thing->class eq "CDD"){$link_url = "http://0-www.ncbi.nlm.nih.gov.library.vu.edu.au:80/Structure/cdd/cddsrv.cgi?uid=$link_id"}
1089        elsif($thing->class eq "PFAM"){$link_url = "http://www.sanger.ac.uk/cgi-bin/Pfam/getacc?$link_id"}
1090        else{$link_url = "NO_URL"}
1091    
1092        $link = {"link_title" => $thing->acc,
1093                 "link" => $link_url};
1094        push(@$links_list,$link);
1095    
1096        my $element_hash = {
1097            "title" => $thing->type,
1098            "start" => $thing->start,
1099            "end" =>  $thing->stop,
1100            "color"=> $color,
1101            "zlayer" => '2',
1102            "links_list" => $links_list,
1103            "description" => $descriptions};
1104    
1105        push(@$line_data,$element_hash);
1106        $gd->add_line($line_data, $line_config);
1107    
1108        return $gd;
1109    
1110    }
1111    
1112    #########################################
1113    #########################################
1114    package Observation::Location;
1115    
1116    use base qw(Observation);
1117    
1118    sub new {
1119    
1120        my ($class,$dataset) = @_;
1121        my $self = $class->SUPER::new($dataset);
1122        $self->{cleavage_prob} = $dataset->{'cleavage_prob'};
1123        $self->{cleavage_loc} = $dataset->{'cleavage_loc'};
1124        $self->{signal_peptide_score} = $dataset->{'signal_peptide_score'};
1125        $self->{cello_location} = $dataset->{'cello_location'};
1126        $self->{cello_score} = $dataset->{'cello_score'};
1127        $self->{tmpred_score} = $dataset->{'tmpred_score'};
1128        $self->{tmpred_locations} = $dataset->{'tmpred_locations'};
1129    
1130        bless($self,$class);
1131        return $self;
1132    }
1133    
1134    sub display {
1135        my ($thing,$gd) = @_;
1136    
1137        my $fid = $thing->fig_id;
1138        my $fig= new FIG;
1139        my $length = length($fig->get_translation($fid));
1140    
1141        my $cleavage_prob;
1142        if($thing->cleavage_prob){$cleavage_prob = $thing->cleavage_prob;}
1143        my ($cleavage_loc_begin,$cleavage_loc_end) = split("-",$thing->cleavage_loc);
1144        my $signal_peptide_score = $thing->signal_peptide_score;
1145        my $cello_location = $thing->cello_location;
1146        my $cello_score = $thing->cello_score;
1147        my $tmpred_score = $thing->tmpred_score;
1148        my @tmpred_locations = split(",",$thing->tmpred_locations);
1149    
1150        my $lines = [];
1151        my $line_config = { 'title' => 'Localization Evidence',
1152                            'short_title' => 'Local',
1153                            'basepair_offset' => '1' };
1154    
1155        #color is
1156        my $color = "5";
1157    
1158        my $line_data = [];
1159    
1160        if($cello_location){
1161            my $cello_descriptions = [];
1162            my $description_cello_location = {"title" => 'Best Cello Location',
1163                                              "value" => $cello_location};
1164    
1165            push(@$cello_descriptions,$description_cello_location);
1166    
1167            my $description_cello_score = {"title" => 'Cello Score',
1168                                           "value" => $cello_score};
1169    
1170            push(@$cello_descriptions,$description_cello_score);
1171    
1172            my $element_hash = {
1173                "title" => "CELLO",
1174                "start" => "1",
1175                "end" =>  $length + 1,
1176                "color"=> $color,
1177                "type" => 'box',
1178                "zlayer" => '2',
1179                "description" => $cello_descriptions};
1180    
1181            push(@$line_data,$element_hash);
1182        }
1183    
1184        my $color = "6";
1185        if($tmpred_score){
1186            foreach my $tmpred (@tmpred_locations){
1187                my $descriptions = [];
1188                my ($begin,$end) =split("-",$tmpred);
1189                my $description_tmpred_score = {"title" => 'TMPRED score',
1190                                 "value" => $tmpred_score};
1191    
1192                push(@$descriptions,$description_tmpred_score);
1193    
1194                my $element_hash = {
1195                "title" => "transmembrane location",
1196                "start" => $begin + 1,
1197                "end" =>  $end + 1,
1198                "color"=> $color,
1199                "zlayer" => '5',
1200                "type" => 'smallbox',
1201                "description" => $descriptions};
1202    
1203                push(@$line_data,$element_hash);
1204            }
1205        }
1206    
1207        my $color = "1";
1208        if($signal_peptide_score){
1209            my $descriptions = [];
1210            my $description_signal_peptide_score = {"title" => 'signal peptide score',
1211                                                    "value" => $signal_peptide_score};
1212    
1213            push(@$descriptions,$description_signal_peptide_score);
1214    
1215            my $description_cleavage_prob = {"title" => 'cleavage site probability',
1216                                             "value" => $cleavage_prob};
1217    
1218            push(@$descriptions,$description_cleavage_prob);
1219    
1220            my $element_hash = {
1221                "title" => "SignalP",
1222                "start" => $cleavage_loc_begin - 2,
1223                "end" =>  $cleavage_loc_end + 3,
1224                "type" => 'bigbox',
1225                "color"=> $color,
1226                "zlayer" => '10',
1227                "description" => $descriptions};
1228    
1229            push(@$line_data,$element_hash);
1230        }
1231    
1232        $gd->add_line($line_data, $line_config);
1233    
1234        return ($gd);
1235    
1236    }
1237    
1238    sub cleavage_loc {
1239    my ($self) = @_;    my ($self) = @_;
1240    
1241    return $self->{feature_id};    return $self->{cleavage_loc};
1242  }  }
1243    
1244    sub cleavage_prob {
1245      my ($self) = @_;
1246    
1247      return $self->{cleavage_prob};
1248    }
1249    
1250    sub signal_peptide_score {
1251      my ($self) = @_;
1252    
1253      return $self->{signal_peptide_score};
1254    }
1255    
1256    sub tmpred_score {
1257      my ($self) = @_;
1258    
1259      return $self->{tmpred_score};
1260    }
1261    
1262    sub tmpred_locations {
1263      my ($self) = @_;
1264    
1265      return $self->{tmpred_locations};
1266    }
1267    
1268    sub cello_location {
1269      my ($self) = @_;
1270    
1271      return $self->{cello_location};
1272    }
1273    
1274    sub cello_score {
1275      my ($self) = @_;
1276    
1277      return $self->{cello_score};
1278    }
1279    
1280    
1281    #########################################
1282    #########################################
1283    package Observation::Sims;
1284    
1285    use base qw(Observation);
1286    
1287    sub new {
1288    
1289        my ($class,$dataset) = @_;
1290        my $self = $class->SUPER::new($dataset);
1291        $self->{identity} = $dataset->{'identity'};
1292        $self->{acc} = $dataset->{'acc'};
1293        $self->{evalue} = $dataset->{'evalue'};
1294        $self->{qstart} = $dataset->{'qstart'};
1295        $self->{qstop} = $dataset->{'qstop'};
1296        $self->{hstart} = $dataset->{'hstart'};
1297        $self->{hstop} = $dataset->{'hstop'};
1298        $self->{database} = $dataset->{'database'};
1299        $self->{organism} = $dataset->{'organism'};
1300        $self->{function} = $dataset->{'function'};
1301        $self->{qlength} = $dataset->{'qlength'};
1302        $self->{hlength} = $dataset->{'hlength'};
1303    
1304        bless($self,$class);
1305        return $self;
1306    }
1307    
1308    =head3 display_table()
1309    
1310    If available use the function specified here to display the "raw" observation.
1311    This code will display a table for the similarities protein
1312    
1313    B<Please note> that URL linked to in display_method() is an external component and needs to added to the code for every class of evidence.
1314    
1315    =cut
1316    
1317    sub display_table {
1318        my ($self,$dataset) = @_;
1319    
1320        my $data = [];
1321        my $count = 0;
1322        my $content;
1323        my $fig = new FIG;
1324        my $cgi = new CGI;
1325        foreach my $thing (@$dataset) {
1326            my $single_domain = [];
1327            next if ($thing->class ne "SIM");
1328            $count++;
1329    
1330            my $id = $thing->acc;
1331    
1332            # add the subsystem information
1333            my @in_sub  = $fig->peg_to_subsystems($id);
1334            my $in_sub;
1335    
1336            if (@in_sub > 0) {
1337                $in_sub = @in_sub;
1338    
1339                # RAE: add a javascript popup with all the subsystems
1340                my $ss_list=join "<br>", map { my $g = $_; $g =~ s/\_/ /g; $_ = $g } sort {$a cmp $b} @in_sub;
1341                $in_sub = $cgi->a( {id=>"subsystems", onMouseover=>"javascript:if(!this.tooltip) this.tooltip=new Popup_Tooltip(this, 'Subsystems', '$ss_list', ''); this.tooltip.addHandler(); return false;"}, $in_sub);
1342            } else {
1343                $in_sub = "&nbsp;";
1344            }
1345    
1346            # add evidence code with tool tip
1347            my $ev_codes=" &nbsp; ";
1348            my @ev_codes = "";
1349            if ($id =~ /^fig\|\d+\.\d+\.peg\.\d+$/) {
1350                my @codes = grep { $_->[1] =~ /^evidence_code/i } $fig->get_attributes($id);
1351                @ev_codes = ();
1352                foreach my $code (@codes) {
1353                    my $pretty_code = $code->[2];
1354                    if ($pretty_code =~ /;/) {
1355                        my ($cd, $ss) = split(";", $code->[2]);
1356                        $ss =~ s/_/ /g;
1357                        $pretty_code = $cd;# . " in " . $ss;
1358                    }
1359                    push(@ev_codes, $pretty_code);
1360                }
1361            }
1362    
1363            if (scalar(@ev_codes) && $ev_codes[0]) {
1364                my $ev_code_help=join("<br />", map {&HTML::evidence_codes_explain($_)} @ev_codes);
1365                $ev_codes = $cgi->a(
1366                                    {
1367                                        id=>"evidence_codes", onMouseover=>"javascript:if(!this.tooltip) this.tooltip=new Popup_Tooltip(this, 'Evidence Codes', '$ev_code_help', ''); this.tooltip.addHandler(); return false;"}, join("<br />", @ev_codes));
1368            }
1369    
1370            # add the aliases
1371            my $aliases = undef;
1372            $aliases = &html_enc( join( ", ", $fig->feature_aliases($id) ) );
1373            $aliases = &HTML::set_prot_links( $cgi, $aliases );
1374            $aliases ||= "&nbsp;";
1375    
1376            my $iden    = $thing->identity;
1377            my $ln1     = $thing->qlength;
1378            my $ln2     = $thing->hlength;
1379            my $b1      = $thing->qstart;
1380            my $e1      = $thing->qstop;
1381            my $b2      = $thing->hstart;
1382            my $e2      = $thing->hstop;
1383            my $d1      = abs($e1 - $b1) + 1;
1384            my $d2      = abs($e2 - $b2) + 1;
1385            my $reg1    = "$b1-$e1 (<b>$d1/$ln1</b>)";
1386            my $reg2    = "$b2-$e2 (<b>$d2/$ln2</b>)";
1387    
1388    
1389            push(@$single_domain,$thing->database);
1390            push(@$single_domain,&HTML::set_prot_links($cgi,$thing->acc));
1391            push(@$single_domain,$thing->evalue);
1392            push(@$single_domain,"$iden\%");
1393            push(@$single_domain,$reg1);
1394            push(@$single_domain,$reg2);
1395            push(@$single_domain,$in_sub);
1396            push(@$single_domain,$ev_codes);
1397            push(@$single_domain,$thing->organism);
1398            push(@$single_domain,$thing->function);
1399            push(@$single_domain,$aliases);
1400            push(@$data,$single_domain);
1401        }
1402    
1403        if ($count >0){
1404            $content = $data;
1405        }
1406        else
1407        {
1408            $content = "<p>This PEG does not have any similarities</p>";
1409        }
1410        return ($content);
1411    }
1412    
1413    sub html_enc { $_ = $_[0]; s/\&/&amp;/g; s/\>/&gt;/g; s/\</&lt;/g; $_ }
1414    
1415    
1416    
1417    ############################
1418    package Observation::Cluster;
1419    
1420    use base qw(Observation);
1421    
1422    sub new {
1423    
1424        my ($class,$dataset) = @_;
1425        my $self = $class->SUPER::new($dataset);
1426        $self->{context} = $dataset->{'context'};
1427        bless($self,$class);
1428        return $self;
1429    }
1430    
1431    sub display {
1432        my ($self,$gd) = @_;
1433    
1434        my $fid = $self->fig_id;
1435        my $compare_or_coupling = $self->context;
1436        my $gd_window_size = $gd->window_size;
1437        my $fig = new FIG;
1438        my $all_regions = [];
1439    
1440        #get the organism genome
1441        my $target_genome = $fig->genome_of($fid);
1442    
1443        # get location of the gene
1444        my $data = $fig->feature_location($fid);
1445        my ($contig, $beg, $end);
1446        my %reverse_flag;
1447    
1448        if ($data =~ /(.*)_(\d+)_(\d+)$/){
1449            $contig = $1;
1450            $beg = $2;
1451            $end = $3;
1452        }
1453    
1454        my $offset;
1455        my ($region_start, $region_end);
1456        if ($beg < $end)
1457        {
1458            $region_start = $beg - 4000;
1459            $region_end = $end+4000;
1460            $offset = ($2+(($3-$2)/2))-($gd_window_size/2);
1461        }
1462        else
1463        {
1464            $region_start = $end-4000;
1465            $region_end = $beg+4000;
1466            $offset = ($3+(($2-$3)/2))-($gd_window_size/2);
1467            $reverse_flag{$target_genome} = 1;
1468        }
1469    
1470        # call genes in region
1471        my ($target_gene_features, $reg_beg, $reg_end) = $fig->genes_in_region($target_genome, $contig, $region_start, $region_end);
1472        push(@$all_regions,$target_gene_features);
1473        my (@start_array_region);
1474        push (@start_array_region, $offset);
1475    
1476        my %all_genes;
1477        my %all_genomes;
1478        foreach my $feature (@$target_gene_features){ $all_genes{$feature} = 1;}
1479    
1480        if ($compare_or_coupling eq "diverse")
1481        {
1482            my @coup = grep { $_->[1]} $fig->coupling_and_evidence($fid,5000,1e-10,4,1);
1483    
1484            my $coup_count = 0;
1485    
1486            foreach my $pair (@{$coup[0]->[2]}) {
1487                #   last if ($coup_count > 10);
1488                my ($peg1,$peg2) = @$pair;
1489    
1490                my ($pair_contig,$pair_beg,$pair_end,$pair_region_start,$pair_region_stop,$pair_genome);
1491                $pair_genome = $fig->genome_of($peg1);
1492    
1493                my $location = $fig->feature_location($peg1);
1494                if($location =~/(.*)_(\d+)_(\d+)$/){
1495                    $pair_contig = $1;
1496                    $pair_beg = $2;
1497                    $pair_end = $3;
1498                    if ($pair_beg < $pair_end)
1499                    {
1500                        $pair_region_start = $pair_beg - 4000;
1501                        $pair_region_stop = $pair_end+4000;
1502                        $offset = ($2+(($3-$2)/2))-($gd_window_size/2);
1503                    }
1504                    else
1505                    {
1506                        $pair_region_start = $pair_end-4000;
1507                        $pair_region_stop = $pair_beg+4000;
1508                        $offset = ($3+(($2-$3)/2))-($gd_window_size/2);
1509                        $reverse_flag{$pair_genome} = 1;
1510                    }
1511    
1512                    push (@start_array_region, $offset);
1513    
1514                    $all_genomes{$pair_genome} = 1;
1515                    my ($pair_features) = $fig->genes_in_region($pair_genome, $pair_contig, $pair_region_start, $pair_region_stop);
1516                    push(@$all_regions,$pair_features);
1517                    foreach my $pair_feature (@$pair_features){ $all_genes{$pair_feature} = 1;}
1518                }
1519                $coup_count++;
1520            }
1521        }
1522    
1523        elsif ($compare_or_coupling eq "close")
1524        {
1525            # make a hash of genomes that are phylogenetically close
1526            #my $close_threshold = ".26";
1527            #my @genomes = $fig->genomes('complete');
1528            #my %close_genomes = ();
1529            #foreach my $compared_genome (@genomes)
1530            #{
1531            #    my $dist = $fig->crude_estimate_of_distance($target_genome,$compared_genome);
1532            #    #$close_genomes{$compared_genome} = $dist;
1533            #    if ($dist <= $close_threshold)
1534            #    {
1535            #       $all_genomes{$compared_genome} = 1;
1536            #    }
1537            #}
1538            $all_genomes{"216592.1"} = 1;
1539            $all_genomes{"79967.1"} = 1;
1540            $all_genomes{"199310.1"} = 1;
1541            $all_genomes{"216593.1"} = 1;
1542            $all_genomes{"155864.1"} = 1;
1543            $all_genomes{"83334.1"} = 1;
1544            $all_genomes{"316407.3"} = 1;
1545    
1546            foreach my $comp_genome (keys %all_genomes){
1547                my $return = $fig->bbh_list($comp_genome,[$fid]);
1548                my $feature_list = $return->{$fid};
1549                foreach my $peg1 (@$feature_list){
1550                    my $location = $fig->feature_location($peg1);
1551                    my ($pair_contig,$pair_beg,$pair_end,$pair_region_start,$pair_region_stop,$pair_genome);
1552                    $pair_genome = $fig->genome_of($peg1);
1553    
1554                    if($location =~/(.*)_(\d+)_(\d+)$/){
1555                        $pair_contig = $1;
1556                        $pair_beg = $2;
1557                        $pair_end = $3;
1558                        if ($pair_beg < $pair_end)
1559                        {
1560                            $pair_region_start = $pair_beg - 4000;
1561                            $pair_region_stop = $pair_end + 4000;
1562                            $offset = ($2+(($3-$2)/2))-($gd_window_size/2);
1563                        }
1564                        else
1565                        {
1566                            $pair_region_start = $pair_end-4000;
1567                            $pair_region_stop = $pair_beg+4000;
1568                            $offset = ($3+(($2-$3)/2))-($gd_window_size/2);
1569                            $reverse_flag{$pair_genome} = 1;
1570                        }
1571    
1572                        push (@start_array_region, $offset);
1573                        $all_genomes{$pair_genome} = 1;
1574                        my ($pair_features) = $fig->genes_in_region($pair_genome, $pair_contig, $pair_region_start, $pair_region_stop);
1575                        push(@$all_regions,$pair_features);
1576                        foreach my $pair_feature (@$pair_features){ $all_genes{$pair_feature} = 1;}
1577                    }
1578                }
1579            }
1580        }
1581    
1582        # get the PCH to each of the genes
1583        my $pch_sets = [];
1584        my %pch_already;
1585        foreach my $gene_peg (keys %all_genes)
1586        {
1587            if ($pch_already{$gene_peg}){next;};
1588            my $gene_set = [$gene_peg];
1589            foreach my $pch_peg ($fig->in_pch_pin_with($gene_peg)) {
1590                $pch_peg =~ s/,.*$//;
1591                my $pch_genome = $fig->genome_of($pch_peg);
1592                if ( ($gene_peg ne $pch_peg) && ($all_genomes{$pch_genome})) {
1593                    push(@$gene_set,$pch_peg);
1594                    $pch_already{$pch_peg}=1;
1595                }
1596                $pch_already{$gene_peg}=1;
1597            }
1598            push(@$pch_sets,$gene_set);
1599        }
1600    
1601        #create a rank of the pch's
1602        my %pch_set_rank;
1603        my $order = 0;
1604        foreach my $set (@$pch_sets){
1605            my $count = scalar(@$set);
1606            $pch_set_rank{$order} = $count;
1607            $order++;
1608        }
1609    
1610        my %peg_rank;
1611        my $counter =  1;
1612        foreach my $pch_order (sort {$pch_set_rank{$b} <=> $pch_set_rank{$a}} keys %pch_set_rank){
1613            my $good_set = @$pch_sets[$pch_order];
1614            my $flag_set = 0;
1615            if (scalar (@$good_set) > 1)
1616            {
1617                foreach my $peg (@$good_set){
1618                    if ((!$peg_rank{$peg})){
1619                        $peg_rank{$peg} = $counter;
1620                        $flag_set = 1;
1621                    }
1622                }
1623                $counter++ if ($flag_set == 1);
1624            }
1625            else
1626            {
1627                foreach my $peg (@$good_set){
1628                    $peg_rank{$peg} = 100;
1629                }
1630            }
1631        }
1632    
1633    
1634    #    my $bbh_sets = [];
1635    #    my %already;
1636    #    foreach my $gene_key (keys(%all_genes)){
1637    #       if($already{$gene_key}){next;}
1638    #       my $gene_set = [$gene_key];
1639    #
1640    #       my $gene_key_genome = $fig->genome_of($gene_key);
1641    #
1642    #       foreach my $genome_key (keys(%all_genomes)){
1643    #           #next if ($gene_key_genome eq $genome_key);
1644    #           my $return = $fig->bbh_list($genome_key,[$gene_key]);
1645    #
1646    #           my $feature_list = $return->{$gene_key};
1647    #           foreach my $fl (@$feature_list){
1648    #               push(@$gene_set,$fl);
1649    #           }
1650    #       }
1651    #       $already{$gene_key} = 1;
1652    #       push(@$bbh_sets,$gene_set);
1653    #    }
1654    #
1655    #    my %bbh_set_rank;
1656    #    my $order = 0;
1657    #    foreach my $set (@$bbh_sets){
1658    #       my $count = scalar(@$set);
1659    #       $bbh_set_rank{$order} = $count;
1660    #       $order++;
1661    #    }
1662    #
1663    #    my %peg_rank;
1664    #    my $counter =  1;
1665    #    foreach my $bbh_order (sort {$bbh_set_rank{$b} <=> $bbh_set_rank{$a}} keys %bbh_set_rank){
1666    #       my $good_set = @$bbh_sets[$bbh_order];
1667    #       my $flag_set = 0;
1668    #       if (scalar (@$good_set) > 1)
1669    #       {
1670    #           foreach my $peg (@$good_set){
1671    #               if ((!$peg_rank{$peg})){
1672    #                   $peg_rank{$peg} = $counter;
1673    #                   $flag_set = 1;
1674    #               }
1675    #           }
1676    #           $counter++ if ($flag_set == 1);
1677    #       }
1678    #       else
1679    #       {
1680    #           foreach my $peg (@$good_set){
1681    #               $peg_rank{$peg} = 100;
1682    #           }
1683    #       }
1684    #    }
1685    
1686        foreach my $region (@$all_regions){
1687            my $sample_peg = @$region[0];
1688            my $region_genome = $fig->genome_of($sample_peg);
1689            my $region_gs = $fig->genus_species($region_genome);
1690            my $abbrev_name = $fig->abbrev($region_gs);
1691            my $line_config = { 'title' => $region_gs,
1692                                'short_title' => $abbrev_name,
1693                                'basepair_offset' => '0'
1694                                };
1695    
1696            my $offsetting = shift @start_array_region;
1697    
1698            my $line_data = [];
1699            foreach my $fid1 (@$region){
1700                my $element_hash;
1701                my $links_list = [];
1702                my $descriptions = [];
1703    
1704                my $color = $peg_rank{$fid1};
1705    
1706                # get subsystem information
1707                my $function = $fig->function_of($fid1);
1708                my $url_link = "http://seed-viewer.theseed.org/index.cgi?action=ShowAnnotation&prot=".$fid1;
1709    
1710                my $link;
1711                $link = {"link_title" => $fid1,
1712                         "link" => $url_link};
1713                push(@$links_list,$link);
1714    
1715                my @subsystems = $fig->peg_to_subsystems($fid1);
1716                foreach my $subsystem (@subsystems){
1717                    my $link;
1718                    $link = {"link" => "http://seed-viewer.theseed.org/index.cgi?action=ShowSubsystem&subsystem_name=$subsystem",
1719                             "link_title" => $subsystem};
1720                    push(@$links_list,$link);
1721                }
1722    
1723                my $description_function;
1724                $description_function = {"title" => "function",
1725                                         "value" => $function};
1726                push(@$descriptions,$description_function);
1727    
1728                my $description_ss;
1729                my $ss_string = join (",", @subsystems);
1730                $description_ss = {"title" => "subsystems",
1731                                   "value" => $ss_string};
1732                push(@$descriptions,$description_ss);
1733    
1734    
1735                my $fid_location = $fig->feature_location($fid1);
1736                if($fid_location =~/(.*)_(\d+)_(\d+)$/){
1737                    my($start,$stop);
1738                    $start = $2 - $offsetting;
1739                    $stop = $3 - $offsetting;
1740    
1741                    if (defined($reverse_flag{$region_genome})){
1742                        $start = $gd_window_size - $start;
1743                        $stop = $gd_window_size - $stop;
1744                    }
1745    
1746                    $element_hash = {
1747                        "title" => $fid1,
1748                        "start" => $start,
1749                        "end" =>  $stop,
1750                        "type"=> 'arrow',
1751                        "color"=> $color,
1752                        "zlayer" => "2",
1753                        "links_list" => $links_list,
1754                        "description" => $descriptions
1755                    };
1756                    push(@$line_data,$element_hash);
1757                }
1758            }
1759            $gd->add_line($line_data, $line_config);
1760        }
1761        return $gd;
1762    }
1763    
1764    

Legend:
Removed from v.1.3  
changed lines
  Added in v.1.24

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3