[Bio] / FigWebServices / homologs_in_clusters.cgi Repository:
ViewVC logotype

Annotation of /FigWebServices/homologs_in_clusters.cgi

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (view) (download)

1 : overbeek 1.1 #### start #####
2 :     use FIG;
3 :     my $fig = new FIG;
4 :    
5 :     use HTML;
6 :     use strict;
7 :     use CGI;
8 :     my $cgi = new CGI;
9 :    
10 :     if (0)
11 :     {
12 :     my $VAR1;
13 :     eval(join("",`cat /tmp/homologs_in_clusters_parms`));
14 :     $cgi = $VAR1;
15 :     # print STDERR &Dumper($cgi);
16 :     }
17 :    
18 :     if (0)
19 :     {
20 :     print $cgi->header;
21 :     my @params = $cgi->param;
22 :     print "<pre>\n";
23 :     foreach $_ (@params)
24 :     {
25 :     print "$_\t:",join(",",$cgi->param($_)),":\n";
26 :     }
27 :    
28 :     if (0)
29 :     {
30 :     if (open(TMP,">/tmp/homologs_in_clusters_parms"))
31 :     {
32 :     print TMP &Dumper($cgi);
33 :     close(TMP);
34 :     }
35 :     }
36 :     exit;
37 :     }
38 :    
39 :     my $html = [];
40 :     unshift @$html, "<TITLE>The SEED: Homologs in Clusters Page</TITLE>\n";
41 :    
42 :     my $prot = $cgi->param('prot');
43 :     if (! $prot)
44 :     {
45 :     push(@$html,"<h1>Sorry, you need to specify a protein</h1>\n");
46 :     &HTML::show_page($cgi,$html);
47 :     exit;
48 :     }
49 :    
50 :    
51 :     if ($prot !~ /^fig\|/)
52 :     {
53 :     my @poss = $fig->by_alias($prot);
54 :     if (@poss > 0)
55 :     {
56 :     $prot = $poss[0];
57 :     }
58 :     else
59 :     {
60 :     push(@$html,"<h1>Sorry, $prot appears not to have a FIG id at this point</h1>\n");
61 :     &HTML::show_page($cgi,$html);
62 :     exit;
63 :     }
64 :     }
65 :    
66 :     &compute_desired_homologs($fig,$cgi,$html,$prot);
67 :     &HTML::show_page($cgi,$html);
68 :     exit;
69 :    
70 :     sub compute_desired_homologs {
71 :     my($fig,$cgi,$html,$peg) = @_;
72 :    
73 :     my @pinned = &relevant_homologs($fig,$cgi,$peg);
74 :     # print STDERR &Dumper(\@pinned);
75 :    
76 : overbeek 1.3 # my @clusters = sort { (@$b <=> @$a) } &sets_of_homologs($fig,$cgi,$peg,\@pinned);
77 : overbeek 1.1 # print STDERR &Dumper(\@clusters);
78 :    
79 : overbeek 1.3 # my @homologs = &extract_homologs($peg,\@pinned,\@clusters);
80 : overbeek 1.1 # print STDERR &Dumper(\@homologs);
81 :    
82 :     my $sc;
83 :     my @tab = map { ($peg,$sc) = @$_; [$sc,
84 :     &HTML::fid_link($cgi,$peg),
85 : redwards 1.2 $fig->genus_species($fig->genome_of($peg)),
86 : overbeek 1.1 scalar $fig->function_of($peg,$cgi->param('user')),
87 :     &HTML::set_prot_links($cgi,join( ', ', $fig->feature_aliases($peg) ))
88 : overbeek 1.3 ] } @pinned;
89 :     push(@$html,&HTML::make_table(["Score","PEG","Genome", "Function","Aliases"],\@tab,"PEGs that Might Be in Clusters"));
90 : overbeek 1.1 }
91 :    
92 :     sub relevant_homologs {
93 :     my($fig,$cgi,$peg) = @_;
94 :     my($maxN,$maxP,$genome1,$sim,$id2,$genome2,%seen);
95 :    
96 :     $maxN = $cgi->param('maxN');
97 :     $maxN = $maxN ? $maxN : 50;
98 :    
99 :     $maxP = $cgi->param('maxP');
100 :     $maxP = $maxP ? $maxP : 1.0e-10;
101 :    
102 :     my @sims = $fig->sims( $peg, $maxN, $maxP, "fig");
103 :    
104 :     my @homologs = ();
105 :     $seen{&FIG::genome_of($peg)} = 1;
106 :     foreach $sim (@sims)
107 :     {
108 :     $id2 = $sim->id2;
109 :     $genome2 = &FIG::genome_of($id2);
110 : overbeek 1.3 my @coup;
111 :     if ((! $seen{$genome2}) && (@coup = $fig->coupled_to($id2)) && (@coup > 0))
112 : overbeek 1.1 {
113 :     $seen{$genome2} = 1;
114 : overbeek 1.3 push(@homologs,[$id2,scalar @coup]);
115 : overbeek 1.1 }
116 :     }
117 : overbeek 1.3 return sort { $b->[1] <=> $a->[1] } @homologs;
118 : overbeek 1.1 }
119 :    
120 :     sub sets_of_homologs {
121 :     my($fig,$cgi,$given_peg,$pinned) = @_;
122 :     my($peg,$mid,$min,$max,$feat,$fid);
123 :    
124 :     my $bound = $cgi->param('bound');
125 :     $bound = $bound ? $bound : 4000;
126 :    
127 :     my @pegs = ();
128 :     foreach $peg (($given_peg,@$pinned))
129 :     {
130 :     my $loc = $fig->feature_location($peg);
131 :     if ($loc)
132 :     {
133 :     my($contig,$beg,$end) = &FIG::boundaries_of($loc);
134 :     if ($contig && $beg && $end)
135 :     {
136 :     $mid = int(($beg + $end) / 2);
137 :     $min = $mid - $bound;
138 :     $max = $mid + $bound;
139 :    
140 :     ($feat,undef,undef) = $fig->genes_in_region(&FIG::genome_of($peg),$contig,$min,$max);
141 :     foreach $fid (@$feat)
142 :     {
143 :     if (&FIG::ftype($fid) eq "peg")
144 :     {
145 :     push(@pegs,$fid);
146 :     }
147 :     }
148 :     }
149 :     }
150 :     }
151 :    
152 :     my %represents;
153 :     foreach $peg (@pegs)
154 :     {
155 :     my $tmp = $fig->maps_to_id($peg);
156 :     push(@{$represents{$tmp}},$peg);
157 :     }
158 :     my($sim,%conn,$x,$y,$i,$j);
159 :     foreach $y (keys(%represents))
160 :     {
161 :     $x = $represents{$y};
162 :     for ($i=0; ($i < @$x); $i++)
163 :     {
164 :     for ($j=$i+1; ($j < @$x); $j++)
165 :     {
166 :     push(@{$conn{$x->[$i]}},$x->[$j]);
167 :     push(@{$conn{$x->[$j]}},$x->[$i]);
168 :     }
169 :     }
170 :     }
171 :    
172 :     my $maxN = $cgi->param('maxN');
173 :     $maxN = $maxN ? $maxN : 50;
174 :    
175 :     my $maxP = $cgi->param('maxP');
176 :     $maxP = $maxP ? $maxP : 1.0e-10;
177 :    
178 :     foreach $peg (@pegs)
179 :     {
180 :     foreach $sim ($fig->sims($peg,$maxN,$maxP,"raw"))
181 :     {
182 :     if (defined($x = $represents{$sim->id2}))
183 :     {
184 :     foreach $y (@$x)
185 :     {
186 :     push(@{$conn{$peg}},$y);
187 :     }
188 :     }
189 :     }
190 :     }
191 :    
192 :     my(%seen,$k,$cluster);
193 :     my @clusters = ();
194 :     for ($i=0; ($i < @pegs); $i++)
195 :     {
196 :     $peg = $pegs[$i];
197 :     if (! $seen{$peg})
198 :     {
199 :     $cluster = [$peg];
200 :     $seen{$peg} = 1;
201 :     for ($j=0; ($j < @$cluster); $j++)
202 :     {
203 :     $x = $conn{$cluster->[$j]};
204 :     foreach $k (@$x)
205 :     {
206 :     if (! $seen{$k})
207 :     {
208 :     push(@$cluster,$k);
209 :     $seen{$k} = 1;
210 :     }
211 :     }
212 :     }
213 :    
214 :     if (@$cluster > 1)
215 :     {
216 :     push(@clusters,$cluster);
217 :     }
218 :     }
219 :     }
220 :     return @clusters;
221 :     }
222 :    
223 :     sub extract_homologs {
224 :     my($given_peg,$pinned,$clusters) = @_;
225 :     my(%main,$cluster,$peg,%counts,@with_counts);
226 :    
227 :     %main = map { $_ => 1 } ($given_peg,@$pinned);
228 :     foreach $cluster (@$clusters)
229 :     {
230 :     foreach $peg (@$cluster)
231 :     {
232 :     if (! $main{$peg})
233 :     {
234 :     $counts{&FIG::genome_of($peg)} += @$cluster - 1;
235 :     }
236 :     }
237 :     }
238 :    
239 :     foreach $peg (($given_peg,@$pinned))
240 :     {
241 :     push(@with_counts,[$peg,$counts{&FIG::genome_of($peg)}]);
242 :     }
243 :    
244 :     return grep { $_->[1] > 2} sort { $b->[1] <=> $a->[1] } @with_counts;
245 :     }

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3