[Bio] / SubsystemExtension / FigKernelScripts / subsystem_extension.pl Repository:
ViewVC logotype

Annotation of /SubsystemExtension/FigKernelScripts/subsystem_extension.pl

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : heiko 1.1 #!/Users/fig/FIGDisk/env/mac/bin/perl
2 :    
3 :     use strict;
4 :     use warnings;
5 :    
6 :     use Getopt::Std;
7 :     use Term::ReadKey;
8 :     use IO::Handle;
9 :     use Data::Dumper;
10 :     use FIG;
11 :     use Subsystem;
12 :     use SubsystemExtension::AEOS;
13 :     use SubsystemExtension::ValidationInterface;
14 :     use SubsystemExtension::SubsystemCandidateFactory::SEED;
15 :     use constant DEBUG => 1;
16 :     no warnings qw(redefine);
17 :    
18 :     #
19 :     # this is necessary if the script is started via rsh(1)
20 :     # otherwise you won't see any output until the first <RETURN>
21 :     #
22 :     STDOUT->autoflush(1);
23 :    
24 :     sub usage {
25 :     print "\t -= subsystem_extension =-\n\n";
26 :     print "This script is used to extend subsystems (either one or a list in a file) to one ore more genomes.\n";
27 :     print "If neither the parameter -t nor -m is defined all genomes not yet present in the subsystem (excluding environmental samples) will be used\n\n";
28 :    
29 :     print "usage:\n subsystem_extension -s <subsystem name> -f <subsystem file> [-t <taxon_id> -m <taxon_id file> -a] -d <validation_directory> -e <e-value cutoff: e.g. 1.0e-10> -r <number of organisms to search in> [-x use AEOS]\n\n";
30 :    
31 :     }
32 :    
33 :     # global variables
34 :     our($opt_s, $opt_t, $opt_f, $opt_d, $opt_m, $opt_a, $opt_x, $opt_e, $opt_r);
35 :    
36 :     getopts('s:t:f:d:m:e:r:ax');
37 :    
38 :    
39 :     if ((!$opt_s) && (!$opt_f)) {
40 :     usage;
41 :     print "ERROR: Can't start script : no subsystem/filename given!\n";
42 :     exit 1;
43 :     }
44 :    
45 :    
46 :    
47 :    
48 :     my $fig = new FIG;
49 :     my @subsystems;
50 :    
51 :     push @subsystems, $opt_s if ($opt_s);
52 :    
53 :    
54 :     ## collect the subsystem names from a flat file
55 :     # each line holds one subsystem name
56 :    
57 :     if ($opt_f && -e $opt_f) {
58 :     open SSF, "<$opt_f";
59 :     while (<SSF>) {
60 :     chomp;
61 :     push @subsystems, $_;
62 :     }
63 :     close SSF;
64 :     }
65 :    
66 :     my $html = [];
67 :    
68 :    
69 :     my @missing_genomes;
70 :    
71 :     # add single or multiple taxon ids (parsed from file) to
72 :     # the list of missing genomes
73 :    
74 :     if ($opt_t) {
75 :     push @missing_genomes, $opt_t;
76 :     } elsif ($opt_m && -e $opt_m) {
77 :     open TF, "<$opt_m";
78 :     while (<TF>) {
79 :     chomp;
80 :     push @missing_genomes, $_ if $_ =~ /\d+\.\d+/;
81 :     }
82 :     close TF;
83 :    
84 :     } elsif ($opt_a) {
85 :     @missing_genomes = grep { $_ !~ /^99999/ } $fig->genomes("complete");
86 :     }
87 :    
88 :    
89 :    
90 :     foreach (@subsystems) {
91 :    
92 :     my $subsystem = $fig->get_subsystem($_);
93 :    
94 :     if (!$subsystem) {
95 :     print STDERR "could not initialize the subsystem $opt_s\n";
96 :     next;
97 :     }
98 :    
99 :     my $extension = SubsystemExtension::SubsystemCandidateFactory::SEED->new($subsystem, "master", $opt_e, $opt_r);
100 :    
101 :    
102 :     if (ref $extension && $extension->isa("SubsystemExtension::SubsystemCandidateFactory")) {
103 :    
104 :     my %in = map { $_ => 1 } $subsystem->get_genomes;
105 :    
106 :    
107 :     unless (scalar @missing_genomes > 0) {
108 :     print STDERR "getting missing genomes for subsystem ".$subsystem->get_name()."\n" if (DEBUG);
109 :     @missing_genomes = grep { ! $in{$_} } grep { $_ !~ /^99999/ } $fig->genomes("complete");
110 :     }
111 :    
112 :    
113 :     foreach my $genome (@missing_genomes) {
114 :     my $subsystem_candidates = $extension->generate([$genome]);
115 :    
116 :     foreach (keys %$subsystem_candidates) {
117 :     my $subsystem_candidate;
118 :     # generate the aeos and extend the subsystem!!
119 :     if ($opt_x) {
120 :     my $aeos = SubsystemExtension::AEOS->new($fig, $subsystem, $genome, $subsystem_candidates->{$_}, {mute=>1});
121 :     $subsystem_candidate = $aeos->extend();
122 :     } else {
123 :     $subsystem_candidate = $subsystem_candidates->{$_};
124 :     }
125 :    
126 :     my $subsystem_name = $subsystem->get_name();
127 :     $subsystem_name =~ s/\(/\\\(/;
128 :     $subsystem_name =~ s/\)/\\\)/;
129 :     print STDERR $subsystem_name."\n\n\n";
130 :    
131 :    
132 :     mkdir $opt_d."/".$subsystem->get_name() unless (-d $opt_d."/".$subsystem->get_name());
133 :    
134 :     $subsystem_candidate->toFile($opt_d."/".$subsystem->get_name()."/".$_.".scs") if (-d $opt_d);
135 :    
136 :     if ($subsystem_candidate->functional_variant() > 0) {
137 :     my $validation_interface = SubsystemExtension::ValidationInterface->new($fig, $subsystem, $_, $subsystem_candidate);
138 :     if (-d $opt_d) {
139 :    
140 :     my $filename = $opt_d."/".$subsystem->get_name()."/".$_;
141 :     $filename .= "_complete" if ($subsystem_candidate->functional_variant_score() >= 0.9);
142 :     $filename .= ".html";
143 :     open VAL, ">$filename";
144 :     print VAL $validation_interface->html_output();
145 :     close VAL;
146 :     # $subsystem_candidate->toFile($opt_d."/".$subsystem->get_name()."/".$_.".scs");
147 :     } else {
148 :     print $validation_interface->html_output();
149 :     }
150 :     }
151 :     }
152 :     }
153 :    
154 :     }
155 :     }
156 :    
157 :    

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3