[Bio] / FigKernelScripts / clustered_hypotheticals.pl Repository:
ViewVC logotype

Annotation of /FigKernelScripts/clustered_hypotheticals.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.2 - (view) (download) (as text)

1 : overbeek 1.1 #!/usr/bin/perl -w
2 :    
3 :     use strict;
4 :     use SeedEnv;
5 :     use Data::Dumper;
6 :    
7 :     my $sapObject = SAPserver->new();
8 :    
9 :     my $genomeHash = $sapObject->all_genomes(-complete => 1);
10 : overbeek 1.2 my @completeG = keys(%$genomeHash);
11 :     my $isP = $sapObject->is_prokaryotic(-ids => \@completeG);
12 :     my @complete_proks = grep { $isP->{$_} } @completeG;
13 :     my %seen;
14 :    
15 :     foreach my $genome (sort { $a <=> $b } @complete_proks)
16 :     {
17 : overbeek 1.1 print STDERR "processing $genome\n";
18 :    
19 :     my $genomeName = $genomeHash->{$genome};
20 :     my $geneHash = $sapObject->feature_assignments(-genome => $genome,
21 :     -type => 'peg');
22 : overbeek 1.2 my @hypotheticalGenes = sort { &SeedUtils::by_fig_id($a,$b) }
23 :     grep { &SeedUtils::hypo($geneHash->{$_}) }
24 :     grep { ! $seen{$_} }
25 :     keys(%$geneHash);
26 :    
27 :     my $couplingHash = $sapObject->conserved_in_neighborhood(-ids => \@hypotheticalGenes);
28 :     my @all_coupled_ids = map { map { $_->[1] } @{$couplingHash->{$_}}} keys(%$couplingHash);
29 :     my $subHash = $sapObject->ids_to_subsystems(-ids => \@all_coupled_ids);
30 : overbeek 1.1
31 : overbeek 1.2 foreach my $gene (@hypotheticalGenes)
32 :     {
33 : overbeek 1.1 my $couplingList = $couplingHash->{$gene};
34 : overbeek 1.2 if (defined $couplingList)
35 :     {
36 :     my ($bestCoupler, $bestScore, $bestFunction) = (undef, 0, '');
37 :     my $best_pairset;
38 : overbeek 1.1
39 : overbeek 1.2 foreach my $coupling (@$couplingList)
40 :     {
41 :     my ($score, $coupler, $function, $pairset) = @$coupling;
42 : overbeek 1.1
43 : overbeek 1.2 if ($subHash->{$coupler} && $score > $bestScore)
44 :     {
45 :     $bestCoupler = $coupler;
46 :     $bestScore = $score;
47 :     $bestFunction = $function ? $function : '';
48 :     $best_pairset = $pairset;
49 : overbeek 1.1 }
50 :     }
51 : overbeek 1.2
52 :     if (defined $bestCoupler)
53 :     {
54 : overbeek 1.1 print join("\t", $gene, $geneHash->{$gene}, $genome, $genomeName,
55 :     $bestCoupler, $bestScore, $bestFunction) . "\n";
56 : overbeek 1.2 &set_seen_for_pairset($gene,$best_pairset,\%seen,$sapObject);
57 : overbeek 1.1 }
58 :     }
59 :     }
60 :     }
61 :    
62 : overbeek 1.2 sub set_seen_for_pairset {
63 :     my($peg,$pairset,$seen,$sapObject) = @_;
64 :    
65 :     my $pairsetsH = $sapObject->pairsets( -ids => [$pairset] );
66 :     if ($_ = $pairsetsH->{$pairset})
67 :     {
68 :     my($sc,$pairs) = @$_;
69 :     my $got = 0;
70 :     my @set = map { if ($_->[0] eq $peg) { $got = 1 }; $_->[0] } @$pairs;
71 :     if (! $got)
72 :     {
73 :     @set = map { $_->[1] } @$pairs;
74 :     }
75 :     foreach $_ (@set)
76 :     {
77 :     $seen->{$_} = 1;
78 :     }
79 :     }
80 :     }
81 :    
82 :    

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3