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

Annotation of /FigWebServices/find_poss_subsys_instances.cgi

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.2 - (view) (download)

1 : overbeek 1.1 # -*- perl -*-
2 :     #
3 :     # Copyright (c) 2003-2006 University of Chicago and Fellowship
4 :     # for Interpretations of Genomes. All Rights Reserved.
5 :     #
6 :     # This file is part of the SEED Toolkit.
7 :     #
8 :     # The SEED Toolkit is free software. You can redistribute
9 :     # it and/or modify it under the terms of the SEED Toolkit
10 :     # Public License.
11 :     #
12 :     # You should have received a copy of the SEED Toolkit Public License
13 :     # along with this program; if not write to the University of Chicago
14 :     # at info@ci.uchicago.edu or the Fellowship for Interpretation of
15 :     # Genomes at veronika@thefig.info or download a copy from
16 :     # http://www.theseed.org/LICENSE.TXT.
17 :     #
18 :    
19 :    
20 :     use FIG;
21 :     my $fig = new FIG;
22 :    
23 :     use Subsystem;
24 :    
25 :     use HTML;
26 :     use strict;
27 :     use Carp;
28 :    
29 :     use CGI;
30 :     my $cgi = new CGI;
31 :    
32 :     if (0)
33 :     {
34 :     my $VAR1;
35 :     eval(join("",`cat /tmp/find_subsys_parms`));
36 :     $cgi = $VAR1;
37 : overbeek 1.2 # print STDERR &Dumper($cgi);
38 : overbeek 1.1 }
39 :    
40 : overbeek 1.2 if (0)
41 : overbeek 1.1 {
42 :     print $cgi->header;
43 :     my @params = $cgi->param;
44 :     print "<pre>\n";
45 :     foreach $_ (@params)
46 :     {
47 :     print "$_\t:",join(",",$cgi->param($_)),":\n";
48 :     }
49 :    
50 : overbeek 1.2 if (0)
51 : overbeek 1.1 {
52 :     if (open(TMP,">/tmp/find_subsys_parms"))
53 :     {
54 :     print TMP &Dumper($cgi);
55 :     close(TMP);
56 :     }
57 :     }
58 :     exit;
59 :     }
60 :    
61 : overbeek 1.2 use FigFams;
62 :     my $figfams = new FigFams($fig);
63 :    
64 :     use ToCall;
65 :    
66 : overbeek 1.1 my $html = [];
67 :     unshift @$html, "<TITLE>Find Instances of a Subsystem</TITLE>\n";
68 :    
69 :     my $subsys = $cgi->param('subsystem');
70 :     my $roles = $cgi->param('roles');
71 :     my $definitions = $cgi->param('definitions');
72 :     my $rules = $cgi->param('rules');
73 :     my @orgs = $cgi->param('korgs');
74 :     @orgs = map { $_ =~ /(\d+\.\d+)/; $1 } @orgs;
75 :     my(@rolesI,@rulesI,@definitionsI);
76 :     if ($roles =~ /\S/)
77 :     {
78 :     $roles =~ tr/\r/\n/;
79 :     @rolesI = grep { $_ } split(/\n/,$roles);
80 :     }
81 :    
82 :     if ($definitions =~ /\S/)
83 :     {
84 :     $definitions =~ tr/\r/\n/;
85 :     @definitionsI = grep { $_ } split(/\n/,$definitions);
86 :     }
87 :    
88 :     if ($rules =~ /\S/)
89 :     {
90 :     $rules =~ tr/\r/\n/;
91 :     @rulesI = grep { $_ } split(/\n/,$rules);
92 :     }
93 :    
94 :     my %org_labels;
95 :    
96 : overbeek 1.2 if ((! $subsys) || (@rolesI < 1) || (@rulesI < 1) || (@orgs < 1))
97 : overbeek 1.1 {
98 :     my @orgs = ();
99 :     foreach my $org ($fig->genomes('complete'))
100 :     {
101 :     my $label = &compute_genome_label($fig, $org);
102 :     $org_labels{$org} = $label;
103 :     push(@orgs, $org);
104 :     }
105 :     @orgs = sort { lc( $org_labels{$a} ) cmp lc( $org_labels{$b} ) } @orgs;
106 :     my @org_names = map { $org_labels{$_} } @orgs;
107 :    
108 :     my @sub = $fig->all_subsystems;
109 :     push(@$html,$cgi->start_form(-action => "find_poss_subsys_instances.cgi", -method => 'get'),
110 :     $cgi->scrolling_list( -name => 'subsystem',
111 :     -values => [ sort @sub ],
112 :     -size => 1
113 :     ),$cgi->br,
114 :     $cgi->h2('Roles'),$cgi->textarea( -name => 'roles', -rows => 20, columns => 200), $cgi->br, $cgi->hr,$cgi->br,
115 :     $cgi->h2('Definitions'),$cgi->textarea( -name => 'definitions', -rows => 20, columns => 200), $cgi->hr,$cgi->br, $cgi->br,
116 :     $cgi->h2('Rules'),$cgi->textarea( -name => 'rules', -rows => 20, columns => 200), $cgi->br, $cgi->hr,$cgi->br,
117 :     $cgi->scrolling_list( -name => 'korgs',
118 :     -values => [ @orgs ],
119 :     -labels => \%org_labels,
120 :     -size => 10,
121 :     ), $cgi->br,
122 :     $cgi->submit( 'Compute Predicted Variant Codes' ),
123 :     $cgi->end_form
124 :     );
125 :     }
126 :     else
127 :     {
128 : overbeek 1.2 my %genomesS = map { $_ => 1 } @orgs;
129 : overbeek 1.1
130 :     my $encoding = [[],0]; # Encoding is a 2-tuple [Memory,NxtAvail]
131 :     my $abbrev_to_loc = {};
132 :    
133 :     my @roles = &load_roles($cgi,$html,$encoding,$abbrev_to_loc,\@rolesI);
134 :     if (@roles < 1)
135 :     {
136 :     push(@$html,$cgi->h1("Roles are invalid"));
137 :     }
138 :     else
139 :     {
140 :     my $rc = &load_definitions($cgi,$html,$encoding,$abbrev_to_loc,\@definitionsI);
141 :     if (! $rc)
142 :     {
143 :     push(@$html,$cgi->h1("Definitions are invalid"));
144 :     }
145 :     else
146 :     {
147 :     my @rules = &parse_rules($cgi,$html,$encoding,$abbrev_to_loc,\@rulesI);
148 :     if (@rules < 1)
149 :     {
150 :     push(@$html,$cgi->h1("Rules are invalid"));
151 :     }
152 :     else
153 :     {
154 :     my $n = @rules;
155 :     my $operational = 0;
156 :     foreach my $genome (map { $_->[0] }
157 :     sort { $a->[1] cmp $b->[1] }
158 :     map { [$_,$fig->genus_species($_)] }
159 : overbeek 1.2 @orgs)
160 : overbeek 1.1 {
161 : overbeek 1.2 my($vcT,$tab) = &find_vc($encoding,\@roles,$figfams,$abbrev_to_loc,\@rules,$fig,$genome,$cgi);
162 :     my $gs = $fig->genus_species($genome);
163 :     push(@$html,$cgi->h1("$vcT for $genome [$gs]"),$tab);
164 : overbeek 1.1 }
165 :     }
166 :     }
167 :     }
168 :     }
169 :     &HTML::show_page($cgi,$html);
170 :    
171 :    
172 :     sub find_vc {
173 : overbeek 1.2 my($encoding,$roles,$figfams,$abbrev_to_loc,$rules,$fig,$genome,$cgi) = @_;
174 : overbeek 1.1
175 :     my $vcT = undef;
176 :     my $rule;
177 :    
178 :     my $role;
179 :     my $relevant_genes = {};
180 : overbeek 1.2 my $tab = [];
181 :     foreach $role (@$roles)
182 : overbeek 1.1 {
183 : overbeek 1.2 $relevant_genes->{$role} = &gather_genes($fig,$genome,$figfams,$role);
184 :     foreach my $x (@{$relevant_genes->{$role}})
185 :     {
186 :     my($peg,$func) = @$x;
187 :     my $link = &HTML::fid_link($cgi,$peg,1);
188 :     push(@$tab,[$link,$func]);
189 :     }
190 : overbeek 1.1 }
191 :    
192 :     my($i,$matched);
193 :     for ($i=0, $matched=undef; (! defined($matched)) && ($i < @$rules); $i++)
194 :     {
195 :     $matched = &is_rule_true($rules->[$i],$relevant_genes);
196 :     }
197 :    
198 :     my $vcT = defined($matched) ? $matched : -1;
199 : overbeek 1.2 my $pegs = &HTML::make_table(["PEG","Proposed Function"],$tab,"PEGs with Relevant Functions");
200 :     return ($vcT,$pegs);
201 : overbeek 1.1 }
202 :    
203 :     sub load_roles {
204 :     my($cgi,$html,$encoding,$abbrev_to_loc,$rolesI) = @_;
205 :    
206 :     my @roles = ();
207 :     foreach $_ (@$rolesI)
208 :     {
209 :     if ($_ =~ /^(\S+)\s+(\S.*\S)/)
210 :     {
211 :     my($abbrev,$role) = ($1,$2);
212 :     my $loc = &add_to_encoding($encoding,['role',$role]);
213 :     $abbrev_to_loc->{$abbrev} = $loc;
214 :     push(@roles,$role);
215 :     }
216 :     elsif ($_ =~ /\S/)
217 :     {
218 :     push(@$html,$cgi->h1("Invalid Role: $_"));
219 :     }
220 :     }
221 :     return @roles;
222 :     }
223 :    
224 :     sub add_to_encoding {
225 :     my($encoding,$val) = @_;
226 :    
227 :     my($mem,$nxt) = @$encoding;
228 :     $mem->[$nxt] = $val;
229 :     $encoding->[1]++;
230 :     return $nxt;
231 :     }
232 :    
233 :     sub load_definitions {
234 :     my($cgi,$html,$encoding,$abbrev_to_loc,$defI) = @_;
235 :    
236 :     my $rc = 1;
237 :     foreach my $def (@$defI)
238 :     {
239 :     if ($def =~ /^(\S+)\s+(\S.*\S)/)
240 :     {
241 :     my($abbrev,$bool) = ($1,$2);
242 :     my $loc = &parse_bool($bool,$encoding,$abbrev_to_loc);
243 :     $abbrev_to_loc->{$abbrev} = $loc;
244 :     }
245 :     elsif ($def =~ /\S/)
246 :     {
247 :     push(@$html,$cgi->h1("Invalid Definition: $def"));
248 :     $rc = 0;
249 :     }
250 :     }
251 :     return $rc;
252 :     }
253 :    
254 :     sub parse_rules {
255 :     my($cgi,$html,$encoding,$abbrev_to_loc,$rulesI) = @_;
256 :    
257 :     my @rules = ();
258 :     foreach $_ (@$rulesI)
259 :     {
260 :     my($boolexp,$variant_code,$loc);
261 :     if (($_ =~ /^\s*(\S+)\s+(\S.*\S)\s*$/) &&
262 :     (($variant_code,$boolexp) = ($1,$2)) &&
263 :     defined($loc = &parse_bool($boolexp,$encoding,$abbrev_to_loc)))
264 :     {
265 :     push(@rules,[$variant_code,[$encoding->[0],$loc]]);
266 :     }
267 :     elsif ($_ =~ /\S/)
268 :     {
269 :     push(@$html,$cgi->h1("Invalid rule: $_"));
270 :     }
271 :     }
272 :     return @rules;
273 :     }
274 :    
275 :     sub parse_bool {
276 :     my($s,$encoding,$abbrev_to_loc) = @_;
277 :    
278 :     my $input = $s;
279 :     my $abbrev;
280 :     foreach $abbrev (sort { length($b) <=> length($a) } keys(%$abbrev_to_loc))
281 :     {
282 :     my $loc = $abbrev_to_loc->{$abbrev};
283 :     my $abbrevQ = quotemeta $abbrev;
284 :     while ($s =~ s/(^|[\s\{,(])($abbrevQ)($|[\s\},)])/$1<$loc>$3/) {}
285 :     }
286 :     my $got = 0;
287 :     while ($s !~ /^\s*<\d+>\s*$/)
288 :     {
289 :     my $nxt = $encoding->[1];
290 :     if ($s =~ s/\(\s*(<\d+>)\s*\)/$1/)
291 :     {
292 :     $got = 1;
293 :     }
294 :     elsif ($s =~ s/not\s+<(\d+)>/<$nxt>/)
295 :     {
296 :     &add_to_encoding($encoding,["not",$1]);
297 :     $got = 1;
298 :     }
299 :     elsif ($s =~ s/<(\d+)>\s+and\s+<(\d+)>/<$nxt>/)
300 :     {
301 :     &add_to_encoding($encoding,["and",$1,$2]);
302 :     $got = 1;
303 :     }
304 :     elsif ($s =~ s/<(\d+)>\s+or\s+<(\d+)>/<$nxt>/)
305 :     {
306 :     &add_to_encoding($encoding,["or",$1,$2]);
307 :     $got = 1;
308 :     }
309 :     elsif ($s =~ s/<(\d+)>\s+->\s+<(\d+)>/<$nxt>/)
310 :     {
311 :     &add_to_encoding($encoding,["->",$1,$2]);
312 :     $got = 1;
313 :     }
314 :    
315 :     elsif ($s =~ s/(\d+)\s+of\s+\{\s*(<\d+>(,\s*<\d+>)*)\s*\}/<$nxt>/)
316 :     {
317 :     my $n = $1;
318 :     my $args = $2;
319 :     my @args = map { $_ =~ /<(\d+)>/; $1 } split(/,\s*/,$args);
320 :     &add_to_encoding($encoding,["of",$n,[@args]]);
321 :     $got = 1;
322 :     }
323 :     last if (! $got);
324 :     }
325 :     return ($s =~ /^\s*<(\d+)>\s*$/) ? $1 : undef;
326 :     }
327 :    
328 : overbeek 1.2 sub gather_genes {
329 :     my($fig,$genome,$figfams,$role) = @_;
330 :     my($fam_id,$peg,$hit,%found);
331 :    
332 :     # print STDERR "Gathering for role=$role\n";
333 :    
334 :     my $genomeO = new ToCall("$FIG_Config::organisms/$genome");
335 :     my @fam_ids = $figfams->families_with_functional_role($role);
336 :     # print STDERR "Families for role: ",join(",",@fam_ids),"\n";
337 : overbeek 1.1
338 : overbeek 1.2 foreach $fam_id (@fam_ids)
339 :     {
340 :     my $figfam = new FigFam($fig,$fam_id);
341 :     my $query;
342 : overbeek 1.1
343 : overbeek 1.2 foreach $query ($figfam->representatives)
344 :     {
345 :     foreach $hit ($genomeO->candidate_orfs(-seq => $query))
346 :     {
347 :     # print STDERR " candidate=",$hit->get_fid,"\n";
348 :     my $hit_seq = $hit->seq;
349 :     my($ok,undef) = $figfam->should_be_member($hit_seq);
350 :     if ($ok)
351 :     {
352 :     # print STDERR " match\n";
353 :     $found{$hit->get_fid} = $figfam->family_function;
354 :     }
355 :     }
356 :     }
357 :     }
358 :     return [sort { &FIG::by_fig_id($a->[0],$b->[0]) } map { [$_,$found{$_}] } keys(%found)];
359 : overbeek 1.1 }
360 :    
361 :     sub is_rule_true {
362 :     my($rule,$relevant_genes) = @_;
363 :    
364 :     my($variant,$exp) = @$rule;
365 :     return &is_true_exp($exp,$relevant_genes) ? $variant : undef;
366 :     }
367 :    
368 :     sub is_true_exp {
369 :     my($bool,$relevant_genes) = @_;
370 :    
371 :     my($nodes,$root) = @$bool;
372 :     my $val = $nodes->[$root];
373 :     if (! ref $val)
374 :     {
375 :     return &is_true_exp([$nodes,$val],$relevant_genes);
376 :     }
377 :     else
378 :     {
379 :     my $op = $val->[0];
380 :    
381 :     if ($op eq 'role')
382 :     {
383 :     my $x;
384 :     return (($x = $relevant_genes->{$val->[1]}) && (@$x > 0)) ? 1 : 0;
385 :     }
386 :     elsif ($op eq "of")
387 :     {
388 :     my $truth_value;
389 :     my $count = 0;
390 :     foreach $truth_value (map { &is_true_exp([$nodes,$_],$relevant_genes) } @{$val->[2]})
391 :     {
392 :     if ($truth_value) { $count++ }
393 :     }
394 :     return $val->[1] <= $count;
395 :     }
396 :     elsif ($op eq "not")
397 :     {
398 :     return &is_true_exp([$nodes,$val->[1]],$relevant_genes) ? 0 : 1;
399 :     }
400 :     else
401 :     {
402 :     my $v1 = &is_true_exp([$nodes,$val->[1]],$relevant_genes);
403 :     my $v2 = &is_true_exp([$nodes,$val->[2]],$relevant_genes);
404 :     if ($op eq "and") { return $v1 && $v2 };
405 :     if ($op eq "or") { return $v1 || $v2 };
406 :     if ($op eq "->") { return ((not $v1) || $v2) }
407 :     else
408 :     {
409 :     print STDERR &Dumper($val);
410 :     die "invalid expression";
411 :     }
412 :     }
413 :     }
414 :     }
415 :    
416 :     sub print_bool {
417 :     my($bool) = @_;
418 :    
419 :     my $s = &printable_bool($bool);
420 :     print $s,"\n";
421 :     }
422 :    
423 :     sub printable_bool {
424 :     my($bool) = @_;
425 :    
426 :     my($nodes,$root) = @$bool;
427 :     my $val = $nodes->[$root];
428 :    
429 :     if (! ref $val)
430 :     {
431 :     return &printable_bool([$nodes,$val]);
432 :     }
433 :     else
434 :     {
435 :     my $op = $val->[0];
436 :    
437 :     if ($op eq 'role')
438 :     {
439 :     return $val->[1];
440 :     }
441 :     elsif ($op eq "of")
442 :     {
443 :     my @expanded_args = map { &printable_bool([$nodes,$_]) } @{$val->[2]};
444 :     my $args = join(',',@expanded_args);
445 :     return "$val->[1] of {$args}";
446 :     }
447 :     elsif ($op eq "not")
448 :     {
449 :     return "($op " . &printable_bool([$nodes,$val->[1]]) . ")";
450 :     }
451 :     else
452 :     {
453 :     return "(" . &printable_bool([$nodes,$val->[1]]) . " $op " . &printable_bool([$nodes,$val->[2]]) . ")";
454 :     }
455 :     }
456 :     }
457 :    
458 :     sub compute_genome_label
459 :     {
460 :     my($fig, $org) = @_;
461 :     my $label;
462 :    
463 :     my $gs = $fig->genus_species($org);
464 :     if ($fig->genome_domain($org) ne "Environmental Sample")
465 :     {
466 :     my $gc=$fig->number_of_contigs($org);
467 :     $label = "$gs ($org) [$gc contigs]";
468 :     }
469 :     else
470 :     {
471 :     $label = "$gs ($org)";
472 :     }
473 :     return $label;
474 :     }
475 :    

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3