[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.14 - (view) (download)

1 : olson 1.9 #
2 :     # Copyright (c) 2003-2006 University of Chicago and Fellowship
3 :     # for Interpretations of Genomes. All Rights Reserved.
4 :     #
5 :     # This file is part of the SEED Toolkit.
6 : parrello 1.10 #
7 : olson 1.9 # The SEED Toolkit is free software. You can redistribute
8 :     # it and/or modify it under the terms of the SEED Toolkit
9 : parrello 1.10 # Public License.
10 : olson 1.9 #
11 :     # You should have received a copy of the SEED Toolkit Public License
12 :     # along with this program; if not write to the University of Chicago
13 :     # at info@ci.uchicago.edu or the Fellowship for Interpretation of
14 :     # Genomes at veronika@thefig.info or download a copy from
15 :     # http://www.theseed.org/LICENSE.TXT.
16 :     #
17 :    
18 : overbeek 1.1 #### start #####
19 : overbeek 1.4 use InterfaceRoutines;
20 :    
21 : overbeek 1.1
22 :     use HTML;
23 :     use strict;
24 :     use CGI;
25 : parrello 1.13 use FIG_CGI;
26 :     use FIG;
27 : overbeek 1.1
28 : overbeek 1.4 my $sproutAvail = eval {
29 :     require SproutFIG;
30 :     require PageBuilder;
31 :     };
32 :    
33 : parrello 1.13 my($fig_or_sprout, $cgi) = FIG_CGI::init();
34 : parrello 1.14 if (ref $fig_or_sprout eq 'SFXlate') {
35 :     my $prot = $cgi->param('prot');
36 :     print $cgi->redirect(-uri => "$FIG_Config::cgi_url/wiki/rest.cgi/NmpdrPlugin/SeedViewer?page=Evidence;feature=$prot",
37 :     -status => 301);
38 :     }
39 : parrello 1.13
40 : overbeek 1.4 my $html = [];
41 :    
42 : parrello 1.13 unshift @$html, "<TITLE>The Homologs in Clusters Page</TITLE>\n";
43 : overbeek 1.4
44 : overbeek 1.1 if (0)
45 :     {
46 :     my $VAR1;
47 :     eval(join("",`cat /tmp/homologs_in_clusters_parms`));
48 :     $cgi = $VAR1;
49 :     # print STDERR &Dumper($cgi);
50 :     }
51 :    
52 :     if (0)
53 :     {
54 :     print $cgi->header;
55 :     my @params = $cgi->param;
56 :     print "<pre>\n";
57 :     foreach $_ (@params)
58 :     {
59 :     print "$_\t:",join(",",$cgi->param($_)),":\n";
60 :     }
61 :    
62 :     if (0)
63 :     {
64 :     if (open(TMP,">/tmp/homologs_in_clusters_parms"))
65 :     {
66 :     print TMP &Dumper($cgi);
67 :     close(TMP);
68 :     }
69 :     }
70 :     exit;
71 :     }
72 :    
73 :     my $prot = $cgi->param('prot');
74 :     if (! $prot)
75 :     {
76 :     push(@$html,"<h1>Sorry, you need to specify a protein</h1>\n");
77 :     &HTML::show_page($cgi,$html);
78 :     exit;
79 :     }
80 :    
81 :    
82 :     if ($prot !~ /^fig\|/)
83 :     {
84 : overbeek 1.4 my @poss = $fig_or_sprout->by_alias($prot);
85 : overbeek 1.1 if (@poss > 0)
86 :     {
87 :     $prot = $poss[0];
88 :     }
89 :     else
90 :     {
91 :     push(@$html,"<h1>Sorry, $prot appears not to have a FIG id at this point</h1>\n");
92 :     &HTML::show_page($cgi,$html);
93 :     exit;
94 :     }
95 :     }
96 :    
97 : overbeek 1.4 &compute_desired_homologs($fig_or_sprout,$cgi,$html,$prot);
98 : olson 1.5
99 : parrello 1.13 if (ref $fig_or_sprout eq 'SFXlate')
100 : olson 1.5 {
101 : parrello 1.12 my $h = { homologs => $html,
102 :     title => "NMPDR Homologs in Clusters Page"};
103 : parrello 1.10
104 : olson 1.5 print "Content-Type: text/html\n";
105 :     print "\n";
106 : parrello 1.12 my $templ = "$FIG_Config::template_url/Homologs_tmpl.php";
107 :     print PageBuilder::Build("$templ", $h,"Html");
108 : olson 1.5 }
109 :     else
110 :     {
111 :     &HTML::show_page($cgi,$html);
112 :     }
113 : overbeek 1.1 exit;
114 :    
115 :     sub compute_desired_homologs {
116 : overbeek 1.4 my($fig_or_sprout,$cgi,$html,$peg) = @_;
117 : overbeek 1.1
118 : overbeek 1.4 my @pinned = &relevant_homologs($fig_or_sprout,$cgi,$peg);
119 : overbeek 1.1 # print STDERR &Dumper(\@pinned);
120 :    
121 : overbeek 1.4 # my @clusters = sort { (@$b <=> @$a) } &sets_of_homologs($fig_or_sprout,$cgi,$peg,\@pinned);
122 : overbeek 1.1 # print STDERR &Dumper(\@clusters);
123 :    
124 : overbeek 1.3 # my @homologs = &extract_homologs($peg,\@pinned,\@clusters);
125 : overbeek 1.1 # print STDERR &Dumper(\@homologs);
126 :    
127 :     my $sc;
128 : overbeek 1.7 my @tab = map { my($peg,$sc,$sim) = @$_; [$sim,$sc,
129 : overbeek 1.1 &HTML::fid_link($cgi,$peg),
130 : overbeek 1.4 $fig_or_sprout->genus_species($fig_or_sprout->genome_of($peg)),
131 :     scalar $fig_or_sprout->function_of($peg,$cgi->param('user')),
132 :     &HTML::set_prot_links($cgi,join( ', ', $fig_or_sprout->feature_aliases($peg) ))
133 : overbeek 1.3 ] } @pinned;
134 : overbeek 1.6 if (@tab > 0)
135 :     {
136 : overbeek 1.8 push(@$html,&HTML::make_table(["Sim. Sc.","Cluster Size","PEG","Genome", "Function","Aliases"],\@tab,"PEGs that Might Be in Clusters"));
137 : overbeek 1.6 }
138 :     else
139 :     {
140 :     push(@$html, $cgi->h1("Sorry, we have no clusters containing homologs of $peg"));
141 :     }
142 : parrello 1.10 }
143 : overbeek 1.1
144 :     sub relevant_homologs {
145 : overbeek 1.4 my($fig_or_sprout,$cgi,$peg) = @_;
146 : overbeek 1.1 my($maxN,$maxP,$genome1,$sim,$id2,$genome2,%seen);
147 :    
148 :     $maxN = $cgi->param('maxN');
149 :     $maxN = $maxN ? $maxN : 50;
150 :    
151 :     $maxP = $cgi->param('maxP');
152 :     $maxP = $maxP ? $maxP : 1.0e-10;
153 :    
154 : overbeek 1.4 my @sims = $fig_or_sprout->sims( $peg, $maxN, $maxP, "fig");
155 : overbeek 1.1
156 :     my @homologs = ();
157 :     $seen{&FIG::genome_of($peg)} = 1;
158 :     foreach $sim (@sims)
159 :     {
160 :     $id2 = $sim->id2;
161 :     $genome2 = &FIG::genome_of($id2);
162 : overbeek 1.3 my @coup;
163 : overbeek 1.4 if ((! $seen{$genome2}) && (@coup = $fig_or_sprout->coupled_to($id2)) && (@coup > 0))
164 : overbeek 1.1 {
165 :     $seen{$genome2} = 1;
166 : overbeek 1.8 push(@homologs,[$id2,@coup+1,$sim->psc]);
167 : overbeek 1.1 }
168 :     }
169 : overbeek 1.3 return sort { $b->[1] <=> $a->[1] } @homologs;
170 : overbeek 1.1 }
171 :    
172 :     sub sets_of_homologs {
173 : overbeek 1.4 my($fig_or_sprout,$cgi,$given_peg,$pinned) = @_;
174 : overbeek 1.1 my($peg,$mid,$min,$max,$feat,$fid);
175 :    
176 :     my $bound = $cgi->param('bound');
177 :     $bound = $bound ? $bound : 4000;
178 :    
179 :     my @pegs = ();
180 :     foreach $peg (($given_peg,@$pinned))
181 :     {
182 : overbeek 1.4 my $loc = $fig_or_sprout->feature_location($peg);
183 : overbeek 1.1 if ($loc)
184 :     {
185 : parrello 1.10 my($contig,$beg,$end) = $fig_or_sprout->boundaries_of($loc);
186 : overbeek 1.1 if ($contig && $beg && $end)
187 :     {
188 :     $mid = int(($beg + $end) / 2);
189 :     $min = $mid - $bound;
190 :     $max = $mid + $bound;
191 :    
192 : overbeek 1.4 ($feat,undef,undef) = &genes_in_region($fig_or_sprout,$cgi,&FIG::genome_of($peg),$contig,$min,$max);
193 : overbeek 1.1 foreach $fid (@$feat)
194 :     {
195 :     if (&FIG::ftype($fid) eq "peg")
196 :     {
197 :     push(@pegs,$fid);
198 :     }
199 :     }
200 :     }
201 :     }
202 :     }
203 :    
204 :     my %represents;
205 :     foreach $peg (@pegs)
206 :     {
207 : overbeek 1.4 my $tmp = $fig_or_sprout->maps_to_id($peg);
208 : overbeek 1.1 push(@{$represents{$tmp}},$peg);
209 : overbeek 1.4 # if ($tmp ne $peg) { push(@{$represents{$peg}},$peg) }
210 : overbeek 1.1 }
211 :     my($sim,%conn,$x,$y,$i,$j);
212 :     foreach $y (keys(%represents))
213 :     {
214 :     $x = $represents{$y};
215 :     for ($i=0; ($i < @$x); $i++)
216 :     {
217 :     for ($j=$i+1; ($j < @$x); $j++)
218 :     {
219 :     push(@{$conn{$x->[$i]}},$x->[$j]);
220 :     push(@{$conn{$x->[$j]}},$x->[$i]);
221 :     }
222 :     }
223 :     }
224 :    
225 :     my $maxN = $cgi->param('maxN');
226 :     $maxN = $maxN ? $maxN : 50;
227 :    
228 :     my $maxP = $cgi->param('maxP');
229 :     $maxP = $maxP ? $maxP : 1.0e-10;
230 :    
231 :     foreach $peg (@pegs)
232 :     {
233 : overbeek 1.4 foreach $sim ($fig_or_sprout->sims( $peg, $maxN, $maxP, "raw"))
234 : overbeek 1.1 {
235 :     if (defined($x = $represents{$sim->id2}))
236 :     {
237 :     foreach $y (@$x)
238 :     {
239 :     push(@{$conn{$peg}},$y);
240 :     }
241 :     }
242 :     }
243 :     }
244 :    
245 :     my(%seen,$k,$cluster);
246 :     my @clusters = ();
247 :     for ($i=0; ($i < @pegs); $i++)
248 :     {
249 :     $peg = $pegs[$i];
250 :     if (! $seen{$peg})
251 :     {
252 :     $cluster = [$peg];
253 :     $seen{$peg} = 1;
254 :     for ($j=0; ($j < @$cluster); $j++)
255 :     {
256 :     $x = $conn{$cluster->[$j]};
257 :     foreach $k (@$x)
258 :     {
259 :     if (! $seen{$k})
260 :     {
261 :     push(@$cluster,$k);
262 :     $seen{$k} = 1;
263 :     }
264 :     }
265 :     }
266 : parrello 1.10
267 : overbeek 1.1 if (@$cluster > 1)
268 :     {
269 :     push(@clusters,$cluster);
270 :     }
271 :     }
272 :     }
273 :     return @clusters;
274 :     }
275 :    
276 :     sub extract_homologs {
277 :     my($given_peg,$pinned,$clusters) = @_;
278 :     my(%main,$cluster,$peg,%counts,@with_counts);
279 :    
280 :     %main = map { $_ => 1 } ($given_peg,@$pinned);
281 :     foreach $cluster (@$clusters)
282 :     {
283 :     foreach $peg (@$cluster)
284 :     {
285 :     if (! $main{$peg})
286 :     {
287 :     $counts{&FIG::genome_of($peg)} += @$cluster - 1;
288 :     }
289 :     }
290 :     }
291 :    
292 :     foreach $peg (($given_peg,@$pinned))
293 :     {
294 :     push(@with_counts,[$peg,$counts{&FIG::genome_of($peg)}]);
295 :     }
296 :    
297 :     return grep { $_->[1] > 2} sort { $b->[1] <=> $a->[1] } @with_counts;
298 :     }

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3