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

Annotation of /FigKernelPackages/FIGMODELmodel.pm

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : chenry 1.2 package FIGMODELmodel;
2 :     use strict;
3 :    
4 :     =head1 FIGMODELmodel object
5 :     =head2 Introduction
6 :     Module for manipulating model objects.
7 :     =head2 Core Object Methods
8 :    
9 :     =head3 new
10 :     Definition:
11 :     FIGMODELmodel = FIGMODELmodel->new();
12 :     Description:
13 :     This is the constructor for the FIGMODELmodel object.
14 :     =cut
15 :     sub new {
16 :     my ($class,$figmodel,$id) = @_;
17 :     #Error checking first
18 :     if (!defined($figmodel)) {
19 :     print STDERR "FIGMODELmodel->new(undef,".$id."):figmodel must be defined to create a model object!\n";
20 :     return undef;
21 :     }
22 :     my $self = {_figmodel => $figmodel};
23 :     bless $self;
24 :     if (!defined($id)) {
25 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,undef):id must be defined to create a model object");
26 :     return undef;
27 :     }
28 :    
29 :     #Checking that the id exists
30 :     my $tbl = $self->figmodel()->database()->GetDBTable("MODELS");
31 :     if (!defined($tbl)) {
32 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not load MODELS table. Check database!");
33 :     return undef;
34 :     }
35 :    
36 :     #If the id is a number, we get the model row by index
37 :     my $index = $id;
38 :     if ($id =~ m/^\d+$/) {
39 :     $self->{_data} = $tbl->get_row($id);
40 :     } else {
41 :     $self->{_data} = $tbl->get_row_by_key($id,"id");
42 :     if (defined($self->{_data})) {
43 :     $index = $tbl->row_index($self->{_data});
44 :     }
45 :     }
46 :     if (!defined($self->{_data})) {
47 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find specified id in database!");
48 :     return undef;
49 :     }
50 :     $self->{_index} = $index;
51 :     $self->figmodel()->{_models}->{$self->id()} = $self;
52 :     $self->figmodel()->{_models}->{$self->index()} = $self;
53 :    
54 :     return $self;
55 :     }
56 :    
57 :     =head3 config
58 :     Definition:
59 :     ref::key value = FIGMODELmodel->config(string::key);
60 :     Description:
61 :     Trying to avoid using calls that assume configuration data is stored in a particular manner.
62 :     Call this function to get file paths etc.
63 :     =cut
64 :     sub config {
65 :     my ($self,$key) = @_;
66 :     return $self->figmodel()->config($key);
67 :     }
68 :    
69 :     =head3 fail
70 :     Definition:
71 :     -1 = FIGMODELmodel->fail();
72 :     Description:
73 :     Standard return for failed functions.
74 :     =cut
75 :     sub fail {
76 :     my ($self) = @_;
77 :     return $self->figmodel()->fail();
78 :     }
79 :    
80 :     =head3 success
81 :     Definition:
82 :     1 = FIGMODELmodel->success();
83 :     Description:
84 :     Standard return for successful functions.
85 :     =cut
86 :     sub success {
87 :     my ($self) = @_;
88 :     return $self->figmodel()->success();
89 :     }
90 :    
91 :     =head3 figmodel
92 :     Definition:
93 :     FIGMODEL = FIGMODELmodel->figmodel();
94 :     Description:
95 :     Returns a FIGMODEL object
96 :     =cut
97 :     sub figmodel {
98 :     my ($self) = @_;
99 :     return $self->{"_figmodel"};
100 :     }
101 :    
102 :     =head3 fig
103 :     Definition:
104 :     FIGMODEL = FIGMODELmodel->fig();
105 :     Description:
106 :     Returns a FIG object
107 :     =cut
108 :     sub fig {
109 :     my ($self) = @_;
110 :    
111 :     if (!defined($self->{_fig}) && $self->source() !~ /^MGRAST/) {
112 :     if ($self->source() =~ /^RAST/) {
113 :     $self->{"_fig"} = $self->figmodel()->fig();#$self->genome());
114 :     } else {
115 :     $self->{"_fig"} = $self->figmodel()->fig();
116 :     }
117 :     }
118 :    
119 :     return $self->{"_fig"};
120 :     }
121 :    
122 :     =head3 mgdata
123 :     Definition:
124 :     FIGMODEL = FIGMODELmodel->mgdata();
125 :     Description:
126 :     Returns a mgrast database object
127 :     =cut
128 :     sub mgdata {
129 :     my ($self) = @_;
130 :    
131 :     if (!defined($self->{_mgdata}) && $self->source() =~ /^MGRAST/) {
132 :     require MGRAST;
133 :     $self->{_mgdata} = $self->figmodel()->mgrast()->Job->get_objects( { 'genome_id' => $self->genome() } )
134 :     }
135 :    
136 :     return $self->{_mgdata};
137 :     }
138 :    
139 :     =head3 id
140 :     Definition:
141 :     string = FIGMODELmodel->id();
142 :     Description:
143 :     Returns model id
144 :     =cut
145 :     sub id {
146 :     my ($self) = @_;
147 :     return $self->{_data}->{id}->[0];
148 :     }
149 :    
150 :     =head3 owner
151 :     Definition:
152 :     string = FIGMODELmodel->owner();
153 :     Description:
154 :     Returns the username for the model owner
155 :     =cut
156 :     sub owner {
157 :     my ($self) = @_;
158 :     return $self->{_data}->{owner}->[0];
159 :     }
160 :    
161 :     =head3 index
162 :     Definition:
163 :     string = FIGMODELmodel->index();
164 :     Description:
165 :     Returns model index
166 :     =cut
167 :     sub index {
168 :     my ($self) = @_;
169 :     return $self->{_index};
170 :     }
171 :    
172 :     =head3 name
173 :     Definition:
174 :     string = FIGMODELmodel->name();
175 :     Description:
176 :     Returns the name of the organism or metagenome sample being modeled
177 :     =cut
178 :     sub name {
179 :     my ($self) = @_;
180 :    
181 :     if (!defined($self->{_name})) {
182 :     my $source = $self->source();
183 :     if ($source =~ /^MGRAST/) {
184 :     $self->{_name} = "NA";
185 :     if (defined($self->mgdata())) {
186 :     $self->{_name} = $self->mgdata()->genome_name;
187 :     }
188 :     } elsif (defined($self->stats())) {
189 :     $self->{_name} = $self->stats()->{'Organism name'}->[0];
190 :     } elsif ($source !~ /^RAST/) {
191 :     $self->{_name} = $self->fig()->orgname_of_orgid($self->genome());
192 :     }
193 :     }
194 :    
195 :     return $self->{_name};
196 :     }
197 :    
198 :     =head3 stats
199 :     Definition:
200 :     string = FIGMODELmodel->stats();
201 :     Description:
202 :     Returns model stats
203 :     =cut
204 :     sub stats {
205 :     my ($self) = @_;
206 :    
207 :     if (!defined($self->{_stats})) {
208 :     $self->{_stats} = $self->figmodel()->database()->GetDBTable("MODEL STATS")->get_row_by_key($self->id(),"Model ID");
209 :     }
210 :     return $self->{_stats};
211 :     }
212 :    
213 :     =head3 get_reaction_class
214 :     Definition:
215 :     string = FIGMODELmodel->get_reaction_class(string::reaction ID);
216 :     Description:
217 :     Returns reaction class
218 :     =cut
219 :     sub get_reaction_class {
220 :     my ($self,$reaction,$nohtml) = @_;
221 :    
222 :     if (!-e $self->directory()."ReactionClassification-".$self->id().".tbl") {
223 :     if (!defined($self->{_reaction_classes})) {
224 :     $self->{_reaction_classes} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id()."-Complete.tbl",";","|",0,["REACTION"]);
225 :     if (!defined($self->{_reaction_classes})) {
226 :     return undef;
227 :     }
228 :     }
229 :    
230 :     my $ClassRow = $self->{_reaction_classes}->get_row_by_key($reaction,"REACTION");
231 :     if (defined($ClassRow) && defined($ClassRow->{CLASS})) {
232 :     my $class;
233 :     if ($ClassRow->{CLASS}->[0] eq "Positive") {
234 :     $class = "Essential =>";
235 :     } elsif ($ClassRow->{CLASS}->[0] eq "Negative") {
236 :     $class = "Essential <=";
237 :     } elsif ($ClassRow->{CLASS}->[0] eq "Positive variable") {
238 :     $class = "Active =>";
239 :     } elsif ($ClassRow->{CLASS}->[0] eq "Negative variable") {
240 :     $class = "Active <=";
241 :     } elsif ($ClassRow->{CLASS}->[0] eq "Variable") {
242 :     $class = "Active <=>";
243 :     } elsif ($ClassRow->{CLASS}->[0] eq "Blocked") {
244 :     $class = "Inactive";
245 :     } elsif ($ClassRow->{CLASS}->[0] eq "Dead") {
246 :     $class = "Disconnected";
247 :     }
248 :     if (!defined($nohtml) || $nohtml ne "1") {
249 :     $class = "<span title=\"Flux:".$ClassRow->{MIN}->[0]." to ".$ClassRow->{MAX}->[0]."\">".$class."</span>";
250 :     }
251 :     return $class;
252 :     }
253 :     return undef;
254 :     }
255 :    
256 :     if (!defined($self->{_reaction_classes})) {
257 :     $self->{_reaction_classes} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id().".tbl",";","|",0,["REACTION"]);
258 :     if (!defined($self->{_reaction_classes})) {
259 :     return undef;
260 :     }
261 :     }
262 :    
263 :     my $ClassRow = $self->{_reaction_classes}->get_row_by_key($reaction,"REACTION");
264 :     my $classstring = "";
265 :     if (defined($ClassRow) && defined($ClassRow->{CLASS})) {
266 :     for (my $i=0; $i < @{$ClassRow->{CLASS}};$i++) {
267 :     if (length($classstring) > 0) {
268 :     $classstring .= "<br>";
269 :     }
270 :     my $NewClass;
271 :     if ($ClassRow->{CLASS}->[$i] eq "Positive") {
272 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Essential =>";
273 :     if (!defined($nohtml) || $nohtml ne "1") {
274 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
275 :     }
276 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Negative") {
277 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Essential <=";
278 :     if (!defined($nohtml) || $nohtml ne "1") {
279 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
280 :     }
281 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Positive variable") {
282 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Active =>";
283 :     if (!defined($nohtml) || $nohtml ne "1") {
284 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
285 :     }
286 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Negative variable") {
287 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=";
288 :     if (!defined($nohtml) || $nohtml ne "1") {
289 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
290 :     }
291 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Variable") {
292 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=>";
293 :     if (!defined($nohtml) || $nohtml ne "1") {
294 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
295 :     }
296 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Blocked") {
297 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Inactive";
298 :     if (!defined($nohtml) || $nohtml ne "1") {
299 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
300 :     }
301 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Dead") {
302 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Disconnected";
303 :     if (!defined($nohtml) || $nohtml ne "1") {
304 :     $NewClass = "<span title=\"Flux:".$ClassRow->{MIN}->[$i]." to ".$ClassRow->{MAX}->[$i]."\">".$NewClass."</span>";
305 :     }
306 :     }
307 :     $classstring .= $NewClass;
308 :     }
309 :     }
310 :     return $classstring;
311 :     }
312 :    
313 :     =head3 get_biomass
314 :     Definition:
315 :     string = FIGMODELmodel->get_biomass();
316 :     Description:
317 :     Returns data for the biomass reaction
318 :     =cut
319 :     sub get_biomass {
320 :     my ($self) = @_;
321 :    
322 :     if (!defined($self->{_biomass})) {
323 :     my $rxntbl = $self->reaction_table();
324 :     for (my $i=0; $i < $rxntbl->size(); $i++) {
325 :     if ($rxntbl->get_row($i)->{"LOAD"}->[0] =~ m/bio\d\d\d\d\d/) {
326 :     $self->{_biomass} = $rxntbl->get_row($i)->{"LOAD"}->[0];
327 :     last;
328 :     }
329 :     }
330 :     }
331 :    
332 :     return $self->get_reaction_data($self->{_biomass});
333 :     }
334 :    
335 :     =head3 get_reaction_data
336 :     Definition:
337 :     string = FIGMODELmodel->get_reaction_data(string::reaction ID);
338 :     Description:
339 :     Returns model reaction data
340 :     =cut
341 :     sub get_reaction_data {
342 :     my ($self,$reaction) = @_;
343 :    
344 :     if (!defined($self->reaction_table())) {
345 :     return undef;
346 :     }
347 :     if ($reaction =~ m/^\d+$/) {
348 :     $self->reaction_table()->get_row($reaction);
349 :     }
350 :     return $self->reaction_table()->get_row_by_key($reaction,"LOAD");
351 :     }
352 :    
353 :     =head3 get_reaction_id
354 :     Definition:
355 :     string = FIGMODELmodel->get_reaction_id(string::reaction ID);
356 :     Description:
357 :     Returns model reaction id
358 :     =cut
359 :     sub get_reaction_id {
360 :     my ($self,$reaction) = @_;
361 :    
362 :     my $Data = $self->get_reaction_data($reaction);
363 :     if (!defined($Data) || !defined($Data->{LOAD}->[0])) {
364 :     return undef;
365 :     }
366 :     return $Data->{LOAD}->[0];
367 :     }
368 :    
369 :     =head3 get_reaction_number
370 :     Definition:
371 :     int = FIGMODELmodel->get_reaction_number();
372 :     Description:
373 :     Returns the number of reactions in the model
374 :     =cut
375 :     sub get_reaction_number {
376 :     my ($self) = @_;
377 :    
378 :     if (!defined($self->reaction_table())) {
379 :     return 0;
380 :     }
381 :    
382 :     return $self->reaction_table()->size();
383 :     }
384 :    
385 :     =head3 reaction_table
386 :     Definition:
387 :     FIGMODELTable = FIGMODELmodel->reaction_table();
388 :     Description:
389 :     Returns FIGMODELTable with the reaction list for the model
390 :     =cut
391 :     sub reaction_table {
392 :     my ($self) = @_;
393 :    
394 :     if (!defined($self->{_reaction_data})) {
395 :     $self->{_reaction_data} = $self->figmodel()->database()->GetDBModel($self->id());
396 :     }
397 :    
398 :     return $self->{_reaction_data};
399 :     }
400 :    
401 :     =head3 reaction_class_table
402 :     Definition:
403 :     FIGMODELTable = FIGMODELmodel->reaction_class_table();
404 :     Description:
405 :     Returns FIGMODELTable with the reaction class data, and creates the table file if it does not exist
406 :     =cut
407 :     sub reaction_class_table {
408 :     my ($self) = @_;
409 :    
410 :     if (!defined($self->{_reaction_class_table})) {
411 :     if (-e $self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl") {
412 :     $self->{_reaction_class_table} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl",";","|",0,["REACTION","CLASS","MEDIA"]);
413 :     } else {
414 :     $self->{_reaction_class_table} = FIGMODELTable->new(["REACTION","MEDIA","CLASS","MIN","MAX"],$self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl",["REACTION","CLASS","MEDIA"],";","|",undef);
415 :     }
416 :     }
417 :    
418 :     return $self->{_reaction_class_table};
419 :     }
420 :    
421 :     =head3 compound_class_table
422 :     Definition:
423 :     FIGMODELTable = FIGMODELmodel->compound_class_table();
424 :     Description:
425 :     Returns FIGMODELTable with the compound class data, and creates the table file if it does not exist
426 :     =cut
427 :     sub compound_class_table {
428 :     my ($self) = @_;
429 :    
430 :     if (!defined($self->{_compound_class_table})) {
431 :     if (-e $self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl") {
432 :     $self->{_compound_class_table} = $self->figmodel()->database()->load_table($self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl",";","|",0,["COMPOUND","CLASS","MEDIA"]);
433 :     } else {
434 :     $self->{_compound_class_table} = FIGMODELTable->new(["COMPOUND","MEDIA","CLASS","MIN","MAX"],$self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl",["COMPOUND","CLASS","MEDIA"],";","|",undef);
435 :     }
436 :     }
437 :    
438 :     return $self->{_compound_class_table};
439 :     }
440 :    
441 :     =head3 compound_table
442 :     Definition:
443 :     FIGMODELTable = FIGMODELmodel->compound_table();
444 :     Description:
445 :     Returns FIGMODELTable with the compound list for the model
446 :     =cut
447 :     sub compound_table {
448 :     my ($self) = @_;
449 :    
450 :     if (!defined($self->{_compound_table})) {
451 :     $self->{_compound_table} = $self->figmodel()->database()->GetDBModelCompounds($self->id());
452 :     }
453 :    
454 :     return $self->{_compound_table};
455 :     }
456 :    
457 :     =head3 get_compound_data
458 :     Definition:
459 :     string = FIGMODELmodel->get_compound_data(string::compound ID);
460 :     Description:
461 :     Returns model compound data
462 :     =cut
463 :     sub get_compound_data {
464 :     my ($self,$compound) = @_;
465 :    
466 :     if (!defined($self->compound_table())) {
467 :     return undef;
468 :     }
469 :     if ($compound =~ m/^\d+$/) {
470 :     $self->compound_table()->get_row($compound);
471 :     }
472 :     return $self->compound_table()->get_row_by_key($compound,"DATABASE");
473 :     }
474 :    
475 :     =head3 genome
476 :     Definition:
477 :     string = FIGMODELmodel->genome();
478 :     Description:
479 :     Returns model genome
480 :     =cut
481 :     sub genome {
482 :     my ($self) = @_;
483 :     return $self->{_data}->{genome}->[0];
484 :     }
485 :    
486 :     =head3 source
487 :     Definition:
488 :     string = FIGMODELmodel->source();
489 :     Description:
490 :     Returns model source
491 :     =cut
492 :     sub source {
493 :     my ($self) = @_;
494 :     return $self->{_data}->{source}->[0];
495 :     }
496 :    
497 :     =head3 rights
498 :     Definition:
499 :     1/0 = FIGMODELmodel->rights(string::username);
500 :     Description:
501 :     Returns 1 if the input user can view the model, and zero otherwise
502 :     =cut
503 :     sub rights {
504 :     my ($self,$username) = @_;
505 :    
506 :     if ($self->public()) {
507 :     return 1;
508 :     }
509 :     if (!defined($username) || $username eq "NONE") {
510 :     return 0;
511 :     }
512 :     if (!defined($self->{_userrights}->{$username})) {
513 :     if (defined($self->{_data}->{master})) {
514 :     for (my $i=0; $i < @{$self->{_data}->{master}};$i++) {
515 :     $self->{_userrights}->{$self->{_data}->{master}->[$i]} = 1;
516 :     }
517 :     }
518 :     if (defined($self->{_data}->{users})) {
519 :     for (my $i=0; $i < @{$self->{_data}->{users}};$i++) {
520 :     $self->{_userrights}->{$self->{_data}->{users}->[$i]} = 1;
521 :     }
522 :     }
523 :     }
524 :     return $self->{_userrights}->{$username};
525 :     }
526 :    
527 :     =head3 public
528 :     Definition:
529 :     1/0 = FIGMODELmodel->public();
530 :     Description:
531 :     Returns 1 if the model is public, and zero otherwise
532 :     =cut
533 :     sub public {
534 :     my ($self) = @_;
535 :    
536 :     if (!defined($self->{_public})) {
537 :     $self->{_public} = 0;
538 :     if (defined($self->{_data}->{users}->[0]) && $self->{_data}->{users}->[0] eq "all") {
539 :     $self->{_public} = 1;
540 :     }
541 :     }
542 :     return $self->{_public};
543 :     }
544 :    
545 :     =head3 directory
546 :     Definition:
547 :     string = FIGMODELmodel->directory();
548 :     Description:
549 :     Returns model directory
550 :     =cut
551 :     sub directory {
552 :     my ($self) = @_;
553 :    
554 :     if (!defined($self->{_directory})) {
555 :     my $userdirectory = "";
556 :     if ($self->owner() ne "master") {
557 :     $userdirectory = $self->owner()."/";
558 :     }
559 :     my $source = $self->source();
560 :     if ($source =~ /^MGRAST/) {
561 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
562 :     } elsif ($source =~ /^RAST/) {
563 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
564 :     } elsif ($source =~ /^SEED/) {
565 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
566 :     } elsif ($source =~ /^PM/) {
567 :     if (length($userdirectory) == 0) {
568 :     $self->{_directory} = $self->figmodel()->config("imported model directory")->[0].$self->id()."/";
569 :     } else {
570 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->id()."/";
571 :     }
572 :     }
573 :     }
574 :    
575 :     return $self->{_directory};
576 :     }
577 :    
578 :     =head3 filename
579 :     Definition:
580 :     string = FIGMODELmodel->filename();
581 :     Description:
582 :     Returns model filename
583 :     =cut
584 :     sub filename {
585 :     my ($self) = @_;
586 :    
587 :     return $self->directory().$self->id().$self->selected_version().".txt";
588 :     }
589 :    
590 :     =head3 set_metagenome_stats
591 :     Definition:
592 :     string = FIGMODELmodel->set_metagenome_stats();
593 :     Description:
594 :     Sets the values of many model stats for a metagenome
595 :     =cut
596 :     sub set_metagenome_stats {
597 :     my ($self) = @_;
598 :    
599 :     $self->{_total_compounds} = 0;
600 :     if (defined($self->compound_table())) {
601 :     $self->{_total_compounds} = $self->compound_table()->size();
602 :     }
603 :     $self->{_gene_reactions} = 0;
604 :     $self->{_gapfilling_reactions} = 0;
605 :     $self->{_model_genes} = 0;
606 :     $self->{_total_reactions} = 0;
607 :     if (defined($self->reaction_table())) {
608 :     $self->{_total_reactions} = $self->reaction_table()->size();
609 :     my $tbl = $self->reaction_table();
610 :     my $spontaneous = 0;
611 :     for (my $i=0; $i < $tbl->size(); $i++) {
612 :     my $row = $tbl->get_row($i);
613 :     if (!defined($row->{"ASSOCIATED PEG"}->[0]) || $row->{"ASSOCIATED PEG"}->[0] !~ m/peg/) {
614 :     if ($row->{"ASSOCIATED PEG"}->[0] =~ m/SPONTANEOUS/) {
615 :     $spontaneous++;
616 :     } else {
617 :     $self->{_gapfilling_reactions}++;
618 :     }
619 :     } else {
620 :     for (my $j=0; $j < @{$row->{"CONFIDENCE"}}; $j++) {
621 :     my @ecores = split(/;/,$row->{"CONFIDENCE"}->[$j]);
622 :     $self->{_model_genes} += @ecores;
623 :     }
624 :     }
625 :     }
626 :     $self->{_gene_reactions} = $tbl->size() - $spontaneous - $self->{_gapfilling_reactions};
627 :     }
628 :     }
629 :    
630 :     =head3 version
631 :     Definition:
632 :     string = FIGMODELmodel->version();
633 :     Description:
634 :     Returns the version of the model
635 :     =cut
636 :     sub version {
637 :     my ($self) = @_;
638 :    
639 :     if (!defined($self->{_version})) {
640 :     if (!defined($self->{_selectedversion})) {
641 :     if (defined($self->stats())) {
642 :     $self->{_version} = "V".$self->stats()->{"Version"}->[0].".".$self->stats()->{"Gap fill version"}->[0];
643 :     }
644 :     } else {
645 :     $self->{_version} = $self->{_selectedversion};
646 :     }
647 :     }
648 :     return $self->{_version};
649 :     }
650 :    
651 :     =head3 selected_version
652 :     Definition:
653 :     string = FIGMODELmodel->selected_version();
654 :     Description:
655 :     Returns the selected version of the model
656 :     =cut
657 :     sub selected_version {
658 :     my ($self) = @_;
659 :    
660 :     if (!defined($self->{_selectedversion})) {
661 :     return "";
662 :     }
663 :     return $self->{_selectedversion};
664 :     }
665 :    
666 :     =head3 modification_time
667 :     Definition:
668 :     string = FIGMODELmodel->modification_time();
669 :     Description:
670 :     Returns the selected version of the model
671 :     =cut
672 :     sub modification_time {
673 :     my ($self) = @_;
674 :     if (!defined($self->{_modification_time})) {
675 :     my $stats = $self->stats();
676 :     if (defined($stats)) {
677 :     $self->{_modification_time} = 0;
678 :     if (defined($stats->{"Build date"}->[0]) && $self->{_modification_time} < $stats->{"Build date"}->[0]) {
679 :     $self->{_modification_time} = $stats->{"Build date"}->[0];
680 :     } elsif (defined($stats->{"Gap fill date"}->[0]) && $self->{_modification_time} < $stats->{"Gap fill date"}->[0]) {
681 :     $self->{_modification_time} = $stats->{"Gap fill date"}->[0];
682 :     }
683 :     } else {
684 :     $self->{_modification_time} = 0;
685 :     }
686 :     }
687 :     return $self->{_modification_time};
688 :     }
689 :    
690 :     =head3 gene_reactions
691 :     Definition:
692 :     string = FIGMODELmodel->gene_reactions();
693 :     Description:
694 :     Returns the number of reactions added by the gap filling
695 :     =cut
696 :     sub gene_reactions {
697 :     my ($self) = @_;
698 :    
699 :     if (!defined($self->{_gene_reactions})) {
700 :     if ($self->source() =~ /^MGRAST/) {
701 :     $self->set_metagenome_stats();
702 :     } elsif (defined($self->stats())) {
703 :     $self->{_gene_reactions} = $self->total_reactions() - $self->gapfilling_reactions() - $self->stats()->{'Spontaneous'}->[0] - $self->stats()->{'Growmatch reactions'}->[0] - $self->stats()->{'Biolog gap filling reactions'}->[0];
704 :     }
705 :     }
706 :     return $self->{_gene_reactions};
707 :     }
708 :    
709 :     =head3 total_compounds
710 :     Definition:
711 :     string = FIGMODELmodel->total_compounds();
712 :     Description:
713 :     Returns the number of compounds in the model
714 :     =cut
715 :     sub total_compounds {
716 :     my ($self) = @_;
717 :    
718 :     if (!defined($self->{_total_compounds})) {
719 :     if ($self->source() =~ /^MGRAST/) {
720 :     $self->set_metagenome_stats();
721 :     } elsif (defined($self->stats())) {
722 :     $self->{_total_compounds} = $self->stats()->{'Metabolites'}->[0];
723 :     }
724 :     }
725 :     return $self->{_total_compounds};
726 :     }
727 :    
728 :     =head3 gapfilling_reactions
729 :     Definition:
730 :     string = FIGMODELmodel->gapfilling_reactions();
731 :     Description:
732 :     Returns the number of reactions added by the gap filling
733 :     =cut
734 :     sub gapfilling_reactions {
735 :     my ($self) = @_;
736 :    
737 :     if (!defined($self->{_gapfilling_reactions})) {
738 :     if ($self->source() =~ /^MGRAST/) {
739 :     $self->set_metagenome_stats();
740 :     } elsif (defined($self->stats())) {
741 :     $self->{_gapfilling_reactions} = $self->stats()->{'Gap filling reactions'}->[0];
742 :     }
743 :     }
744 :     return $self->{_gapfilling_reactions};
745 :     }
746 :    
747 :     =head3 total_reactions
748 :     Definition:
749 :     string = FIGMODELmodel->total_reactions();
750 :     Description:
751 :     Returns the total number of reactions in the model
752 :     =cut
753 :     sub total_reactions {
754 :     my ($self) = @_;
755 :    
756 :     if (!defined($self->{_total_reactions})) {
757 :     if ($self->source() =~ /^MGRAST/) {
758 :     $self->set_metagenome_stats();
759 :     } elsif (defined($self->stats())) {
760 :     $self->{_total_reactions} = $self->stats()->{'Number of reactions'}->[0];
761 :     }
762 :     }
763 :     return $self->{_total_reactions};
764 :     }
765 :    
766 :     =head3 model_genes
767 :     Definition:
768 :     string = FIGMODELmodel->model_genes();
769 :     Description:
770 :     Returns the number of genes mapped to one or more reactions in the model
771 :     =cut
772 :     sub model_genes {
773 :     my ($self) = @_;
774 :    
775 :     if (!defined($self->{_model_genes})) {
776 :     if ($self->source() =~ /^MGRAST/) {
777 :     $self->set_metagenome_stats();
778 :     } elsif (defined($self->stats())) {
779 :     $self->{_model_genes} = $self->stats()->{'Genes with reactions'}->[0];
780 :     }
781 :     }
782 :     return $self->{_model_genes};
783 :     }
784 :    
785 :     =head3 class
786 :     Definition:
787 :     string = FIGMODELmodel->class();
788 :     Description:
789 :     Returns the class of the model: gram positive, gram negative, other
790 :     =cut
791 :     sub class {
792 :     my ($self) = @_;
793 :    
794 :     if (!defined($self->{_class})) {
795 :     if ($self->source() =~ /^MGRAST/) {
796 :     $self->{_class} = "Other";
797 :     } elsif (defined($self->stats())) {
798 :     $self->{_class} = $self->stats()->{Class}->[0];
799 :     }
800 :     }
801 :     return $self->{_class};
802 :     }
803 :    
804 :     =head3 taxonomy
805 :     Definition:
806 :     string = FIGMODELmodel->taxonomy();
807 :     Description:
808 :     Returns model taxonomy or biome if this is an metagenome model
809 :     =cut
810 :     sub taxonomy {
811 :     my ($self) = @_;
812 :    
813 :     if (!defined($self->{_taxonomy})) {
814 :     my $source = $self->source();
815 :     if ($source =~ /^MGRAST/) {
816 :     $self->{_taxonomy} = "NA";
817 :     my $mgmodeldata = $self->figmodel()->database()->mg_model_data($self->genome());
818 :     if (defined($mgmodeldata)) {
819 :     $self->{_taxonomy} = $mgmodeldata->biome();
820 :     }
821 :     } else {
822 :     $self->{_taxonomy} = $self->fig()->taxonomy_of($self->genome());
823 :     }
824 :     }
825 :    
826 :     return $self->{_taxonomy};
827 :     }
828 :    
829 :     =head3 genome_size
830 :     Definition:
831 :     string = FIGMODELmodel->genome_size();
832 :     Description:
833 :     Returns size of the modeled genome in KB
834 :     =cut
835 :     sub genome_size {
836 :     my ($self) = @_;
837 :    
838 :     if (!defined($self->{_genome_size})) {
839 :     my $source = $self->source();
840 :     if ($source =~ /^MGRAST/) {
841 :     $self->{_genome_size} = "NA";
842 :     if (defined($self->mgdata())) {
843 :     $self->{_genome_size} = $self->mgdata()->size;
844 :     }
845 :     } else {
846 :     $self->{_genome_size} = $self->fig()->genome_szdna($self->genome());
847 :     }
848 :     }
849 :    
850 :     return $self->{_genome_size};
851 :     }
852 :    
853 :     =head3 genome_genes
854 :     Definition:
855 :     string = FIGMODELmodel->genome_genes();
856 :     Description:
857 :     Returns the number of genes in the modeled genome
858 :     =cut
859 :     sub genome_genes {
860 :     my ($self) = @_;
861 :    
862 :     if (!defined($self->{_genome_genes})) {
863 :     my $source = $self->source();
864 :     if ($source =~ /^MGRAST/) {
865 :     $self->{_genome_genes} = "NA";
866 :     if (defined($self->mgdata())) {
867 :     $self->{_genome_genes} = $self->mgdata()->genome_contig_count;
868 :     }
869 :     } elsif (defined($self->stats())) {
870 :     $self->{_genome_genes} = $self->stats()->{'Total genes'}->[0];
871 :     }
872 :     }
873 :    
874 :     return $self->{_genome_genes};
875 :     }
876 :    
877 :     =head3 run_default_model_predictions
878 :     Definition:
879 :     0/1::status = FIGMODELmodel->run_default_model_predictions(string::media ID);
880 :     Description:
881 :     =cut
882 :     sub run_default_model_predictions {
883 :     my ($self,$Media) = @_;
884 :    
885 :     #Assuming complete media if none is provided
886 :     if (!defined($Media)) {
887 :     $Media = "Complete";
888 :     }
889 :    
890 :     #Predicting essentiality
891 :     my $result = $self->figmodel()->RunFBASimulation($self->id(),"SINGLEKO",undef,undef,[$self->id()],[$Media]);
892 :     #Checking that the table is defined and the output file exists
893 :     if (defined($result) && defined($result->get_row(0)->{"ESSENTIALGENES"})) {
894 :     $self->figmodel()->database()->print_array_to_file($self->directory()."EssentialGenes-".$self->id()."-".$Media.".tbl",[join("\n",@{$result->get_row(0)->{"ESSENTIALGENES"}})]);
895 :     }
896 :    
897 :     #Classifying reactions and compounds
898 :     my $tbl = $self->classify_model_reactions($Media);
899 :     $tbl->save();
900 :    
901 :     return $self->figmodel()->success();
902 :     }
903 :    
904 :     =head3 save_obsolete_stats
905 :     Definition:
906 :     FIGMODELmodel->save_obsolete_stats();
907 :     Description:
908 :     =cut
909 :     sub save_obsolete_stats {
910 :     my ($self) = @_;
911 :    
912 :     #checking if stats exists
913 :     my $stats = $self->stats();
914 :     if (defined($stats)) {
915 :     $stats->{"Model ID"}->[0] = $self->id()."V".$stats->{"Version"}->[0].".".$stats->{"Gap fill version"}->[0];
916 :     $self->figmodel()->database()->update_row("OBSOLETE MODEL STATS",$stats,"Model ID");
917 :     $stats->{"Model ID"}->[0] = $self->id();
918 :     }
919 :     }
920 :    
921 :     =head3 update_stats_for_gap_filling
922 :     Definition:
923 :     {string => [string]} = FIGMODELmodel->update_stats_for_gap_filling(int::gapfill time);
924 :     Description:
925 :     =cut
926 :     sub update_stats_for_gap_filling {
927 :     my ($self,$gapfilltime) = @_;
928 :    
929 :     #preserving the stats for the now obselete model
930 :     $self->save_obsolete_stats();
931 :     my $stats = $self->update_model_stats(0);
932 :     $stats->{"Gap filling time"}->[0] = $gapfilltime;
933 :     $stats->{"Gap fill date"}->[0] = time();
934 :     if (!defined($stats->{"Gap fill version"}->[0])) {
935 :     $stats->{"Gap fill version"}->[0] = 0;
936 :     }
937 :     $stats->{"Gap fill version"}->[0]++;
938 :    
939 :     #Updating the stats stored in the table
940 :     $self->figmodel()->database()->update_row("MODEL STATS",$stats,"Model ID");
941 :     return $stats;
942 :     }
943 :    
944 :     =head3 update_stats_for_build
945 :     Definition:
946 :     {string => [string]} = FIGMODELmodel->update_stats_for_build();
947 :     Description:
948 :     =cut
949 :     sub update_stats_for_build {
950 :     my ($self) = @_;
951 :    
952 :     #preserving the stats for the now obselete model
953 :     $self->save_obsolete_stats();
954 :     my $stats = $self->update_model_stats(0);
955 :     $stats->{"Build date"}->[0] = time();
956 :     if (!defined($stats->{"Version"}->[0])) {
957 :     $stats->{"Version"}->[0] = 0;
958 :     }
959 :     $stats->{"Version"}->[0]++;
960 :    
961 :     #Updating the stats stored in the table
962 :     $self->figmodel()->database()->update_row("MODEL STATS",$stats,"Model ID");
963 :     return $stats;
964 :     }
965 :    
966 :     =head3 update_model_stats
967 :     Definition:
968 :     FIGMODELmodel->update_model_stats();
969 :     Description:
970 :     =cut
971 :     sub update_model_stats {
972 :     my ($self) = @_;
973 :    
974 :     #Getting reaction table
975 :     my $rxntbl = $self->reaction_table();
976 :     if (!defined($rxntbl)) {
977 :     $self->figmodel()->error_message("FIGMODELmodel:update_model_stats:Could not load reaction list for ".$self->id());
978 :     return undef;
979 :     }
980 :     my $cpdtbl = $self->compound_table();
981 :    
982 :     #Creating empty status row
983 :     my $CurrentStats = {"Genes with reactions" => [0],
984 :     "Metabolites" => [$cpdtbl->size()],
985 :     "Growmatch reactions" => [0],
986 :     "Spontaneous" => [0],
987 :     "Biolog gap filling reactions" => [0],
988 :     "Number of reactions" => [$rxntbl->size()],
989 :     "Transport reaction"=>[0],
990 :     "Gap filling reactions" => [0],
991 :     "Model ID" => [$self->id()],
992 :     "Subsystem genes with reactions" => [0],
993 :     "Nonsubsystem genes with reactions" => [0],
994 :     Version => [0],
995 :     "Gap fill version" => [0],
996 :     Source => [$self->source()],
997 :     "Genes with one reaction" => [0],
998 :     "Genome ID" => [$self->genome()]};
999 :    
1000 :     my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
1001 :     if (defined($genomestats)) {
1002 :     $CurrentStats->{SOURCE}->[0] = $genomestats->{SOURCE}->[0];
1003 :     $CurrentStats->{"Organism name"}->[0] = $genomestats->{NAME}->[0];
1004 :     $CurrentStats->{"Total genes"}->[0] = $genomestats->{"TOTAL GENES"}->[0];
1005 :     $CurrentStats->{"Gram positive genes"}->[0] = $genomestats->{"GRAM POSITIVE GENES"}->[0];
1006 :     $CurrentStats->{"Gram negative genes"}->[0] = $genomestats->{"GRAM NEGATIVE GENES"}->[0];
1007 :     $CurrentStats->{"Class"}->[0] = $genomestats->{CLASS}->[0];
1008 :     $CurrentStats->{"Genes with functions"}->[0] = $genomestats->{"GENES WITH FUNCTIONS"}->[0];
1009 :     $CurrentStats->{"Subsystem genes"}->[0] = $genomestats->{"SUBSYSTEM GENES"}->[0];
1010 :     $CurrentStats->{"Nonsubsystem genes"}->[0] = $genomestats->{"NON SUBSYSTEM GENES"}->[0];
1011 :     if (defined($genomestats->{TAXONOMY})) {
1012 :     $CurrentStats->{"Taxonomy 0"}->[0] = $genomestats->{TAXONOMY}->[0];
1013 :     $CurrentStats->{"Taxonomy 1"}->[0] = $genomestats->{TAXONOMY}->[1];
1014 :     $CurrentStats->{"Taxonomy 2"}->[0] = $genomestats->{TAXONOMY}->[2];
1015 :     $CurrentStats->{"Taxonomy 3"}->[0] = $genomestats->{TAXONOMY}->[3];
1016 :     $CurrentStats->{"Taxonomy 4"}->[0] = $genomestats->{TAXONOMY}->[4];
1017 :     $CurrentStats->{"Taxonomy 5"}->[0] = $genomestats->{TAXONOMY}->[5];
1018 :     }
1019 :     }
1020 :    
1021 :     #Transfering build, version, and gap fill data from existing stats
1022 :     if (defined($self->stats())) {
1023 :     $CurrentStats->{Version}->[0] = $self->stats()->{Version}->[0];
1024 :     $CurrentStats->{"Gap fill version"}->[0] = $self->stats()->{"Gap fill version"}->[0];
1025 :     $CurrentStats->{"Gap filling time"}->[0] = $self->stats()->{"Gap filling time"}->[0];
1026 :     $CurrentStats->{"Gap fill date"}->[0] = $self->stats()->{"Gap fill date"}->[0];
1027 :     $CurrentStats->{"Build date"}->[0] = $self->stats()->{"Build date"}->[0];
1028 :     }
1029 :    
1030 :     my %GeneHash;
1031 :     my %NonpegHash;
1032 :     my %CompoundHash;
1033 :     for (my $i=0; $i < $rxntbl->size(); $i++) {
1034 :     my $Row = $rxntbl->get_row($i);
1035 :     if (defined($Row) && defined($Row->{"ASSOCIATED PEG"})) {
1036 :     my $ReactionRow = $self->figmodel()->get_reaction($Row->{"LOAD"}->[0]);
1037 :     if (defined($ReactionRow->{"EQUATION"}->[0])) {
1038 :     #Checking for extracellular metabolites which indicate that this is a transporter
1039 :     if ($ReactionRow->{"EQUATION"}->[0] =~ m/\[e\]/) {
1040 :     $CurrentStats->{"Transport reaction"}->[0]++;
1041 :     }
1042 :     }
1043 :     #Identifying spontaneous/biolog/gapfilling/gene associated reactions
1044 :     if ($Row->{"ASSOCIATED PEG"}->[0] =~ m/BIOLOG/i) {
1045 :     $CurrentStats->{"Biolog gap filling reactions"}->[0]++;
1046 :     } elsif ($Row->{"ASSOCIATED PEG"}->[0] =~ m/SPONTANEOUS/i) {
1047 :     $CurrentStats->{"Spontaneous"}->[0]++;
1048 :     } elsif ($Row->{"ASSOCIATED PEG"}->[0] =~ m/GAP/ || $Row->{"ASSOCIATED PEG"}->[0] =~ m/UNIVERSAL/i || $Row->{"ASSOCIATED PEG"}->[0] =~ m/UNKNOWN/i) {
1049 :     $CurrentStats->{"Gap filling reactions"}->[0]++;
1050 :     } else {
1051 :     foreach my $GeneSet (@{$Row->{"ASSOCIATED PEG"}}) {
1052 :     my @GeneList = split(/\+/,$GeneSet);
1053 :     foreach my $Gene (@GeneList) {
1054 :     if ($Gene =~ m/(peg\.\d+)/) {
1055 :     $GeneHash{$1} = 1;
1056 :     } else {
1057 :     $NonpegHash{$Gene} = 1;
1058 :     }
1059 :     }
1060 :     }
1061 :     }
1062 :     }
1063 :     }
1064 :     my @genes = keys(%GeneHash);
1065 :     my @othergenes = keys(%NonpegHash);
1066 :     $CurrentStats->{"Genes with reactions"}->[0] = @genes + @othergenes;
1067 :    
1068 :     #Updating the stats stored in the table
1069 :     $self->figmodel()->database()->update_row("MODEL STATS",$CurrentStats,"Model ID");
1070 :     $self->{_stats} = $CurrentStats;
1071 :     return $CurrentStats;
1072 :     }
1073 :    
1074 :     =head3 GapFillModel
1075 :     Definition:
1076 :     (success/fail) = FIGMODELmodel->GapFillModel();
1077 :     Description:
1078 :     This function performs an optimization to identify the minimal set of reactions that must be added to a model in order for biomass to be produced by the biomass reaction in the model.
1079 :     Before running the gap filling, the existing model is backup in the same directory with the current version numbers appended.
1080 :     If the model has been gap filled previously, the previous gap filling reactions are removed prior to running the gap filling again.
1081 :     =cut
1082 :     sub GapFillModel {
1083 :     my ($self,$donotclear) = @_;
1084 :    
1085 :     #Setting status of model to gap filling
1086 :     my $OrganismID = $self->genome();
1087 :     $self->set_status(1,"Auto completion running");
1088 :     my $UniqueFilename = $self->figmodel()->filename();
1089 :     my $StartTime = time();
1090 :    
1091 :     #Archiving the existing model
1092 :     $self->ArchiveModel();
1093 :    
1094 :     #Removing any gapfilling reactions that may be currently present in the model
1095 :     if (!defined($donotclear) || $donotclear != 1) {
1096 :     my $ModelTable = $self->reaction_table();
1097 :     for (my $i=0; $i < $ModelTable->size(); $i++) {
1098 :     if (!defined($ModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0]) || $ModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] eq "GAP FILLING" || $ModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] =~ m/BIOLOG/ || $ModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] =~ m/GROWMATCH/) {
1099 :     $ModelTable->delete_row($i);
1100 :     $i--;
1101 :     }
1102 :     }
1103 :     $ModelTable->save();
1104 :     }
1105 :    
1106 :     #Calling the MFAToolkit to run the gap filling optimization
1107 :     my $MinimalMediaTable = $self->figmodel()->database()->GetDBTable("MINIMAL MEDIA TABLE");
1108 :     my $Media = "Complete";
1109 :     if (defined($MinimalMediaTable->get_row_by_key($self->genome(),"Organism"))) {
1110 :     $Media = $MinimalMediaTable->get_row_by_key($self->genome(),"Organism")->{"Minimal media"}->[0];
1111 :     #Loading media, changing bounds, saving media as a test media
1112 :     my $MediaTable = FIGMODELTable::load_table($self->config("Media directory")->[0].$Media.".txt",";","",0,["VarName"]);
1113 :     for (my $i=0; $i < $MediaTable->size(); $i++) {
1114 :     if ($MediaTable->get_row($i)->{"Min"}->[0] < 0) {
1115 :     $MediaTable->get_row($i)->{"Min"}->[0] = -10000;
1116 :     }
1117 :     if ($MediaTable->get_row($i)->{"Max"}->[0] > 0) {
1118 :     $MediaTable->get_row($i)->{"Max"}->[0] = 10000;
1119 :     }
1120 :     }
1121 :     $MediaTable->save($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");
1122 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$UniqueFilename."TestMedia",["GapFilling"],{"Default max drain flux" => 0},"GapFill".$self->id().".log",undef));
1123 :     unlink($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");
1124 :     } else {
1125 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),undef,["GapFilling"],undef,"GapFill".$self->id().".log",undef));
1126 :     }
1127 :    
1128 :     #Parse the solutions file for the model and get the reaction list from it
1129 :     my $SolutionData = $self->figmodel()->LoadProblemReport($UniqueFilename);
1130 :     #Clearing the mfatoolkit output and log file
1131 :     $self->figmodel()->clearing_output($UniqueFilename,"GapFill".$self->id().".log");
1132 :     if (!defined($SolutionData) || $SolutionData->size() == 0) {
1133 :     $self->set_status(1,"Auto completion failed: auto completion solution not found");
1134 :     $self->figmodel()->error_message("FIGMODEL:GapFillModel: Gap filling report not found!");
1135 :     return $self->fail();
1136 :     }
1137 :     $SolutionData->save("/home/chenry/Solution.txt");
1138 :    
1139 :     #Looking for the last printed recursive MILP solution
1140 :     for (my $i=($SolutionData->size()-1); $i >=0; $i--) {
1141 :     if (defined($SolutionData->get_row($i)->{"Notes"}) && $SolutionData->get_row($i)->{"Notes"}->[0] =~ m/^Recursive/) {
1142 :     my $AllSolutions = substr($SolutionData->get_row($i)->{"Notes"}->[0],15);
1143 :     print "Solution:".$AllSolutions."\n";
1144 :     my @TempThree = split(/\|/,$AllSolutions);
1145 :     if (@TempThree > 0 && $TempThree[0] =~ m/.+:(.+)/) {
1146 :     my @TempFour = split(/,/,$1);
1147 :     my $DirectionList;
1148 :     my $ReactionList;
1149 :     for (my $j=0; $j < @TempFour; $j++) {
1150 :     my $ID = "";
1151 :     my $Sign = "=>";
1152 :     if ($TempFour[$j] =~ m/\-(rxn\d\d\d\d\d)/) {
1153 :     $ID = $1;
1154 :     $Sign = "<=";
1155 :     } elsif ($TempFour[$j] =~ m/\+(rxn\d\d\d\d\d)/) {
1156 :     $ID = $1;
1157 :     }
1158 :     if ($ID ne "") {
1159 :     if ($self->figmodel()->reversibility_of_reaction($ID) ne $Sign) {
1160 :     $Sign = "<=>";
1161 :     }
1162 :     push(@{$DirectionList},$Sign);
1163 :     push(@{$ReactionList},$ID);
1164 :     }
1165 :     }
1166 :     print "Solution:".$TempThree[0]."\n";
1167 :     $self->figmodel()->IntegrateGrowMatchSolution($self->id(),undef,$ReactionList,$DirectionList,"GAP FILLING",0,1);
1168 :     }
1169 :     last;
1170 :     }
1171 :     }
1172 :     #Updating model stats with gap filling results
1173 :     my $ElapsedTime = time() - $StartTime;
1174 :     $self->figmodel()->ClearDBModel($self->id(),1);
1175 :     #Determining why each gap filling reaction was added
1176 :     $self->figmodel()->IdentifyDependancyOfGapFillingReactions($self->id(),$Media);
1177 :     if (!defined($donotclear) || $donotclear != 1) {
1178 :     $self->figmodel()->AddBiologTransporters($self->id());
1179 :     if ($self->id() !~ m/MGRast/) {
1180 :     $self->update_stats_for_gap_filling($ElapsedTime);
1181 :     }
1182 :     } else {
1183 :     $self->update_model_stats();
1184 :     }
1185 :     #Printing the updated SBML file
1186 :     $self->PrintSBMLFile();
1187 :     $self->PrintModelLPFile($self->id());
1188 :     $self->set_status(1,"Auto completion successfully finished");
1189 :     $self->run_default_model_predictions();
1190 :    
1191 :     return $self->success();
1192 :     }
1193 :    
1194 :     =head3 datagapfill
1195 :     Definition:
1196 :     success()/fail() = FIGMODELmodel->datagapfill();
1197 :     Description:
1198 :     Run gapfilling on the input run specifications
1199 :     =cut
1200 :     sub datagapfill {
1201 :     my ($self,$GapFillingRunSpecs,$TansferFileSuffix) = @_;
1202 :     my $UniqueFilename = $self->figmodel()->filename();
1203 :     if (defined($GapFillingRunSpecs) && @{$GapFillingRunSpecs} > 0) {
1204 :     print "Gap filling specs:\n".join("\n",@{$GapFillingRunSpecs})."\n\n";
1205 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id().$self->selected_version(),"NoBounds",["DataGapFilling"],{"Reactions to knockout" => $self->config("permanently knocked out reactions")->[0],"Gap filling runs" => join(";",@{$GapFillingRunSpecs})},"GapFilling-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,undef));
1206 :    
1207 :     #Checking that the solution exists
1208 :     if (!-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingSolutionTable.txt") {
1209 :     $self->figmodel()->error_message("FIGMODEL:GapFillingAlgorithm: Could not find MFA output file!");
1210 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-GFS.txt",["Experiment;Solution index;Solution cost;Solution reactions"]);
1211 :     return undef;
1212 :     }
1213 :     my $GapFillResultTable = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingSolutionTable.txt",";","",0,undef);
1214 :     if (defined($TansferFileSuffix)) {
1215 :     system("cp ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingSolutionTable.txt ".$self->directory().$self->id().$self->selected_version()."-".$TansferFileSuffix.".txt");
1216 :     system("cp ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingReport.txt ".$self->directory().$self->id().$self->selected_version()."-".$TansferFileSuffix.".txt");
1217 :     }
1218 :    
1219 :     #If the system is not configured to preserve all logfiles, then the mfatoolkit output folder is deleted
1220 :     $self->figmodel()->clearing_output($UniqueFilename,"GapFilling-".$self->id().$self->selected_version()."-".$UniqueFilename.".log");
1221 :     return $GapFillingRunSpecs;
1222 :     }
1223 :     if (defined($TansferFileSuffix)) {
1224 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-".$TansferFileSuffix.".txt",["Experiment;Solution index;Solution cost;Solution reactions"]);
1225 :     }
1226 :     return undef;
1227 :     }
1228 :    
1229 :     =head3 TestSolutions
1230 :     Definition:
1231 :     $model->TestSolutions($ModelID,$NumProcessors,$ProcessorIndex,$GapFill);
1232 :     Description:
1233 :     Example:
1234 :     =cut
1235 :    
1236 :     sub TestSolutions {
1237 :     my ($self,$OriginalErrorFilename,$GapFillResultTable) = @_;
1238 :     #Getting the filename
1239 :     my $UniqueFilename = $self->figmodel()->filename();
1240 :     #Reading in the original error matrix which has the headings for the original model simulation
1241 :     my $OriginalErrorData;
1242 :     if (!defined($OriginalErrorFilename) || !-e $self->directory().$OriginalErrorFilename) {
1243 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector) = $self->RunAllStudiesWithDataFast("All");
1244 :     $OriginalErrorData = [$HeadingVector,$Errorvector];
1245 :     } else {
1246 :     $OriginalErrorData = $self->figmodel()->database()->load_single_column_file($self->directory().$OriginalErrorFilename,"");
1247 :     }
1248 :     my $HeadingHash;
1249 :     my @HeadingArray = split(/;/,$OriginalErrorData->[0]);
1250 :     my @OrigErrorArray = split(/;/,$OriginalErrorData->[1]);
1251 :     for (my $i=0; $i < @HeadingArray; $i++) {
1252 :     my @SubArray = split(/:/,$HeadingArray[$i]);
1253 :     $HeadingHash->{$SubArray[0].":".$SubArray[1].":".$SubArray[2]} = $i;
1254 :     }
1255 :     #Scanning through the gap filling solutions
1256 :     my $TempVersion = "V".$UniqueFilename;
1257 :     my $ErrorMatrixLines;
1258 :     for (my $i=0; $i < $GapFillResultTable->size(); $i++) {
1259 :     print "Starting problem solving ".$i."\n";
1260 :     my $ErrorLine = $GapFillResultTable->get_row($i)->{"Experiment"}->[0].";".$GapFillResultTable->get_row($i)->{"Solution index"}->[0].";".$GapFillResultTable->get_row($i)->{"Solution cost"}->[0].";".$GapFillResultTable->get_row($i)->{"Solution reactions"}->[0];
1261 :     #Integrating solution into test model
1262 :     my $ReactionArray;
1263 :     my $DirectionArray;
1264 :     my @ReactionList = split(/,/,$GapFillResultTable->get_row($i)->{"Solution reactions"}->[0]);
1265 :     my %SolutionHash;
1266 :     for (my $k=0; $k < @ReactionList; $k++) {
1267 :     if ($ReactionList[$k] =~ m/(.+)(rxn\d\d\d\d\d)/) {
1268 :     my $Reaction = $2;
1269 :     my $Sign = $1;
1270 :     if (defined($SolutionHash{$Reaction})) {
1271 :     $SolutionHash{$Reaction} = "<=>";
1272 :     } elsif ($Sign eq "-") {
1273 :     $SolutionHash{$Reaction} = "<=";
1274 :     } elsif ($Sign eq "+") {
1275 :     $SolutionHash{$Reaction} = "=>";
1276 :     } else {
1277 :     $SolutionHash{$Reaction} = $Sign;
1278 :     }
1279 :     }
1280 :     }
1281 :     @ReactionList = keys(%SolutionHash);
1282 :     for (my $k=0; $k < @ReactionList; $k++) {
1283 :     push(@{$ReactionArray},$ReactionList[$k]);
1284 :     push(@{$DirectionArray},$SolutionHash{$ReactionList[$k]});
1285 :     }
1286 :     print "Integrating solution!\n";
1287 :     $self->IntegrateGrowMatchSolution($self->id().$self->selected_version(),$self->directory().$self->id().$TempVersion.".txt",$ReactionArray,$DirectionArray,"Gapfilling ".$GapFillResultTable->get_row($i)->{"Experiment"}->[0],1,1);
1288 :     my $model = $self->get_model($self->id().$TempVersion);
1289 :     $model->PrintModelLPFile();
1290 :     #Running the model against all available experimental data
1291 :     print "Running test model!\n";
1292 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector) = $self->RunAllStudiesWithDataFast("All");
1293 :    
1294 :     @HeadingArray = split(/;/,$HeadingVector);
1295 :     my @ErrorArray = @OrigErrorArray;
1296 :     my @TempArray = split(/;/,$Errorvector);
1297 :     for (my $j=0; $j < @HeadingArray; $j++) {
1298 :     my @SubArray = split(/:/,$HeadingArray[$j]);
1299 :     $ErrorArray[$HeadingHash->{$SubArray[0].":".$SubArray[1].":".$SubArray[2]}] = $TempArray[$j];
1300 :     }
1301 :     $ErrorLine .= ";".$FalsePostives."/".$FalseNegatives.";".join(";",@ErrorArray);
1302 :     push(@{$ErrorMatrixLines},$ErrorLine);
1303 :     print "Finishing problem solving ".$i."\n";
1304 :     }
1305 :     #Clearing out the test model
1306 :     if (-e $self->directory().$self->id().$TempVersion.".txt") {
1307 :     unlink($self->directory().$self->id().$TempVersion.".txt");
1308 :     unlink($self->directory()."SimulationOutput".$self->id().$TempVersion.".txt");
1309 :     }
1310 :     return $ErrorMatrixLines;
1311 :     }
1312 :    
1313 :     =head3 status
1314 :     Definition:
1315 :     int::model status = FIGMODELmodel->status();
1316 :     Description:
1317 :     Returns the current status of the SEED model associated with the input genome ID.
1318 :     model status = 1: model exists
1319 :     model status = 0: model is being built
1320 :     model status = -1: model does not exist
1321 :     model status = -2: model build failed
1322 :     =cut
1323 :     sub status {
1324 :     my ($self) = @_;
1325 :     return $self->{_data}->{status}->[0];
1326 :     }
1327 :    
1328 :     =head3 message
1329 :     Definition:
1330 :     string::model message = FIGMODELmodel->message();
1331 :     Description:
1332 :     Returns a message associated with the models current status
1333 :     =cut
1334 :     sub message {
1335 :     my ($self) = @_;
1336 :     return $self->{_data}->{message}->[0];
1337 :     }
1338 :    
1339 :     =head3 set_status
1340 :     Definition:
1341 :     (success/fail) = FIGMODELmodel->set_status(int::new status,string::status message);
1342 :     Description:
1343 :     Changes the current status of the SEED model
1344 :     new status = 1: model exists
1345 :     new status = 0: model is being built
1346 :     new status = -1: model does not exist
1347 :     new status = -2: model build failed
1348 :     =cut
1349 :     sub set_status {
1350 :     my ($self,$NewStatus,$Message) = @_;
1351 :    
1352 :     #Getting the model row from the MODELS table
1353 :     $self->{_data}->{status}->[0] = $NewStatus;
1354 :     $self->{_data}->{message}->[0] = $Message;
1355 :     $self->figmodel()->database()->update_row("MODELS",$self->{_data},"id");
1356 :     return $self->config("SUCCESS")->[0];
1357 :     }
1358 :    
1359 :     =head3 CreateSingleGenomeReactionList
1360 :     Definition:
1361 :     FIGMODEL->CreateSingleGenomeReactionList();
1362 :     Description:
1363 :     This function uses fig calls to obtain a list of genes and functions for a genome, and it uses a file mapping reactions and functional roles to produce a reaction list.
1364 :     Example:
1365 :     =cut
1366 :    
1367 :     sub CreateSingleGenomeReactionList {
1368 :     my ($self,$RunGapFilling) = @_;
1369 :    
1370 :     #Getting genome stats
1371 :     my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
1372 :     my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
1373 :     if (!defined($FeatureTable)) {
1374 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." genome features could not be accessed!");
1375 :     return $self->fail();
1376 :     }
1377 :     #Checking that the number of genes exceeds the minimum size
1378 :     if ($FeatureTable->size() < $self->config("minimum genome size for modeling")->[0]) {
1379 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." genome rejected as too small for modeling!");
1380 :     return $self->fail();
1381 :     }
1382 :     #Setting up needed variables
1383 :     my $OriginalModelTable = undef;
1384 :     #Locking model status table
1385 :     my $ModelTable = $self->figmodel()->database()->LockDBTable("MODELS");
1386 :     my $Row = $ModelTable->get_row_by_key($self->id(),"id");
1387 :     if (!defined($Row) || !defined($Row->{status}->[0])) {
1388 :     $self->figmodel()->database()->UnlockDBTable("MODELS");
1389 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." has no row in models table!");
1390 :     return $self->fail();
1391 :     } elsif ($Row->{status}->[0] == 0) {
1392 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList:Model is already being built. Canceling current build.");
1393 :     $self->figmodel()->database()->UnlockDBTable("MODELS");
1394 :     return $self->fail();
1395 :     }elsif ($Row->{status}->[0] == 1) {
1396 :     $OriginalModelTable = $self->reaction_table();
1397 :     $self->ArchiveModel();
1398 :     $Row->{status}->[0] = 0;
1399 :     $Row->{message}->[0] = "Rebuilding preliminary reconstruction";
1400 :     } else {
1401 :     $Row->{status}->[0] = 0;
1402 :     $Row->{message}->[0] = "Preliminary reconstruction";
1403 :     }
1404 :     #Updating the status table
1405 :     $self->figmodel()->database()->save_table($ModelTable);
1406 :     $self->figmodel()->database()->UnlockDBTable("MODELS");
1407 :     #Sorting GenomeData by gene location on the chromosome
1408 :     $FeatureTable->sort_rows("MIN LOCATION");
1409 :     my ($ComplexHash,$SuggestedMappings,$UnrecognizedReactions,$ReactionHash);
1410 :     my %LastGenePosition;
1411 :     my $GeneRoles;
1412 :     for (my $j=0; $j < $FeatureTable->size(); $j++) {
1413 :     my $CurrentRow = $FeatureTable->get_row($j);
1414 :     #"ID","ALIASES","MIN LOCATION","MAX LOCATION","ROLES","SUBSYSTEMS","SUBSYSTEM CLASS"
1415 :     if (defined($CurrentRow)) {
1416 :     my $GeneID = $CurrentRow->{"ID"}->[0];
1417 :     if ($GeneID =~ m/(peg\.\d+)/) {
1418 :     $GeneID = $1;
1419 :     }
1420 :     foreach my $Role (@{$CurrentRow->{"ROLES"}}) {
1421 :     if ($self->figmodel()->role_is_valid($Role) != 0) {
1422 :     push(@{$GeneRoles->{$GeneID}},$Role);
1423 :     my $ReactionHashArrayRef = $self->figmodel()->reactions_of_role($Role);
1424 :     if ($ReactionHashArrayRef != 0) {
1425 :     foreach my $Mapping (@{$ReactionHashArrayRef}) {
1426 :     if (defined($Mapping->{"REACTION"}) && defined($Mapping->{"MASTER"}) && defined($Mapping->{"SUBSYSTEM"}) && defined($Mapping->{"SOURCE"})) {
1427 :     if ($Mapping->{"REACTION"}->[0] =~ m/rxn\d\d\d\d\d/) {
1428 :     if ($Mapping->{"MASTER"}->[0] eq 1) {
1429 :     #Creating a complex if consecutive genes have been assigned to the same reaction
1430 :     $ComplexHash->{$Mapping->{"REACTION"}->[0]}->{$Mapping->{"COMPLEX"}->[0]}->{$Role}->{$GeneID} = 1;
1431 :     if (!defined($LastGenePosition{$Mapping->{"REACTION"}->[0]})) {
1432 :     $LastGenePosition{$Mapping->{"REACTION"}->[0]} = $j;
1433 :     push(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}},$GeneID);
1434 :     } elsif (($j-$LastGenePosition{$Mapping->{"REACTION"}->[0]}) < 3 && $LastGenePosition{$Mapping->{"REACTION"}->[0]} != $j) {
1435 :     my $CurrentComplex = pop(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}});
1436 :     push(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}},$CurrentComplex."+".$GeneID);
1437 :     } elsif ($LastGenePosition{$Mapping->{"REACTION"}->[0]} != $j) {
1438 :     push(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}},$GeneID);
1439 :     }
1440 :     $LastGenePosition{$Mapping->{"REACTION"}->[0]} = $j;
1441 :     #Adding a subsystem for the reaction
1442 :     if ($self->figmodel()->subsystem_is_valid($Mapping->{"SUBSYSTEM"}->[0]) == 1) {
1443 :     ($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SUBSYSTEMS"},my $NumMatches) = $self->figmodel()->add_elements_unique($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SUBSYSTEMS"},$Mapping->{"SUBSYSTEM"}->[0]);
1444 :     if (!defined($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}) || $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] > 1) {
1445 :     if ($Mapping->{"SOURCE"}->[0] =~ m/Hope\sFiles/) {
1446 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 1;
1447 :     } elsif ($Mapping->{"SOURCE"}->[0] =~ m/SEED/) {
1448 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 2;
1449 :     } elsif (!defined($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}) || $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] > 2) {
1450 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 3;
1451 :     }
1452 :     }
1453 :     }
1454 :     #Handling confidence
1455 :     if (!defined($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}) || $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] > 2) {
1456 :     if ($Mapping->{"SOURCE"}->[0] =~ m/MATT/) {
1457 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 3;
1458 :     } elsif ($Mapping->{"SOURCE"}->[0] =~ m/CHRIS/) {
1459 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 4;
1460 :     } else {
1461 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 5;
1462 :     }
1463 :     }
1464 :     #Parsing sources
1465 :     ($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SOURCE"},my $NumMatches) = $self->figmodel()->add_elements_unique($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SOURCE"},split(/\|/,$Mapping->{"SOURCE"}->[0]));
1466 :     } else {
1467 :     push(@{$SuggestedMappings},$GeneID."\t".$Mapping->{"REACTION"}->[0]."\t".$Role);
1468 :     }
1469 :     } else {
1470 :     push(@{$UnrecognizedReactions},$GeneID."\t".$Mapping->{"REACTION"}->[0]."\t".$Role);
1471 :     }
1472 :     }
1473 :     }
1474 :     }
1475 :     }
1476 :     }
1477 :     }
1478 :     }
1479 :    
1480 :     #Creating nonadjacent complexes
1481 :     my @ReactionList = keys(%{$ReactionHash});
1482 :     foreach my $Reaction (@ReactionList) {
1483 :     #If multiple genes are assigned to the reaction, we check if they should should be in a complex
1484 :     if (@{$ReactionHash->{$Reaction}->{"GENES"}} > 0 && defined($ComplexHash->{$Reaction})) {
1485 :     my $GeneArray;
1486 :     foreach my $Complex (keys(%{$ComplexHash->{$Reaction}})) {
1487 :     my %ComplexComponents;
1488 :     foreach my $CurrentGeneSet (@{$ReactionHash->{$Reaction}->{"GENES"}}) {
1489 :     my @GeneList = split(/\+/,$CurrentGeneSet);
1490 :     my %RoleHash;
1491 :     foreach my $Gene (@GeneList) {
1492 :     foreach my $Role (@{$GeneRoles->{$Gene}}) {
1493 :     if (defined($ComplexHash->{$Reaction}->{$Complex}->{$Role})) {
1494 :     $RoleHash{$Role} = 1;
1495 :     }
1496 :     }
1497 :     }
1498 :     if (keys(%RoleHash) > 0) {
1499 :     if (!defined($ComplexComponents{join("|",sort(keys(%RoleHash)))})) {
1500 :     my @RoleList = keys(%RoleHash);
1501 :     my @ComplexList = keys(%ComplexComponents);
1502 :     foreach my $ComplexSet (@ComplexList) {
1503 :     my @RoleList = split(/\|/,$ComplexSet);
1504 :     my $Match = 0;
1505 :     foreach my $SingleRole (@RoleList) {
1506 :     if (defined($RoleHash{$SingleRole})) {
1507 :     $Match = 1;
1508 :     last;
1509 :     }
1510 :     }
1511 :     if ($Match == 1) {
1512 :     foreach my $SingleRole (@RoleList) {
1513 :     $RoleHash{$SingleRole} = 1
1514 :     }
1515 :     push(@{$ComplexComponents{join("|",sort(keys(%RoleHash)))}},@{$ComplexComponents{$ComplexSet}});
1516 :     delete $ComplexComponents{$ComplexSet};
1517 :     }
1518 :     }
1519 :     }
1520 :     push(@{$ComplexComponents{join("|",sort(keys(%RoleHash)))}},$CurrentGeneSet);
1521 :     }
1522 :     }
1523 :     my @Position;
1524 :     my @Options;
1525 :     my $Count = 0;
1526 :     foreach my $RoleSet (keys(%ComplexComponents)) {
1527 :     push(@Position,0);
1528 :     push(@{$Options[$Count]},@{$ComplexComponents{$RoleSet}});
1529 :     $Count++;
1530 :     }
1531 :     my $Done = 0;
1532 :     $Count = 0;
1533 :     my $NewRelationship;
1534 :     while($Done == 0) {
1535 :     #Creating complex with current indecies
1536 :     $NewRelationship->[$Count] = $Options[0]->[$Position[0]];
1537 :     for (my $i=1; $i < @Position; $i++) {
1538 :     $NewRelationship->[$Count] .= "+".$Options[$i]->[$Position[$i]];
1539 :     }
1540 :     $NewRelationship->[$Count] = join("+",$self->figmodel()->remove_duplicates(split(/\+/,$NewRelationship->[$Count])));
1541 :     $Count++;
1542 :     #Iterating indecies
1543 :     my $CurrentIndex = 0;
1544 :     while($CurrentIndex >= 0) {
1545 :     if ($CurrentIndex >= @Position) {
1546 :     $CurrentIndex = -1000;
1547 :     } elsif ($Position[$CurrentIndex]+1 == @{$Options[$CurrentIndex]}) {
1548 :     $Position[$CurrentIndex] = -1;
1549 :     $CurrentIndex++;
1550 :     } else {
1551 :     $Position[$CurrentIndex]++;
1552 :     $CurrentIndex--;
1553 :     }
1554 :     }
1555 :     if ($CurrentIndex == -1000) {
1556 :     $Done = 1;
1557 :     }
1558 :     }
1559 :     push(@{$GeneArray},@{$NewRelationship});
1560 :     }
1561 :     @{$ReactionHash->{$Reaction}->{"GENES"}} = $self->figmodel()->remove_duplicates(@{$GeneArray});
1562 :     }
1563 :     }
1564 :    
1565 :     #Getting the reaction table
1566 :     my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");
1567 :    
1568 :     #Creating the model reaction table
1569 :     my $NewModelTable = FIGMODELTable->new(["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG","SUBSYSTEM","CONFIDENCE","REFERENCE","NOTES"],$self->directory().$self->id().".txt",["LOAD"],";","|","REACTIONS\n");
1570 :     @ReactionList = keys(%{$ReactionHash});
1571 :     foreach my $ReactionID (@ReactionList) {
1572 :     #Getting the thermodynamic reversibility from the database
1573 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1574 :     my $Subsystem = "NONE";
1575 :     if (defined($ReactionHash->{$ReactionID}->{"SUBSYSTEMS"})) {
1576 :     $Subsystem = join("|",@{$ReactionHash->{$ReactionID}->{"SUBSYSTEMS"}});
1577 :     }
1578 :     my $Source = "NONE";
1579 :     if (defined($ReactionHash->{$ReactionID}->{"SOURCE"})) {
1580 :     $Source = join("|",@{$ReactionHash->{$ReactionID}->{"SOURCE"}});
1581 :     }
1582 :     $NewModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => [join("|",@{$ReactionHash->{$ReactionID}->{"GENES"}})],"SUBSYSTEM" => [$Subsystem],"CONFIDENCE" => [$ReactionHash->{$ReactionID}->{"CONFIDENCE"}->[0]],"REFERENCE" => [$Source],"NOTES" => ["NONE"]});
1583 :     }
1584 :    
1585 :     #Adding the spontaneous and universal reactions
1586 :     foreach my $ReactionID (@{$self->config("spontaneous reactions")}) {
1587 :     #Getting the thermodynamic reversibility from the database
1588 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1589 :     #Checking if the reaction is already in the model
1590 :     if (!defined($NewModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1591 :     $NewModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["SPONTANEOUS"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["SPONTANEOUS"],"NOTES" => ["NONE"]});
1592 :     }
1593 :     }
1594 :     foreach my $ReactionID (@{$self->config("universal reactions")}) {
1595 :     #Getting the thermodynamic reversibility from the database
1596 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1597 :     #Checking if the reaction is already in the model
1598 :     if (!defined($NewModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1599 :     $NewModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["UNIVERSAL"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["UNIVERSAL"],"NOTES" => ["NONE"]});
1600 :     }
1601 :     }
1602 :    
1603 :     #Checking if a biomass reaction already exists
1604 :     my $BiomassReactionRow = $self->figmodel()->database()->get_row_by_key("BIOMASS TABLE",$self->id(),"MODELS");
1605 :     if (!defined($BiomassReactionRow)) {
1606 :     $BiomassReactionRow = $self->BuildSpecificBiomassReaction();
1607 :     if (!defined($BiomassReactionRow)) {
1608 :     $self->set_status(-2,"Preliminary reconstruction failed: could not generate biomass reaction");
1609 :     $self->figmodel()->error_message("FIGMODEL:CreateModelReactionList: Could not generate biomass function for ".$self->id()."!");
1610 :     return $self->fail();
1611 :     }
1612 :     }
1613 :     my $ReactionList = $BiomassReactionRow->{"ESSENTIAL REACTIONS"};
1614 :     push(@{$ReactionList},$BiomassReactionRow->{DATABASE}->[0]);
1615 :    
1616 :     #Adding biomass reactions to the model table
1617 :     foreach my $BOFReaction (@{$ReactionList}) {
1618 :     #Getting the thermodynamic reversibility from the database
1619 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($BOFReaction);
1620 :     #Checking if the reaction is already in the model
1621 :     if (!defined($NewModelTable->get_row_by_key($BOFReaction,"LOAD"))) {
1622 :     if ($BOFReaction =~ m/bio/) {
1623 :     $NewModelTable->add_row({"LOAD" => [$BOFReaction],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["BOF"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [1],"REFERENCE" => ["Biomass objective function"],"NOTES" => ["NONE"]});
1624 :     } else {
1625 :     $NewModelTable->add_row({"LOAD" => [$BOFReaction],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["INITIAL GAP FILLING"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [5],"REFERENCE" => ["Initial gap filling"],"NOTES" => ["NONE"]});
1626 :     }
1627 :     }
1628 :     }
1629 :    
1630 :     #Completing any incomplete reactions sets
1631 :     my $ReactionSetTable = $self->figmodel()->database()->GetDBTable("REACTION SETS");
1632 :     for (my $i=0; $i < $ReactionSetTable->size(); $i++) {
1633 :     if (defined($NewModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1634 :     foreach my $Reaction (@{$ReactionSetTable->get_row($i)->{"Dependant reactions"}}) {
1635 :     if (!defined($NewModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1636 :     #Getting the thermodynamic reversibility from the database
1637 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($Reaction);
1638 :     $NewModelTable->add_row({"LOAD" => [$Reaction],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["REACTION SET GAP FILLING"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [5],"REFERENCE" => ["Added due to presence of ".$ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0]],"NOTES" => ["NONE"]});
1639 :     }
1640 :     }
1641 :     }
1642 :     }
1643 :    
1644 :     #Now we compare the model to the previous model to determine if any differences exist that aren't gap filling reactions
1645 :     if (defined($OriginalModelTable)) {
1646 :     my $PerfectMatch = 1;
1647 :     my $ReactionCount = 0;
1648 :     for (my $i=0; $i < $OriginalModelTable->size(); $i++) {
1649 :     #We only check that nongapfilling reactions exist in the new model
1650 :     if ($OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] !~ m/GAP/ || $OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] eq "INITIAL GAP FILLING") {
1651 :     $ReactionCount++;
1652 :     my $Row = $NewModelTable->get_row_by_key($OriginalModelTable->get_row($i)->{"LOAD"}->[0],"LOAD");
1653 :     if (defined($Row)) {
1654 :     #We check that the reaction directionality is identical
1655 :     if ($Row->{"DIRECTIONALITY"}->[0] ne $OriginalModelTable->get_row($i)->{"DIRECTIONALITY"}->[0]) {
1656 :     if (defined($OriginalModelTable->get_row($i)->{"NOTES"}->[0]) && $OriginalModelTable->get_row($i)->{"NOTES"}->[0] =~ m/Directionality\sswitched\sfrom\s([^\s])/) {
1657 :     if ($1 ne $Row->{"DIRECTIONALITY"}->[0]) {
1658 :     print "Directionality mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".$1." vs ".$Row->{"DIRECTIONALITY"}->[0]."\n";
1659 :     $PerfectMatch = 0;
1660 :     last;
1661 :     }
1662 :     } else {
1663 :     print "Directionality mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".$OriginalModelTable->get_row($i)->{"DIRECTIONALITY"}->[0]." vs ".$Row->{"DIRECTIONALITY"}->[0]."\n";
1664 :     $PerfectMatch = 0;
1665 :     last;
1666 :     }
1667 :     }
1668 :     #We check that the genes assigned to the reaction are identical
1669 :     if ($PerfectMatch == 1 && @{$OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}} != @{$Row->{"ASSOCIATED PEG"}}) {
1670 :     print "Gene associatation mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".@{$OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}}." vs ".@{$Row->{"ASSOCIATED PEG"}}."\n";
1671 :     $PerfectMatch = 0;
1672 :     last;
1673 :     }
1674 :     if ($PerfectMatch == 1) {
1675 :     my @GeneSetOne = sort(@{$OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}});
1676 :     my @GeneSetTwo = sort(@{$Row->{"ASSOCIATED PEG"}});
1677 :     for (my $j=0; $j < @GeneSetOne; $j++) {
1678 :     if ($GeneSetOne[$j] ne $GeneSetTwo[$j]) {
1679 :     print "Gene mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".$GeneSetOne[$j]." vs ".$GeneSetTwo[$j]."\n";
1680 :     $PerfectMatch = 0;
1681 :     $i = $OriginalModelTable->size();
1682 :     last;
1683 :     }
1684 :     }
1685 :     }
1686 :     } else {
1687 :     print "Original model contains an extra reaction:".$OriginalModelTable->get_row($i)->{"LOAD"}->[0]."\n";
1688 :     $PerfectMatch = 0;
1689 :     last;
1690 :     }
1691 :     }
1692 :     }
1693 :     if ($PerfectMatch == 1 && $ReactionCount == $NewModelTable->size()) {
1694 :     #Bailing out of function as the model has not changed
1695 :     $self->set_status(1,"rebuild canceled because model has not changed");
1696 :     return $self->success();
1697 :     }
1698 :     }
1699 :    
1700 :     #Saving the new model to file
1701 :     $self->set_status(1,"Preliminary reconstruction complete");
1702 :     $self->figmodel()->database()->save_table($NewModelTable);
1703 :     $self->{_reaction_data} = $NewModelTable;
1704 :     #Clearing the previous model from the cache
1705 :     $self->figmodel()->ClearDBModel($self->id(),1);
1706 :     #Updating the model stats table
1707 :     $self->update_stats_for_build();
1708 :     $self->PrintSBMLFile();
1709 :    
1710 :     #Adding model to gapfilling queue
1711 :     if (defined($RunGapFilling) && $RunGapFilling == 1) {
1712 :     $self->set_status(1,"Autocompletion queued");
1713 :     $self->figmodel()->add_job_to_queue("gapfillmodel?".$self->id(),"QSUB","cplex",$self->owner(),"BACK");
1714 :     }
1715 :     return $self->success();
1716 :     }
1717 :    
1718 :     =head3 CreateMetaGenomeReactionList
1719 :     Definition:
1720 :     (success/fail) = FIGMODELmodel->CreateMetaGenomeReactionList();
1721 :     Description:
1722 :     This is the code called to create or update the reaction list for a metgenome model
1723 :     =cut
1724 :    
1725 :     sub CreateMetaGenomeReactionList {
1726 :     my ($self) = @_;
1727 :    
1728 :     #Checking if the metagenome file exists
1729 :     if (!-e $self->config("raw MGRAST directory")->[0].$self->genome().".summary") {
1730 :     $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());
1731 :     }
1732 :     #Loading metagenome data
1733 :     my $MGRASTData = $self->figmodel()->database()->load_multiple_column_file($self->config("raw MGRAST directory")->[0].$self->genome().".summary","\t");
1734 :     if (!defined($MGRASTData)) {
1735 :     $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());
1736 :     }
1737 :    
1738 :     #Setting up needed variables
1739 :     my $OriginalModelTable = undef;
1740 :    
1741 :     #Checking status
1742 :     if ($self->status() < 0) {
1743 :     $self->set_status(0,"Preliminary reconstruction");
1744 :     } elsif ($self->status() == 0) {
1745 :     $self->error_message("FIGMODEL->CreateModelReactionList:Model is already being built. Canceling current build.");
1746 :     return $self->fail();
1747 :     } else {
1748 :     $OriginalModelTable = $self->reaction_table();
1749 :     $self->ArchiveModel();
1750 :     $self->set_status(0,"Rebuilding preliminary reconstruction");
1751 :     }
1752 :    
1753 :     #Getting the reaction table
1754 :     my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");
1755 :     #Creating model table
1756 :     my $ModelTable = FIGMODELTable->new(["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG","SUBSYSTEM","CONFIDENCE","REFERENCE","NOTES"],$self->directory().$self->id().".txt",["LOAD"],";","|","REACTIONS\n");
1757 :     for (my $i=0; $i < @{$MGRASTData};$i++) {
1758 :     #MD5,PEG,number of sims,role,sim e-scores
1759 :     my $Role = $MGRASTData->[$i]->[3];
1760 :     my $MD5 = $MGRASTData->[$i]->[0];
1761 :     my $peg = $MGRASTData->[$i]->[1];
1762 :     my $sims = $MGRASTData->[$i]->[4];
1763 :     $sims =~ s/;/,/g;
1764 :    
1765 :     #Checking for subsystems
1766 :     my $GeneSubsystems = $self->figmodel()->subsystems_of_role($Role);
1767 :     #Checking if there are reactions associated with this role
1768 :     my $ReactionHashArrayRef = $self->figmodel()->reactions_of_role($Role);
1769 :     if ($ReactionHashArrayRef != 0) {
1770 :     foreach my $Mapping (@{$ReactionHashArrayRef}) {
1771 :     if (defined($Mapping->{"REACTION"}) && defined($Mapping->{"MASTER"}) && defined($Mapping->{"SUBSYSTEM"}) && defined($Mapping->{"SOURCE"})) {
1772 :     if ($Mapping->{"REACTION"}->[0] =~ m/rxn\d\d\d\d\d/) {
1773 :     if ($Mapping->{"MASTER"}->[0] eq 1) {
1774 :     #Checking if the reaction is already in the model
1775 :     my $ReactionRow = $ModelTable->get_row_by_key($Mapping->{"REACTION"}->[0],"LOAD");
1776 :     if (!defined($ReactionRow)) {
1777 :     $ReactionRow = {"LOAD" => [$Mapping->{"REACTION"}->[0]],"DIRECTIONALITY" => [$self->figmodel()->reversibility_of_reaction($Mapping->{"REACTION"}->[0])],"COMPARTMENT" => ["c"]};
1778 :     $ModelTable->add_row($ReactionRow);
1779 :     }
1780 :     push(@{$ReactionRow->{"ASSOCIATED PEG"}},substr($peg,4));
1781 :     push(@{$ReactionRow->{"REFERENCE"}},$MD5.":".$Role);
1782 :     push(@{$ReactionRow->{"CONFIDENCE"}},$sims);
1783 :     if (defined($GeneSubsystems)) {
1784 :     push(@{$ReactionRow->{"SUBSYSTEM"}},@{$GeneSubsystems});
1785 :     }
1786 :     }
1787 :     }
1788 :     }
1789 :     }
1790 :     }
1791 :     }
1792 :    
1793 :     #Adding the spontaneous and universal reactions
1794 :     foreach my $ReactionID (@{$self->config("spontaneous reactions")}) {
1795 :     #Getting the thermodynamic reversibility from the database
1796 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1797 :     #Checking if the reaction is already in the model
1798 :     if (!defined($ModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1799 :     $ModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["SPONTANEOUS"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["SPONTANEOUS"],"NOTES" => ["NONE"]});
1800 :     }
1801 :     }
1802 :     foreach my $ReactionID (@{$self->config("universal reactions")}) {
1803 :     #Getting the thermodynamic reversibility from the database
1804 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1805 :     #Checking if the reaction is already in the model
1806 :     if (!defined($ModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1807 :     $ModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["UNIVERSAL"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["UNIVERSAL"],"NOTES" => ["NONE"]});
1808 :     }
1809 :     }
1810 :    
1811 :     #Completing any incomplete reactions sets
1812 :     my $ReactionSetTable = $self->figmodel()->database()->GetDBTable("REACTION SETS");
1813 :     for (my $i=0; $i < $ReactionSetTable->size(); $i++) {
1814 :     if (defined($ModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1815 :     foreach my $Reaction (@{$ReactionSetTable->get_row($i)->{"Dependant reactions"}}) {
1816 :     if (!defined($ModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1817 :     #Getting the thermodynamic reversibility from the database
1818 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($Reaction);
1819 :     $ModelTable->add_row({"LOAD" => [$Reaction],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["REACTION SET GAP FILLING"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [5],"REFERENCE" => ["Added due to presence of ".$ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0]],"NOTES" => ["NONE"]});
1820 :     }
1821 :     }
1822 :     }
1823 :     }
1824 :    
1825 :     #Clearing the previous model from the cache
1826 :     $self->figmodel()->ClearDBModel($self->id(),1);
1827 :     $ModelTable->save();
1828 :    
1829 :     return $self->success();
1830 :     }
1831 :    
1832 :     =head3 ArchiveModel
1833 :     Definition:
1834 :     (success/fail) = FIGMODELmodel->ArchiveModel();
1835 :     Description:
1836 :     This function archives the specified model in the model directory with the current version numbers appended.
1837 :     This function is used to preserve old versions of models prior to overwriting so new versions may be compared with old versions.
1838 :     =cut
1839 :     sub ArchiveModel {
1840 :     my ($self) = @_;
1841 :    
1842 :     #Making sure the model exists
1843 :     if (!defined($self->stats())) {
1844 :     $self->figmodel()->error_message("FIGMODEL:ArchiveModel: Could not find specified ".$self->id()." in database!");
1845 :     return $self->fail();
1846 :     }
1847 :    
1848 :     #Checking that the model file exists
1849 :     if (!(-e $self->filename())) {
1850 :     $self->figmodel()->error_message("FIGMODEL:ArchiveModel: Model file ".$self->filename()." not found!");
1851 :     return $self->fail();
1852 :     }
1853 :    
1854 :     #Copying the model file
1855 :     system("cp ".$self->filename()." ".$self->directory().$self->id().$self->version().".txt");
1856 :     }
1857 :    
1858 :     =head3 PrintModelDataToFile
1859 :     Definition:
1860 :     (success/fail) = FIGMODELmodel->PrintModelDataToFile();
1861 :     Description:
1862 :     This function uses the MFAToolkit to print out all of the compound and reaction data for the input model.
1863 :     Some of the data printed by the toolkit is calculated internally in the toolkit and not stored in any files, so this data can only be retrieved through this
1864 :     function. The LoadModel function for example would not give you this data.
1865 :     =cut
1866 :     sub PrintModelDataToFile {
1867 :     my($self) = @_;
1868 :    
1869 :     #Running the MFAToolkit on the model file
1870 :     my $OutputIndex = $self->figmodel()->filename();
1871 :     my $Command = $self->config("MFAToolkit executable")->[0]." parameterfile ../Parameters/Printing.txt resetparameter output_folder ".$OutputIndex.'/ LoadCentralSystem "'.$self->filename().'"';
1872 :     system($Command);
1873 :    
1874 :     #Copying the model file printed by the toolkit out of the output directory and into the model directory
1875 :     if (!-e $self->config("MFAToolkit output directory")->[0].$OutputIndex."/".$self->id().$self->selected_version().".txt") {
1876 :     $self->figmodel()->error_message("New model file not created due to an error. Check that the input modelfile exists.");
1877 :     $self->figmodel()->cleardirectory($OutputIndex);
1878 :     return $self->fail();
1879 :     }
1880 :    
1881 :     $Command = 'cp "'.$self->config("MFAToolkit output directory")->[0].$OutputIndex."/".$self->id().$self->selected_version().'.txt" "'.$self->directory().$self->id().$self->selected_version().'Data.txt"';
1882 :     system($Command);
1883 :     $Command = 'cp "'.$self->config("MFAToolkit output directory")->[0].$OutputIndex.'/ErrorLog0.txt" "'.$self->directory().'ModelErrors.txt"';
1884 :     system($Command);
1885 :     $self->figmodel()->cleardirectory($OutputIndex);
1886 :     return $self->success();
1887 :     }
1888 :    
1889 :     =head2 Analysis Functions
1890 :    
1891 :     =head3 run_microarray_analysis
1892 :     Definition:
1893 :     int::status = FIGMODEL->run_microarray_analysis(string::media,string::job id,string::gene calls);
1894 :     Description:
1895 :     Runs microarray analysis attempting to turn off genes that are inactive in the microarray
1896 :     =cut
1897 :     sub run_microarray_analysis {
1898 :     my ($self,$media,$jobid,$index,$genecall) = @_;
1899 :     $genecall =~ s/_/:/g;
1900 :     $genecall =~ s/\//;/g;
1901 :     #print "\n\n".$genecall."\n\n";
1902 :     my $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($jobid,$self->id(),$media,["MFA","MicroarrayAssertions"],{"Microarray assertions" => $self->id().";".$index.";".$genecall,"MFASolver" => "CPLEX","Network output location" => "/scratch/"},"MicroarrayAnalysis-".$jobid.".txt",undef,$self->selected_version());
1903 :     #print $command."\n";
1904 :     system($command);
1905 :     #system("/home/chenry/Software/MFAToolkitRepository/Linux/mfatoolkit resetparameter \"user bounds filename\" \"Carbon-D-Glucose.txt\" resetparameter output_folder \"32749.354149.0/\" resetparameter \"Microarray assertions\" \"Seed83333.1;0;peg.3130:-1;peg.4035:-1\" resetparameter \"MFASolver\" \"CPLEX\" resetparameter \"Network output location\" \"/scratch/\" LoadCentralSystem \"/vol/model-dev/MODEL_DEV_DB/Models/83333.1/Seed83333.1.txt\" > \"/vol/model-dev/MODEL_DEV_DB/ReactionDB/log/MicroarrayAnalysis-32749.354149.0.txt\"");
1906 :     }
1907 :    
1908 :     =head3 find_minimal_pathways
1909 :     Definition:
1910 :     int::status = FIGMODEL->find_minimal_pathways(string::media,string::objective);
1911 :     Description:
1912 :     Runs microarray analysis attempting to turn off genes that are inactive in the microarray
1913 :     =cut
1914 :     sub find_minimal_pathways {
1915 :     my ($self,$media,$objective,$solutionnum) = @_;
1916 :    
1917 :     #Setting default media
1918 :     if (!defined($media)) {
1919 :     $media = "Complete";
1920 :     }
1921 :    
1922 :     #Setting default solution number
1923 :     if (!defined($solutionnum)) {
1924 :     $solutionnum = "5";
1925 :     }
1926 :    
1927 :     #Translating objective
1928 :     my $objectivestring;
1929 :     if ($objective eq "ALL") {
1930 :     #Getting the list of universal building blocks
1931 :     my $buildingblocks = $self->config("universal building blocks");
1932 :     my @objectives = keys(%{$buildingblocks});
1933 :     #Getting the nonuniversal building blocks
1934 :     my $otherbuildingblocks = $self->config("nonuniversal building blocks");
1935 :     my @array = keys(%{$otherbuildingblocks});
1936 :     if (defined($self->get_biomass()) && defined($self->figmodel()->get_reaction($self->get_biomass()->{"LOAD"}->[0]))) {
1937 :     my $equation = $self->figmodel()->get_reaction($self->get_biomass()->{"LOAD"}->[0])->{"EQUATION"}->[0];
1938 :     if (defined($equation)) {
1939 :     for (my $i=0; $i < @array; $i++) {
1940 :     if (CORE::index($equation,$array[$i]) > 0) {
1941 :     push(@objectives,$array[$i]);
1942 :     }
1943 :     }
1944 :     }
1945 :     }
1946 :     for (my $i=0; $i < @objectives; $i++) {
1947 :     $self->find_minimal_pathways($media,$objectives[$i])
1948 :     }
1949 :     } elsif ($objective eq "ENERGY") {
1950 :     $objectivestring = "MAX;FLUX;rxn00062;c;1";
1951 :     } elsif ($objective =~ m/cpd\d\d\d\d\d/) {
1952 :     $objectivestring = "MIN;DRAIN_FLUX;".$objective.";c;1";
1953 :     }
1954 :    
1955 :     #Setting additional exchange fluxes
1956 :     my $additionalexchange = $self->config("default exchange fluxes")->[0];
1957 :     if ($objective eq "cpd15665" || $objective eq "cpd15667" || $objective eq "cpd15668" || $objective eq "cpd15669") {
1958 :     $additionalexchange .= ";cpd15666[c]:0:100";
1959 :     } elsif ($objective eq "cpd11493") {
1960 :     $additionalexchange .= ";cpd12370[c]:0:100";
1961 :     } elsif ($objective eq "cpd00166") {
1962 :     $additionalexchange .= ";cpd01997[c]:0:100;cpd03422[c]:0:100";
1963 :     }
1964 :    
1965 :     #Running MFAToolkit
1966 :     my $filename = $self->figmodel()->filename();
1967 :     my $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($filename,$self->id(),$media,["MFA"],{"Recursive MILP solution limit" => $solutionnum,"Generate pathways to objective" => 1,"MFASolver" => "CPLEX","objective" => $objectivestring,"exchange species" => $additionalexchange},"MinimalPathways-".$media."-".$self->id().$self->selected_version().".txt",undef,$self->selected_version());
1968 :     system($command);
1969 :    
1970 :     #Loading problem report
1971 :     my $results = $self->figmodel()->LoadProblemReport($filename);
1972 :     #Clearing output
1973 :     $self->figmodel()->clearing_output($filename,"MinimalPathways-".$media."-".$self->id()."-".$objective.".txt");
1974 :     if (!defined($results)) {
1975 :     return;
1976 :     }
1977 :    
1978 :     #Parsing output
1979 :     my @Array;
1980 :     my $row = $results->get_row(1);
1981 :     if (defined($row->{"Notes"}->[0])) {
1982 :     $_ = $row->{"Notes"}->[0];
1983 :     @Array = /\d+:([^\|]+)\|/g;
1984 :     }
1985 :    
1986 :     # Storing data in a figmodel table
1987 :     my $TableObject;
1988 :     if (-e $self->directory()."MinimalPathways-".$media."-".$self->id().$self->selected_version().".txt") {
1989 :     $TableObject->load_table($self->directory()."MinimalPathways-".$media."-".$self->id().$self->selected_version().".txt",";","|",0,["OBJECTIVE"]);
1990 :     } else {
1991 :     $TableObject = FIGMODELTable->new(["OBJECTIVE","REACTIONS"],$self->directory()."MinimalPathways-".$media."-".$self->id().$self->selected_version().".txt",["OBJECTIVE"],";","|",undef);
1992 :     }
1993 :     my $tablerow = $TableObject->get_row_by_key($objective,"OBJECTIVE",1);
1994 :     push(@{$tablerow->{"REACTIONS"}},@Array);
1995 :     $TableObject->save();
1996 :     }
1997 :    
1998 :     =head3 calculate_growth
1999 :     Definition:
2000 :     string::growth = FIGMODELmodel->calculate_growth(string:media);
2001 :     Description:
2002 :     Calculating growth in the input media
2003 :     =cut
2004 :     sub calculate_growth {
2005 :     my ($self,$Media) = @_;
2006 :     my $UniqueFilename = $self->figmodel()->filename();
2007 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$Media,["ProductionMFA"],{"optimize metabolite production if objective is zero" => 1},$self->id()."-".$Media."-GrowthTest.txt",undef,$self->selected_version()));
2008 :     my $ProblemReport = $self->figmodel()->LoadProblemReport($UniqueFilename);
2009 :     my $Result;
2010 :     if (defined($ProblemReport)) {
2011 :     my $Row = $ProblemReport->get_row(0);
2012 :     if (defined($Row) && defined($Row->{"Objective"}->[0])) {
2013 :     if ($Row->{"Objective"}->[0] < 0.00000001) {
2014 :     $Result = "NOGROWTH:".$Row->{"Individual metabolites with zero production"}->[0];
2015 :     } else {
2016 :     $Result = $Row->{"Objective"}->[0];
2017 :     }
2018 :     }
2019 :     }
2020 :     return $Result;
2021 :     }
2022 :    
2023 :     =head3 classify_model_reactions
2024 :     Definition:
2025 :     (FIGMODELTable:Reaction classes,FIGMODELTable:Compound classes) = FIGMODELmodel->classify_model_reactions(string:media);
2026 :     Description:
2027 :     This function uses the MFAToolkit to minimize and maximize the flux through every reaction in the input model during minimal growth on the input media.
2028 :     The results are returned in a hash of strings where the keys are the reaction IDs and the strings are structured as follows: "Class;Min flux;Max flux".
2029 :     Possible values for "Class" include:
2030 :     1.) Positive: these reactions are essential in the forward direction.
2031 :     2.) Negative: these reactions are essential in the reverse direction.
2032 :     3.) Positive variable: these reactions are nonessential, but they only ever proceed in the forward direction.
2033 :     4.) Negative variable: these reactions are nonessential, but they only ever proceed in the reverse direction.
2034 :     5.) Variable: these reactions are nonessential and proceed in the forward or reverse direction.
2035 :     6.) Blocked: these reactions never carry any flux at all in the media condition tested.
2036 :     7.) Dead: these reactions are disconnected from the network.
2037 :     =cut
2038 :     sub classify_model_reactions {
2039 :     my ($self,$Media) = @_;
2040 :    
2041 :     #Getting unique file for printing model output
2042 :     my $UniqueFilename = $self->figmodel()->filename();
2043 :     #Running the MFAToolkit
2044 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$Media,["ProductionMFA"],{"identify dead ends" => 1,"find tight bounds" => 1,"MFASolver" => "GLPK"},"Classify-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
2045 :     #Reading in the output bounds file
2046 :     my ($ReactionTB,$CompoundTB,$DeadCompounds,$DeadEndCompounds,$DeadReactions);
2047 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsReactionData0.txt") {
2048 :     $ReactionTB = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsReactionData0.txt",";","|",1,["DATABASE ID"]);
2049 :     }
2050 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsCompoundData0.txt") {
2051 :     $CompoundTB = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsCompoundData0.txt",";","|",1,["DATABASE ID"]);
2052 :     }
2053 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadReactions.txt") {
2054 :     $DeadReactions = $self->figmodel()->put_array_in_hash($self->figmodel()->database()->load_single_column_file($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadReactions.txt",""));
2055 :     }
2056 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadMetabolites.txt") {
2057 :     $DeadCompounds = $self->figmodel()->put_array_in_hash($self->figmodel()->database()->load_single_column_file($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadMetabolites.txt",""));
2058 :     }
2059 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadEndMetabolites.txt") {
2060 :     $DeadEndCompounds = $self->figmodel()->put_array_in_hash($self->figmodel()->database()->load_single_column_file($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadEndMetabolites.txt",""));
2061 :     }
2062 :     if (!defined($ReactionTB) && !defined($CompoundTB)) {
2063 :     print STDERR "FIGMODEL:ClassifyModelReactions: Classification file not found when classifying reactions in ".$self->id().$self->selected_version()." with ".$Media." media. Most likely the model did not grow.\n";
2064 :     return (undef,undef);
2065 :     }
2066 :    
2067 :     #Clearing output
2068 :     $self->figmodel()->clearing_output($UniqueFilename,"Classify-".$self->id().$self->selected_version()."-".$UniqueFilename.".log");
2069 :     #Creating the table objects that will hold the results of the reaction classification
2070 :     my $rxnclasstable = $self->reaction_class_table();
2071 :     my $cpdclasstable = $self->compound_class_table();
2072 :     #Loading the compound table
2073 :     if (defined($CompoundTB)) {
2074 :     for (my $i=0; $i < $CompoundTB->size(); $i++) {
2075 :     my $Row = $CompoundTB->get_row($i);
2076 :     if (defined($Row->{"DATABASE ID"})) {
2077 :     #Getting the compound row
2078 :     my $CpdRow = $cpdclasstable->get_row_by_key($Row->{"DATABASE ID"}->[0].$Row->{COMPARTMENT}->[0],"COMPOUND",1);
2079 :     #Setting row values
2080 :     my $Max = 0;
2081 :     my $Min = 0;
2082 :     my $Class = "Unknown";
2083 :     if (defined($DeadCompounds) && defined($DeadCompounds->{$Row->{"DATABASE ID"}->[0]})) {
2084 :     $Class = "Dead";
2085 :     } elsif (defined($DeadEndCompounds) && defined($DeadEndCompounds->{$Row->{"DATABASE ID"}->[0]})) {
2086 :     $Class = "Deadend";
2087 :     } elsif (defined($Row->{"Min DRAIN_FLUX"}) && defined($Row->{"Max DRAIN_FLUX"}) && $Row->{"Min DRAIN_FLUX"}->[0] ne "1e+07") {
2088 :     $Max = $Row->{"Max DRAIN_FLUX"}->[0];
2089 :     $Min = $Row->{"Min DRAIN_FLUX"}->[0];
2090 :     if ($Row->{"Min DRAIN_FLUX"}->[0] > 0.00000001) {
2091 :     $Class = "Positive";
2092 :     } elsif ($Row->{"Max DRAIN_FLUX"}->[0] < -0.00000001) {
2093 :     $Class = "Negative";
2094 :     } elsif ($Row->{"Min DRAIN_FLUX"}->[0] < -0.00000001) {
2095 :     if ($Row->{"Max DRAIN_FLUX"}->[0] > 0.00000001) {
2096 :     $Class = "Variable";
2097 :     } else {
2098 :     $Max = 0;
2099 :     $Class = "Negative variable";
2100 :     }
2101 :     } elsif ($Row->{"Max DRAIN_FLUX"}->[0] > 0.00000001) {
2102 :     $Min = 0;
2103 :     $Class = "Positive variable";
2104 :     } else {
2105 :     $Min = 0;
2106 :     $Max = 0;
2107 :     $Class = "Blocked";
2108 :     }
2109 :     }
2110 :     my $index = 0;
2111 :     if (defined($CpdRow->{MEDIA})) {
2112 :     for (my $i=0; $i < @{$CpdRow->{MEDIA}};$i++) {
2113 :     $index++;
2114 :     if ($CpdRow->{MEDIA}->[$i] eq $Media) {
2115 :     $index = $i;
2116 :     last;
2117 :     }
2118 :     }
2119 :     }
2120 :     $CpdRow->{MIN}->[$index] = $Min;
2121 :     $CpdRow->{MAX}->[$index] = $Max;
2122 :     $CpdRow->{CLASS}->[$index] = $Class;
2123 :     $CpdRow->{MEDIA}->[$index] = $Media;
2124 :     }
2125 :     }
2126 :     $cpdclasstable->save();
2127 :     }
2128 :     if (defined($ReactionTB)) {
2129 :     for (my $i=0; $i < $ReactionTB->size(); $i++) {
2130 :     my $Row = $ReactionTB->get_row($i);
2131 :     if (defined($Row->{"DATABASE ID"})) {
2132 :     #Getting the compound row
2133 :     my $Compartment = "c";
2134 :     if (defined($Row->{COMPARTMENT}->[0])) {
2135 :     $Compartment = $Row->{COMPARTMENT}->[0];
2136 :     }
2137 :     my $RxnRow = $rxnclasstable->get_row_by_key($Row->{"DATABASE ID"}->[0],"REACTION",1);
2138 :     my $Max = 0;
2139 :     my $Min = 0;
2140 :     my $Class = "Unknown";
2141 :     if (defined($DeadReactions) && defined($DeadReactions->{$Row->{"DATABASE ID"}->[0]})) {
2142 :     $Class = "Dead";
2143 :     } elsif (defined($Row->{"Min FLUX"}) && defined($Row->{"Max FLUX"})) {
2144 :     $Max = $Row->{"Max FLUX"}->[0];
2145 :     $Min = $Row->{"Min FLUX"}->[0];
2146 :     if ($Row->{"Min FLUX"}->[0] > 0.00000001) {
2147 :     $Class = "Positive";
2148 :     } elsif ($Row->{"Max FLUX"}->[0] < -0.00000001) {
2149 :     $Class = "Negative";
2150 :     } elsif ($Row->{"Min FLUX"}->[0] < -0.00000001) {
2151 :     if ($Row->{"Max FLUX"}->[0] > 0.00000001) {
2152 :     $Class = "Variable";
2153 :     } else {
2154 :     $Max = 0;
2155 :     $Class = "Negative variable";
2156 :     }
2157 :     } elsif ($Row->{"Max FLUX"}->[0] > 0.00000001) {
2158 :     $Min = 0;
2159 :     $Class = "Positive variable";
2160 :     } else {
2161 :     $Min = 0;
2162 :     $Max = 0;
2163 :     $Class = "Blocked";
2164 :     }
2165 :     }
2166 :     my $index = 0;
2167 :     if (defined($RxnRow->{MEDIA})) {
2168 :     for (my $i=0; $i < @{$RxnRow->{MEDIA}};$i++) {
2169 :     $index++;
2170 :     if ($RxnRow->{MEDIA}->[$i] eq $Media) {
2171 :     $index = $i;
2172 :     last;
2173 :     }
2174 :     }
2175 :     }
2176 :     $RxnRow->{MIN}->[$index] = $Min;
2177 :     $RxnRow->{MAX}->[$index] = $Max;
2178 :     $RxnRow->{CLASS}->[$index] = $Class;
2179 :     $RxnRow->{MEDIA}->[$index] = $Media;
2180 :     }
2181 :     }
2182 :     $rxnclasstable->save();
2183 :     }
2184 :     return ($rxnclasstable,$cpdclasstable);
2185 :     }
2186 :    
2187 :     =head3 RunAllStudiesWithDataFast
2188 :     Definition:
2189 :     (integer::false positives,integer::false negatives,integer::correct negatives,integer::correct positives,string::error vector,string heading vector) = FIGMODELmodel->RunAllStudiesWithDataFast(string::experiment,0/1::print result);
2190 :     Description:
2191 :     Simulates every experimental condition currently available for the model.
2192 :     =cut
2193 :    
2194 :     sub RunAllStudiesWithDataFast {
2195 :     my ($self,$Experiment,$PrintResults) = @_;
2196 :    
2197 :     #Printing lp and key file for model
2198 :     if (!-e $self->directory()."FBA-".$self->id().$self->selected_version().".lp") {
2199 :     $self->PrintModelLPFile();
2200 :     }
2201 :     my $UniqueFilename = $self->figmodel()->filename();
2202 :    
2203 :     #Determing the simulations that need to be run
2204 :     my $ExperimentalDataTable = $self->figmodel()->GetExperimentalDataTable($self->genome(),$Experiment);
2205 :     #Creating the table of jobs to submit
2206 :     my $JobArray = $self->GetSimulationJobTable($ExperimentalDataTable,$Experiment,$UniqueFilename);
2207 :     #Printing the job file
2208 :     if (!-d $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/") {
2209 :     system("mkdir ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/");
2210 :     }
2211 :     $JobArray->save();
2212 :    
2213 :     #Running simulations
2214 :     system($self->config("mfalite executable")->[0]." ".$self->config("Reaction database directory")->[0]."masterfiles/MediaTable.txt ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Jobfile.txt ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt");
2215 :     print "simulation complete for ".$self->id().$self->selected_version()."\n";
2216 :     #Parsing the results
2217 :     my $Results = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt",";","\\|",0,undef);
2218 :     if (!defined($Results)) {
2219 :     $self->figmodel()->error_message("FIGMODELmodel:RunAllStudiesWithDataFast:Could not find simulation results: ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt");
2220 :     return undef;
2221 :     }
2222 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector,$SimulationResults) = $self->EvaluateSimulationResults($Results,$ExperimentalDataTable);
2223 :     #Printing results to file
2224 :     $self->figmodel()->database()->save_table($SimulationResults,undef,undef,undef,"False negatives\tFalse positives\tCorrect negatives\tCorrect positives\n".$FalseNegatives."\t".$FalsePostives."\t".$CorrectNegatives."\t".$CorrectPositives."\n");
2225 :     $self->figmodel()->clearing_output($UniqueFilename);
2226 :    
2227 :     return ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector);
2228 :     }
2229 :    
2230 :     =head3 GetSimulationJobTable
2231 :     Definition:
2232 :     my $JobTable = $model->GetSimulationJobTable($Experiment,$PrintResults,$Version);
2233 :     Description:
2234 :     =cut
2235 :    
2236 :     sub GetSimulationJobTable {
2237 :     my ($self,$SimulationTable,$Experiment,$Folder) = @_;
2238 :    
2239 :     #Determing the simulations that need to be run
2240 :     if (!defined($SimulationTable)) {
2241 :     $SimulationTable = $self->figmodel()->GetExperimentalDataTable($self->genome(),$Experiment);
2242 :     if (!defined($SimulationTable)) {
2243 :     return undef;
2244 :     }
2245 :     }
2246 :    
2247 :     #Creating the job table
2248 :     my $JobTable = $self->figmodel()->CreateJobTable($Folder);
2249 :     for (my $i=0; $i < $SimulationTable->size(); $i++) {
2250 :     if ($SimulationTable->get_row($i)->{"Heading"}->[0] =~ m/Gene\sKO/) {
2251 :     my $Row = $JobTable->get_row_by_key("Gene KO","LABEL",1);
2252 :     $JobTable->add_data($Row,"MEDIA",$SimulationTable->get_row($i)->{"Media"}->[0],1);
2253 :     } elsif ($SimulationTable->get_row($i)->{"Heading"}->[0] =~ m/Media\sgrowth/) {
2254 :     my $Row = $JobTable->get_row_by_key("Growth phenotype","LABEL",1);
2255 :     $JobTable->add_data($Row,"MEDIA",$SimulationTable->get_row($i)->{"Media"}->[0],1);
2256 :     } elsif ($SimulationTable->get_row($i)->{"Heading"}->[0] =~ m/Interval\sKO/) {
2257 :     my $Row = $JobTable->get_row_by_key($SimulationTable->get_row($i)->{"Heading"}->[0],"LABEL",1);
2258 :     $JobTable->add_data($Row,"MEDIA",$SimulationTable->get_row($i)->{"Media"}->[0],1);
2259 :     $JobTable->add_data($Row,"GENE KO",$SimulationTable->get_row($i)->{"Experiment type"}->[0],1);
2260 :     }
2261 :     }
2262 :    
2263 :     #Filling in model specific elements of the job table
2264 :     for (my $i=0; $i < $JobTable->size(); $i++) {
2265 :     if ($JobTable->get_row($i)->{"LABEL"}->[0] =~ m/Gene\sKO/) {
2266 :     $JobTable->get_row($i)->{"RUNTYPE"}->[0] = "SINGLEKO";
2267 :     $JobTable->get_row($i)->{"SAVE NONESSENTIALS"}->[0] = 1;
2268 :     } else {
2269 :     $JobTable->get_row($i)->{"RUNTYPE"}->[0] = "GROWTH";
2270 :     $JobTable->get_row($i)->{"SAVE NONESSENTIALS"}->[0] = 0;
2271 :     }
2272 :     $JobTable->get_row($i)->{"LP FILE"}->[0] = $self->directory()."FBA-".$self->id().$self->selected_version();
2273 :     $JobTable->get_row($i)->{"MODEL"}->[0] = $self->directory().$self->id().$self->selected_version().".txt";
2274 :     $JobTable->get_row($i)->{"SAVE FLUXES"}->[0] = 0;
2275 :     }
2276 :    
2277 :     return $JobTable;
2278 :     }
2279 :    
2280 :     =head3 EvaluateSimulationResults
2281 :     Definition:
2282 :     (integer::false positives,integer::false negatives,integer::correct negatives,integer::correct positives,string::error vector,string heading vector,FIGMODELtable::simulation results) = FIGMODELmodel->EvaluateSimulationResults(FIGMODELtable::raw simulation results,FIGMODELtable::experimental data);
2283 :     Description:
2284 :     Compares simulation results with experimental data to produce a table indicating where predictions are incorrect.
2285 :     =cut
2286 :    
2287 :     sub EvaluateSimulationResults {
2288 :     my ($self,$Results,$ExperimentalDataTable) = @_;
2289 :    
2290 :     #Comparing experimental results with simulation results
2291 :     my $SimulationResults = FIGMODELTable->new(["Run result","Experiment type","Media","Experiment ID","Reactions knocked out"],$self->directory()."SimulationOutput".$self->id().$self->selected_version().".txt",["Experiment ID","Media"],"\t",",",undef);
2292 :     my $FalsePostives = 0;
2293 :     my $FalseNegatives = 0;
2294 :     my $CorrectNegatives = 0;
2295 :     my $CorrectPositives = 0;
2296 :     my @Errorvector;
2297 :     my @HeadingVector;
2298 :     my $ReactionKOWithGeneHash;
2299 :     for (my $i=0; $i < $Results->size(); $i++) {
2300 :     if ($Results->get_row($i)->{"LABEL"}->[0] eq "Gene KO") {
2301 :     if (defined($Results->get_row($i)->{"REACTION KO WITH GENES"})) {
2302 :     for (my $j=0; $j < @{$Results->get_row($i)->{"REACTION KO WITH GENES"}}; $j++) {
2303 :     my @Temp = split(/:/,$Results->get_row($i)->{"REACTION KO WITH GENES"}->[$j]);
2304 :     if (defined($Temp[1]) && length($Temp[1]) > 0) {
2305 :     $ReactionKOWithGeneHash->{$Temp[0]} = $Temp[1];
2306 :     }
2307 :     }
2308 :     }
2309 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] == 0) {
2310 :     for (my $j=0; $j < @{$Results->get_row($i)->{"NONESSENTIALGENES"}}; $j++) {
2311 :     my $Row = $ExperimentalDataTable->get_row_by_key("Gene KO:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j],"Heading");
2312 :     if (defined($Row)) {
2313 :     my $KOReactions = "none";
2314 :     if (defined($ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]})) {
2315 :     $KOReactions = $ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]};
2316 :     }
2317 :     push(@HeadingVector,$Row->{"Heading"}->[0].":".$KOReactions);
2318 :     my $Status = "Unknown";
2319 :     if ($Row->{"Growth"}->[0] > 0) {
2320 :     $Status = "False negative";
2321 :     $FalseNegatives++;
2322 :     push(@Errorvector,3);
2323 :     } else {
2324 :     $Status = "False positive";
2325 :     $FalsePostives++;
2326 :     push(@Errorvector,2);
2327 :     }
2328 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Gene KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => [$KOReactions]});
2329 :     }
2330 :     }
2331 :     } else {
2332 :     for (my $j=0; $j < @{$Results->get_row($i)->{"ESSENTIALGENES"}}; $j++) {
2333 :     #print $j."\t".$Results->get_row($i)->{"ESSENTIALGENES"}->[$j]."\n";
2334 :     my $Row = $ExperimentalDataTable->get_row_by_key("Gene KO:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"ESSENTIALGENES"}->[$j],"Heading");
2335 :     if (defined($Row)) {
2336 :     my $KOReactions = "none";
2337 :     if (defined($ReactionKOWithGeneHash->{$Results->get_row($i)->{"ESSENTIALGENES"}->[$j]})) {
2338 :     $KOReactions = $ReactionKOWithGeneHash->{$Results->get_row($i)->{"ESSENTIALGENES"}->[$j]};
2339 :     }
2340 :     push(@HeadingVector,$Row->{"Heading"}->[0].":".$KOReactions);
2341 :     my $Status = "Unknown";
2342 :     if ($Row->{"Growth"}->[0] > 0) {
2343 :     $Status = "False negative";
2344 :     $FalseNegatives++;
2345 :     push(@Errorvector,3);
2346 :     } else {
2347 :     $Status = "Correct negative";
2348 :     $CorrectNegatives++;
2349 :     push(@Errorvector,1);
2350 :     }
2351 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Gene KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => [$KOReactions]});
2352 :     }
2353 :     }
2354 :     for (my $j=0; $j < @{$Results->get_row($i)->{"NONESSENTIALGENES"}}; $j++) {
2355 :     my $Row = $ExperimentalDataTable->get_row_by_key("Gene KO:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j],"Heading");
2356 :     if (defined($Row)) {
2357 :     my $KOReactions = "none";
2358 :     if (defined($ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]})) {
2359 :     $KOReactions = $ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]};
2360 :     }
2361 :     push(@HeadingVector,$Row->{"Heading"}->[0].":".$KOReactions);
2362 :     my $Status = "Unknown";
2363 :     if ($Row->{"Growth"}->[0] > 0) {
2364 :     $Status = "Correct positive";
2365 :     $CorrectPositives++;
2366 :     push(@Errorvector,0);
2367 :     } else {
2368 :     $Status = "False positive";
2369 :     $FalsePostives++;
2370 :     push(@Errorvector,2);
2371 :     }
2372 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Gene KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => [$KOReactions]});
2373 :     }
2374 :     }
2375 :     }
2376 :     } elsif ($Results->get_row($i)->{"LABEL"}->[0] eq "Growth phenotype") {
2377 :     my $Row = $ExperimentalDataTable->get_row_by_key("Media growth:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"MEDIA"}->[0],"Heading");
2378 :     if (defined($Row)) {
2379 :     push(@HeadingVector,$Row->{"Heading"}->[0].":none");
2380 :     my $Status = "Unknown";
2381 :     if ($Row->{"Growth"}->[0] > 0) {
2382 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2383 :     $Status = "Correct positive";
2384 :     $CorrectPositives++;
2385 :     push(@Errorvector,0);
2386 :     } else {
2387 :     $Status = "False negative";
2388 :     $FalseNegatives++;
2389 :     push(@Errorvector,3);
2390 :     }
2391 :     } else {
2392 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2393 :     $Status = "False positive";
2394 :     $FalsePostives++;
2395 :     push(@Errorvector,2);
2396 :     } else {
2397 :     $Status = "Correct negative";
2398 :     $CorrectNegatives++;
2399 :     push(@Errorvector,1);
2400 :     }
2401 :     }
2402 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Media growth"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Media"}->[0]],"Reactions knocked out" => ["none"]});
2403 :     }
2404 :     } elsif ($Results->get_row($i)->{"LABEL"}->[0] =~ m/Interval\sKO/ && defined($Results->get_row($i)->{"KOGENES"}->[0])) {
2405 :     my $Row = $ExperimentalDataTable->get_row_by_key($Results->get_row($i)->{"LABEL"}->[0],"Heading");
2406 :     if (defined($Row)) {
2407 :     my $Status = "Unknown";
2408 :     if ($Row->{"Growth"}->[0] > 0) {
2409 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2410 :     $Status = "Correct positive";
2411 :     $CorrectPositives++;
2412 :     push(@Errorvector,0);
2413 :     } else {
2414 :     $Status = "False negative";
2415 :     $FalseNegatives++;
2416 :     push(@Errorvector,3);
2417 :     }
2418 :     } else {
2419 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2420 :     $Status = "False positive";
2421 :     $FalsePostives++;
2422 :     push(@Errorvector,2);
2423 :     } else {
2424 :     $Status = "Correct negative";
2425 :     $CorrectNegatives++;
2426 :     push(@Errorvector,1);
2427 :     }
2428 :     }
2429 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Interval KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => ["none"]});
2430 :     }
2431 :     }
2432 :     }
2433 :    
2434 :     return ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,join(";",@Errorvector),join(";",@HeadingVector),$SimulationResults);
2435 :     }
2436 :    
2437 :     =head3 InspectSolution
2438 :     Definition:
2439 :     $model->InspectSolution(string::gene knocked out,string::media condition,[string]::list of reactions);
2440 :     Description:
2441 :     =cut
2442 :    
2443 :     sub InspectSolution {
2444 :     my ($self,$GeneKO,$Media,$ReactionList) = @_;
2445 :    
2446 :     #Getting a directory for the results
2447 :     my $UniqueFilename = $self->figmodel()->filename();
2448 :     system("mkdir ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/");
2449 :     my $TempVersion = "V".$UniqueFilename;
2450 :    
2451 :     #Setting gene ko to none if no genes are to be knocked out
2452 :     if ($GeneKO !~ m/^peg\./) {
2453 :     $GeneKO = "none";
2454 :     }
2455 :    
2456 :     #Implementing the input solution in the test model
2457 :     my $ReactionArray;
2458 :     my $DirectionArray;
2459 :     my %SolutionHash;
2460 :     for (my $k=0; $k < @{$ReactionList}; $k++) {
2461 :     if ($ReactionList->[$k] =~ m/(.+)(rxn\d\d\d\d\d)/) {
2462 :     my $Reaction = $2;
2463 :     my $Sign = $1;
2464 :     if (defined($SolutionHash{$Reaction})) {
2465 :     $SolutionHash{$Reaction} = "<=>";
2466 :     } elsif ($Sign eq "-") {
2467 :     $SolutionHash{$Reaction} = "<=";
2468 :     } elsif ($Sign eq "+") {
2469 :     $SolutionHash{$Reaction} = "=>";
2470 :     } else {
2471 :     $SolutionHash{$Reaction} = $Sign;
2472 :     }
2473 :     }
2474 :     }
2475 :     my @TempList = keys(%SolutionHash);
2476 :     for (my $k=0; $k < @TempList; $k++) {
2477 :     push(@{$ReactionArray},$TempList[$k]);
2478 :     push(@{$DirectionArray},$SolutionHash{$TempList[$k]});
2479 :     }
2480 :    
2481 :     print "Integrating solution!\n";
2482 :     $self->figmodel()->IntegrateGrowMatchSolution($self->id().$self->selected_version(),$self->directory().$self->id().$TempVersion.".txt",$ReactionArray,$DirectionArray,"SolutionInspection",1,1);
2483 :    
2484 :     #Printing lp and key file for model
2485 :     $self->PrintModelLPFile();
2486 :    
2487 :     #Running FBA on the test model
2488 :     my $JobTable = $self->figmodel()->CreateJobTable($UniqueFilename);
2489 :     $JobTable->add_row({"LABEL" => ["TEST"],"RUNTYPE" => ["GROWTH"],"LP FILE" => [$self->directory()."FBA-".$self->id().$TempVersion],"MODEL" => [$self->directory().$self->id().$TempVersion.".txt"],"MEDIA" => [$Media],"REACTION KO" => ["none|".join("|",@{$ReactionList})],"GENE KO" => [$GeneKO],"SAVE FLUXES" => [0],"SAVE NONESSENTIALS" => [0]});
2490 :     $JobTable->save();
2491 :    
2492 :     #Running simulations
2493 :     system($self->config("mfalite executable")->[0]." ".$self->config("Reaction database directory")->[0]."masterfiles/MediaTable.txt ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Jobfile.txt ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt");
2494 :    
2495 :     #Parsing the results
2496 :     my $Results = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt",";","\\|",0,undef);
2497 :     if (!defined($Results)) {
2498 :     $self->figmodel()->error_message("FIGMODELmodel:InspectSolution:Could not load problem report ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt");
2499 :     return undef;
2500 :     }
2501 :    
2502 :     #Making sure that the model grew with all reactions present
2503 :     my $Found = 0;
2504 :     for (my $i=0; $i < $Results->size(); $i++) {
2505 :     if (defined($Results->get_row($i)->{"KOGENES"}->[0]) && defined($Results->get_row($i)->{"KOREACTIONS"}->[0]) && $Results->get_row($i)->{"KOREACTIONS"}->[0] eq "none" && $Results->get_row($i)->{"KOGENES"}->[0] eq $GeneKO && $Results->get_row($i)->{"OBJECTIVE"}->[0] > 0.00001) {
2506 :     $Found = 1;
2507 :     }
2508 :     }
2509 :     if ($Found == 0) {
2510 :     print "Solution no longer valid\n";
2511 :     return undef;
2512 :     }
2513 :    
2514 :     #Making sure all of the reactions added are still necessary
2515 :     my $FinalReactionList;
2516 :     for (my $k=0; $k < $Results->size(); $k++) {
2517 :     if (defined($Results->get_row($k)->{"KOGENES"}->[0]) && $Results->get_row($k)->{"KOGENES"}->[0] eq $GeneKO) {
2518 :     if (defined($Results->get_row($k)->{"KOREACTIONS"}->[0]) && $Results->get_row($k)->{"KOREACTIONS"}->[0] =~ m/rxn\d\d\d\d\d/ && $Results->get_row($k)->{"OBJECTIVE"}->[0] < 0.000001) {
2519 :     push(@{$FinalReactionList},$Results->get_row($k)->{"KOREACTIONS"}->[0]);
2520 :     }
2521 :     }
2522 :     }
2523 :    
2524 :     #Deleting extra files created
2525 :     unlink($self->directory()."FBA-".$self->id().$TempVersion.".lp");
2526 :     unlink($self->directory()."FBA-".$self->id().$TempVersion.".key");
2527 :     unlink($self->directory().$self->id().$TempVersion.".txt");
2528 :    
2529 :     #Deleting the test model and the MFA folder
2530 :     $self->figmodel()->clearing_output($UniqueFilename);
2531 :    
2532 :     return $FinalReactionList;
2533 :     }
2534 :    
2535 :     =head3 GapFillingAlgorithm
2536 :    
2537 :     Definition:
2538 :     FIGMODELmodel->GapFillingAlgorithm();
2539 :    
2540 :     Description:
2541 :     This is a wrapper for running the gap filling algorithm on any model in the database.
2542 :     The algorithm performs a gap filling for any false negative prediction of the avialable experimental data.
2543 :     This function is threaded to improve efficiency: one thread does nothing but using the MFAToolkit to fill gaps for every false negative prediction.
2544 :     The other thread reads in the gap filling solutions, builds a test model for each solution, and runs the test model against all available experimental data.
2545 :     This function prints two important output files in the Model directory:
2546 :     1.) GapFillingOutput.txt: this is a summary of the results of the gap filling analysis
2547 :     2.) GapFillingErrorMatrix.txt: this lists the correct and incorrect predictions for each gapfilling solution implemented in a test model.
2548 :     =cut
2549 :    
2550 :     sub GapFillingAlgorithm {
2551 :     my ($self) = @_;
2552 :    
2553 :     #First the input model version and model filename should be simulated and the false negatives identified
2554 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector) = $self->RunAllStudiesWithDataFast("All");
2555 :    
2556 :     #Getting the filename
2557 :     my $UniqueFilename = $self->figmodel()->filename();
2558 :    
2559 :     #Printing the original performance vector
2560 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-OPEM".".txt",[$HeadingVector,$Errorvector]);
2561 :    
2562 :     my $PreviousGapFilling;
2563 :     if (-e $self->directory().$self->id().$self->selected_version()."-GFS.txt") {
2564 :     #Backing up the old solution file
2565 :     system("cp ".$self->directory().$self->id().$self->selected_version()."-GFS.txt ".$self->directory().$self->id().$self->selected_version()."-OldGFS.txt");
2566 :     unlink($self->directory().$self->id().$self->selected_version()."-GFS.txt");
2567 :     }
2568 :     if (-e $self->directory().$self->id().$self->selected_version()."-OldGFS.txt") {
2569 :     #Reading in the solution file from the previous gap filling if it exists
2570 :     $PreviousGapFilling = $self->figmodel()->database()->load_table($self->directory().$self->id().$self->selected_version()."-OldGFS.txt",";",",",0,["Experiment"]);
2571 :     }
2572 :    
2573 :     #Now we use the simulation output to make the gap filling run data
2574 :     my @Errors = split(/;/,$Errorvector);
2575 :     my @Headings = split(/;/,$HeadingVector);
2576 :     my $GapFillingRunSpecs = "";
2577 :     my $Count = 0;
2578 :     my $RescuedPreviousResults;
2579 :     my $RunCount = 0;
2580 :     my $SolutionExistedCount = 0;
2581 :     my $AcceptedSolutions = 0;
2582 :     my $RejectedSolutions = 0;
2583 :     my $NoExistingSolutions = 0;
2584 :     for (my $i=0; $i < @Errors; $i++) {
2585 :     if ($Errors[$i] == 3) {
2586 :     my @HeadingDataArray = split(/:/,$Headings[$i]);
2587 :     if ($HeadingDataArray[2] !~ m/^peg\./ || $HeadingDataArray[3] ne "none") {
2588 :     my $SolutionFound = 0;
2589 :     if (defined($PreviousGapFilling) && defined($PreviousGapFilling->get_row_by_key($HeadingDataArray[2],"Experiment"))) {
2590 :     my @Rows = $PreviousGapFilling->get_rows_by_key($HeadingDataArray[2],"Experiment");
2591 :     for (my $j=0; $j < @Rows; $j++) {
2592 :     if ($HeadingDataArray[2] =~ m/^peg\./) {
2593 :     my $ReactionList = $self->InspectSolution($HeadingDataArray[2],$HeadingDataArray[1],$Rows[$j]->{"Solution reactions"});
2594 :     if (defined($ReactionList)) {
2595 :     print join(",",@{$Rows[$j]->{"Solution reactions"}})."\t".join(",",@{$ReactionList})."\n";
2596 :     $SolutionFound++;
2597 :     push(@{$RescuedPreviousResults},$Rows[$j]->{"Experiment"}->[0].";".$Rows[$j]->{"Solution index"}->[0].";".$Rows[$j]->{"Solution cost"}->[0].";".join(",",@{$ReactionList}));
2598 :     $AcceptedSolutions++;
2599 :     } else {
2600 :     $RejectedSolutions++;
2601 :     }
2602 :     } else {
2603 :     my $ReactionList = $self->InspectSolution($HeadingDataArray[2],$HeadingDataArray[1],$Rows[$j]->{"Solution reactions"});
2604 :     if (defined($ReactionList)) {
2605 :     print join(",",@{$Rows[$j]->{"Solution reactions"}})."\t".join(",",@{$ReactionList})."\n";
2606 :     $SolutionFound++;
2607 :     push(@{$RescuedPreviousResults},$Rows[$j]->{"Experiment"}->[0].";".$Rows[$j]->{"Solution index"}->[0].";".$Rows[$j]->{"Solution cost"}->[0].";".join(",",@{$ReactionList}));
2608 :     $AcceptedSolutions++;
2609 :     } else {
2610 :     $RejectedSolutions++;
2611 :     }
2612 :     }
2613 :     }
2614 :     } else {
2615 :     $NoExistingSolutions++;
2616 :     }
2617 :     if ($SolutionFound == 0) {
2618 :     $RunCount++;
2619 :     if (length($GapFillingRunSpecs) > 0) {
2620 :     $GapFillingRunSpecs .= ";";
2621 :     }
2622 :     $GapFillingRunSpecs .= $HeadingDataArray[2].":".$HeadingDataArray[1].".txt:".$HeadingDataArray[3];
2623 :     } else {
2624 :     $SolutionExistedCount++;
2625 :     }
2626 :     }
2627 :     $Count++;
2628 :     }
2629 :     }
2630 :    
2631 :     #Updating the growmatch progress table
2632 :     my $Row = $self->figmodel()->database()->get_row_by_key("GROWMATCH TABLE",$self->genome(),"ORGANISM",1);
2633 :     $Row->{"INITIAL FP"}->[0] = $FalsePostives;
2634 :     $Row->{"INITIAL FN"}->[0] = $FalseNegatives;
2635 :     $Row->{"GF TIMING"}->[0] = time()."-";
2636 :     $Row->{"FN WITH SOL"}->[0] = $FalseNegatives-$NoExistingSolutions;
2637 :     $Row->{"FN WITH ACCEPTED SOL"}->[0] = $SolutionExistedCount;
2638 :     $Row->{"TOTAL ACCEPTED GF SOL"}->[0] = $AcceptedSolutions;
2639 :     $Row->{"TOTAL REJECTED GF SOL"}->[0] = $RejectedSolutions;
2640 :     $Row->{"FN WITH NO SOL"}->[0] = $NoExistingSolutions+$RejectedSolutions;
2641 :     $self->figmodel()->database()->update_row("GROWMATCH TABLE",$Row,"ORGANISM");
2642 :    
2643 :     #Running the gap filling once to correct all false negative errors
2644 :     my $SolutionsFound = 0;
2645 :     my $GapFillingArray;
2646 :     push(@{$GapFillingArray},split(/;/,$GapFillingRunSpecs));
2647 :     $self->datagapfill($GapFillingArray);
2648 :    
2649 :     if (defined($RescuedPreviousResults) && @{$RescuedPreviousResults} > 0) {
2650 :     #Printing previous solutions to GFS file
2651 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-GFS.txt",$RescuedPreviousResults,1);
2652 :     $SolutionsFound = 1;
2653 :     }
2654 :    
2655 :     #Recording the finishing of the gapfilling
2656 :     $Row = $self->figmodel()->database()->get_row_by_key("GROWMATCH TABLE",$self->genome(),"ORGANISM",1);
2657 :     $Row->{"GF TIMING"}->[0] .= time();
2658 :     $self->figmodel()->database()->update_row("GROWMATCH TABLE",$Row,"ORGANISM");
2659 :    
2660 :     if ($SolutionsFound == 1) {
2661 :     #Scheduling solution testing
2662 :     $self->figmodel()->add_job_to_queue("testsolutions?".$self->id().$self->selected_version()."?-1?GF","QSUB","fast","master","BACK");
2663 :     } else {
2664 :     $self->figmodel()->error_message("No false negative predictions found. Data gap filling not necessary!");
2665 :     }
2666 :    
2667 :     return $self->success();
2668 :     }
2669 :    
2670 :     =head3 BuildSpecificBiomassReaction
2671 :     Definition:
2672 :     FIGMODELmodel->BuildSpecificBiomassReaction();
2673 :     Description:
2674 :     =cut
2675 :     sub BuildSpecificBiomassReaction {
2676 :     my ($self) = @_;
2677 :    
2678 :     my $biomassrxn;
2679 :     my $OrganismID = $self->genome();
2680 :     #Checking for a biomass override
2681 :     if (defined($self->config("biomass reaction override")->{$OrganismID})) {
2682 :     $biomassrxn = $self->config("biomass reaction override")->{$OrganismID};
2683 :     print "Overriding biomass template and selecting ".$biomassrxn." for ".$OrganismID.".\n";
2684 :     } else {#Creating biomass reaction from the template
2685 :     #Getting the genome stats
2686 :     my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
2687 :     my $Class = $genomestats->{CLASS}->[0];
2688 :     my $Name = $genomestats->{NAME}->[0];
2689 :    
2690 :     #Checking for phoenix variants
2691 :     my $PhoenixVariantTable = $self->figmodel()->database()->GetDBTable("Phoenix variants table");
2692 :     my $Phoenix = 0;
2693 :     my @Rows = $PhoenixVariantTable->get_rows_by_key($OrganismID,"GENOME");
2694 :     my $VariantHash;
2695 :     for (my $i=0; $i < @Rows; $i++) {
2696 :     $Phoenix = 1;
2697 :     if (defined($Rows[$i]->{"SUBSYSTEM"}) && defined($Rows[$i]->{"VARIANT"})) {
2698 :     $VariantHash->{$Rows[$i]->{"SUBSYSTEM"}->[0]} = $Rows[$i]->{"VARIANT"}->[0];
2699 :     }
2700 :     }
2701 :    
2702 :     #Collecting genome data
2703 :     my $RoleHash;
2704 :     my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
2705 :     for (my $i=0; $i < $FeatureTable->size(); $i++) {
2706 :     for (my $j=0; $j < @{$FeatureTable->get_row($i)->{"ROLES"}}; $j++) {
2707 :     $RoleHash->{$FeatureTable->get_row($i)->{"ROLES"}->[$j]} = 1;
2708 :     my $Subsystems = $self->figmodel()->subsystems_of_role($FeatureTable->get_row($i)->{"ROLES"}->[$j]);
2709 :     if (defined($Subsystems)) {
2710 :     for (my $k=0; $k < @{$Subsystems}; $k++) {
2711 :     if ($Phoenix == 0) {
2712 :     $VariantHash->{$Subsystems->[$k]} = 1;
2713 :     }
2714 :     }
2715 :     }
2716 :     }
2717 :     }
2718 :    
2719 :     #Scanning through the template item by item and determinine which biomass components should be added
2720 :     my $ComponentTypes;
2721 :     my $EquationData;
2722 :     my $BiomassReactionTemplateTable = $self->figmodel()->database()->GetDBTable("BIOMASS TEMPLATE");
2723 :     for (my $i=0; $i < $BiomassReactionTemplateTable->size(); $i++) {
2724 :     my $Row = $BiomassReactionTemplateTable->get_row($i);
2725 :     if (defined($Row->{"INCLUSION CRITERIA"}->[0]) && $Row->{"INCLUSION CRITERIA"}->[0] eq "UNIVERSAL") {
2726 :     push(@{$EquationData},$Row);
2727 :     $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;
2728 :     } else {
2729 :     my $Criteria = $Row->{"INCLUSION CRITERIA"}->[0];
2730 :     my $End = 0;
2731 :     while ($End == 0) {
2732 :     if ($Criteria =~ m/^(.+)(AND)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(AND)\{([^{^}]+)\}$/ || $Criteria =~ m/^(.+)(OR)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(OR)\{([^{^}]+)\}$/) {
2733 :     print $Criteria."\n";
2734 :     my $Start = "";
2735 :     my $End = "";
2736 :     my $Condition = $1;
2737 :     my $Data = $2;
2738 :     if ($1 ne "AND" && $1 ne "OR") {
2739 :     $Start = $1;
2740 :     $End = $4;
2741 :     $Condition = $2;
2742 :     $Data = $3;
2743 :     }
2744 :     my $Result = "YES";
2745 :     if ($Condition eq "OR") {
2746 :     $Result = "NO";
2747 :     }
2748 :     my @Array = split(/\|/,$Data);
2749 :     for (my $j=0; $j < @Array; $j++) {
2750 :     if ($Array[$j] eq "YES" && $Condition eq "OR") {
2751 :     $Result = "YES";
2752 :     last;
2753 :     } elsif ($Array[$j] eq "NO" && $Condition eq "AND") {
2754 :     $Result = "NO";
2755 :     last;
2756 :     } elsif ($Array[$j] =~ m/^COMPOUND:(.+)/) {
2757 :     my $Match = 0;
2758 :     for (my $k=0; $k < @{$EquationData}; $k++) {
2759 :     if ($EquationData->[$k]->{"ID"}->[0] eq $1) {
2760 :     $Match = 1;
2761 :     last;
2762 :     }
2763 :     }
2764 :     if ($Match == 1 && $Condition eq "OR") {
2765 :     $Result = "YES";
2766 :     last;
2767 :     } elsif ($Match != 1 && $Condition eq "AND") {
2768 :     $Result = "NO";
2769 :     last;
2770 :     }
2771 :     } elsif ($Array[$j] =~ m/^!COMPOUND:(.+)/) {
2772 :     my $Match = 0;
2773 :     for (my $k=0; $k < @{$EquationData}; $k++) {
2774 :     if ($EquationData->[$k]->{"ID"}->[0] eq $1) {
2775 :     $Match = 1;
2776 :     last;
2777 :     }
2778 :     }
2779 :     if ($Match != 1 && $Condition eq "OR") {
2780 :     $Result = "YES";
2781 :     last;
2782 :     } elsif ($Match == 1 && $Condition eq "AND") {
2783 :     $Result = "NO";
2784 :     last;
2785 :     }
2786 :     } elsif ($Array[$j] =~ m/^NAME:(.+)/) {
2787 :     my $Comparison = $1;
2788 :     if ((!defined($Comparison) || !defined($Name) || $Name =~ m/$Comparison/) && $Condition eq "OR") {
2789 :     $Result = "YES";
2790 :     last;
2791 :     } elsif (defined($Comparison) && defined($Name) && $Name !~ m/$Comparison/ && $Condition eq "AND") {
2792 :     $Result = "NO";
2793 :     last;
2794 :     }
2795 :     } elsif ($Array[$j] =~ m/^!NAME:(.+)/) {
2796 :     my $Comparison = $1;
2797 :     if ((!defined($Comparison) || !defined($Name) || $Name !~ m/$Comparison/) && $Condition eq "OR") {
2798 :     $Result = "YES";
2799 :     last;
2800 :     } elsif (defined($Comparison) && defined($Name) && $Name =~ m/$Comparison/ && $Condition eq "AND") {
2801 :     $Result = "NO";
2802 :     last;
2803 :     }
2804 :     } elsif ($Array[$j] =~ m/^SUBSYSTEM:(.+)/) {
2805 :     my @SubsystemArray = split(/`/,$1);
2806 :     if (@SubsystemArray == 1) {
2807 :     if (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} ne -1 && $Condition eq "OR") {
2808 :     $Result = "YES";
2809 :     last;
2810 :     } elsif ((!defined($VariantHash->{$SubsystemArray[0]}) || $VariantHash->{$SubsystemArray[0]} eq -1) && $Condition eq "AND") {
2811 :     $Result = "NO";
2812 :     last;
2813 :     }
2814 :     } else {
2815 :     my $Match = 0;
2816 :     for (my $k=1; $k < @SubsystemArray; $k++) {
2817 :     if (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} eq $SubsystemArray[$k]) {
2818 :     $Match = 1;
2819 :     last;
2820 :     }
2821 :     }
2822 :     if ($Match == 1 && $Condition eq "OR") {
2823 :     $Result = "YES";
2824 :     last;
2825 :     } elsif ($Match != 1 && $Condition eq "AND") {
2826 :     $Result = "NO";
2827 :     last;
2828 :     }
2829 :     }
2830 :     } elsif ($Array[$j] =~ m/^!SUBSYSTEM:(.+)/) {
2831 :     my @SubsystemArray = split(/`/,$1);
2832 :     if (@SubsystemArray == 1) {
2833 :     if ((!defined($VariantHash->{$SubsystemArray[0]}) || $VariantHash->{$SubsystemArray[0]} eq -1) && $Condition eq "OR") {
2834 :     $Result = "YES";
2835 :     last;
2836 :     } elsif (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} ne -1 && $Condition eq "AND") {
2837 :     $Result = "NO";
2838 :     last;
2839 :     }
2840 :     } else {
2841 :     my $Match = 0;
2842 :     for (my $k=1; $k < @SubsystemArray; $k++) {
2843 :     if (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} eq $SubsystemArray[$k]) {
2844 :     $Match = 1;
2845 :     last;
2846 :     }
2847 :     }
2848 :     if ($Match != 1 && $Condition eq "OR") {
2849 :     $Result = "YES";
2850 :     last;
2851 :     } elsif ($Match == 1 && $Condition eq "AND") {
2852 :     $Result = "NO";
2853 :     last;
2854 :     }
2855 :     }
2856 :     } elsif ($Array[$j] =~ m/^ROLE:(.+)/) {
2857 :     if (defined($RoleHash->{$1}) && $Condition eq "OR") {
2858 :     $Result = "YES";
2859 :     last;
2860 :     } elsif (!defined($RoleHash->{$1}) && $Condition eq "AND") {
2861 :     $Result = "NO";
2862 :     last;
2863 :     }
2864 :     } elsif ($Array[$j] =~ m/^!ROLE:(.+)/) {
2865 :     if (!defined($RoleHash->{$1}) && $Condition eq "OR") {
2866 :     $Result = "YES";
2867 :     last;
2868 :     } elsif (defined($RoleHash->{$1}) && $Condition eq "AND") {
2869 :     $Result = "NO";
2870 :     last;
2871 :     }
2872 :     } elsif ($Array[$j] =~ m/^CLASS:(.+)/) {
2873 :     if ($Class eq $1 && $Condition eq "OR") {
2874 :     $Result = "YES";
2875 :     last;
2876 :     } elsif ($Class ne $1 && $Condition eq "AND") {
2877 :     $Result = "NO";
2878 :     last;
2879 :     }
2880 :     } elsif ($Array[$j] =~ m/^!CLASS:(.+)/) {
2881 :     if ($Class ne $1 && $Condition eq "OR") {
2882 :     $Result = "YES";
2883 :     last;
2884 :     } elsif ($Class eq $1 && $Condition eq "AND") {
2885 :     $Result = "NO";
2886 :     last;
2887 :     }
2888 :     }
2889 :     }
2890 :     $Criteria = $Start.$Result.$End;
2891 :     print $Criteria."\n";
2892 :     } else {
2893 :     $End = 1;
2894 :     last;
2895 :     }
2896 :     }
2897 :     if ($Criteria eq "YES") {
2898 :     push(@{$EquationData},$Row);
2899 :     $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;
2900 :     }
2901 :     }
2902 :     }
2903 :    
2904 :     #Building biomass equation
2905 :     my %Reactants;
2906 :     my %Products;
2907 :     foreach my $EquationRow (@{$EquationData}) {
2908 :     #First determine what the coefficient should be
2909 :     my $Coefficient;
2910 :     if (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^[0-9\.]+$/) {
2911 :     $Coefficient = $EquationRow->{"COEFFICIENT"}->[0];
2912 :     } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/cpd\d\d\d\d\d/) {
2913 :     $Coefficient = 0;
2914 :     my @CompoundList = split(/,/,$EquationRow->{"COEFFICIENT"}->[0]);
2915 :     foreach my $Compound (@CompoundList) {
2916 :     if (defined($Reactants{$Compound})) {
2917 :     $Coefficient += $Reactants{$Compound};
2918 :     }
2919 :     }
2920 :     } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^([0-9\.]+)\/(.+)$/) {
2921 :     my @Keys = keys(%{$ComponentTypes->{$2}});
2922 :     my $MW = 1;
2923 :     my $CompoundData = $self->figmodel()->LoadObject($EquationRow->{"ID"}->[0]);
2924 :     if (defined($CompoundData->{"MASS"})) {
2925 :     $MW = $CompoundData->{"MASS"}->[0];
2926 :     }
2927 :     $Coefficient = $1/@Keys/$MW;
2928 :     }
2929 :     if (defined($EquationRow->{"REACTANT"}) && $EquationRow->{"REACTANT"}->[0] eq "YES") {
2930 :     if (defined($Reactants{$EquationRow->{"ID"}->[0]})) {
2931 :     $Reactants{$EquationRow->{"ID"}->[0]} += $Coefficient;
2932 :     } elsif (defined($Products{$EquationRow->{"ID"}->[0]}) && $Products{$EquationRow->{"ID"}->[0]} > $Coefficient) {
2933 :     $Products{$EquationRow->{"ID"}->[0]} -= $Coefficient;
2934 :     } elsif (defined($Products{$EquationRow->{"ID"}->[0]}) && $Products{$EquationRow->{"ID"}->[0]} < $Coefficient) {
2935 :     $Reactants{$EquationRow->{"ID"}->[0]} = $Coefficient - $Products{$EquationRow->{"ID"}->[0]};
2936 :     delete $Products{$EquationRow->{"ID"}->[0]};
2937 :     } else {
2938 :     $Reactants{$EquationRow->{"ID"}->[0]} = $Coefficient;
2939 :     }
2940 :     } else {
2941 :     if (defined($Products{$EquationRow->{"ID"}->[0]})) {
2942 :     $Products{$EquationRow->{"ID"}->[0]} += $Coefficient;
2943 :     } elsif (defined($Reactants{$EquationRow->{"ID"}->[0]}) && $Reactants{$EquationRow->{"ID"}->[0]} > $Coefficient) {
2944 :     $Reactants{$EquationRow->{"ID"}->[0]} -= $Coefficient;
2945 :     } elsif (defined($Reactants{$EquationRow->{"ID"}->[0]}) && $Reactants{$EquationRow->{"ID"}->[0]} < $Coefficient) {
2946 :     $Products{$EquationRow->{"ID"}->[0]} = $Coefficient - $Reactants{$EquationRow->{"ID"}->[0]};
2947 :     delete $Reactants{$EquationRow->{"ID"}->[0]};
2948 :     } else {
2949 :     $Products{$EquationRow->{"ID"}->[0]} = $Coefficient;
2950 :     }
2951 :     }
2952 :     }
2953 :     my $Equation = "";
2954 :     my @ReactantList = sort(keys(%Reactants));
2955 :     for (my $i=0; $i < @ReactantList; $i++) {
2956 :     if (length($Equation) > 0) {
2957 :     $Equation .= " + ";
2958 :     }
2959 :     $Equation .= $self->figmodel()->format_coefficient($Reactants{$ReactantList[$i]})." ".$ReactantList[$i];
2960 :     }
2961 :     $Equation .= " => ";
2962 :     my $First = 1;
2963 :     @ReactantList = sort(keys(%Products));
2964 :     for (my $i=0; $i < @ReactantList; $i++) {
2965 :     if ($First == 0) {
2966 :     $Equation .= " + ";
2967 :     }
2968 :     $First = 0;
2969 :     $Equation .= $self->figmodel()->format_coefficient($Products{$ReactantList[$i]})." ".$ReactantList[$i];
2970 :     }
2971 :    
2972 :     #Adding the biomass equation to the biomass table
2973 :     my $NewRow = $self->figmodel()->add_biomass_reaction($Equation,$self->id(),"Template:".$self->genome());
2974 :     $biomassrxn = $NewRow->{DATABASE}->[0];
2975 :     print $biomassrxn."\n";
2976 :     }
2977 :     print $biomassrxn."\n";
2978 :     my $BiomassRow = $self->figmodel()->add_model_to_biomass_reaction($biomassrxn,$self->id());
2979 :     return $BiomassRow;
2980 :     }
2981 :    
2982 :     =head3 PrintSBMLFile
2983 :     Definition:
2984 :     FIGMODELmodel->PrintSBMLFile();
2985 :     Description:
2986 :     Printing file with model data in SBML format
2987 :     =cut
2988 :     sub PrintSBMLFile {
2989 :     my($self) = @_;
2990 :    
2991 :     #Opening the SBML file for printing
2992 :     my $Filename = $self->config("SBML files")->[0].$self->id().".xml";
2993 :     if ($self->owner() ne "master") {
2994 :     if (!-d $self->config("SBML files")->[0].$self->owner()."/") {
2995 :     system("mkdir ".$self->config("SBML files")->[0].$self->owner()."/");
2996 :     }
2997 :     $Filename = $self->config("SBML files")->[0].$self->owner()."/".$self->id().".xml";
2998 :     }
2999 :     if (!open (SBMLOUTPUT, ">$Filename")) {
3000 :     return;
3001 :     }
3002 :    
3003 :     #Loading and parsing the model data
3004 :     my $ModelTable = $self->reaction_table();
3005 :     if (!defined($ModelTable) || !defined($ModelTable->{"array"})) {
3006 :     print "Failed to load ".$self->id()."\n";
3007 :     return;
3008 :     }
3009 :    
3010 :     #Adding intracellular metabolites that also need exchange fluxes to the exchange hash
3011 :     my $ExchangeHash = {"cpd11416" => "c"};
3012 :    
3013 :     my %CompartmentsPresent;
3014 :     $CompartmentsPresent{"c"} = 1;
3015 :     my %CompoundList;
3016 :     my @ReactionList;
3017 :     my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");
3018 :     for (my $i=0; $i < $ModelTable->size(); $i++) {
3019 :     my $Reaction = $ModelTable->get_row($i)->{"LOAD"}->[0];
3020 :     if (defined($ReactionTable->get_row_by_key($Reaction,"DATABASE")) && defined($ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0])) {
3021 :     push(@ReactionList,$Reaction);
3022 :     $_ = $ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0];
3023 :     my @MatchArray = /(cpd\d\d\d\d\d)/g;
3024 :     for (my $j=0; $j < @MatchArray; $j++) {
3025 :     $CompoundList{$MatchArray[$j]}->{"c"} = 1;
3026 :     }
3027 :     $_ = $ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0];
3028 :     @MatchArray = /(cpd\d\d\d\d\d\[\D\])/g;
3029 :     for (my $j=0; $j < @MatchArray; $j++) {
3030 :     if ($MatchArray[$j] =~ m/(cpd\d\d\d\d\d)\[(\D)\]/) {
3031 :     $CompartmentsPresent{lc($2)} = 1;
3032 :     $CompoundList{$1}->{lc($2)} = 1;
3033 :     }
3034 :     }
3035 :     }
3036 :     }
3037 :    
3038 :     #Printing header to SBML file
3039 :     my $ModelName = $self->id();
3040 :     $ModelName =~ s/\./_/;
3041 :     print SBMLOUTPUT '<?xml version="1.0" encoding="UTF-8"?>'."\n";
3042 :     print SBMLOUTPUT '<sbml xmlns="http://www.sbml.org/sbml/level2" level="2" version="1" xmlns:html="http://www.w3.org/1999/xhtml">' . "\n";
3043 :     if (defined($self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Organism name"}->[0])) {
3044 :     print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Organism name"}->[0].' SEED model">'."\n";
3045 :     } else {
3046 :     print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->id().' SEED model">'."\n";
3047 :     }
3048 :    
3049 :     #Printing the unit data
3050 :     print SBMLOUTPUT "<listOfUnitDefinitions>\n";
3051 :     print SBMLOUTPUT "\t<unitDefinition id=\"mmol_per_gDW_per_hr\">\n";
3052 :     print SBMLOUTPUT "\t\t<listOfUnits>\n";
3053 :     print SBMLOUTPUT "\t\t\t<unit kind=\"mole\" scale=\"-3\"/>\n";
3054 :     print SBMLOUTPUT "\t\t\t<unit kind=\"gram\" exponent=\"-1\"/>\n";
3055 :     print SBMLOUTPUT "\t\t\t<unit kind=\"second\" multiplier=\".00027777\" exponent=\"-1\"/>\n";
3056 :     print SBMLOUTPUT "\t\t</listOfUnits>\n";
3057 :     print SBMLOUTPUT "\t</unitDefinition>\n";
3058 :     print SBMLOUTPUT "</listOfUnitDefinitions>\n";
3059 :    
3060 :     #Printing compartments for SBML file
3061 :     print SBMLOUTPUT '<listOfCompartments>'."\n";
3062 :     foreach my $Compartment (keys(%CompartmentsPresent)) {
3063 :     if (defined($self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0])) {
3064 :     my @OutsideList = split(/\//,$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Outside"}->[0]);
3065 :     my $Printed = 0;
3066 :     foreach my $Outside (@OutsideList) {
3067 :     if (defined($CompartmentsPresent{$Outside}) && defined($self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Outside}->[0]->{"Name"}->[0])) {
3068 :     print SBMLOUTPUT '<compartment id="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'" outside="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Outside}->[0]->{"Name"}->[0].'"/>'."\n";
3069 :     $Printed = 1;
3070 :     last;
3071 :     }
3072 :     }
3073 :     if ($Printed eq 0) {
3074 :     print SBMLOUTPUT '<compartment id="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'"/>'."\n";
3075 :     }
3076 :     }
3077 :     }
3078 :     print SBMLOUTPUT '</listOfCompartments>'."\n";
3079 :    
3080 :     #Printing the list of metabolites involved in the model
3081 :     print SBMLOUTPUT '<listOfSpecies>'."\n";
3082 :     foreach my $Compound (keys(%CompoundList)) {
3083 :     my $Name = $Compound;
3084 :     my $Formula = "";
3085 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0])) {
3086 :     $Formula = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0];
3087 :     }
3088 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0])) {
3089 :     $Name = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0];
3090 :     $Name =~ s/\s/_/;
3091 :     $Name .= "_".$Formula;
3092 :     }
3093 :     $Name =~ s/[<>;&\*]//;
3094 :     my $Charge = 0;
3095 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0])) {
3096 :     $Charge = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0];
3097 :     }
3098 :     foreach my $Compartment (keys(%{$CompoundList{$Compound}})) {
3099 :     if ($Compartment eq "e") {
3100 :     $ExchangeHash->{$Compound} = "e";
3101 :     }
3102 :     print SBMLOUTPUT '<species id="'.$Compound.'_'.$Compartment.'" name="'.$Name.'" compartment="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'" charge="'.$Charge.'" boundaryCondition="false"/>'."\n";
3103 :     }
3104 :     }
3105 :     #Printing the boundary species
3106 :     foreach my $Compound (keys(%{$ExchangeHash})) {
3107 :     my $Name = $Compound;
3108 :     my $Formula = "";
3109 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0])) {
3110 :     $Formula = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0];
3111 :     }
3112 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0])) {
3113 :     $Name = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0];
3114 :     $Name =~ s/\s/_/;
3115 :     $Name .= "_".$Formula;
3116 :     }
3117 :     $Name =~ s/[<>;&\*]//;
3118 :     my $Charge = 0;
3119 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0])) {
3120 :     $Charge = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0];
3121 :     }
3122 :     print SBMLOUTPUT '<species id="'.$Compound.'_b" name="'.$Name.'" compartment="Extracellular" charge="'.$Charge.'" boundaryCondition="true"/>'."\n";
3123 :     }
3124 :     print SBMLOUTPUT '</listOfSpecies>'."\n";
3125 :    
3126 :     #Printing the list of reactions involved in the model
3127 :     my $ObjectiveCoef;
3128 :     print SBMLOUTPUT '<listOfReactions>'."\n";
3129 :     foreach my $Reaction (@ReactionList) {
3130 :     $ObjectiveCoef = "0.0";
3131 :     if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"EQUATION"}->[0])) {
3132 :     if ($Reaction =~ m/^bio/) {
3133 :     $ObjectiveCoef = "1.0";
3134 :     }
3135 :     my $LowerBound = -10000;
3136 :     my $UpperBound = 10000;
3137 :     my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateData($Reaction);
3138 :     my $Name = $Reaction;
3139 :     if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"NAME"}->[0])) {
3140 :     $Name = $self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"NAME"}->[0];
3141 :     $Name =~ s/[<>;&]//g;
3142 :     }
3143 :     my $Reversibility = "true";
3144 :     if (defined($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0])) {
3145 :     if ($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0] ne "<=>") {
3146 :     $LowerBound = 0;
3147 :     $Reversibility = "false";
3148 :     }
3149 :     if ($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0] eq "<=") {
3150 :     my $Temp = $Products;
3151 :     $Products = $Reactants;
3152 :     $Reactants = $Temp;
3153 :     }
3154 :     }
3155 :     print SBMLOUTPUT '<reaction id="'.$Reaction.'" name="'.$Name.'" reversible="'.$Reversibility.'">'."\n";
3156 :     print SBMLOUTPUT "<notes>\n";
3157 :     my $ECData = "";
3158 :     if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"ENZYME"}->[0])) {
3159 :     $ECData = $self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"ENZYME"}->[0];
3160 :     }
3161 :     my $SubsystemData = "";
3162 :     if (defined($ModelTable->{$Reaction}->[0]->{"SUBSYSTEM"}->[0])) {
3163 :     $SubsystemData = $ModelTable->{$Reaction}->[0]->{"SUBSYSTEM"}->[0];
3164 :     }
3165 :     my $GeneAssociation = "";
3166 :     my $ProteinAssociation = "";
3167 :     if (defined($ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0])) {
3168 :     if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} == 1 && $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0] !~ m/\+/) {
3169 :     $GeneAssociation = $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0];
3170 :     } else {
3171 :     if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} > 1) {
3172 :     $GeneAssociation = "( ";
3173 :     }
3174 :     for (my $i=0; $i < @{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}}; $i++) {
3175 :     if ($i > 0) {
3176 :     $GeneAssociation .= " ) or ( ";
3177 :     }
3178 :     $GeneAssociation .= $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[$i];
3179 :     }
3180 :     if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} > 1) {
3181 :     $GeneAssociation .= " )";
3182 :     }
3183 :     }
3184 :     $GeneAssociation =~ s/\+/ and /g;
3185 :     if ($GeneAssociation =~ m/\sor\s/ || $GeneAssociation =~ m/\sand\s/) {
3186 :     $GeneAssociation = "( ".$GeneAssociation." )";
3187 :     }
3188 :     $ProteinAssociation = $GeneAssociation;
3189 :     if (defined($self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Genome ID"}->[0])) {
3190 :     $ProteinAssociation = $self->figmodel()->translate_gene_to_protein($ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0],$self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Genome ID"}->[0]);
3191 :     }
3192 :     }
3193 :     print SBMLOUTPUT "<html:p>GENE_ASSOCIATION:".$GeneAssociation."</html:p>\n";
3194 :     print SBMLOUTPUT "<html:p>PROTEIN_ASSOCIATION:".$ProteinAssociation."</html:p>\n";
3195 :     print SBMLOUTPUT "<html:p>SUBSYSTEM: ".$SubsystemData."</html:p>\n";
3196 :     print SBMLOUTPUT "<html:p>PROTEIN_CLASS: ".$ECData."</html:p>\n";
3197 :     print SBMLOUTPUT "</notes>\n";
3198 :     print SBMLOUTPUT "<listOfReactants>\n";
3199 :     foreach my $Reactant (@{$Reactants}) {
3200 :     print SBMLOUTPUT '<speciesReference species="'.$Reactant->{"DATABASE"}->[0]."_".$Reactant->{"COMPARTMENT"}->[0].'" stoichiometry="'.$Reactant->{"COEFFICIENT"}->[0].'"/>'."\n";
3201 :     }
3202 :     print SBMLOUTPUT "</listOfReactants>\n";
3203 :     print SBMLOUTPUT "<listOfProducts>\n";
3204 :     foreach my $Product (@{$Products}) {
3205 :     print SBMLOUTPUT '<speciesReference species="'.$Product->{"DATABASE"}->[0]."_".$Product->{"COMPARTMENT"}->[0].'" stoichiometry="'.$Product->{"COEFFICIENT"}->[0].'"/>'."\n";
3206 :     }
3207 :     print SBMLOUTPUT "</listOfProducts>\n";
3208 :     print SBMLOUTPUT "<kineticLaw>\n";
3209 :     print SBMLOUTPUT "\t<math xmlns=\"http://www.w3.org/1998/Math/MathML\">\n";
3210 :     print SBMLOUTPUT "\t\t\t<ci> FLUX_VALUE </ci>\n";
3211 :     print SBMLOUTPUT "\t</math>\n";
3212 :     print SBMLOUTPUT "\t<listOfParameters>\n";
3213 :     print SBMLOUTPUT "\t\t<parameter id=\"LOWER_BOUND\" value=\"".$LowerBound."\" units=\"mmol_per_gDW_per_hr\"/>\n";
3214 :     print SBMLOUTPUT "\t\t<parameter id=\"UPPER_BOUND\" value=\"".$UpperBound."\" units=\"mmol_per_gDW_per_hr\"/>\n";
3215 :     print SBMLOUTPUT "\t\t<parameter id=\"OBJECTIVE_COEFFICIENT\" value=\"".$ObjectiveCoef."\"/>\n";
3216 :     print SBMLOUTPUT "\t\t<parameter id=\"FLUX_VALUE\" value=\"0.0\" units=\"mmol_per_gDW_per_hr\"/>\n";
3217 :     print SBMLOUTPUT "\t</listOfParameters>\n";
3218 :     print SBMLOUTPUT "</kineticLaw>\n";
3219 :     print SBMLOUTPUT '</reaction>'."\n";
3220 :     }
3221 :     }
3222 :    
3223 :     my @ExchangeList = keys(%{$ExchangeHash});
3224 :     foreach my $ExCompound (@ExchangeList) {
3225 :     my $ExCompoundName = $ExCompound;
3226 :     my $Row = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->get_row_by_key($ExCompound,"DATABASE");
3227 :     if (defined($Row) && defined($Row->{"NAME"}->[0])) {
3228 :     $ExCompoundName = $Row->{"NAME"}->[0];
3229 :     $ExCompoundName =~ s/[<>;&]//g;
3230 :     }
3231 :     $ObjectiveCoef = "0.0";
3232 :     print SBMLOUTPUT '<reaction id="EX_'.$ExCompound.'_'.$ExchangeHash->{$ExCompound}.'" name="EX_'.$ExCompoundName.'_'.$ExchangeHash->{$ExCompound}.'" reversible="true">'."\n";
3233 :     print SBMLOUTPUT "\t".'<notes>'."\n";
3234 :     print SBMLOUTPUT "\t\t".'<html:p>GENE_ASSOCIATION: </html:p>'."\n";
3235 :     print SBMLOUTPUT "\t\t".'<html:p>PROTEIN_ASSOCIATION: </html:p>'."\n";
3236 :     print SBMLOUTPUT "\t\t".'<html:p>SUBSYSTEM: S_</html:p>'."\n";
3237 :     print SBMLOUTPUT "\t\t".'<html:p>PROTEIN_CLASS: </html:p>'."\n";
3238 :     print SBMLOUTPUT "\t".'</notes>'."\n";
3239 :     print SBMLOUTPUT "\t".'<listOfReactants>'."\n";
3240 :     print SBMLOUTPUT "\t\t".'<speciesReference species="'.$ExCompound.'_'.$ExchangeHash->{$ExCompound}.'" stoichiometry="1.000000"/>'."\n";
3241 :     print SBMLOUTPUT "\t".'</listOfReactants>'."\n";
3242 :     print SBMLOUTPUT "\t".'<listOfProducts>'."\n";
3243 :     print SBMLOUTPUT "\t\t".'<speciesReference species="'.$ExCompound.'_b" stoichiometry="1.000000"/>'."\n";
3244 :     print SBMLOUTPUT "\t".'</listOfProducts>'."\n";
3245 :     print SBMLOUTPUT "\t".'<kineticLaw>'."\n";
3246 :     print SBMLOUTPUT "\t\t".'<math xmlns="http://www.w3.org/1998/Math/MathML">'."\n";
3247 :     print SBMLOUTPUT "\t\t\t\t".'<ci> FLUX_VALUE </ci>'."\n";
3248 :     print SBMLOUTPUT "\t\t".'</math>'."\n";
3249 :     print SBMLOUTPUT "\t\t".'<listOfParameters>'."\n";
3250 :     print SBMLOUTPUT "\t\t\t".'<parameter id="LOWER_BOUND" value="-10000.000000" units="mmol_per_gDW_per_hr"/>'."\n";
3251 :     print SBMLOUTPUT "\t\t\t".'<parameter id="UPPER_BOUND" value="10000.000000" units="mmol_per_gDW_per_hr"/>'."\n";
3252 :     print SBMLOUTPUT "\t\t\t".'<parameter id="OBJECTIVE_COEFFICIENT" value="'.$ObjectiveCoef.'"/>'."\n";
3253 :     print SBMLOUTPUT "\t\t\t".'<parameter id="FLUX_VALUE" value="0.000000" units="mmol_per_gDW_per_hr"/>'."\n";
3254 :     print SBMLOUTPUT "\t\t".'</listOfParameters>'."\n";
3255 :     print SBMLOUTPUT "\t".'</kineticLaw>'."\n";
3256 :     print SBMLOUTPUT '</reaction>'."\n";
3257 :     }
3258 :    
3259 :     #Closing out the file
3260 :     print SBMLOUTPUT '</listOfReactions>'."\n";
3261 :     print SBMLOUTPUT '</model>'."\n";
3262 :     print SBMLOUTPUT "</sbml>\n";
3263 :     close(SBMLOUTPUT);
3264 :     }
3265 :    
3266 :     =head3 PrintModelLPFile
3267 :     Definition:
3268 :     success()/fail() FIGMODELmodel->PrintModelLPFile();
3269 :     Description:
3270 :     Prints the lp file needed to run the model using the mpifba program
3271 :     =cut
3272 :     sub PrintModelLPFile {
3273 :     my ($self) = @_;
3274 :     #Printing lp and key file for model
3275 :     my $UniqueFilename = $self->figmodel()->filename();
3276 :     #Printing the standard FBA file
3277 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),"NoBounds",["ProdFullFBALP"],undef,$self->id().$self->selected_version()."-LPPrint.log",undef,$self->selected_version()));
3278 :     system("cp ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/CurrentProblem.lp ".$self->directory()."FBA-".$self->id().$self->selected_version().".lp");
3279 :     my $KeyTable = FIGMODELTable::load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/VariableKey.txt",";","|",0,undef);
3280 :     if (!defined($KeyTable)) {
3281 :     print STDERR "FIGMODEL:RunAllStudiesWithDataFast: ".$self->id()." LP file could not be printed.\n";
3282 :     return 0;
3283 :     }
3284 :     $KeyTable->headings(["Variable type","Variable ID"]);
3285 :     $KeyTable->save($self->directory()."FBA-".$self->id().$self->selected_version().".key");
3286 :     unlink($self->config("database message file directory")->[0].$self->id().$self->selected_version()."-LPPrint.log");
3287 :     $self->figmodel()->clearing_output($UniqueFilename,"FBA-".$self->id().$self->selected_version().".lp");
3288 :     }
3289 :    
3290 :     1;

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3