[Bio] / FigKernelScripts / check_model.pl Repository:
ViewVC logotype

Annotation of /FigKernelScripts/check_model.pl

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : formsma 1.1 # _*_ Perl _*_
2 :     #
3 :     #
4 :     #
5 :     use strict;
6 : formsma 1.2 use FIGV;
7 : formsma 1.1 use Scenario;
8 :    
9 : formsma 1.2 my ($genome_id,$orgdir) = @ARGV;
10 : formsma 1.1
11 :     my $debug = 0;
12 :    
13 : formsma 1.2 unless (scalar (@ARGV) >= 1)
14 : formsma 1.1 {
15 :     print STDERR "usage: check_model <genome_id>\n";
16 :     exit(-1);
17 :     }
18 :     chomp $genome_id;
19 : formsma 1.2 my $fig;
20 :     if($orgdir)
21 :     {
22 :     $fig = new FIGV($orgdir);
23 :     }
24 :     else
25 :     {
26 :     $fig = new FIGV;
27 :     }
28 : formsma 1.3 Scenario::set_fig($fig,$genome_id);
29 : formsma 1.1 print STDERR "Running check_model on genome $genome_id\n" if($debug);
30 :     print STDERR "\tStage 1 - Loading Scenarios, model inputs and biomass\n" if($debug);
31 :     my @scenarios = @{Scenario->get_genome_scenarios($genome_id,0)};
32 : formsma 1.2 my $filebase = $fig->model_directory($genome_id)."/Analysis/";
33 : formsma 1.1
34 :     my %scenarios_used;
35 : formsma 1.4 my %scenarios_id;
36 : formsma 1.1 my %all_substrates;
37 :     foreach (@scenarios)
38 :     {
39 :     map {$all_substrates{$_} = 1} @{$_->get_substrates_all()};
40 :     $scenarios_used{$_->get_scenario_name()} = 0;
41 : formsma 1.4 $scenarios_id{$_->get_id()} = 0;
42 : formsma 1.1 }
43 :    
44 :     my %input_cpds;
45 :     open(INPUTS,$filebase."inputs.txt")
46 :     or die ("Failed to open ".$filebase."inputs.txt. Run find_model_information $genome_id inputs");
47 :     while(<INPUTS>)
48 :     {
49 :     chomp;
50 :     my @line = split "\t", $_;
51 :     chomp @line;
52 :     $input_cpds{$line[0]} = 1;
53 :     }
54 :     close(INPUTS);
55 :    
56 : formsma 1.4 #ATP Snythase HACK, lets add ATP as a input
57 :     $input_cpds{"C00002"} = 1;
58 :    
59 : formsma 1.1 my %cofactor_cpds;
60 : formsma 1.4 open(INPUTCO,$filebase."secondary_compounds.txt")
61 :     or die("Failed to open ".$filebase."secondary_compounds.txt. Run find_model_cofactors $genome_id\n");
62 : formsma 1.1 while(<INPUTCO>)
63 :     {
64 :     chomp;
65 :     my @line = split "\t", $_;
66 :     chomp @line;
67 :     $cofactor_cpds{$line[0]} = 1;
68 :     }
69 :     close(INPUTCO);
70 :    
71 :     my %biomass_cpds;
72 :     open(BIOMASS,$filebase."biomass.txt")
73 :     or die ("Failed to open ".$filebase."biomass.txt Run find_model_information $genome_id biomass");
74 :     while(<BIOMASS>)
75 :     {
76 :     chomp;
77 :     my @temp = split "\t", $_;
78 :     chomp @temp;
79 :     $biomass_cpds{$temp[0]} = 1;
80 :     }
81 :     close(BIOMASS);
82 : formsma 1.4
83 :     my %compound_id_to_name;
84 :     open(IN,"$FIG_Config::global/Models/compound_id_to_name.txt") or die("Failed to open compound_id_to_name.txt");
85 :     while($_ = <IN>){
86 :     chomp($_);
87 :     my($id,$name) = split("\t",$_);
88 :     $compound_id_to_name{$id} = $name;
89 :     }
90 :     close(IN);
91 :    
92 : formsma 1.1 print STDERR "\tStage 1 - Complete\n" if($debug);
93 :    
94 :     print STDERR "\tStage 2 - Run Analysis Tools (fueled by caffine)\n" if($debug);
95 :     #Total list of what compounds we reach.
96 :     my @list_generated = keys %input_cpds;
97 : formsma 1.4 push @list_generated , keys %cofactor_cpds;
98 : formsma 1.1 my %possible_problems;
99 :     my %good_compounds;
100 :    
101 :     # Here we are running all the scenarios we have inputs for
102 :     # to analyze what biomass we can reach. We are also identifying
103 :     # any Scenarios that are not being utilized.
104 :    
105 :     my %current_inputs;
106 :     foreach (keys %input_cpds)
107 :     {
108 :     $current_inputs{$_} = 1;
109 :     }
110 :     my %current_cofactors;
111 :     foreach (keys %cofactor_cpds)
112 :     {
113 :     $current_cofactors{$_} = 1;
114 :     }
115 :    
116 :     while(scalar(keys %current_inputs))
117 :     {
118 : formsma 1.4 my ($temp_in,$temp_co,$break) = generate_new_inputs(\%current_inputs,\%current_cofactors);
119 : formsma 1.1
120 : formsma 1.4 map {$current_inputs{$_} = 1;} keys %$temp_in;
121 :     map {$current_cofactors{$_} = 1;} keys %$temp_co;
122 : formsma 1.1
123 :     push @list_generated, keys %$temp_in;
124 :     push @list_generated, keys %$temp_co;
125 : formsma 1.4
126 :     print STDERR "Rechecking with inputs\n" if($debug);
127 :     if($break)
128 :     {
129 :     last;
130 :     }
131 : formsma 1.1 }
132 : formsma 1.4 print STDERR "Finished checking scenarios\n" if($debug);
133 : formsma 1.1 foreach my $found (@list_generated)
134 :     {
135 :     if(defined $biomass_cpds{$found})
136 :     {
137 :     $biomass_cpds{$found} = 2;
138 :     }
139 :     }
140 :    
141 :     # Analyze input compounds
142 :     # and also identify dead-end Scenarios
143 :     my %dead_end_scenarios;
144 :    
145 :     foreach my $scenario (@scenarios)
146 :     {
147 :     my $count = 0;
148 :     my %used_inputs;
149 :     foreach my $substrate (@{$scenario->get_substrates_all()})
150 :     {
151 :     if($input_cpds{$substrate} == 1)
152 :     {
153 :     $used_inputs{$substrate} = 1;
154 :     $count++;
155 :     }
156 :     }
157 :    
158 :     if($count !=0 && $count != scalar(@{$scenario->get_substrates_all()}))
159 :     {
160 :     map {$possible_problems{$_} = 1;} keys %used_inputs;
161 :     }
162 :     else
163 :     {
164 :     map {$good_compounds{$_} = 0;} keys %used_inputs;
165 :     }
166 :    
167 :     #dead-end analysis
168 :     my $bad_end = 1;
169 :     foreach my $product (@{$scenario->get_products_all()})
170 :     {
171 :     if(defined $biomass_cpds{$product} || $all_substrates{$product} == 1)
172 :     {
173 :     $bad_end = 0;
174 :     last;
175 :     }
176 :     }
177 :     #Iff we reach this point, then this scenario produces no biomass
178 :     # or compounds used by another scenario. It's all waste, could be a error.
179 :     if($bad_end)
180 :     {
181 :     $dead_end_scenarios{$scenario->get_scenario_name()} = 1;
182 :     }
183 :     }
184 :    
185 :     open(RESULTS,">".$filebase."info_biomass.txt");
186 :     open(OUT, ">".$filebase."biomass_valid.txt");
187 :     foreach (sort keys %biomass_cpds)
188 :     {
189 :     if($biomass_cpds{$_} != 2)
190 :     {
191 : formsma 1.4 print RESULTS "Biomass $_ $compound_id_to_name{$_} failed to be generated using inputs.txt\n";
192 : formsma 1.1 }
193 :     else
194 :     {
195 : formsma 1.4 print OUT "$_\t$compound_id_to_name{$_}\n";
196 : formsma 1.1 }
197 :     }
198 :     close(RESULTS);
199 :     close(OUT);
200 :    
201 :     open(OUT,">$filebase\scenarios_valid.txt");
202 :     open(SCEN_RESULTS,">".$filebase."info_scenarios.txt");
203 :     foreach my $scen_id (sort keys %scenarios_used)
204 :     {
205 :     if($scenarios_used{$scen_id} == 0)
206 :     {
207 :     print SCEN_RESULTS "Scenario $scen_id unable to run with current inputs\n";
208 :     }
209 :     else
210 :     {
211 :     print OUT "$scen_id\n";
212 :     }
213 :     }
214 :     print SCEN_RESULTS "\n\n#################\n\n\n";
215 :     foreach (sort keys %dead_end_scenarios)
216 :     {
217 :     print SCEN_RESULTS "Scenario $_ is possibly a dead end(no outputs used by other scenarios or biomass.)\n";
218 :     }
219 :     close(SCEN_RESULTS);
220 :     close(OUT);
221 :    
222 : formsma 1.4 #open(RESULTS_INPUT,">".$filebase."info_inputs.txt");
223 :     #foreach (sort keys %possible_problems)
224 :     #{
225 :     # if($possible_problems{$_} == 1 && !defined $good_compounds{$_})
226 :     # {
227 :     # print RESULTS_INPUT "Input $_ may be a false-positive. Perhaps a Scenario should generate this.\n";
228 :     # }
229 :     #}
230 :     #close(RESULTS_INPUT);
231 : formsma 1.1
232 :     print STDERR "\tStage 2 - Complete! Results in info_inputs.txt, info_biomass.txt and info_scenarios.txt\n" if($debug);
233 :    
234 :     sub generate_new_inputs
235 :     {
236 :     my ($inputs,$cofactors) = @_;
237 :     my %new_inputs;
238 :     my %new_cofactors;
239 : formsma 1.4 my $break = 1;
240 : formsma 1.1 foreach my $scenario (@scenarios)
241 :     {
242 : formsma 1.4 if($scenarios_id{$scenario->get_id()} == 1)
243 :     {
244 :     next;
245 :     }
246 : formsma 1.1 my $count = 0;
247 :     foreach my $substrate (@{$scenario->get_substrates_all()})
248 :     {
249 :     if($inputs->{$substrate})
250 :     {
251 :     $count++;
252 :     }
253 : formsma 1.4 elsif($cofactors->{$substrate})
254 : formsma 1.1 {
255 :     $count++;
256 :     }
257 :     }
258 :    
259 :     if($count != scalar(@{$scenario->get_substrates_all()}))
260 :     {
261 :     next;
262 :     }
263 :     else
264 :     {
265 :     print STDERR "Scenario ".$scenario->get_scenario_name()." ran.\n" if($debug);
266 : formsma 1.4 $break = 0;
267 :     $scenarios_id{$scenario->get_id()} = 1;
268 : formsma 1.1 $scenarios_used{$scenario->get_scenario_name()} = 1;
269 :     foreach my $product (@{$scenario->get_products()})
270 :     {
271 :     $new_inputs{$product} = 1;
272 :     }
273 :     foreach my $product (@{$scenario->get_product_cofactors()})
274 :     {
275 :     $new_cofactors{$product} = 1;
276 :     }
277 :     }
278 :     }
279 : formsma 1.4 return (\%new_inputs,\%new_cofactors,$break);
280 : formsma 1.1 }

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3