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

Annotation of /FigKernelPackages/UnvSubsys.pm

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : overbeek 1.1 package UnvSubsys;
2 :    
3 : overbeek 1.2 use Subsystem;
4 : overbeek 1.1 use Carp;
5 :     use FIG;
6 :    
7 :     use Data::Dumper;
8 :     use strict;
9 :    
10 :     sub new
11 :     {
12 : overbeek 1.7 my($class, $ssa, $fig, $active_subsetR, $focus, $show_clusters, $aliases) = @_;
13 : overbeek 1.1
14 :     $ssa =~ s/ /_/g;
15 :    
16 :    
17 : parrello 1.5 ### { Roles =>Roles,
18 : overbeek 1.2 ### RoleIndex => ToRoleIndexHash,
19 :     ### RoleSubsets => ColSubsets,
20 :     ### Genomes => Genomes,
21 :     ### GenomeIndex => ToGenomeIndexHash,
22 :     ### PegHash => PegHash,
23 :     ### Colors => ColorHash,
24 :     ### Aliases => AliasHash,
25 :     ### Curator => Curator,
26 :     ### Notes => Notes,
27 :     ### Reactions => ReactionHash
28 :     ### }
29 : parrello 1.5 ###
30 :     ### Roles = pointer to a list of [Role,Abbrev,[ReactionURLs]]
31 :     ###
32 :     ### ToRoleIndexHash = a pointer to a hash: key=Role Value=RoleIndex
33 :     ###
34 :     ### ColSubsets = pointer to a list of [SubsetName,[RoleIndexesFrom0]]
35 :     ###
36 : overbeek 1.6 ### RowSubSets = pointer to a list of [SubsetName,[GenomeIndexesFrom0]]
37 :     ###
38 : parrello 1.5 ### Genomes is a pointer to a list of [Genome,Variant]
39 :     ###
40 :     ### ToGenomeIndexHash = a pointer to a hash: key=Genome value=GenomeIndex
41 :     ###
42 :     ### PegHash = a pointer to a hash of hashes such that $peg_hash->{$genome_index}->{$role_index} = a
43 :     ### pointer to a list of PEGs
44 :     ###
45 :     ### ColorHash is a hash: key=PEG value=color
46 :     ###
47 :     ### AliasHash is a hash: key=PEG value=aliases
48 :     ###
49 : overbeek 1.2 ### ReactionHash is a hash: key=Role value=[reaction-ids]
50 : overbeek 1.1
51 : overbeek 1.6 if ((ref($fig) eq "FIG") || ((ref($fig) eq 'SFXlate') && ($fig = $fig->{'fig'})))
52 : overbeek 1.1 {
53 : parrello 1.5 my $subsystem = new Subsystem($ssa,$fig,0);
54 :     my $curator = $subsystem->get_curator;
55 :     my $notes = $subsystem->get_notes;
56 :     $notes =~ s/ /\n/g;
57 :     my @roles = $subsystem->get_roles;
58 :     my $reactions = $subsystem->get_reactions;
59 :     my @genomes = $subsystem->get_genomes;
60 : overbeek 1.1 my @col_subsets = $subsystem->get_subset_namesC;
61 :    
62 : parrello 1.5 my $role_info = [];
63 :     my $roleH = {};
64 : overbeek 1.1
65 : parrello 1.5 my($i,$j,$subset,$peg);
66 :     for ($i=0; ($i < @roles); $i++)
67 :     {
68 :     my $role = $roles[$i];
69 :     my $abbrev = $subsystem->get_role_abbr( $subsystem->get_role_index( $role ) );
70 :     my $react = $reactions ? join(",", map { &HTML::reaction_link($_) } @{$reactions->{$role}}) : [];
71 :     push(@$role_info,[$role,$abbrev,$react]);
72 :     $roleH->{$role} = $i;
73 :     }
74 :    
75 :     my $subset_info = [];
76 :     foreach $subset (@col_subsets)
77 :     {
78 :     if ($subset ne 'All')
79 :     {
80 :     push(@$subset_info,[$subset,[map { $roleH->{$_} } $subsystem->get_subsetC_roles($subset)]]);
81 :     }
82 :     }
83 :    
84 :     my $genomes_info = [];
85 :     my $genomeH = {};
86 :     for ($i=0; ($i < @genomes); $i++)
87 :     {
88 :     my $genome = $genomes[$i];
89 :     my $variant = $subsystem->get_variant_code( $subsystem->get_genome_index( $genome ) );
90 :     push(@$genomes_info,[$genome,$variant]);
91 :     $genomeH->{$genome} = $i;
92 :     }
93 :    
94 :     my $pegH = {};
95 :     for ($i=0; ($i < @genomes); $i++)
96 :     {
97 :     my $genome = $genomes[$i];
98 :     for ($j=0; ($j < @roles); $j++)
99 :     {
100 :     my $role = $roles[$j];
101 :     my @pegs = $subsystem->get_pegs_from_cell($genome,$role);
102 :     $pegH->{$i}->{$j} = [@pegs];
103 :     }
104 :     }
105 : overbeek 1.7
106 :     my $row_subsets = &row_subsets($fig,$genomeH,$genomes_info);
107 :     my $active_genomes = &active_genomes($row_subsets,$active_subsetR,$focus,$genomeH,$genomes_info);
108 :     my $colorsH = $show_clusters ? &set_colors($fig,$pegH,$active_genomes) : {};
109 :     my $aliasesH = $aliases ? &set_aliases($fig,$pegH,$active_genomes) : {};
110 : overbeek 1.6 my $reactions = $subsystem->get_reactions;
111 : parrello 1.5 my $self = { Roles => $role_info,
112 :     RoleIndex => $roleH,
113 :     RoleSubsets => $subset_info,
114 :     Genomes => $genomes_info,
115 :     GenomeIndex => $genomeH,
116 : overbeek 1.6 GenomeSubsets => $row_subsets,
117 : parrello 1.5 PegHash => $pegH,
118 :     Colors => $colorsH,
119 :     Aliases => $aliasesH,
120 :     Curator => $curator,
121 :     Notes => $notes,
122 :     Reactions => $reactions
123 :     };
124 :     bless($self, $class);
125 :     return $self;
126 : overbeek 1.1 }
127 :     else
128 :     {
129 : parrello 1.5 return undef;
130 : overbeek 1.1 }
131 :     }
132 :    
133 : overbeek 1.6 sub get_subset_namesR {
134 :     my($self) = @_;
135 :    
136 :     return map { $_->[0] } @{$self->{GenomeSubsets}};
137 :     }
138 :    
139 :     sub get_subsetR {
140 :     my($self,$set) = @_;
141 :     my($i);
142 :    
143 :     my $sets = $self->{GenomeSubsets};
144 :     for ($i=0; ($i < @$sets) && ($sets->[$i]->[0] ne $set); $i++) {}
145 :     if ($i < @$sets)
146 :     {
147 :     return map { $self->{Genomes}->[$_]->[0] } @{$sets->[$i]->[1]}
148 :     }
149 :     return undef;
150 :     }
151 :    
152 :     sub get_subsetsR {
153 :     my($self) = @_;
154 :    
155 :     my $sets = $self->{GenomeSubsets};
156 :     my @pairs = ();
157 :     my $pair;
158 :     foreach $pair (@$sets)
159 :     {
160 :     my($id,$members) = @$pair;
161 :     push(@pairs,[$id,[map { $self->{Genomes}->[$_]->[0] } @$members]]);
162 :     }
163 :     return @pairs;
164 :     }
165 :    
166 :     sub row_subsets {
167 :     my ($fig,$genomeH,$genomes_info) = @_;
168 :    
169 :     my $subsets = [];
170 :     my $taxonomic_groups = $fig->taxonomic_groups_of_complete(10);
171 :    
172 :     my($pair,$id,$members);
173 :     foreach $pair (@$taxonomic_groups)
174 :     {
175 :     ($id,$members) = @$pair;
176 :     my @mem = grep { defined($_) } map { $genomeH->{$_} } @$members;
177 :     if (@mem > 10)
178 :     {
179 :     push(@$subsets,[$id,[@mem]]);
180 :     }
181 :     }
182 :     return $subsets;
183 :     }
184 :    
185 : overbeek 1.4 sub set_aliases {
186 : overbeek 1.7 my($fig,$pegH,$active_genomes) = @_;
187 : overbeek 1.4 my($genomeI,$roleI,$pegs,$peg,$roleH);
188 :    
189 :     my $aliasesH = {};
190 :    
191 : overbeek 1.7 foreach $genomeI (grep { $active_genomes->{$_} } keys(%$pegH))
192 : overbeek 1.4 {
193 : parrello 1.5 $roleH = $pegH->{$genomeI};
194 :     foreach $roleI (keys(%$roleH))
195 :     {
196 :     $pegs = $roleH->{$roleI};
197 :     foreach $peg (@$pegs)
198 :     {
199 :     if (! $aliasesH->{$peg})
200 :     {
201 :     $aliasesH->{$peg} = scalar &ext_id($fig,$peg);
202 :     }
203 :     }
204 :     }
205 : overbeek 1.4 }
206 :     return $aliasesH;
207 :     }
208 :    
209 :     sub set_colors {
210 : overbeek 1.7 my($fig,$pegH,$active_genomes) = @_;
211 : overbeek 1.4 my($genomeI,$roleI,$pegs,$peg,$roleH,$peg,%pegs_in_genome);
212 :    
213 :     my $colorsH = {};
214 :    
215 : overbeek 1.7 foreach $genomeI (grep { $active_genomes->{$_} } keys(%$pegH))
216 : overbeek 1.4 {
217 : parrello 1.5 undef %pegs_in_genome;
218 :     $roleH = $pegH->{$genomeI};
219 :     foreach $roleI (keys(%$roleH))
220 :     {
221 :     $pegs = $roleH->{$roleI};
222 :     foreach $peg (@$pegs)
223 :     {
224 :     $pegs_in_genome{$peg} = 1;
225 :     }
226 :     }
227 :    
228 :     my @pegs = keys(%pegs_in_genome);
229 :     my($tuple,$peg,$color);
230 :     my $colors_for_one_genome = &set_colors_for_genome($fig,\@pegs);
231 :    
232 :     while (($peg,$color) = each %$colors_for_one_genome)
233 :     {
234 :     $colorsH->{$peg} = $colors_for_one_genome->{$peg};
235 :     }
236 : overbeek 1.4 }
237 :     return $colorsH;
238 :     }
239 :    
240 :     sub set_colors_for_genome {
241 :     my($fig,$pegs) = @_;
242 :     my($peg,@clusters,$cluster,@colors,$color,%seen,%conn,$x,$peg1,@pegs,$i);
243 :    
244 :     my $color_of = {};
245 :     foreach $peg (@$pegs) { $color_of->{$peg} = '#FFFFFF' }
246 :    
247 :     @pegs = keys(%$color_of); # Use of keys makes @pegs entries unique
248 :    
249 :     foreach $peg (@pegs)
250 :     {
251 : parrello 1.5 $conn{$peg} = [grep { $color_of->{$_} && ($_ ne $peg) } $fig->close_genes($peg,5000)];
252 : overbeek 1.4 }
253 :    
254 :     @clusters = ();
255 :     while ($peg = shift @pegs)
256 :     {
257 : parrello 1.5 if (! $seen{$peg})
258 :     {
259 :     $cluster = [$peg];
260 :     $seen{$peg} = 1;
261 :     for ($i=0; ($i < @$cluster); $i++)
262 :     {
263 :     my @tmp = grep { ! $seen{$_} } @{$conn{$cluster->[$i]}};
264 :     if (@tmp > 0)
265 :     {
266 :     foreach my $peg1 (@tmp) { $seen{$peg1} = 1 }
267 :     push(@$cluster,@tmp);
268 :     }
269 :     }
270 :     push(@clusters,$cluster);
271 :     }
272 : overbeek 1.4 }
273 :    
274 :     @colors = &cool_colors();
275 :    
276 :     @clusters = grep { @$_ > 1 } sort { @$a <=> @$b } @clusters;
277 :    
278 :     if (@clusters > @colors) { splice(@clusters,0,(@clusters - @colors)) } # make sure we have enough colors
279 :    
280 :     my($cluster);
281 :     foreach $cluster (@clusters)
282 :     {
283 : parrello 1.5 $color = shift @colors;
284 :     foreach $peg (@$cluster)
285 :     {
286 :     $color_of->{$peg} = $color;
287 :     }
288 : overbeek 1.4 }
289 :     return $color_of;
290 :     }
291 :    
292 :     sub cool_colors {
293 :     # just an array of "websafe" colors or whatever colors we want to use. Feel free to remove bad colors (hence the lines not being equal length!)
294 :     return (
295 :     '#C0C0C0', '#FF40C0', '#FF8040', '#FF0080', '#FFC040', '#40C0FF', '#40FFC0', '#C08080', '#C0FF00', '#00FF80', '#00C040',
296 :     "#6B8E23", "#483D8B", "#2E8B57", "#008000", "#006400", "#800000", "#00FF00", "#7FFFD4",
297 :     "#87CEEB", "#A9A9A9", "#90EE90", "#D2B48C", "#8DBC8F", "#D2691E", "#87CEFA", "#E9967A", "#FFE4C4", "#FFB6C1",
298 :     "#E0FFFF", "#FFA07A", "#DB7093", "#9370DB", "#008B8B", "#FFDEAD", "#DA70D6", "#DCDCDC", "#FF00FF", "#6A5ACD",
299 :     "#00FA9A", "#228B22", "#1E90FF", "#FA8072", "#CD853F", "#DC143C", "#FF6347", "#98FB98", "#4682B4",
300 :     "#D3D3D3", "#7B68EE", "#2F4F4F", "#FF7F50", "#FF69B4", "#BC8F8F", "#A0522D", "#DEB887", "#00DED1",
301 :     "#6495ED", "#800080", "#FFD700", "#F5DEB3", "#66CDAA", "#FF4500", "#4B0082", "#CD5C5C",
302 :     "#EE82EE", "#7CFC00", "#FFFF00", "#191970", "#FFFFE0", "#DDA0DD", "#00BFFF", "#DAA520", "#008080",
303 :     "#00FF7F", "#9400D3", "#BA55D3", "#D8BFD8", "#8B4513", "#3CB371", "#00008B", "#5F9EA0",
304 :     "#4169E1", "#20B2AA", "#8A2BE2", "#ADFF2F", "#556B2F",
305 :     "#F0FFFF", "#B0E0E6", "#FF1493", "#B8860B", "#FF0000", "#F08080", "#7FFF00", "#8B0000",
306 :     "#40E0D0", "#0000CD", "#48D1CC", "#8B008B", "#696969", "#AFEEEE", "#FF8C00", "#EEE8AA", "#A52A2A",
307 :     "#FFE4B5", "#B0C4DE", "#FAF0E6", "#9ACD32", "#B22222", "#FAFAD2", "#808080", "#0000FF",
308 :     "#000080", "#32CD32", "#FFFACD", "#9932CC", "#FFA500", "#F0E68C", "#E6E6FA", "#F4A460", "#C71585",
309 :     "#BDB76B", "#00FFFF", "#FFDAB9", "#ADD8E6", "#778899",
310 :     );
311 :     }
312 :    
313 :     sub ext_id {
314 :     my($fig,$peg) = @_;
315 :    
316 :     my @tmp;
317 :     my @aliases = $fig->feature_aliases($peg);
318 :     if ((@tmp = grep { $_ =~ /^uni\|/ } @aliases) > 0)
319 :     {
320 : parrello 1.5 @aliases = @tmp;
321 : overbeek 1.4 }
322 :     elsif ((@tmp = grep { $_ =~ /^sp\|/ } @aliases) > 0)
323 :     {
324 : parrello 1.5 @aliases = @tmp;
325 : overbeek 1.4 }
326 :     elsif ((@tmp = grep { $_ =~ /^gi\|/ } @aliases) > 0)
327 :     {
328 : parrello 1.5 @aliases = @tmp;
329 : overbeek 1.4 }
330 :     elsif ((@tmp = grep { $_ =~ /^kegg\|/ } @aliases) > 0)
331 :     {
332 : parrello 1.5 @aliases = @tmp;
333 : overbeek 1.4 }
334 :     else
335 :     {
336 : parrello 1.5 @aliases = ();
337 : overbeek 1.4 }
338 :    
339 :     if (wantarray())
340 :     {
341 : parrello 1.5 return @aliases;
342 : overbeek 1.4 }
343 :     else
344 :     {
345 : parrello 1.5 return $aliases[0];
346 : overbeek 1.4 }
347 :     }
348 :    
349 :    
350 : overbeek 1.2 sub subsystem_curator {
351 :     my($self) = @_;
352 :    
353 :     my $curator = $self->{Curator};
354 :     $curator =~ s/master://;
355 :     return $curator;
356 :     }
357 :    
358 :     sub get_roles {
359 :     my($self) = @_;
360 :    
361 :     return map { $_->[0] } @{$self->{Roles}};
362 :     }
363 :    
364 : overbeek 1.3 sub get_genome_index {
365 :     my($self,$genome) = @_;
366 :    
367 :     return $self->{GenomeIndex}->{$genome};
368 :     }
369 :    
370 :     sub get_genomes {
371 :     my($self) = @_;
372 :    
373 :     return map { $_->[0] } @{$self->{Genomes}};
374 :     }
375 :    
376 :     sub get_variant_code {
377 :     my($self,$genome) = @_;
378 :    
379 :     if ($genome =~ /^\d+$/)
380 :     {
381 : parrello 1.5 return $self->{Genomes}->[$genome]->[1];
382 : overbeek 1.3 }
383 :     else
384 :     {
385 : parrello 1.5 my $genomeI = $self->{GenomeIndex}->{$genome};
386 :     return $self->{Genomes}->[$genomeI]->[1];
387 : overbeek 1.3 }
388 :     }
389 :    
390 :     sub get_pegs_from_cell {
391 :     my($self,$genome,$role) = @_;
392 :    
393 :     my $genomeI = $self->{GenomeIndex}->{$genome};
394 :     my $roleI = $self->{RoleIndex}->{$role};
395 :    
396 :     my $pegs = $self->{PegHash}->{$genomeI}->{$roleI};
397 :     return $pegs ? @$pegs : ();
398 :     }
399 :    
400 :     sub get_notes {
401 :     my($self) = @_;
402 :    
403 :     return $self->{Notes};
404 :     }
405 :    
406 : overbeek 1.2 sub get_role_index {
407 :     my($self,$role) = @_;
408 :    
409 :     return $self->{RoleIndex}->{$role};
410 :     }
411 :    
412 :     sub get_role_abbr {
413 :     my($self,$roleI) = @_;
414 :    
415 :     if ($roleI !~ /^\d+$/)
416 :     {
417 : parrello 1.5 $roleI = $self->{RoleIndex}->{$roleI};
418 : overbeek 1.2 }
419 :     my $roles = $self->{Roles};
420 :     return $roles->[$roleI]->[1];
421 :     }
422 :    
423 :     sub get_reactions {
424 :     my($self) = @_;
425 :    
426 :     return $self->{Reactions};
427 :     }
428 :    
429 :     sub get_subset_namesC {
430 :     my($self) = @_;
431 :    
432 :     return map { $_->[0] } @{$self->{RoleSubsets}};
433 :     }
434 :    
435 :     sub get_subsetC_roles {
436 :     my($self,$subset) = @_;
437 :     my($i,$j);
438 :    
439 :     my $subset_info = $self->{RoleSubsets};
440 :     for ($i=0; ($i < @$subset_info) && ($subset_info->[$i]->[0] ne $subset); $i++) {}
441 :     if ($i < @$subset_info)
442 :     {
443 : parrello 1.5 my @roles = ();
444 :     foreach $j (@{$subset_info->[$i]->[1]})
445 :     {
446 :     push(@roles,$self->{Roles}->[$j]->[0]);
447 :     }
448 :     return @roles;
449 : overbeek 1.2 }
450 :     return undef;
451 :     }
452 :    
453 : overbeek 1.4 sub get_color_of {
454 :     my($self,$peg) = @_;
455 :    
456 :     return $self->{Colors}->{$peg};
457 :     }
458 :    
459 : overbeek 1.7 sub active_genomes {
460 :     my($row_subsets,$active_subsetR,$focus,$genomeH,$genomes_info) = @_;
461 :     my($i,@bestL);
462 :    
463 :     my $active_genomes = {};
464 :     if ($active_subsetR)
465 :     {
466 :     for ($i=0; ($i < @$row_subsets) && ($row_subsets->[$i]->[0] ne $active_subsetR); $i++) {}
467 :     if ($i < @$row_subsets)
468 :     {
469 :     @bestL = @{$row_subsets->[$i]->[1]};
470 :     }
471 :     else
472 :     {
473 :     $active_subsetR = 'All';
474 :     @bestL = keys(%$genomeH);
475 :     }
476 :     }
477 :     else
478 :     {
479 :     if (! $focus)
480 :     {
481 :     $active_subsetR = 'All';
482 :     @bestL = keys(%$genomeH);
483 :     }
484 :     else
485 :     {
486 :     my $bestN = undef;
487 :     @bestL = ();
488 :    
489 :     my $tuple;
490 :     foreach $tuple (@$row_subsets)
491 :     {
492 :     my($id,$genomeIs) = @$tuple;
493 :     for ($i=0; ($i < @$genomeIs) && ($genomes_info->[$genomeIs->[$i]]->[0] ne $focus); $i++) {}
494 :     if ($i < @$genomeIs)
495 :     {
496 :     if ((! $bestN) || (@$genomeIs < @bestL))
497 :     {
498 :     $bestN = $id;
499 :     @bestL = @$genomeIs;
500 :     }
501 :     }
502 :     }
503 :     $active_subsetR = $bestN;
504 :     }
505 :     }
506 :    
507 :     my %active_genomes = map { $_ => 1 } @bestL;
508 :     return \%active_genomes;
509 :     }
510 :    
511 : overbeek 1.1 1;
512 :    
513 :    

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3