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

Annotation of /FigKernelPackages/FIGMODELmodel.pm

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : chenry 1.18 use strict;
2 :     package FIGMODELmodel;
3 : chenry 1.12
4 : chenry 1.18 =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 :    
18 :     #Error checking first
19 :     if (!defined($figmodel)) {
20 :     print STDERR "FIGMODELmodel->new(undef,".$id."):figmodel must be defined to create a model object!\n";
21 :     return undef;
22 :     }
23 :     my $self = {_figmodel => $figmodel};
24 :     bless $self;
25 :     if (!defined($id)) {
26 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,undef):id must be defined to create a model object");
27 :     return undef;
28 :     }
29 :    
30 :     #Checking that the id exists
31 :     my $modelHandle = $self->figmodel()->database()->get_object_manager("model");
32 :     if (!defined($modelHandle)) {
33 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not load MODELS table. Check database!");
34 :     return undef;
35 :     }
36 :    
37 :     #If the id is a number, we get the model row by index
38 :     if ($id =~ m/^\d+$/) {
39 :     my $objects = $modelHandle->get_objects();
40 :     $self->{_data} = $objects->[$id];
41 :     $self->figmodel()->{_models}->{$id} = $self;
42 :     } else {
43 :     my $objects = $modelHandle->get_objects({id => $id});
44 :     if (!defined($objects->[0]) && $id =~ m/(.+)(V[^V]+)/) {
45 :     $objects = $modelHandle->get_objects({id => $1});
46 :     if (!defined($objects->[0])) {
47 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find model ".$id." in database!");
48 :     return undef;
49 :     }
50 :     $self->{_selectedversion} = $2;
51 :     } elsif (!defined($objects->[0])) {
52 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find model ".$id." in database!");
53 :     return undef;
54 :     }
55 :     $self->{_data} = $objects->[0];
56 :     }
57 :     if (!defined($self->{_data})) {
58 :     $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find specified id in database!");
59 :     return undef;
60 :     }
61 :     $self->figmodel()->{_models}->{$self->id()} = $self;
62 :    
63 :     return $self;
64 :     }
65 :    
66 :     =head3 config
67 :     Definition:
68 :     ref::key value = FIGMODELmodel->config(string::key);
69 :     Description:
70 :     Trying to avoid using calls that assume configuration data is stored in a particular manner.
71 :     Call this function to get file paths etc.
72 :     =cut
73 :     sub config {
74 :     my ($self,$key) = @_;
75 :     return $self->figmodel()->config($key);
76 :     }
77 :    
78 :     =head3 fail
79 :     Definition:
80 :     -1 = FIGMODELmodel->fail();
81 :     Description:
82 :     Standard return for failed functions.
83 :     =cut
84 :     sub fail {
85 :     my ($self) = @_;
86 :     return $self->figmodel()->fail();
87 :     }
88 :    
89 :     =head3 success
90 :     Definition:
91 :     1 = FIGMODELmodel->success();
92 :     Description:
93 :     Standard return for successful functions.
94 :     =cut
95 :     sub success {
96 :     my ($self) = @_;
97 :     return $self->figmodel()->success();
98 :     }
99 :    
100 :     =head3 figmodel
101 :     Definition:
102 :     FIGMODEL = FIGMODELmodel->figmodel();
103 :     Description:
104 :     Returns a FIGMODEL object
105 :     =cut
106 :     sub figmodel {
107 :     my ($self) = @_;
108 :     return $self->{"_figmodel"};
109 :     }
110 :    
111 :     =head3 fig
112 :     Definition:
113 :     FIGMODEL = FIGMODELmodel->fig();
114 :     Description:
115 :     Returns a FIG object
116 :     =cut
117 :     sub fig {
118 :     my ($self) = @_;
119 :     if (!defined($self->{_fig}) && $self->source() !~ /^MGRAST/) {
120 :     if ($self->source() =~ /^RAST/) {
121 :     $self->{"_fig"} = $self->figmodel()->fig();#$self->genome());
122 :     } else {
123 :     $self->{"_fig"} = $self->figmodel()->fig();
124 :     }
125 :     }
126 :     return $self->{"_fig"};
127 :     }
128 :    
129 :     =head3 aquireModelLock
130 :    
131 :     Definition:
132 :    
133 :     FIGMODELmodel->aquireModelLock();
134 :    
135 :     Description:
136 :    
137 :     Locks the database for alterations relating to the current model object
138 :    
139 :     =cut
140 :     sub aquireModelLock {
141 :     my ($self) = @_;
142 :     $self->figmodel()->database()->genericLock($self->id());
143 :     }
144 :    
145 :     =head3 releaseModelLock
146 :    
147 :     Definition:
148 :    
149 :     FIGMODELmodel->releaseModelLock();
150 :    
151 :     Description:
152 :    
153 :     Unlocks the database for alterations relating to the current model object
154 :    
155 :     =cut
156 :     sub releaseModelLock {
157 :     my ($self) = @_;
158 :     $self->figmodel()->database()->genericUnlock($self->id());
159 :     }
160 :    
161 :     =head3 mgdata
162 :     Definition:
163 :     FIGMODEL = FIGMODELmodel->mgdata();
164 :     Description:
165 :     Returns a mgrast database object
166 :     =cut
167 :     sub mgdata {
168 :     my ($self) = @_;
169 :     if (!defined($self->{_mgdata}) && $self->source() =~ /^MGRAST/) {
170 :     require MGRAST;
171 :     $self->{_mgdata} = $self->figmodel()->mgrast()->Job->get_objects( { 'genome_id' => $self->genome() } )
172 :     }
173 :     return $self->{_mgdata};
174 :     }
175 :    
176 :     =head3 id
177 :     Definition:
178 :     string = FIGMODELmodel->id();
179 :     Description:
180 :     Returns model id
181 :     =cut
182 :     sub id {
183 :     my ($self) = @_;
184 :     return $self->{_data}->id();
185 :     }
186 :    
187 :     =head3 owner
188 :     Definition:
189 :     string = FIGMODELmodel->owner();
190 :     Description:
191 :     Returns the username for the model owner
192 :     =cut
193 :     sub owner {
194 :     my ($self) = @_;
195 :     return $self->{_data}->owner();
196 :     }
197 :    
198 :     =head3 name
199 :     Definition:
200 :     string = FIGMODELmodel->name();
201 :     Description:
202 :     Returns the name of the organism or metagenome sample being modeled
203 :     =cut
204 :     sub name {
205 :     my ($self) = @_;
206 :    
207 :     if (!defined($self->{_name})) {
208 :     my $source = $self->source();
209 :     if ($source =~ /^MGRAST/) {
210 :     $self->{_name} = "NA";
211 :     if (defined($self->mgdata())) {
212 :     $self->{_name} = $self->mgdata()->genome_name;
213 :     }
214 :     } elsif ($source !~ /^RAST/) {
215 :     $self->{_name} = $self->fig()->orgname_of_orgid($self->genome());
216 :     } else {
217 :     $self->{_name} = $self->figmodel()->get_genome_stats($self->genome())->{NAME}->[0];
218 :     }
219 :     }
220 :    
221 :     return $self->{_name};
222 :     }
223 :    
224 :     =head3 get_reaction_class
225 :     Definition:
226 :     string = FIGMODELmodel->get_reaction_class(string::reaction ID);
227 :     Description:
228 :     Returns reaction class
229 :     =cut
230 :     sub get_reaction_class {
231 :     my ($self,$reaction,$nohtml,$brief_flux) = @_;
232 :    
233 :     if (!-e $self->directory()."ReactionClassification-".$self->id().".tbl") {
234 :     if (!defined($self->{_reaction_classes})) {
235 :     $self->{_reaction_classes} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id()."-Complete.tbl",";","|",0,["REACTION"]);
236 :     if (!defined($self->{_reaction_classes})) {
237 :     return undef;
238 :     }
239 :     }
240 :    
241 :     my $ClassRow = $self->{_reaction_classes}->get_row_by_key($reaction,"REACTION");
242 :     if (defined($ClassRow) && defined($ClassRow->{CLASS})) {
243 :     my $class;
244 :     my $min = $ClassRow->{MIN}->[0];
245 :     my $max = $ClassRow->{MAX}->[0];
246 :     if ($ClassRow->{CLASS}->[0] eq "Positive") {
247 :     $class = "Essential =>";
248 :     $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
249 :     } elsif ($ClassRow->{CLASS}->[0] eq "Negative") {
250 :     $class = "Essential <=";
251 :     $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
252 :     } elsif ($ClassRow->{CLASS}->[0] eq "Positive variable") {
253 :     $class = "Active =>";
254 :     $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
255 :     } elsif ($ClassRow->{CLASS}->[0] eq "Negative variable") {
256 :     $class = "Active <=";
257 :     $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
258 :     } elsif ($ClassRow->{CLASS}->[0] eq "Variable") {
259 :     $class = "Active <=>";
260 :     $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
261 :     } elsif ($ClassRow->{CLASS}->[0] eq "Blocked") {
262 :     $class = "Inactive";
263 :     } elsif ($ClassRow->{CLASS}->[0] eq "Dead") {
264 :     $class = "Disconnected";
265 :     }
266 :    
267 :     if (!defined($nohtml) || $nohtml ne "1") {
268 :     $class = "<span title=\"Flux:".$min." to ".$max."\">".$class."</span>";
269 :     }
270 :    
271 :     return $class;
272 :     }
273 :     return undef;
274 :     }
275 :    
276 :     if (!defined($self->{_reaction_classes})) {
277 :     $self->{_reaction_classes} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id().".tbl",";","|",0,["REACTION"]);
278 :     if (!defined($self->{_reaction_classes})) {
279 :     return undef;
280 :     }
281 :     }
282 :    
283 :     my $ClassRow = $self->{_reaction_classes}->get_row_by_key($reaction,"REACTION");
284 :     my $classstring = "";
285 :     if (defined($ClassRow) && defined($ClassRow->{CLASS})) {
286 :     for (my $i=0; $i < @{$ClassRow->{CLASS}};$i++) {
287 :     if (length($classstring) > 0) {
288 :     $classstring .= "<br>";
289 :     }
290 :     my $NewClass;
291 :     my $min = $ClassRow->{MIN}->[$i];
292 :     my $max = $ClassRow->{MAX}->[$i];
293 :     if ($ClassRow->{CLASS}->[$i] eq "Positive") {
294 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Essential =>";
295 :     $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
296 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Negative") {
297 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Essential <=";
298 :     $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
299 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Positive variable") {
300 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Active =>";
301 :     $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
302 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Negative variable") {
303 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=";
304 :     $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
305 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Variable") {
306 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=>";
307 :     $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
308 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Blocked") {
309 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Inactive";
310 :     } elsif ($ClassRow->{CLASS}->[$i] eq "Dead") {
311 :     $NewClass = $ClassRow->{MEDIA}->[$i].":Disconnected";
312 :     }
313 :    
314 :     if (!defined($nohtml) || $nohtml ne "1") {
315 :     $NewClass = "<span title=\"Flux:".$min." to ".$max."\">".$NewClass."</span>";
316 :     }
317 :     $classstring .= $NewClass;
318 :     }
319 :     }
320 :     return $classstring;
321 :     }
322 :    
323 :     =head3 get_biomass
324 :     Definition:
325 :     string = FIGMODELmodel->get_biomass();
326 :     Description:
327 :     Returns data for the biomass reaction
328 :     =cut
329 :     sub get_biomass {
330 :     my ($self) = @_;
331 :    
332 :     if (!defined($self->{_biomass})) {
333 :     my $rxntbl = $self->reaction_table();
334 :     if (defined($rxntbl)) {
335 :     for (my $i=0; $i < $rxntbl->size(); $i++) {
336 :     if ($rxntbl->get_row($i)->{"LOAD"}->[0] =~ m/bio\d\d\d\d\d/) {
337 :     $self->{_biomass} = $rxntbl->get_row($i)->{"LOAD"}->[0];
338 :     last;
339 :     }
340 :     }
341 :     }
342 :     }
343 :    
344 :     return $self->get_reaction_data($self->{_biomass});
345 :     }
346 :    
347 :     =head3 get_reaction_data
348 :     Definition:
349 :     string = FIGMODELmodel->get_reaction_data(string::reaction ID);
350 :     Description:
351 :     Returns model reaction data
352 :     =cut
353 :     sub get_reaction_data {
354 :     my ($self,$reaction) = @_;
355 :    
356 :     if (!defined($self->reaction_table())) {
357 :     return undef;
358 :     }
359 :     if ($reaction =~ m/^\d+$/) {
360 :     $self->reaction_table()->get_row($reaction);
361 :     }
362 :     return $self->reaction_table()->get_row_by_key($reaction,"LOAD");
363 :     }
364 :    
365 :     =head3 get_reaction_id
366 :     Definition:
367 :     string = FIGMODELmodel->get_reaction_id(string::reaction ID);
368 :     Description:
369 :     Returns model reaction id
370 :     =cut
371 :     sub get_reaction_id {
372 :     my ($self,$reaction) = @_;
373 :    
374 :     my $Data = $self->get_reaction_data($reaction);
375 :     if (!defined($Data) || !defined($Data->{LOAD}->[0])) {
376 :     return undef;
377 :     }
378 :     return $Data->{LOAD}->[0];
379 :     }
380 :    
381 :     =head3 load_model_table
382 :    
383 :     Definition: FIGMODELTable = FIGMODELmodel->load_model_table(string:table name,0/1:refresh the table));
384 :    
385 :     Description: Returns the table specified by the input filename. Table will be stored in a file in the model directory.
386 :    
387 :     =cut
388 :     sub load_model_table {
389 :     my ($self,$name,$refresh) = @_;
390 :     if (!defined($refresh)) {
391 :     $refresh = 1;
392 :     }
393 :     if ($refresh == 1) {
394 :     delete $self->{"_".$name};
395 :     }
396 :     if (!defined($self->{"_".$name})) {
397 :     my $tbldef = $self->figmodel()->config("$name");
398 :     if (!defined($tbldef)) {
399 :     return undef;
400 :     }
401 :     my $filename = $self->directory().$name."-".$self->id().$self->selected_version().".tbl";
402 :     $self->{"_".$name} = $self->figmodel()->database()->load_table($filename,"\t","|",$tbldef->{headingline}->[0],$tbldef->{hashcolumns});
403 :     if (!defined($self->{"_".$name})) {
404 :     if (defined($tbldef->{prefix})) {
405 :     $self->{"_".$name} = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},"\t","|",join(@{$tbldef->{prefix}},"\n"));
406 :     } else {
407 :     $self->{"_".$name} = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},"\t","|");
408 :     }
409 :     }
410 :     }
411 :     return $self->{"_".$name};
412 :     }
413 :    
414 :     =head3 create_table_prototype
415 :    
416 :     Definition:
417 :     FIGMODELTable::table = FIGMODELmodel->create_table_prototype(string::table);
418 :     Description:
419 :     Returns a empty FIGMODELTable with all the metadata associated with the input table name
420 :    
421 :     =cut
422 :     sub create_table_prototype {
423 :     my ($self,$TableName) = @_;
424 :    
425 :     #Checking if the table definition exists in the FIGMODELconfig file
426 :     if (!defined($self->figmodel()->config($TableName))) {
427 :     $self->figmodel()->error_message("FIGMODELdatabase:create_table_prototype:Definition not found for ".$TableName);
428 :     return undef;
429 :     }
430 :     #Checking that this is a database table
431 :     if (!defined($self->config($TableName)->{tabletype}) || $self->config($TableName)->{tabletype}->[0] ne "ModelTable") {
432 :     $self->figmodel()->error_message("FIGMODELdatabase:create_table_prototype:".$TableName." is not a model table!");
433 :     return undef;
434 :     }
435 :     if (!defined($self->config($TableName)->{delimiter})) {
436 :     $self->config($TableName)->{delimiter}->[0] = ";";
437 :     }
438 :     if (!defined($self->config($TableName)->{itemdelimiter})) {
439 :     $self->config($TableName)->{itemdelimiter}->[0] = "|";
440 :     }
441 :     my $prefix;
442 :     if (defined($self->config($TableName)->{prefix})) {
443 :     $prefix = join("\n",@{$self->config($TableName)->{prefix}});
444 :     }
445 :     my $tbl = FIGMODELTable->new($self->config($TableName)->{columns},$self->directory().$self->config($TableName)->{filename_prefix}->[0]."-".$self->id().$self->selected_version().".txt",$self->config($TableName)->{hashcolumns},$self->config($TableName)->{delimiter}->[0],$self->config($TableName)->{itemdelimiter}->[0],$prefix);
446 :     return $tbl;
447 :     }
448 :    
449 :     =head3 get_reaction_number
450 :     Definition:
451 :     int = FIGMODELmodel->get_reaction_number();
452 :     Description:
453 :     Returns the number of reactions in the model
454 :     =cut
455 :     sub get_reaction_number {
456 :     my ($self) = @_;
457 :     if (!defined($self->reaction_table())) {
458 :     return 0;
459 :     }
460 :     return $self->reaction_table()->size();
461 :     }
462 :    
463 :     =head3 reaction_table
464 :     Definition:
465 :     FIGMODELTable = FIGMODELmodel->reaction_table();
466 :     Description:
467 :     Returns FIGMODELTable with the reaction list for the model
468 :     =cut
469 :     sub reaction_table {
470 :     my ($self) = @_;
471 :    
472 :     if (!defined($self->{_reaction_data})) {
473 :     $self->{_reaction_data} = $self->figmodel()->database()->GetDBModel($self->id());
474 :     my $classTbl = $self->reaction_class_table();
475 :     for (my $i=0; $i < $classTbl->size(); $i++) {
476 :     my $row = $classTbl->get_row($i);
477 :     my $rxnRow = $self->{_reaction_data}->get_row_by_key($row->{"REACTION"}->[0],"LOAD");
478 :     for (my $j=0; $j < @{$row->{MEDIA}};$j++) {
479 :     my $class = "Active <=>";
480 :     if ($row->{CLASS}->[$j] eq "Positive") {
481 :     $class = "Essential =>";
482 :     } elsif ($row->{CLASS}->[$j] eq "Negative") {
483 :     $class = "Essential <=";
484 :     } elsif ($row->{CLASS}->[$j] eq "Blocked") {
485 :     $class = "Inactive";
486 :     } elsif ($row->{CLASS}->[$j] eq "Positive variable") {
487 :     $class = "Active =>";
488 :     } elsif ($row->{CLASS}->[$j] eq "Negative variable") {
489 :     $class = "Active <=";
490 :     } elsif ($row->{CLASS}->[$j] eq "Variable") {
491 :     $class = "Active <=>";
492 :     } elsif ($row->{CLASS}->[$j] eq "Dead") {
493 :     $class = "Dead";
494 :     }
495 :     push(@{$rxnRow->{PREDICTIONS}},$row->{MEDIA}->[$j].":".$class);
496 :     }
497 :     }
498 :     }
499 :    
500 :     return $self->{_reaction_data};
501 :     }
502 :    
503 :     =head3 feature_table
504 :     Definition:
505 :     FIGMODELTable = FIGMODELmodel->feature_table();
506 :     Description:
507 :     Returns FIGMODELTable with the feature list for the model
508 :     =cut
509 :     sub feature_table {
510 :     my ($self) = @_;
511 :    
512 :     if (!defined($self->{_feature_data})) {
513 :     #Getting the genome feature list
514 :     my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
515 :     if (!defined($FeatureTable)) {
516 :     print STDERR "FIGMODELmodel:feature_table:Could not get features for genome ".$self->genome()." in database!";
517 :     return undef;
518 :     }
519 :     #Getting the reaction table for the model
520 :     my $rxnTable = $self->reaction_table();
521 :     if (!defined($rxnTable)) {
522 :     print STDERR "FIGMODELmodel:feature_table:Could not get reaction table for model ".$self->id()." in database!";
523 :     return undef;
524 :     }
525 :     #Cloning the feature table
526 :     $self->{_feature_data} = $FeatureTable->clone_table_def();
527 :     $self->{_feature_data}->add_headings(($self->id()."REACTIONS",$self->id()."PREDICTIONS"));
528 :     for (my $i=0; $i < $rxnTable->size(); $i++) {
529 :     my $Row = $rxnTable->get_row($i);
530 :     if (defined($Row) && defined($Row->{"ASSOCIATED PEG"})) {
531 :     foreach my $GeneSet (@{$Row->{"ASSOCIATED PEG"}}) {
532 :     my $temp = $GeneSet;
533 :     $temp =~ s/\+/|/g;
534 :     $temp =~ s/\sAND\s/|/gi;
535 :     $temp =~ s/\sOR\s/|/gi;
536 :     $temp =~ s/[\(\)\s]//g;
537 :     my @GeneList = split(/\|/,$temp);
538 :     foreach my $Gene (@GeneList) {
539 :     my $FeatureRow = $self->{_feature_data}->get_row_by_key("fig|".$self->genome().".".$Gene,"ID");
540 :     if (!defined($FeatureRow)) {
541 :     $FeatureRow = $FeatureTable->get_row_by_key("fig|".$self->genome().".".$Gene,"ID");
542 :     if (defined($FeatureRow)) {
543 :     $self->{_feature_data}->add_row($FeatureRow);
544 :     }
545 :     }
546 :     if (defined($FeatureRow)) {
547 :     $self->{_feature_data}->add_data($FeatureRow,$self->id()."REACTIONS",$Row->{"LOAD"}->[0],1);
548 :     }
549 :     }
550 :     }
551 :     }
552 :     }
553 :     #Loading predictions
554 :     my @files = glob($self->directory()."EssentialGenes-".$self->id()."-*");
555 :     foreach my $file (@files) {
556 :     if ($file =~ m/\-([^\-]+)\.tbl/) {
557 :     my $media = $1;
558 :     my $list = $self->figmodel()->database()->load_single_column_file($file,"");
559 :     my $hash;
560 :     foreach my $gene (@{$list}) {
561 :     $hash->{$gene} = 1;
562 :     }
563 :     for (my $i=0; $i < $self->{_feature_data}->size(); $i++) {
564 :     my $Row = $self->{_feature_data}->get_row($i);
565 :     if ($Row->{ID}->[0] =~ m/(peg\.\d+)/) {
566 :     my $gene = $1;
567 :     if (defined($hash->{$gene})) {
568 :     push(@{$Row->{$self->id()."PREDICTIONS"}},$media.":essential");
569 :     } else {
570 :     push(@{$Row->{$self->id()."PREDICTIONS"}},$media.":nonessential");
571 :     }
572 :     }
573 :     }
574 :     }
575 :     }
576 :     }
577 :     return $self->{_feature_data};
578 :     }
579 :    
580 :     =head3 reaction_class_table
581 :     Definition:
582 :     FIGMODELTable = FIGMODELmodel->reaction_class_table();
583 :     Description:
584 :     Returns FIGMODELTable with the reaction class data, and creates the table file if it does not exist
585 :     =cut
586 :     sub reaction_class_table {
587 :     my ($self) = @_;
588 :    
589 :     if (!defined($self->{_reaction_class_table})) {
590 :     if (-e $self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl") {
591 :     $self->{_reaction_class_table} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl",";","|",0,["REACTION","CLASS","MEDIA"]);
592 :     } else {
593 :     $self->{_reaction_class_table} = FIGMODELTable->new(["REACTION","MEDIA","CLASS","MIN","MAX"],$self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl",["REACTION","CLASS","MEDIA"],";","|",undef);
594 :     }
595 :     }
596 :    
597 :     return $self->{_reaction_class_table};
598 :     }
599 :    
600 :     =head3 compound_class_table
601 :     Definition:
602 :     FIGMODELTable = FIGMODELmodel->compound_class_table();
603 :     Description:
604 :     Returns FIGMODELTable with the compound class data, and creates the table file if it does not exist
605 :     =cut
606 :     sub compound_class_table {
607 :     my ($self) = @_;
608 :    
609 :     if (!defined($self->{_compound_class_table})) {
610 :     if (-e $self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl") {
611 :     $self->{_compound_class_table} = $self->figmodel()->database()->load_table($self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl",";","|",0,["COMPOUND","CLASS","MEDIA"]);
612 :     } else {
613 :     $self->{_compound_class_table} = FIGMODELTable->new(["COMPOUND","MEDIA","CLASS","MIN","MAX"],$self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl",["COMPOUND","CLASS","MEDIA"],";","|",undef);
614 :     }
615 :     }
616 :    
617 :     return $self->{_compound_class_table};
618 :     }
619 :    
620 :     =head3 get_essential_genes
621 :     Definition:
622 :     [string::peg ID] = FIGMODELmodel->get_essential_genes(string::media condition);
623 :     Description:
624 :     Returns an reference to an array of the predicted essential genes during growth in the input media condition
625 :     =cut
626 :     sub get_essential_genes {
627 :     my ($self,$media) = @_;
628 :    
629 :     if (!defined($media)) {
630 :     $media = "Complete";
631 :     }
632 :     if (!defined($self->{_essential_genes}->{$media})) {
633 :     $self->{_essential_genes}->{$media} = undef;
634 :     if (-e $self->directory()."EssentialGenes-".$self->id().$self->selected_version()."-".$media.".tbl") {
635 :     $self->{_essential_genes}->{$media} = $self->figmodel()->database()->load_single_column_file($self->directory()."EssentialGenes-".$self->id().$self->selected_version()."-".$media.".tbl","");
636 :     }
637 :     }
638 :    
639 :     return $self->{_essential_genes}->{$media};
640 :     }
641 :    
642 :     =head3 compound_table
643 :     Definition:
644 :     FIGMODELTable = FIGMODELmodel->compound_table();
645 :     Description:
646 :     Returns FIGMODELTable with the compound list for the model
647 :     =cut
648 :     sub compound_table {
649 :     my ($self) = @_;
650 :    
651 :     if (!defined($self->{_compound_table})) {
652 :     $self->{_compound_table} = $self->figmodel()->database()->GetDBModelCompounds($self->id());
653 :     }
654 :    
655 :     return $self->{_compound_table};
656 :     }
657 :    
658 :     =head3 get_compound_data
659 :     Definition:
660 :     {string:key=>[string]:values} = FIGMODELmodel->get_compound_data(string::compound ID);
661 :     Description:
662 :     Returns model compound data
663 :     =cut
664 :     sub get_compound_data {
665 :     my ($self,$compound) = @_;
666 :    
667 :     if (!defined($self->compound_table())) {
668 :     return undef;
669 :     }
670 :     if ($compound =~ m/^\d+$/) {
671 :     return $self->compound_table()->get_row($compound);
672 :     }
673 :     return $self->compound_table()->get_row_by_key($compound,"DATABASE");
674 :     }
675 :    
676 :     =head3 get_feature_data
677 :     Definition:
678 :     {string:key=>[string]:values} = FIGMODELmodel->get_feature_data(string::feature ID);
679 :     Description:
680 :     Returns model feature data
681 :     =cut
682 :     sub get_feature_data {
683 :     my ($self,$feature) = @_;
684 :     if (!defined($self->feature_table())) {
685 :     return undef;
686 :     }
687 :     if ($feature =~ m/^\d+$/) {
688 :     return $self->feature_table()->get_row($feature);
689 :     }
690 :     if ($feature =~ m/(peg\.\d+)/) {
691 :     $feature = $1;
692 :     }
693 :     return $self->feature_table()->get_row_by_key("fig|".$self->genome().".".$feature,"ID");
694 :     }
695 :    
696 :     =head3 genome
697 :     Definition:
698 :     string = FIGMODELmodel->genome();
699 :     Description:
700 :     Returns model genome
701 :     =cut
702 :     sub genome {
703 :     my ($self) = @_;
704 :     return $self->{_data}->genome();
705 :     }
706 :    
707 :     =head3 source
708 :     Definition:
709 :     string = FIGMODELmodel->source();
710 :     Description:
711 :     Returns model source
712 :     =cut
713 :     sub source {
714 :     my ($self) = @_;
715 :     return $self->{_data}->source();
716 :     }
717 :    
718 :     =head3 rights
719 :     Definition:
720 :     1/0 = FIGMODELmodel->rights(string::username);
721 :     Description:
722 :     Returns 1 if the input user can view the model, and zero otherwise
723 :     =cut
724 :     sub rights {
725 :     my ($self,$username) = @_;
726 :     if ($self->public()) {
727 :     return 1;
728 :     }
729 :     if (!defined($username) || $username eq "NONE") {
730 :     return 0;
731 :     }
732 :     if (!defined($self->{_userrights}->{$username})) {
733 :     $self->{_userrights}->{$self->{_data}->owner()} = 1;
734 :     my @users = split(/\|/,$self->{_data}->users());
735 :     for (my $i=0; $i < @users;$i++) {
736 :     $self->{_userrights}->{$users[$i]} = 1;
737 :     }
738 :     }
739 :     return $self->{_userrights}->{$username};
740 :     }
741 :    
742 :     =head3 public
743 :     Definition:
744 :     1/0 = FIGMODELmodel->public();
745 :     Description:
746 :     Returns 1 if the model is public, and zero otherwise
747 :     =cut
748 :     sub public {
749 :     my ($self) = @_;
750 :     if ($self->{_data}->users() eq "all") {
751 :     return 1;
752 :     }
753 :     return 0;
754 :     }
755 :    
756 :     =head3 directory
757 :     Definition:
758 :     string = FIGMODELmodel->directory();
759 :     Description:
760 :     Returns model directory
761 :     =cut
762 :     sub directory {
763 :     my ($self) = @_;
764 :    
765 :     if (!defined($self->{_directory})) {
766 :     my $userdirectory = "";
767 :     if ($self->owner() ne "master") {
768 :     $userdirectory = $self->owner()."/";
769 :     }
770 :     my $source = $self->source();
771 :     if ($source =~ /^MGRAST/) {
772 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
773 :     } elsif ($source =~ /^RAST/) {
774 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
775 :     } elsif ($source =~ /^SEED/) {
776 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
777 :     } elsif ($source =~ /^PM/) {
778 :     if (length($userdirectory) == 0) {
779 :     $self->{_directory} = $self->figmodel()->config("imported model directory")->[0].$self->id()."/";
780 :     } else {
781 :     $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->id()."/";
782 :     }
783 :     }
784 :     }
785 :    
786 :     return $self->{_directory};
787 :     }
788 :    
789 :     =head3 filename
790 :     Definition:
791 :     string = FIGMODELmodel->filename();
792 :     Description:
793 :     Returns model filename
794 :     =cut
795 :     sub filename {
796 :     my ($self) = @_;
797 :    
798 :     return $self->directory().$self->id().$self->selected_version().".txt";
799 :     }
800 :    
801 :     =head3 version
802 :     Definition:
803 :     string = FIGMODELmodel->version();
804 :     Description:
805 :     Returns the version of the model
806 :     =cut
807 :     sub version {
808 :     my ($self) = @_;
809 :    
810 :     if (!defined($self->{_version})) {
811 :     if (!defined($self->{_selectedversion})) {
812 :     $self->{_version} = "V".$self->{_data}->version().".".$self->{_data}->autocompleteVersion();
813 :     } else {
814 :     $self->{_version} = $self->{_selectedversion};
815 :     }
816 :     }
817 :     return $self->{_version};
818 :     }
819 :    
820 :     =head3 selected_version
821 :     Definition:
822 :     string = FIGMODELmodel->selected_version();
823 :     Description:
824 :     Returns the selected version of the model
825 :     =cut
826 :     sub selected_version {
827 :     my ($self) = @_;
828 :    
829 :     if (!defined($self->{_selectedversion})) {
830 :     return "";
831 :     }
832 :     return $self->{_selectedversion};
833 :     }
834 :    
835 :     =head3 modification_time
836 :     Definition:
837 :     string = FIGMODELmodel->modification_time();
838 :     Description:
839 :     Returns the selected version of the model
840 :     =cut
841 :     sub modification_time {
842 :     my ($self) = @_;
843 :     return $self->{_data}->modificationDate();
844 :     }
845 :    
846 :     =head3 gene_reactions
847 :     Definition:
848 :     string = FIGMODELmodel->gene_reactions();
849 :     Description:
850 :     Returns the number of reactions added by the gap filling
851 :     =cut
852 :     sub gene_reactions {
853 :     my ($self) = @_;
854 :     return ($self->{_data}->reactions() - $self->{_data}->autoCompleteReactions() - $self->{_data}->spontaneousReactions() - $self->{_data}->gapFillReactions());
855 :     }
856 :    
857 :     =head3 total_compounds
858 :     Definition:
859 :     string = FIGMODELmodel->total_compounds();
860 :     Description:
861 :     Returns the number of compounds in the model
862 :     =cut
863 :     sub total_compounds {
864 :     my ($self) = @_;
865 :     return $self->{_data}->compounds();
866 :     }
867 :    
868 :     =head3 gapfilling_reactions
869 :     Definition:
870 :     string = FIGMODELmodel->gapfilling_reactions();
871 :     Description:
872 :     Returns the number of reactions added by the gap filling
873 :     =cut
874 :     sub gapfilling_reactions {
875 :     my ($self) = @_;
876 :     return ($self->{_data}->autoCompleteReactions()+$self->{_data}->gapFillReactions());
877 :     }
878 :    
879 :     =head3 total_reactions
880 :     Definition:
881 :     string = FIGMODELmodel->total_reactions();
882 :     Description:
883 :     Returns the total number of reactions in the model
884 :     =cut
885 :     sub total_reactions {
886 :     my ($self) = @_;
887 :     return $self->{_data}->reactions();
888 :     }
889 :    
890 :     =head3 model_genes
891 :     Definition:
892 :     string = FIGMODELmodel->model_genes();
893 :     Description:
894 :     Returns the number of genes mapped to one or more reactions in the model
895 :     =cut
896 :     sub model_genes {
897 :     my ($self) = @_;
898 :     return $self->{_data}->associatedGenes();
899 :     }
900 :    
901 :     =head3 class
902 :     Definition:
903 :     string = FIGMODELmodel->class();
904 :     Description:
905 :     Returns the class of the model: gram positive, gram negative, other
906 :     =cut
907 :     sub class {
908 :     my ($self) = @_;
909 :     return $self->{_data}->cellwalltype();
910 :     }
911 :    
912 :     =head3 taxonomy
913 :     Definition:
914 :     string = FIGMODELmodel->taxonomy();
915 :     Description:
916 :     Returns model taxonomy or biome if this is an metagenome model
917 :     =cut
918 :     sub taxonomy {
919 :     my ($self) = @_;
920 :    
921 :     if (!defined($self->{_taxonomy})) {
922 :     my $source = $self->source();
923 :     if ($source =~ /^MGRAST/) {
924 :     $self->{_taxonomy} = "NA";
925 :     my $mgmodeldata = $self->figmodel()->database()->mg_model_data($self->genome());
926 :     if (defined($mgmodeldata)) {
927 :     $self->{_taxonomy} = $mgmodeldata->biome();
928 :     }
929 :     } else {
930 :     $self->{_taxonomy} = $self->fig()->taxonomy_of($self->genome());
931 :     }
932 :     }
933 :    
934 :     return $self->{_taxonomy};
935 :     }
936 :    
937 :     =head3 genome_size
938 :     Definition:
939 :     string = FIGMODELmodel->genome_size();
940 :     Description:
941 :     Returns size of the modeled genome in KB
942 :     =cut
943 :     sub genome_size {
944 :     my ($self) = @_;
945 :    
946 :     if (!defined($self->{_genome_size})) {
947 :     my $source = $self->source();
948 :     if ($source =~ /^MGRAST/) {
949 :     $self->{_genome_size} = "NA";
950 :     if (defined($self->mgdata())) {
951 :     $self->{_genome_size} = $self->mgdata()->size;
952 :     }
953 :     } else {
954 :     $self->{_genome_size} = $self->fig()->genome_szdna($self->genome());
955 :     }
956 :     }
957 :    
958 :     return $self->{_genome_size};
959 :     }
960 :    
961 :     =head3 genome_genes
962 :     Definition:
963 :     string = FIGMODELmodel->genome_genes();
964 :     Description:
965 :     Returns the number of genes in the modeled genome
966 :     =cut
967 :     sub genome_genes {
968 :     my ($self) = @_;
969 :    
970 :     if (!defined($self->{_genome_genes})) {
971 :     my $source = $self->source();
972 :     if ($source =~ /^MGRAST/) {
973 :     $self->{_genome_genes} = "NA";
974 :     if (defined($self->mgdata())) {
975 :     $self->{_genome_genes} = $self->mgdata()->genome_contig_count;
976 :     }
977 :     } else {
978 :     $self->{_genome_genes} = $self->figmodel()->get_genome_stats($self->genome())->{"TOTAL GENES"}->[0];
979 :     }
980 :     }
981 :    
982 :     return $self->{_genome_genes};
983 :     }
984 :    
985 :     =head3 run_default_model_predictions
986 :     Definition:
987 :     0/1::status = FIGMODELmodel->run_default_model_predictions(string::media ID);
988 :     Description:
989 :     =cut
990 :     sub run_default_model_predictions {
991 :     my ($self,$Media) = @_;
992 :    
993 :     #Assuming complete media if none is provided
994 :     if (!defined($Media)) {
995 :     $Media = "Complete";
996 :     }
997 :    
998 :     #Predicting essentiality
999 :     my $result = $self->figmodel()->RunFBASimulation($self->id(),"SINGLEKO",undef,undef,[$self->id()],[$Media]);
1000 :     #Checking that the table is defined and the output file exists
1001 :     if (defined($result) && defined($result->get_row(0)->{"ESSENTIALGENES"})) {
1002 :     $self->figmodel()->database()->print_array_to_file($self->directory()."EssentialGenes-".$self->id()."-".$Media.".tbl",[join("\n",@{$result->get_row(0)->{"ESSENTIALGENES"}})]);
1003 :     } else {
1004 :     $self->figmodel()->error_message("FIGMODELmodel:run_default_model_predictions:could not identify essential reactions for model ".$self->id().$self->selected_version().".");
1005 :     return $self->figmodel()->fail();
1006 :     }
1007 :    
1008 :     #Classifying reactions and compounds
1009 :     my $tbl = $self->classify_model_reactions($Media);
1010 :     if (!defined($tbl)) {
1011 :     $self->figmodel()->error_message("FIGMODELmodel:run_default_model_predictions:could not classify reactions for model ".$self->id().$self->selected_version().".");
1012 :     return $self->figmodel()->fail();
1013 :     }
1014 :     $tbl->save();
1015 :    
1016 :     return $self->figmodel()->success();
1017 :     }
1018 :    
1019 :     =head3 update_stats_for_gap_filling
1020 :     Definition:
1021 :     {string => [string]} = FIGMODELmodel->update_stats_for_gap_filling(int::gapfill time);
1022 :     Description:
1023 :     =cut
1024 :     sub update_stats_for_gap_filling {
1025 :     my ($self,$gapfilltime) = @_;
1026 :     $self->{_data}->autoCompleteTime($gapfilltime);
1027 :     $self->{_data}->autocompleteDate(time());
1028 :     $self->{_data}->modificationDate(time());
1029 :     my $version = $self->{_data}->autocompleteVersion();
1030 :     $self->{_data}->autocompleteVersion($version+1);
1031 :     }
1032 :    
1033 :     =head3 update_stats_for_build
1034 :     Definition:
1035 :     {string => [string]} = FIGMODELmodel->update_stats_for_build();
1036 :     Description:
1037 :     =cut
1038 :     sub update_stats_for_build {
1039 :     my ($self) = @_;
1040 :     $self->{_data}->builtDate(time());
1041 :     $self->{_data}->modificationDate(time());
1042 :     my $version = $self->{_data}->version();
1043 :     $self->{_data}->version($version+1);
1044 :     }
1045 :    
1046 :     =head3 update_model_stats
1047 :     Definition:
1048 :     FIGMODELmodel->update_model_stats();
1049 :     Description:
1050 :     =cut
1051 :     sub update_model_stats {
1052 :     my ($self) = @_;
1053 :    
1054 :     #Getting reaction table
1055 :     my $rxntbl = $self->reaction_table();
1056 :     if (!defined($rxntbl)) {
1057 :     $self->figmodel()->error_message("FIGMODELmodel:update_model_stats:Could not load reaction list for ".$self->id());
1058 :     return undef;
1059 :     }
1060 :     my $cpdtbl = $self->compound_table();
1061 :    
1062 :     #Calculating all necessary stats
1063 :     my %GeneHash;
1064 :     my %NonpegHash;
1065 :     my %CompoundHash;
1066 :     my $spontaneousReactions = 0;
1067 :     my $gapFillReactions = 0;
1068 :     my $biologReactions = 0;
1069 :     my $transporters = 0;
1070 :     my $autoCompleteReactions = 0;
1071 :     my $associatedSubsystemGenes = 0;
1072 :     for (my $i=0; $i < $rxntbl->size(); $i++) {
1073 :     my $Row = $rxntbl->get_row($i);
1074 :     if (defined($Row) && defined($Row->{"ASSOCIATED PEG"})) {
1075 :     my $ReactionRow = $self->figmodel()->get_reaction($Row->{"LOAD"}->[0]);
1076 :     if (defined($ReactionRow->{"EQUATION"}->[0])) {
1077 :     #Checking for extracellular metabolites which indicate that this is a transporter
1078 :     if ($ReactionRow->{"EQUATION"}->[0] =~ m/\[e\]/) {
1079 :     $transporters++;
1080 :     }
1081 :     }
1082 :     #Identifying spontaneous/biolog/gapfilling/gene associated reactions
1083 :     if ($Row->{"ASSOCIATED PEG"}->[0] =~ m/BIOLOG/i) {
1084 :     $biologReactions++;
1085 :     } elsif ($Row->{"ASSOCIATED PEG"}->[0] =~ m/GROW/i) {
1086 :     $gapFillReactions++;
1087 :     } elsif ($Row->{"ASSOCIATED PEG"}->[0] =~ m/SPONTANEOUS/i) {
1088 :     $spontaneousReactions++;
1089 :     } elsif ($Row->{"ASSOCIATED PEG"}->[0] =~ m/GAP/ || $Row->{"ASSOCIATED PEG"}->[0] =~ m/UNIVERSAL/i || $Row->{"ASSOCIATED PEG"}->[0] =~ m/UNKNOWN/i) {
1090 :     $autoCompleteReactions++;
1091 :     } else {
1092 :     foreach my $GeneSet (@{$Row->{"ASSOCIATED PEG"}}) {
1093 :     $_ = $GeneSet;
1094 :     my @GeneList = /(peg\.\d+)/g;
1095 :     foreach my $Gene (@GeneList) {
1096 :     if ($Gene =~ m/(peg\.\d+)/) {
1097 :     $GeneHash{$1} = 1;
1098 :     } else {
1099 :     $NonpegHash{$Gene} = 1;
1100 :     }
1101 :     }
1102 :     }
1103 :     }
1104 :     }
1105 :     }
1106 :     my @genes = keys(%GeneHash);
1107 :     my @othergenes = keys(%NonpegHash);
1108 :    
1109 :     #Setting the reaction count
1110 :     $self->{_data}->reactions($rxntbl->size());
1111 :     #Setting the metabolite count
1112 :     $self->{_data}->compounds($rxntbl->size());
1113 :     #Setting the gene count
1114 :     my $geneCount = @genes + @othergenes;
1115 :     $self->{_data}->associatedGenes($geneCount);
1116 :     #Setting remaining stats
1117 :     $self->{_data}->spontaneousReactions($spontaneousReactions);
1118 :     $self->{_data}->gapFillReactions($gapFillReactions);
1119 :     $self->{_data}->biologReactions($biologReactions);
1120 :     $self->{_data}->transporters($transporters);
1121 :     $self->{_data}->autoCompleteReactions($autoCompleteReactions);
1122 :     $self->{_data}->associatedSubsystemGenes($associatedSubsystemGenes);
1123 :     #Setting the model class
1124 :     my $class = "";
1125 :     for (my $i=0; $i < @{$self->figmodel()->config("class list")}; $i++) {
1126 :     if (defined($self->figmodel()->config($self->figmodel()->config("class list")->[$i]))) {
1127 :     if (defined($self->figmodel()->config($self->figmodel()->config("class list")->[$i])->{$self->id()})) {
1128 :     $class = $self->figmodel()->config("class list")->[$i];
1129 :     last;
1130 :     }
1131 :     if ($class eq "" && defined($self->figmodel()->config($self->figmodel()->config("class list")->[$i])->{$self->genome()})) {
1132 :     $class = $self->figmodel()->config("class list")->[$i];
1133 :     }
1134 :     }
1135 :     }
1136 :     if ($class eq "") {
1137 :     $class = $self->figmodel()->get_genome_stats($self->genome())->{CLASS}->[0];
1138 :     }
1139 :     if ($class eq "") {
1140 :     $class = "unknown";
1141 :     }
1142 :     $self->{_data}->cellwalltype($class);
1143 :     }
1144 :    
1145 :     =head3 GapFillModel
1146 :     Definition:
1147 :     (success/fail) = FIGMODELmodel->GapFillModel();
1148 :     Description:
1149 :     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.
1150 :     Before running the gap filling, the existing model is backup in the same directory with the current version numbers appended.
1151 :     If the model has been gap filled previously, the previous gap filling reactions are removed prior to running the gap filling again.
1152 :     =cut
1153 :     sub GapFillModel {
1154 :     my ($self,$donotclear) = @_;
1155 :    
1156 :     #Setting status of model to gap filling
1157 :     my $OrganismID = $self->genome();
1158 :     $self->set_status(1,"Auto completion running");
1159 :     my $UniqueFilename = $self->figmodel()->filename();
1160 :     my $StartTime = time();
1161 :    
1162 :     #Archiving the existing model
1163 :     $self->ArchiveModel();
1164 :    
1165 :     #Removing any gapfilling reactions that may be currently present in the model
1166 :     if (!defined($donotclear) || $donotclear != 1) {
1167 :     my $ModelTable = $self->reaction_table();
1168 :     for (my $i=0; $i < $ModelTable->size(); $i++) {
1169 :     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/) {
1170 :     $ModelTable->delete_row($i);
1171 :     $i--;
1172 :     }
1173 :     }
1174 :     $ModelTable->save();
1175 :     }
1176 :    
1177 :     #Calling the MFAToolkit to run the gap filling optimization
1178 :     my $MinimalMediaTable = $self->figmodel()->database()->GetDBTable("MINIMAL MEDIA TABLE");
1179 :     my $Media = "Complete";
1180 :     if (defined($MinimalMediaTable->get_row_by_key($self->genome(),"Organism"))) {
1181 :     $Media = $MinimalMediaTable->get_row_by_key($self->genome(),"Organism")->{"Minimal media"}->[0];
1182 :     #Loading media, changing bounds, saving media as a test media
1183 :     my $MediaTable = FIGMODELTable::load_table($self->config("Media directory")->[0].$Media.".txt",";","",0,["VarName"]);
1184 :     for (my $i=0; $i < $MediaTable->size(); $i++) {
1185 :     if ($MediaTable->get_row($i)->{"Min"}->[0] < 0) {
1186 :     $MediaTable->get_row($i)->{"Min"}->[0] = -10000;
1187 :     }
1188 :     if ($MediaTable->get_row($i)->{"Max"}->[0] > 0) {
1189 :     $MediaTable->get_row($i)->{"Max"}->[0] = 10000;
1190 :     }
1191 :     }
1192 :     $MediaTable->save($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");
1193 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$UniqueFilename."TestMedia",["GapFilling"],{"Default max drain flux" => 0,"Reactions to knockout" => $self->config("permanently knocked out reactions")->[0]},"GapFill".$self->id().".log",undef));
1194 :     unlink($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");
1195 :     } else {
1196 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),undef,["GapFilling"],{"Reactions to knockout" => $self->config("permanently knocked out reactions")->[0]},"GapFill".$self->id().".log",undef));
1197 :     }
1198 :    
1199 :     #Parse the solutions file for the model and get the reaction list from it
1200 :     my $SolutionData = $self->figmodel()->LoadProblemReport($UniqueFilename);
1201 :     #Clearing the mfatoolkit output and log file
1202 :     $self->figmodel()->clearing_output($UniqueFilename,"GapFill".$self->id().".log");
1203 :     if (!defined($SolutionData) || $SolutionData->size() == 0) {
1204 :     $self->set_status(1,"Auto completion failed: auto completion solution not found");
1205 :     $self->figmodel()->error_message("FIGMODEL:GapFillModel: Gap filling report not found!");
1206 :     return $self->fail();
1207 :     }
1208 :     $SolutionData->save("/home/chenry/Solution.txt");
1209 :    
1210 :     #Looking for the last printed recursive MILP solution
1211 :     for (my $i=($SolutionData->size()-1); $i >=0; $i--) {
1212 :     if (defined($SolutionData->get_row($i)->{"Notes"}) && $SolutionData->get_row($i)->{"Notes"}->[0] =~ m/^Recursive/) {
1213 :     my $AllSolutions = substr($SolutionData->get_row($i)->{"Notes"}->[0],15);
1214 :     my @TempThree = split(/\|/,$AllSolutions);
1215 :     if (@TempThree > 0 && $TempThree[0] =~ m/.+:(.+)/) {
1216 :     my @TempFour = split(/,/,$1);
1217 :     my $DirectionList;
1218 :     my $ReactionList;
1219 :     for (my $j=0; $j < @TempFour; $j++) {
1220 :     my $ID = "";
1221 :     my $Sign = "=>";
1222 :     if ($TempFour[$j] =~ m/\-(rxn\d\d\d\d\d)/) {
1223 :     $ID = $1;
1224 :     $Sign = "<=";
1225 :     } elsif ($TempFour[$j] =~ m/\+(rxn\d\d\d\d\d)/) {
1226 :     $ID = $1;
1227 :     }
1228 :     if ($ID ne "") {
1229 :     if ($self->figmodel()->reversibility_of_reaction($ID) ne $Sign) {
1230 :     $Sign = "<=>";
1231 :     }
1232 :     push(@{$DirectionList},$Sign);
1233 :     push(@{$ReactionList},$ID);
1234 :     }
1235 :     }
1236 :     $self->figmodel()->IntegrateGrowMatchSolution($self->id(),undef,$ReactionList,$DirectionList,"GAP FILLING",0,1);
1237 :     }
1238 :     last;
1239 :     }
1240 :     }
1241 :     #Updating model stats with gap filling results
1242 :     my $ElapsedTime = time() - $StartTime;
1243 :     $self->figmodel()->database()->ClearDBModel($self->id(),1);
1244 :     #Determining why each gap filling reaction was added
1245 :     $self->figmodel()->IdentifyDependancyOfGapFillingReactions($self->id(),$Media);
1246 :     if (!defined($donotclear) || $donotclear != 1) {
1247 :     $self->figmodel()->AddBiologTransporters($self->id());
1248 :     if ($self->id() !~ m/MGRast/) {
1249 :     $self->update_stats_for_gap_filling($ElapsedTime);
1250 :     }
1251 :     } else {
1252 :     $self->update_model_stats();
1253 :     }
1254 :     #Printing the updated SBML file
1255 :     $self->PrintSBMLFile();
1256 :     $self->PrintModelLPFile($self->id());
1257 :     $self->set_status(1,"Auto completion successfully finished");
1258 :     $self->run_default_model_predictions();
1259 :    
1260 :     return $self->success();
1261 :     }
1262 :    
1263 :     =head3 GapGenModel
1264 :     Definition:
1265 :     FIGMODELmodel->GapGenModel();
1266 :     Description:
1267 :     Runs the gap generation algorithm to correct a single false positive prediction. Results are loaded into a table.
1268 :     =cut
1269 :    
1270 :     sub GapGenModel {
1271 :     my ($self,$Media,$KOList,$NoKOList,$Experiment,$SolutionLimit) = @_;
1272 :    
1273 :     #Enforcing nonoptional arguments
1274 :     if (!defined($Media)) {
1275 :     return undef;
1276 :     }
1277 :     if (!defined($KOList)) {
1278 :     $KOList->[0] = "none";
1279 :     }
1280 :     if (!defined($NoKOList)) {
1281 :     $NoKOList->[0] = "none";
1282 :     }
1283 :     if (!defined($Experiment)) {
1284 :     $Experiment= "ReactionKO";
1285 :     }
1286 :     if (!defined($SolutionLimit)) {
1287 :     $SolutionLimit = "10";
1288 :     }
1289 :    
1290 :     #Translating the KO lists into arrays
1291 :     if (ref($KOList) ne "ARRAY") {
1292 :     my $temp = $KOList;
1293 :     $KOList = ();
1294 :     push(@{$KOList},split(/[,;]/,$temp));
1295 :     }
1296 :     my $noKOHash;
1297 :     if (defined($NoKOList) && ref($NoKOList) ne "ARRAY") {
1298 :     my $temp = $NoKOList;
1299 :     $NoKOList = ();
1300 :     push(@{$NoKOList},split(/[,;]/,$temp));
1301 :     foreach my $rxn (@{$NoKOList}) {
1302 :     $noKOHash->{$rxn} = 1;
1303 :     }
1304 :     }
1305 :    
1306 :     #Checking if solutions exist for the input parameters
1307 :     $self->aquireModelLock();
1308 :     my $tbl = $self->load_model_table("GapGenSolutions");
1309 :     my $solutionRow = $tbl->get_table_by_key($Experiment,"Experiment")->get_table_by_key($Media,"Media")->get_row_by_key(join(",",@{$KOList}),"KOlist");
1310 :     my $solutions;
1311 :     if (defined($solutionRow)) {
1312 :     #Checking if any solutions conform to the no KO list
1313 :     foreach my $solution (@{$solutionRow->{Solutions}}) {
1314 :     my @reactions = split(/,/,$solution);
1315 :     my $include = 1;
1316 :     foreach my $rxn (@reactions) {
1317 :     if ($rxn =~ m/(rxn\d\d\d\d\d)/) {
1318 :     if (defined($noKOHash->{$1})) {
1319 :     $include = 0;
1320 :     }
1321 :     }
1322 :     }
1323 :     if ($include == 1) {
1324 :     push(@{$solutions},$solution);
1325 :     }
1326 :     }
1327 :     } else {
1328 :     $solutionRow = {Media => [$Media],Experiment => [$Experiment],KOlist => [join(",",@{$KOList})]};
1329 :     $tbl->add_row($solutionRow);
1330 :     $self->figmodel()->database()->save_table($tbl);
1331 :     }
1332 :     $self->releaseModelLock();
1333 :    
1334 :     #Returning solution list of solutions were found
1335 :     if (defined($solutions) && @{$solutions} > 0) {
1336 :     return $solutions;
1337 :     }
1338 :    
1339 :     #Getting unique filename
1340 :     my $Filename = $self->figmodel()->filename();
1341 :    
1342 :     #Running the gap generation
1343 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($Filename,$self->id().$self->selected_version(),$Media,["GapGeneration"],{"Recursive MILP solution limit" => $SolutionLimit ,"Reactions that should always be active" => join(";",@{$NoKOList}),"Reactions to knockout" => join(";",@{$KOList}),"Reactions that are always blocked" => "none"},"Gapgeneration-".$self->id().$self->selected_version()."-".$Filename.".log",undef,undef));
1344 :     my $ProblemReport = $self->figmodel()->LoadProblemReport($Filename);
1345 :     if (!defined($ProblemReport)) {
1346 :     $self->figmodel()->error_message("FIGMODEL:GapGenerationAlgorithm;No problem report;".$Filename.";".$self->id().$self->selected_version().";".$Media.";".$KOList.";".$NoKOList);
1347 :     return undef;
1348 :     }
1349 :    
1350 :     #Clearing the output folder and log file
1351 :     $self->figmodel()->clearing_output($Filename,"Gapgeneration-".$self->id().$self->selected_version()."-".$Filename.".log");
1352 :    
1353 :     #Saving the solution
1354 :     $self->aquireModelLock();
1355 :     $tbl = $self->load_model_table("GapGenSolutions");
1356 :     $solutionRow = $tbl->get_table_by_key($Experiment,"Experiment")->get_table_by_key($Media,"Media")->get_row_by_key(join(",",@{$KOList}),"KOlist");
1357 :     for (my $j=0; $j < $ProblemReport->size(); $j++) {
1358 :     if ($ProblemReport->get_row($j)->{"Notes"}->[0] =~ m/^Recursive\sMILP\s([^)]+)/) {
1359 :     my @SolutionList = split(/\|/,$1);
1360 :     for (my $k=0; $k < @SolutionList; $k++) {
1361 :     if ($SolutionList[$k] =~ m/(\d+):(.+)/) {
1362 :     push(@{$solutionRow->{Solutions}},$2);
1363 :     push(@{$solutions},$2);
1364 :     }
1365 :     }
1366 :     }
1367 :     }
1368 :     $self->figmodel()->database()->save_table($tbl);
1369 :     $self->releaseModelLock();
1370 :    
1371 :     return $solutions;
1372 :     }
1373 :    
1374 :     =head3 datagapfill
1375 :     Definition:
1376 :     success()/fail() = FIGMODELmodel->datagapfill();
1377 :     Description:
1378 :     Run gapfilling on the input run specifications
1379 :     =cut
1380 :     sub datagapfill {
1381 :     my ($self,$GapFillingRunSpecs,$TansferFileSuffix) = @_;
1382 :     my $UniqueFilename = $self->figmodel()->filename();
1383 :     if (defined($GapFillingRunSpecs) && @{$GapFillingRunSpecs} > 0) {
1384 :     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));
1385 :     #Checking that the solution exists
1386 :     if (!-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingSolutionTable.txt") {
1387 :     $self->figmodel()->error_message("FIGMODEL:GapFillingAlgorithm: Could not find MFA output file!");
1388 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-GFS.txt",["Experiment;Solution index;Solution cost;Solution reactions"]);
1389 :     return undef;
1390 :     }
1391 :     my $GapFillResultTable = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingSolutionTable.txt",";","",0,undef);
1392 :     if (defined($TansferFileSuffix)) {
1393 :     system("cp ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/GapFillingSolutionTable.txt ".$self->directory().$self->id().$self->selected_version()."-".$TansferFileSuffix.".txt");
1394 :     }
1395 :     #If the system is not configured to preserve all logfiles, then the mfatoolkit output folder is deleted
1396 :     $self->figmodel()->clearing_output($UniqueFilename,"GapFilling-".$self->id().$self->selected_version()."-".$UniqueFilename.".log");
1397 :     return $GapFillResultTable;
1398 :     }
1399 :     if (defined($TansferFileSuffix)) {
1400 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-".$TansferFileSuffix.".txt",["Experiment;Solution index;Solution cost;Solution reactions"]);
1401 :     }
1402 :     return undef;
1403 :     }
1404 :    
1405 :     =head3 TestSolutions
1406 :     Definition:
1407 :     $model->TestSolutions($ModelID,$NumProcessors,$ProcessorIndex,$GapFill);
1408 :     Description:
1409 :     Example:
1410 :     =cut
1411 :    
1412 :     sub TestSolutions {
1413 :     my ($self,$OriginalErrorFilename,$GapFillResultTable) = @_;
1414 :     #Getting the filename
1415 :     my $UniqueFilename = $self->figmodel()->filename();
1416 :     #Reading in the original error matrix which has the headings for the original model simulation
1417 :     my $OriginalErrorData;
1418 :     if (!defined($OriginalErrorFilename) || !-e $self->directory().$OriginalErrorFilename) {
1419 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector) = $self->RunAllStudiesWithDataFast("All");
1420 :     $OriginalErrorData = [$HeadingVector,$Errorvector];
1421 :     } else {
1422 :     $OriginalErrorData = $self->figmodel()->database()->load_single_column_file($self->directory().$OriginalErrorFilename,"");
1423 :     }
1424 :     my $HeadingHash;
1425 :     my @HeadingArray = split(/;/,$OriginalErrorData->[0]);
1426 :     my @OrigErrorArray = split(/;/,$OriginalErrorData->[1]);
1427 :     for (my $i=0; $i < @HeadingArray; $i++) {
1428 :     my @SubArray = split(/:/,$HeadingArray[$i]);
1429 :     $HeadingHash->{$SubArray[0].":".$SubArray[1].":".$SubArray[2]} = $i;
1430 :     }
1431 :     #Scanning through the gap filling solutions
1432 :     my $TempVersion = "V".$UniqueFilename;
1433 :     my $ErrorMatrixLines;
1434 :     for (my $i=0; $i < $GapFillResultTable->size(); $i++) {
1435 :     print "Starting problem solving ".$i."\n";
1436 :     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];
1437 :     #Integrating solution into test model
1438 :     my $ReactionArray;
1439 :     my $DirectionArray;
1440 :     my @ReactionList = split(/,/,$GapFillResultTable->get_row($i)->{"Solution reactions"}->[0]);
1441 :     my %SolutionHash;
1442 :     for (my $k=0; $k < @ReactionList; $k++) {
1443 :     if ($ReactionList[$k] =~ m/(.+)(rxn\d\d\d\d\d)/) {
1444 :     my $Reaction = $2;
1445 :     my $Sign = $1;
1446 :     if (defined($SolutionHash{$Reaction})) {
1447 :     $SolutionHash{$Reaction} = "<=>";
1448 :     } elsif ($Sign eq "-") {
1449 :     $SolutionHash{$Reaction} = "<=";
1450 :     } elsif ($Sign eq "+") {
1451 :     $SolutionHash{$Reaction} = "=>";
1452 :     } else {
1453 :     $SolutionHash{$Reaction} = $Sign;
1454 :     }
1455 :     }
1456 :     }
1457 :     @ReactionList = keys(%SolutionHash);
1458 :     for (my $k=0; $k < @ReactionList; $k++) {
1459 :     push(@{$ReactionArray},$ReactionList[$k]);
1460 :     push(@{$DirectionArray},$SolutionHash{$ReactionList[$k]});
1461 :     }
1462 :     print "Integrating solution!\n";
1463 :     $self->figmodel()->IntegrateGrowMatchSolution($self->id().$self->selected_version(),$self->directory().$self->id().$TempVersion.".txt",$ReactionArray,$DirectionArray,"Gapfilling ".$GapFillResultTable->get_row($i)->{"Experiment"}->[0],1,1);
1464 :     $self->PrintModelLPFile();
1465 :     #Running the model against all available experimental data
1466 :     print "Running test model!\n";
1467 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector) = $self->RunAllStudiesWithDataFast("All");
1468 :    
1469 :     @HeadingArray = split(/;/,$HeadingVector);
1470 :     my @ErrorArray = @OrigErrorArray;
1471 :     my @TempArray = split(/;/,$Errorvector);
1472 :     for (my $j=0; $j < @HeadingArray; $j++) {
1473 :     my @SubArray = split(/:/,$HeadingArray[$j]);
1474 :     $ErrorArray[$HeadingHash->{$SubArray[0].":".$SubArray[1].":".$SubArray[2]}] = $TempArray[$j];
1475 :     }
1476 :     $ErrorLine .= ";".$FalsePostives."/".$FalseNegatives.";".join(";",@ErrorArray);
1477 :     push(@{$ErrorMatrixLines},$ErrorLine);
1478 :     print "Finishing problem solving ".$i."\n";
1479 :     }
1480 :     #Clearing out the test model
1481 :     if (-e $self->directory().$self->id().$TempVersion.".txt") {
1482 :     unlink($self->directory().$self->id().$TempVersion.".txt");
1483 :     unlink($self->directory()."SimulationOutput".$self->id().$TempVersion.".txt");
1484 :     }
1485 :     return $ErrorMatrixLines;
1486 :     }
1487 :    
1488 :     =head3 status
1489 :     Definition:
1490 :     int::model status = FIGMODELmodel->status();
1491 :     Description:
1492 :     Returns the current status of the SEED model associated with the input genome ID.
1493 :     model status = 1: model exists
1494 :     model status = 0: model is being built
1495 :     model status = -1: model does not exist
1496 :     model status = -2: model build failed
1497 :     =cut
1498 :     sub status {
1499 :     my ($self) = @_;
1500 :     return $self->{_data}->{status}->[0];
1501 :     }
1502 :    
1503 :     =head3 message
1504 :     Definition:
1505 :     string::model message = FIGMODELmodel->message();
1506 :     Description:
1507 :     Returns a message associated with the models current status
1508 :     =cut
1509 :     sub message {
1510 :     my ($self) = @_;
1511 :     return $self->{_data}->{message}->[0];
1512 :     }
1513 :    
1514 :     =head3 set_status
1515 :     Definition:
1516 :     (success/fail) = FIGMODELmodel->set_status(int::new status,string::status message);
1517 :     Description:
1518 :     Changes the current status of the SEED model
1519 :     new status = 1: model exists
1520 :     new status = 0: model is being built
1521 :     new status = -1: model does not exist
1522 :     new status = -2: model build failed
1523 :     =cut
1524 :     sub set_status {
1525 :     my ($self,$NewStatus,$Message) = @_;
1526 :    
1527 :     #Getting the model row from the MODELS table
1528 :     $self->{_data}->{status}->[0] = $NewStatus;
1529 :     $self->{_data}->{message}->[0] = $Message;
1530 :     $self->figmodel()->database()->update_row("MODELS",$self->{_data},"id");
1531 :     return $self->config("SUCCESS")->[0];
1532 :     }
1533 :    
1534 :     =head3 CreateSingleGenomeReactionList
1535 :     Definition:
1536 :     FIGMODEL->CreateSingleGenomeReactionList();
1537 :     Description:
1538 :     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.
1539 :     Example:
1540 :     =cut
1541 :    
1542 :     sub CreateSingleGenomeReactionList {
1543 :     my ($self,$RunGapFilling) = @_;
1544 :    
1545 :     #Creating directory
1546 :     if ($self->owner() ne "master" && !-d $self->figmodel()->config("organism directory")->[0].$self->owner()."/") {
1547 :     system("mkdir ".$self->figmodel()->config("organism directory")->[0].$self->owner()."/");
1548 :     } elsif ($self->owner() eq "master" && !-d $self->figmodel()->config("organism directory")->[0].$self->genome()."/") {
1549 :     system("mkdir ".$self->figmodel()->config("organism directory")->[0].$self->genome()."/");
1550 :     }
1551 :     if ($self->owner() ne "master" && !-d $self->figmodel()->config("organism directory")->[0].$self->owner()."/".$self->genome()."/") {
1552 :     system("mkdir ".$self->figmodel()->config("organism directory")->[0].$self->owner()."/".$self->genome()."/");
1553 :     }
1554 :    
1555 :     #Getting genome stats
1556 :     my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
1557 :     my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
1558 :     if (!defined($FeatureTable)) {
1559 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." genome features could not be accessed!");
1560 :     return $self->fail();
1561 :     }
1562 :     #Checking that the number of genes exceeds the minimum size
1563 :     if ($FeatureTable->size() < $self->config("minimum genome size for modeling")->[0]) {
1564 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." genome rejected as too small for modeling!");
1565 :     return $self->fail();
1566 :     }
1567 :     #Setting up needed variables
1568 :     my $OriginalModelTable = undef;
1569 :     #Locking model status table
1570 :     my $ModelTable = $self->figmodel()->database()->LockDBTable("MODELS");
1571 :     my $Row = $ModelTable->get_row_by_key($self->id(),"id");
1572 :     if (!defined($Row) || !defined($Row->{status}->[0])) {
1573 :     $self->figmodel()->database()->UnlockDBTable("MODELS");
1574 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." has no row in models table!");
1575 :     return $self->fail();
1576 :     } elsif ($Row->{status}->[0] == 0) {
1577 :     $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList:Model is already being built. Canceling current build.");
1578 :     $self->figmodel()->database()->UnlockDBTable("MODELS");
1579 :     return $self->fail();
1580 :     }elsif ($Row->{status}->[0] == 1) {
1581 :     $OriginalModelTable = $self->reaction_table();
1582 :     $self->ArchiveModel();
1583 :     $Row->{status}->[0] = 0;
1584 :     $Row->{message}->[0] = "Rebuilding preliminary reconstruction";
1585 :     } else {
1586 :     $Row->{status}->[0] = 0;
1587 :     $Row->{message}->[0] = "Preliminary reconstruction";
1588 :     }
1589 :     #Updating the status table
1590 :     $self->figmodel()->database()->save_table($ModelTable);
1591 :     $self->figmodel()->database()->UnlockDBTable("MODELS");
1592 :     #Sorting GenomeData by gene location on the chromosome
1593 :     $FeatureTable->sort_rows("MIN LOCATION");
1594 :     my ($ComplexHash,$SuggestedMappings,$UnrecognizedReactions,$ReactionHash);
1595 :     my %LastGenePosition;
1596 :     my $GeneRoles;
1597 :     for (my $j=0; $j < $FeatureTable->size(); $j++) {
1598 :     my $CurrentRow = $FeatureTable->get_row($j);
1599 :     #"ID","ALIASES","MIN LOCATION","MAX LOCATION","ROLES","SUBSYSTEMS","SUBSYSTEM CLASS"
1600 :     if (defined($CurrentRow)) {
1601 :     my $GeneID = $CurrentRow->{"ID"}->[0];
1602 :     if ($GeneID =~ m/(peg\.\d+)/) {
1603 :     $GeneID = $1;
1604 :     }
1605 :     foreach my $Role (@{$CurrentRow->{"ROLES"}}) {
1606 :     if ($self->figmodel()->role_is_valid($Role) != 0) {
1607 :     push(@{$GeneRoles->{$GeneID}},$Role);
1608 :     my $ReactionHashArrayRef = $self->figmodel()->reactions_of_role($Role);
1609 :     if ($ReactionHashArrayRef != 0) {
1610 :     foreach my $Mapping (@{$ReactionHashArrayRef}) {
1611 :     if (defined($Mapping->{"REACTION"}) && defined($Mapping->{"MASTER"}) && defined($Mapping->{"SUBSYSTEM"}) && defined($Mapping->{"SOURCE"})) {
1612 :     if ($Mapping->{"REACTION"}->[0] =~ m/rxn\d\d\d\d\d/) {
1613 :     if ($Mapping->{"MASTER"}->[0] eq 1) {
1614 :     #Creating a complex if consecutive genes have been assigned to the same reaction
1615 :     $ComplexHash->{$Mapping->{"REACTION"}->[0]}->{$Mapping->{"COMPLEX"}->[0]}->{$Role}->{$GeneID} = 1;
1616 :     if (!defined($LastGenePosition{$Mapping->{"REACTION"}->[0]})) {
1617 :     $LastGenePosition{$Mapping->{"REACTION"}->[0]} = $j;
1618 :     push(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}},$GeneID);
1619 :     } elsif (($j-$LastGenePosition{$Mapping->{"REACTION"}->[0]}) < 3 && $LastGenePosition{$Mapping->{"REACTION"}->[0]} != $j) {
1620 :     my $CurrentComplex = pop(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}});
1621 :     push(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}},$CurrentComplex."+".$GeneID);
1622 :     } elsif ($LastGenePosition{$Mapping->{"REACTION"}->[0]} != $j) {
1623 :     push(@{$ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"GENES"}},$GeneID);
1624 :     }
1625 :     $LastGenePosition{$Mapping->{"REACTION"}->[0]} = $j;
1626 :     #Adding a subsystem for the reaction
1627 :     if ($self->figmodel()->subsystem_is_valid($Mapping->{"SUBSYSTEM"}->[0]) == 1) {
1628 :     ($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SUBSYSTEMS"},my $NumMatches) = $self->figmodel()->add_elements_unique($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SUBSYSTEMS"},$Mapping->{"SUBSYSTEM"}->[0]);
1629 :     if (!defined($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}) || $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] > 1) {
1630 :     if ($Mapping->{"SOURCE"}->[0] =~ m/Hope\sFiles/) {
1631 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 1;
1632 :     } elsif ($Mapping->{"SOURCE"}->[0] =~ m/SEED/) {
1633 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 2;
1634 :     } elsif (!defined($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}) || $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] > 2) {
1635 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 3;
1636 :     }
1637 :     }
1638 :     }
1639 :     #Handling confidence
1640 :     if (!defined($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}) || $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] > 2) {
1641 :     if ($Mapping->{"SOURCE"}->[0] =~ m/MATT/) {
1642 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 3;
1643 :     } elsif ($Mapping->{"SOURCE"}->[0] =~ m/CHRIS/) {
1644 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 4;
1645 :     } else {
1646 :     $ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"CONFIDENCE"}->[0] = 5;
1647 :     }
1648 :     }
1649 :     #Parsing sources
1650 :     ($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SOURCE"},my $NumMatches) = $self->figmodel()->add_elements_unique($ReactionHash->{$Mapping->{"REACTION"}->[0]}->{"SOURCE"},split(/\|/,$Mapping->{"SOURCE"}->[0]));
1651 :     } else {
1652 :     push(@{$SuggestedMappings},$GeneID."\t".$Mapping->{"REACTION"}->[0]."\t".$Role);
1653 :     }
1654 :     } else {
1655 :     push(@{$UnrecognizedReactions},$GeneID."\t".$Mapping->{"REACTION"}->[0]."\t".$Role);
1656 :     }
1657 :     }
1658 :     }
1659 :     }
1660 :     }
1661 :     }
1662 :     }
1663 :     }
1664 :    
1665 :     #Creating nonadjacent complexes
1666 :     my @ReactionList = keys(%{$ReactionHash});
1667 :     foreach my $Reaction (@ReactionList) {
1668 :     #If multiple genes are assigned to the reaction, we check if they should should be in a complex
1669 :     if (@{$ReactionHash->{$Reaction}->{"GENES"}} > 0 && defined($ComplexHash->{$Reaction})) {
1670 :     my $GeneArray;
1671 :     foreach my $Complex (keys(%{$ComplexHash->{$Reaction}})) {
1672 :     my %ComplexComponents;
1673 :     foreach my $CurrentGeneSet (@{$ReactionHash->{$Reaction}->{"GENES"}}) {
1674 :     my @GeneList = split(/\+/,$CurrentGeneSet);
1675 :     my %RoleHash;
1676 :     foreach my $Gene (@GeneList) {
1677 :     foreach my $Role (@{$GeneRoles->{$Gene}}) {
1678 :     if (defined($ComplexHash->{$Reaction}->{$Complex}->{$Role})) {
1679 :     $RoleHash{$Role} = 1;
1680 :     }
1681 :     }
1682 :     }
1683 :     if (keys(%RoleHash) > 0) {
1684 :     if (!defined($ComplexComponents{join("|",sort(keys(%RoleHash)))})) {
1685 :     my @RoleList = keys(%RoleHash);
1686 :     my @ComplexList = keys(%ComplexComponents);
1687 :     foreach my $ComplexSet (@ComplexList) {
1688 :     my @RoleList = split(/\|/,$ComplexSet);
1689 :     my $Match = 0;
1690 :     foreach my $SingleRole (@RoleList) {
1691 :     if (defined($RoleHash{$SingleRole})) {
1692 :     $Match = 1;
1693 :     last;
1694 :     }
1695 :     }
1696 :     if ($Match == 1) {
1697 :     foreach my $SingleRole (@RoleList) {
1698 :     $RoleHash{$SingleRole} = 1
1699 :     }
1700 :     push(@{$ComplexComponents{join("|",sort(keys(%RoleHash)))}},@{$ComplexComponents{$ComplexSet}});
1701 :     delete $ComplexComponents{$ComplexSet};
1702 :     }
1703 :     }
1704 :     }
1705 :     push(@{$ComplexComponents{join("|",sort(keys(%RoleHash)))}},$CurrentGeneSet);
1706 :     }
1707 :     }
1708 :     my @Position;
1709 :     my @Options;
1710 :     my $Count = 0;
1711 :     foreach my $RoleSet (keys(%ComplexComponents)) {
1712 :     push(@Position,0);
1713 :     push(@{$Options[$Count]},@{$ComplexComponents{$RoleSet}});
1714 :     $Count++;
1715 :     }
1716 :     my $Done = 0;
1717 :     $Count = 0;
1718 :     my $NewRelationship;
1719 :     while($Done == 0) {
1720 :     #Creating complex with current indecies
1721 :     $NewRelationship->[$Count] = $Options[0]->[$Position[0]];
1722 :     for (my $i=1; $i < @Position; $i++) {
1723 :     $NewRelationship->[$Count] .= "+".$Options[$i]->[$Position[$i]];
1724 :     }
1725 :     $NewRelationship->[$Count] = join("+",$self->figmodel()->remove_duplicates(split(/\+/,$NewRelationship->[$Count])));
1726 :     $Count++;
1727 :     #Iterating indecies
1728 :     my $CurrentIndex = 0;
1729 :     while($CurrentIndex >= 0) {
1730 :     if ($CurrentIndex >= @Position) {
1731 :     $CurrentIndex = -1000;
1732 :     } elsif ($Position[$CurrentIndex]+1 == @{$Options[$CurrentIndex]}) {
1733 :     $Position[$CurrentIndex] = -1;
1734 :     $CurrentIndex++;
1735 :     } else {
1736 :     $Position[$CurrentIndex]++;
1737 :     $CurrentIndex--;
1738 :     }
1739 :     }
1740 :     if ($CurrentIndex == -1000) {
1741 :     $Done = 1;
1742 :     }
1743 :     }
1744 :     push(@{$GeneArray},@{$NewRelationship});
1745 :     }
1746 :     @{$ReactionHash->{$Reaction}->{"GENES"}} = $self->figmodel()->remove_duplicates(@{$GeneArray});
1747 :     }
1748 :     }
1749 :    
1750 :     #Getting the reaction table
1751 :     my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");
1752 :    
1753 :     #Creating the model reaction table
1754 :     my $NewModelTable = FIGMODELTable->new(["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG","SUBSYSTEM","CONFIDENCE","REFERENCE","NOTES"],$self->directory().$self->id().".txt",["LOAD"],";","|","REACTIONS\n");
1755 :     @ReactionList = keys(%{$ReactionHash});
1756 :     foreach my $ReactionID (@ReactionList) {
1757 :     #Getting the thermodynamic reversibility from the database
1758 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1759 :     my $Subsystem = "NONE";
1760 :     if (defined($ReactionHash->{$ReactionID}->{"SUBSYSTEMS"})) {
1761 :     $Subsystem = join("|",@{$ReactionHash->{$ReactionID}->{"SUBSYSTEMS"}});
1762 :     }
1763 :     my $Source = "NONE";
1764 :     if (defined($ReactionHash->{$ReactionID}->{"SOURCE"})) {
1765 :     $Source = join("|",@{$ReactionHash->{$ReactionID}->{"SOURCE"}});
1766 :     }
1767 :     $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"]});
1768 :     }
1769 :    
1770 :     #Adding the spontaneous and universal reactions
1771 :     foreach my $ReactionID (@{$self->config("spontaneous reactions")}) {
1772 :     #Getting the thermodynamic reversibility from the database
1773 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1774 :     #Checking if the reaction is already in the model
1775 :     if (!defined($NewModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1776 :     $NewModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["SPONTANEOUS"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["SPONTANEOUS"],"NOTES" => ["NONE"]});
1777 :     }
1778 :     }
1779 :     foreach my $ReactionID (@{$self->config("universal reactions")}) {
1780 :     #Getting the thermodynamic reversibility from the database
1781 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1782 :     #Checking if the reaction is already in the model
1783 :     if (!defined($NewModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1784 :     $NewModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["UNIVERSAL"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["UNIVERSAL"],"NOTES" => ["NONE"]});
1785 :     }
1786 :     }
1787 :    
1788 :     #Checking if a biomass reaction already exists
1789 :     my $BiomassReactionRow = $self->BuildSpecificBiomassReaction();
1790 :     if (!defined($BiomassReactionRow)) {
1791 :     $self->set_status(-2,"Preliminary reconstruction failed: could not generate biomass reaction");
1792 :     $self->figmodel()->error_message("FIGMODEL:CreateModelReactionList: Could not generate biomass function for ".$self->id()."!");
1793 :     return $self->fail();
1794 :     }
1795 :     my $ReactionList = $BiomassReactionRow->{"ESSENTIAL REACTIONS"};
1796 :     push(@{$ReactionList},$BiomassReactionRow->{DATABASE}->[0]);
1797 :    
1798 :     #Adding biomass reactions to the model table
1799 :     foreach my $BOFReaction (@{$ReactionList}) {
1800 :     #Getting the thermodynamic reversibility from the database
1801 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($BOFReaction);
1802 :     #Checking if the reaction is already in the model
1803 :     if (!defined($NewModelTable->get_row_by_key($BOFReaction,"LOAD"))) {
1804 :     if ($BOFReaction =~ m/bio/) {
1805 :     $NewModelTable->add_row({"LOAD" => [$BOFReaction],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["BOF"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [1],"REFERENCE" => ["Biomass objective function"],"NOTES" => ["NONE"]});
1806 :     } else {
1807 :     $NewModelTable->add_row({"LOAD" => [$BOFReaction],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["INITIAL GAP FILLING"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [5],"REFERENCE" => ["Initial gap filling"],"NOTES" => ["NONE"]});
1808 :     }
1809 :     }
1810 :     }
1811 :    
1812 :     #Completing any incomplete reactions sets
1813 :     my $ReactionSetTable = $self->figmodel()->database()->GetDBTable("REACTION SETS");
1814 :     for (my $i=0; $i < $ReactionSetTable->size(); $i++) {
1815 :     if (defined($NewModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1816 :     foreach my $Reaction (@{$ReactionSetTable->get_row($i)->{"Dependant reactions"}}) {
1817 :     if (!defined($NewModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1818 :     #Getting the thermodynamic reversibility from the database
1819 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($Reaction);
1820 :     $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"]});
1821 :     }
1822 :     }
1823 :     }
1824 :     }
1825 :    
1826 :     #Now we compare the model to the previous model to determine if any differences exist that aren't gap filling reactions
1827 :     if (defined($OriginalModelTable)) {
1828 :     my $PerfectMatch = 1;
1829 :     my $ReactionCount = 0;
1830 :     for (my $i=0; $i < $OriginalModelTable->size(); $i++) {
1831 :     #We only check that nongapfilling reactions exist in the new model
1832 :     if ($OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] !~ m/GAP/ || $OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] eq "INITIAL GAP FILLING") {
1833 :     $ReactionCount++;
1834 :     my $Row = $NewModelTable->get_row_by_key($OriginalModelTable->get_row($i)->{"LOAD"}->[0],"LOAD");
1835 :     if (defined($Row)) {
1836 :     #We check that the reaction directionality is identical
1837 :     if ($Row->{"DIRECTIONALITY"}->[0] ne $OriginalModelTable->get_row($i)->{"DIRECTIONALITY"}->[0]) {
1838 :     if (defined($OriginalModelTable->get_row($i)->{"NOTES"}->[0]) && $OriginalModelTable->get_row($i)->{"NOTES"}->[0] =~ m/Directionality\sswitched\sfrom\s([^\s])/) {
1839 :     if ($1 ne $Row->{"DIRECTIONALITY"}->[0]) {
1840 :     print "Directionality mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".$1." vs ".$Row->{"DIRECTIONALITY"}->[0]."\n";
1841 :     $PerfectMatch = 0;
1842 :     last;
1843 :     }
1844 :     } else {
1845 :     print "Directionality mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".$OriginalModelTable->get_row($i)->{"DIRECTIONALITY"}->[0]." vs ".$Row->{"DIRECTIONALITY"}->[0]."\n";
1846 :     $PerfectMatch = 0;
1847 :     last;
1848 :     }
1849 :     }
1850 :     #We check that the genes assigned to the reaction are identical
1851 :     if ($PerfectMatch == 1 && @{$OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}} != @{$Row->{"ASSOCIATED PEG"}}) {
1852 :     print "Gene associatation mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".@{$OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}}." vs ".@{$Row->{"ASSOCIATED PEG"}}."\n";
1853 :     $PerfectMatch = 0;
1854 :     last;
1855 :     }
1856 :     if ($PerfectMatch == 1) {
1857 :     my @GeneSetOne = sort(@{$OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}});
1858 :     my @GeneSetTwo = sort(@{$Row->{"ASSOCIATED PEG"}});
1859 :     for (my $j=0; $j < @GeneSetOne; $j++) {
1860 :     if ($GeneSetOne[$j] ne $GeneSetTwo[$j]) {
1861 :     print "Gene mismatch for reaction ".$OriginalModelTable->get_row($i)->{"LOAD"}->[0].": ".$GeneSetOne[$j]." vs ".$GeneSetTwo[$j]."\n";
1862 :     $PerfectMatch = 0;
1863 :     $i = $OriginalModelTable->size();
1864 :     last;
1865 :     }
1866 :     }
1867 :     }
1868 :     } else {
1869 :     print "Original model contains an extra reaction:".$OriginalModelTable->get_row($i)->{"LOAD"}->[0]."\n";
1870 :     $PerfectMatch = 0;
1871 :     last;
1872 :     }
1873 :     }
1874 :     }
1875 :     if ($PerfectMatch == 1 && $ReactionCount == $NewModelTable->size()) {
1876 :     #Bailing out of function as the model has not changed
1877 :     $self->set_status(1,"rebuild canceled because model has not changed");
1878 :     return $self->success();
1879 :     }
1880 :     }
1881 :    
1882 :     #Saving the new model to file
1883 :     $self->set_status(1,"Preliminary reconstruction complete");
1884 :     $self->figmodel()->database()->save_table($NewModelTable);
1885 :     $self->{_reaction_data} = $NewModelTable;
1886 :     #Clearing the previous model from the cache
1887 :     $self->figmodel()->database()->ClearDBModel($self->id(),1);
1888 :     #Updating the model stats table
1889 :     $self->update_stats_for_build();
1890 :     $self->PrintSBMLFile();
1891 :    
1892 :     #Adding model to gapfilling queue
1893 :     if (defined($RunGapFilling) && $RunGapFilling == 1) {
1894 :     $self->set_status(1,"Autocompletion queued");
1895 :     $self->figmodel()->add_job_to_queue("gapfillmodel?".$self->id(),"QSUB","cplex",$self->owner(),"BACK");
1896 :     }
1897 :     return $self->success();
1898 :     }
1899 :    
1900 :     =head3 CreateMetaGenomeReactionList
1901 :     Definition:
1902 :     (success/fail) = FIGMODELmodel->CreateMetaGenomeReactionList();
1903 :     Description:
1904 :     This is the code called to create or update the reaction list for a metgenome model
1905 :     =cut
1906 :    
1907 :     sub CreateMetaGenomeReactionList {
1908 :     my ($self) = @_;
1909 :    
1910 :     #Checking if the metagenome file exists
1911 :     if (!-e $self->config("raw MGRAST directory")->[0].$self->genome().".summary") {
1912 :     $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());
1913 :     }
1914 :     #Loading metagenome data
1915 :     my $MGRASTData = $self->figmodel()->database()->load_multiple_column_file($self->config("raw MGRAST directory")->[0].$self->genome().".summary","\t");
1916 :     if (!defined($MGRASTData)) {
1917 :     $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());
1918 :     }
1919 :    
1920 :     #Setting up needed variables
1921 :     my $OriginalModelTable = undef;
1922 :    
1923 :     #Checking status
1924 :     if ($self->status() < 0) {
1925 :     $self->set_status(0,"Preliminary reconstruction");
1926 :     } elsif ($self->status() == 0) {
1927 :     $self->error_message("FIGMODEL->CreateModelReactionList:Model is already being built. Canceling current build.");
1928 :     return $self->fail();
1929 :     } else {
1930 :     $OriginalModelTable = $self->reaction_table();
1931 :     $self->ArchiveModel();
1932 :     $self->set_status(0,"Rebuilding preliminary reconstruction");
1933 :     }
1934 :    
1935 :     #Getting the reaction table
1936 :     my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");
1937 :     #Creating model table
1938 :     my $ModelTable = FIGMODELTable->new(["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG","SUBSYSTEM","CONFIDENCE","REFERENCE","NOTES"],$self->directory().$self->id().".txt",["LOAD"],";","|","REACTIONS\n");
1939 :     for (my $i=0; $i < @{$MGRASTData};$i++) {
1940 :     #MD5,PEG,number of sims,role,sim e-scores
1941 :     my $Role = $MGRASTData->[$i]->[3];
1942 :     my $MD5 = $MGRASTData->[$i]->[0];
1943 :     my $peg = $MGRASTData->[$i]->[1];
1944 :     my $sims = $MGRASTData->[$i]->[4];
1945 :     $sims =~ s/;/,/g;
1946 :    
1947 :     #Checking for subsystems
1948 :     my $GeneSubsystems = $self->figmodel()->subsystems_of_role($Role);
1949 :     #Checking if there are reactions associated with this role
1950 :     my $ReactionHashArrayRef = $self->figmodel()->reactions_of_role($Role);
1951 :     if ($ReactionHashArrayRef != 0) {
1952 :     foreach my $Mapping (@{$ReactionHashArrayRef}) {
1953 :     if (defined($Mapping->{"REACTION"}) && defined($Mapping->{"MASTER"}) && defined($Mapping->{"SUBSYSTEM"}) && defined($Mapping->{"SOURCE"})) {
1954 :     if ($Mapping->{"REACTION"}->[0] =~ m/rxn\d\d\d\d\d/) {
1955 :     if ($Mapping->{"MASTER"}->[0] eq 1) {
1956 :     #Checking if the reaction is already in the model
1957 :     my $ReactionRow = $ModelTable->get_row_by_key($Mapping->{"REACTION"}->[0],"LOAD");
1958 :     if (!defined($ReactionRow)) {
1959 :     $ReactionRow = {"LOAD" => [$Mapping->{"REACTION"}->[0]],"DIRECTIONALITY" => [$self->figmodel()->reversibility_of_reaction($Mapping->{"REACTION"}->[0])],"COMPARTMENT" => ["c"]};
1960 :     $ModelTable->add_row($ReactionRow);
1961 :     }
1962 :     push(@{$ReactionRow->{"ASSOCIATED PEG"}},substr($peg,4));
1963 :     push(@{$ReactionRow->{"REFERENCE"}},$MD5.":".$Role);
1964 :     push(@{$ReactionRow->{"CONFIDENCE"}},$sims);
1965 :     if (defined($GeneSubsystems)) {
1966 :     push(@{$ReactionRow->{"SUBSYSTEM"}},@{$GeneSubsystems});
1967 :     }
1968 :     }
1969 :     }
1970 :     }
1971 :     }
1972 :     }
1973 :     }
1974 :    
1975 :     #Adding the spontaneous and universal reactions
1976 :     foreach my $ReactionID (@{$self->config("spontaneous reactions")}) {
1977 :     #Getting the thermodynamic reversibility from the database
1978 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1979 :     #Checking if the reaction is already in the model
1980 :     if (!defined($ModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1981 :     $ModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["SPONTANEOUS"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["SPONTANEOUS"],"NOTES" => ["NONE"]});
1982 :     }
1983 :     }
1984 :     foreach my $ReactionID (@{$self->config("universal reactions")}) {
1985 :     #Getting the thermodynamic reversibility from the database
1986 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($ReactionID);
1987 :     #Checking if the reaction is already in the model
1988 :     if (!defined($ModelTable->get_row_by_key($ReactionID,"LOAD"))) {
1989 :     $ModelTable->add_row({"LOAD" => [$ReactionID],"DIRECTIONALITY" => [$Directionality],"COMPARTMENT" => ["c"],"ASSOCIATED PEG" => ["UNIVERSAL"],"SUBSYSTEM" => ["NONE"],"CONFIDENCE" => [4],"REFERENCE" => ["UNIVERSAL"],"NOTES" => ["NONE"]});
1990 :     }
1991 :     }
1992 :    
1993 :     #Completing any incomplete reactions sets
1994 :     my $ReactionSetTable = $self->figmodel()->database()->GetDBTable("REACTION SETS");
1995 :     for (my $i=0; $i < $ReactionSetTable->size(); $i++) {
1996 :     if (defined($ModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1997 :     foreach my $Reaction (@{$ReactionSetTable->get_row($i)->{"Dependant reactions"}}) {
1998 :     if (!defined($ModelTable->get_row_by_key($ReactionSetTable->get_row($i)->{"Trigger reaction"}->[0],"LOAD"))) {
1999 :     #Getting the thermodynamic reversibility from the database
2000 :     my $Directionality = $self->figmodel()->reversibility_of_reaction($Reaction);
2001 :     $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"]});
2002 :     }
2003 :     }
2004 :     }
2005 :     }
2006 :    
2007 :     #Clearing the previous model from the cache
2008 :     $self->figmodel()->database()->ClearDBModel($self->id(),1);
2009 :     $ModelTable->save();
2010 :    
2011 :     return $self->success();
2012 :     }
2013 :    
2014 :     =head3 ArchiveModel
2015 :     Definition:
2016 :     (success/fail) = FIGMODELmodel->ArchiveModel();
2017 :     Description:
2018 :     This function archives the specified model in the model directory with the current version numbers appended.
2019 :     This function is used to preserve old versions of models prior to overwriting so new versions may be compared with old versions.
2020 :     =cut
2021 :     sub ArchiveModel {
2022 :     my ($self) = @_;
2023 :    
2024 :     #Checking that the model file exists
2025 :     if (!(-e $self->filename())) {
2026 :     $self->figmodel()->error_message("FIGMODEL:ArchiveModel: Model file ".$self->filename()." not found!");
2027 :     return $self->fail();
2028 :     }
2029 :    
2030 :     #Copying the model file
2031 :     system("cp ".$self->filename()." ".$self->directory().$self->id().$self->version().".txt");
2032 :     }
2033 :    
2034 :     =head3 PrintModelDataToFile
2035 :     Definition:
2036 :     (success/fail) = FIGMODELmodel->PrintModelDataToFile();
2037 :     Description:
2038 :     This function uses the MFAToolkit to print out all of the compound and reaction data for the input model.
2039 :     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
2040 :     function. The LoadModel function for example would not give you this data.
2041 :     =cut
2042 :     sub PrintModelDataToFile {
2043 :     my($self) = @_;
2044 :    
2045 :     #Running the MFAToolkit on the model file
2046 :     my $OutputIndex = $self->figmodel()->filename();
2047 :     my $Command = $self->config("MFAToolkit executable")->[0]." parameterfile ../Parameters/Printing.txt resetparameter output_folder ".$OutputIndex.'/ LoadCentralSystem "'.$self->filename().'"';
2048 :     system($Command);
2049 :    
2050 :     #Copying the model file printed by the toolkit out of the output directory and into the model directory
2051 :     if (!-e $self->config("MFAToolkit output directory")->[0].$OutputIndex."/".$self->id().$self->selected_version().".txt") {
2052 :     $self->figmodel()->error_message("New model file not created due to an error. Check that the input modelfile exists.");
2053 :     $self->figmodel()->cleardirectory($OutputIndex);
2054 :     return $self->fail();
2055 :     }
2056 :    
2057 :     $Command = 'cp "'.$self->config("MFAToolkit output directory")->[0].$OutputIndex."/".$self->id().$self->selected_version().'.txt" "'.$self->directory().$self->id().$self->selected_version().'Data.txt"';
2058 :     system($Command);
2059 :     $Command = 'cp "'.$self->config("MFAToolkit output directory")->[0].$OutputIndex.'/ErrorLog0.txt" "'.$self->directory().'ModelErrors.txt"';
2060 :     system($Command);
2061 :     $self->figmodel()->cleardirectory($OutputIndex);
2062 :     return $self->success();
2063 :     }
2064 :    
2065 :     =head2 Analysis Functions
2066 :    
2067 :     =head3 run_microarray_analysis
2068 :     Definition:
2069 :     int::status = FIGMODEL->run_microarray_analysis(string::media,string::job id,string::gene calls);
2070 :     Description:
2071 :     Runs microarray analysis attempting to turn off genes that are inactive in the microarray
2072 :     =cut
2073 :     sub run_microarray_analysis {
2074 :     my ($self,$media,$label,$index,$genecall) = @_;
2075 :     $genecall =~ s/_/:/g;
2076 :     $genecall =~ s/\//;/g;
2077 :     my $uniqueFilename = $self->figmodel()->filename();
2078 :     my $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($uniqueFilename,$self->id(),$media,["ProductionMFA","ShewenellaExperiment"],{"Microarray assertions" => $label.";".$index.";".$genecall,"MFASolver" => "CPLEX","Network output location" => "/scratch/"},"MicroarrayAnalysis-".$uniqueFilename.".txt",undef,$self->selected_version());
2079 :     system($command);
2080 :     my $filename = $self->figmodel()->config("MFAToolkit output directory")->[0].$uniqueFilename."/MicroarrayOutput-".$index.".txt";
2081 :     if (-e $filename) {
2082 :     my $output = $self->figmodel()->database()->load_single_column_file($filename);
2083 :     if (defined($output->[0])) {
2084 :     my @array = split(/;/,$output->[0]);
2085 :     $self->figmodel()->clearing_output($uniqueFilename,"MicroarrayAnalysis-".$uniqueFilename.".txt");
2086 :     return ($array[0],$array[1],$array[8].":".$array[2],$array[9].":".$array[3],$array[10].":".$array[4],$array[11].":".$array[5],$array[12].":".$array[6],$array[13].":".$array[7]);
2087 :     }
2088 :     print STDERR $filename." is empty!";
2089 :     }
2090 :     print STDERR $filename." not found!";
2091 :     $self->figmodel()->clearing_output($uniqueFilename,"MicroarrayAnalysis-".$uniqueFilename.".txt");
2092 :    
2093 :     return undef;
2094 :     }
2095 :    
2096 :     =head3 find_minimal_pathways
2097 :     Definition:
2098 :     int::status = FIGMODEL->find_minimal_pathways(string::media,string::objective);
2099 :     Description:
2100 :     Runs microarray analysis attempting to turn off genes that are inactive in the microarray
2101 :     =cut
2102 :     sub find_minimal_pathways {
2103 :     my ($self,$media,$objective,$solutionnum,$AllReversible,$additionalexchange) = @_;
2104 :    
2105 :     #Setting default media
2106 :     if (!defined($media)) {
2107 :     $media = "Complete";
2108 :     }
2109 :    
2110 :     #Setting default solution number
2111 :     if (!defined($solutionnum)) {
2112 :     $solutionnum = "5";
2113 :     }
2114 :    
2115 :     #Setting additional exchange fluxes
2116 :     if (!defined($additionalexchange) || length($additionalexchange) == 0) {
2117 :     if ($self->id() eq "iAF1260") {
2118 :     $additionalexchange = "cpd03422[c]:-100:100;cpd01997[c]:-100:100;cpd11416[c]:-100:0;cpd15378[c]:-100:0;cpd15486[c]:-100:0";
2119 :     } else {
2120 :     $additionalexchange = $self->figmodel()->config("default exchange fluxes")->[0];
2121 :     }
2122 :     }
2123 :    
2124 :     #Translating objective
2125 :     my $objectivestring;
2126 :     if ($objective eq "ALL") {
2127 :     #Getting the list of universal building blocks
2128 :     my $buildingblocks = $self->config("universal building blocks");
2129 :     my @objectives = keys(%{$buildingblocks});
2130 :     #Getting the nonuniversal building blocks
2131 :     my $otherbuildingblocks = $self->config("nonuniversal building blocks");
2132 :     my @array = keys(%{$otherbuildingblocks});
2133 :     if (defined($self->get_biomass()) && defined($self->figmodel()->get_reaction($self->get_biomass()->{"LOAD"}->[0]))) {
2134 :     my $equation = $self->figmodel()->get_reaction($self->get_biomass()->{"LOAD"}->[0])->{"EQUATION"}->[0];
2135 :     if (defined($equation)) {
2136 :     for (my $i=0; $i < @array; $i++) {
2137 :     if (CORE::index($equation,$array[$i]) > 0) {
2138 :     push(@objectives,$array[$i]);
2139 :     }
2140 :     }
2141 :     }
2142 :     }
2143 :     for (my $i=0; $i < @objectives; $i++) {
2144 :     $self->find_minimal_pathways($media,$objectives[$i]);
2145 :     }
2146 :     return;
2147 :     } elsif ($objective eq "ENERGY") {
2148 :     $objectivestring = "MAX;FLUX;rxn00062;c;1";
2149 :     } elsif ($objective =~ m/cpd\d\d\d\d\d/) {
2150 :     if ($objective =~ m/\[(\w)\]/) {
2151 :     $objectivestring = "MIN;DRAIN_FLUX;".$objective.";".$1.";1";
2152 :     $additionalexchange .= ";".$objective."[".$1."]:-100:0";
2153 :     } else {
2154 :     $objectivestring = "MIN;DRAIN_FLUX;".$objective.";c;1";
2155 :     $additionalexchange .= ";".$objective."[c]:-100:0";
2156 :     }
2157 :     } elsif ($objective =~ m/(rxn\d\d\d\d\d)/) {
2158 :     my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateData($objective);
2159 :     for (my $i=0; $i < @{$Products};$i++) {
2160 :     my $temp = $Products->[$i]->{"DATABASE"}->[0];
2161 :     if ($additionalexchange !~ m/$temp/) {
2162 :     #$additionalexchange .= ";".$temp."[c]:-100:0";
2163 :     }
2164 :     }
2165 :     for (my $i=0; $i < @{$Reactants};$i++) {
2166 :     print $Reactants->[$i]->{"DATABASE"}->[0]." started\n";
2167 :     $self->find_minimal_pathways($media,$Reactants->[$i]->{"DATABASE"}->[0],$additionalexchange);
2168 :     print $Reactants->[$i]->{"DATABASE"}->[0]." done\n";
2169 :     }
2170 :     return;
2171 :     }
2172 :    
2173 :     #Adding additional drains
2174 :     if (($objective eq "cpd15665" || $objective eq "cpd15667" || $objective eq "cpd15668" || $objective eq "cpd15669") && $additionalexchange !~ m/cpd15666/) {
2175 :     $additionalexchange .= ";cpd15666[c]:0:100";
2176 :     } elsif ($objective eq "cpd11493" && $additionalexchange !~ m/cpd12370/) {
2177 :     $additionalexchange .= ";cpd12370[c]:0:100";
2178 :     } elsif ($objective eq "cpd00166" && $additionalexchange !~ m/cpd01997/) {
2179 :     $additionalexchange .= ";cpd01997[c]:0:100;cpd03422[c]:0:100";
2180 :     }
2181 :    
2182 :     #Running MFAToolkit
2183 :     my $filename = $self->figmodel()->filename();
2184 :     my $command;
2185 :     if (defined($AllReversible) && $AllReversible == 1) {
2186 :     $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($filename,$self->id(),$media,["ProductionMFA"],{"Make all reactions reversible in MFA"=>1, "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());
2187 :     } else {
2188 :     $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($filename,$self->id(),$media,["ProductionMFA"],{"Make all reactions reversible in MFA"=>0, "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());
2189 :     }
2190 :     system($command);
2191 :    
2192 :     #Loading problem report
2193 :     my $results = $self->figmodel()->LoadProblemReport($filename);
2194 :     #Clearing output
2195 :     $self->figmodel()->clearing_output($filename,"MinimalPathways-".$media."-".$self->id()."-".$objective.".txt");
2196 :     if (!defined($results)) {
2197 :     print STDERR $objective." pathway results not found!\n";
2198 :     return;
2199 :     }
2200 :    
2201 :     #Parsing output
2202 :     my @Array;
2203 :     my $row = $results->get_row(1);
2204 :     if (defined($row->{"Notes"}->[0])) {
2205 :     $_ = $row->{"Notes"}->[0];
2206 :     @Array = /\d+:([^\|]+)\|/g;
2207 :     }
2208 :    
2209 :     #Writing output to file
2210 :     $self->figmodel()->database()->print_array_to_file($self->directory()."MinimalPathways-".$media."-".$objective."-".$self->id()."-".$AllReversible."-".$self->selected_version().".txt",[join("|",@Array)]);
2211 :     }
2212 :    
2213 :     =head3 find_minimal_pathways
2214 :     Definition:
2215 :     int::status = FIGMODEL->find_minimal_pathways(string::media,string::objective);
2216 :     Description:
2217 :     Runs microarray analysis attempting to turn off genes that are inactive in the microarray
2218 :     =cut
2219 :     sub find_minimal_pathways_two {
2220 :     my ($self,$media,$objective,$solutionnum,$AllReversible,$additionalexchange) = @_;
2221 :    
2222 :     #Setting default media
2223 :     if (!defined($media)) {
2224 :     $media = "Complete";
2225 :     }
2226 :    
2227 :     #Setting default solution number
2228 :     if (!defined($solutionnum)) {
2229 :     $solutionnum = "5";
2230 :     }
2231 :    
2232 :     #Setting additional exchange fluxes
2233 :     if (!defined($additionalexchange) || length($additionalexchange) == 0) {
2234 :     if ($self->id() eq "iAF1260") {
2235 :     $additionalexchange = "cpd03422[c]:-100:100;cpd01997[c]:-100:100;cpd11416[c]:-100:0;cpd15378[c]:-100:0;cpd15486[c]:-100:0";
2236 :     } else {
2237 :     $additionalexchange = $self->figmodel()->config("default exchange fluxes")->[0];
2238 :     }
2239 :     }
2240 :    
2241 :     #Translating objective
2242 :     my $objectivestring;
2243 :     if ($objective eq "ALL") {
2244 :     #Getting the list of universal building blocks
2245 :     my $buildingblocks = $self->config("universal building blocks");
2246 :     my @objectives = keys(%{$buildingblocks});
2247 :     #Getting the nonuniversal building blocks
2248 :     my $otherbuildingblocks = $self->config("nonuniversal building blocks");
2249 :     my @array = keys(%{$otherbuildingblocks});
2250 :     if (defined($self->get_biomass()) && defined($self->figmodel()->get_reaction($self->get_biomass()->{"LOAD"}->[0]))) {
2251 :     my $equation = $self->figmodel()->get_reaction($self->get_biomass()->{"LOAD"}->[0])->{"EQUATION"}->[0];
2252 :     if (defined($equation)) {
2253 :     for (my $i=0; $i < @array; $i++) {
2254 :     if (CORE::index($equation,$array[$i]) > 0) {
2255 :     push(@objectives,$array[$i]);
2256 :     }
2257 :     }
2258 :     }
2259 :     }
2260 :     for (my $i=0; $i < @objectives; $i++) {
2261 :     $self->find_minimal_pathways($media,$objectives[$i]);
2262 :     }
2263 :     return;
2264 :     } elsif ($objective eq "ENERGY") {
2265 :     $objectivestring = "MAX;FLUX;rxn00062;c;1";
2266 :     } elsif ($objective =~ m/cpd\d\d\d\d\d/) {
2267 :     if ($objective =~ m/\[(\w)\]/) {
2268 :     $objectivestring = "MIN;DRAIN_FLUX;".$objective.";".$1.";1";
2269 :     $additionalexchange .= ";".$objective."[".$1."]:-100:0";
2270 :     } else {
2271 :     $objectivestring = "MIN;DRAIN_FLUX;".$objective.";c;1";
2272 :     $additionalexchange .= ";".$objective."[c]:-100:0";
2273 :     }
2274 :     } elsif ($objective =~ m/(rxn\d\d\d\d\d)/) {
2275 :     my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateData($objective);
2276 :     for (my $i=0; $i < @{$Products};$i++) {
2277 :     my $temp = $Products->[$i]->{"DATABASE"}->[0];
2278 :     if ($additionalexchange !~ m/$temp/) {
2279 :     #$additionalexchange .= ";".$temp."[c]:-100:0";
2280 :     }
2281 :     }
2282 :     for (my $i=0; $i < @{$Reactants};$i++) {
2283 :     print $Reactants->[$i]->{"DATABASE"}->[0]." started\n";
2284 :     $self->find_minimal_pathways($media,$Reactants->[$i]->{"DATABASE"}->[0],$additionalexchange);
2285 :     print $Reactants->[$i]->{"DATABASE"}->[0]." done\n";
2286 :     }
2287 :     return;
2288 :     }
2289 :    
2290 :     #Adding additional drains
2291 :     if (($objective eq "cpd15665" || $objective eq "cpd15667" || $objective eq "cpd15668" || $objective eq "cpd15669") && $additionalexchange !~ m/cpd15666/) {
2292 :     $additionalexchange .= ";cpd15666[c]:0:100";
2293 :     } elsif ($objective eq "cpd11493" && $additionalexchange !~ m/cpd12370/) {
2294 :     $additionalexchange .= ";cpd12370[c]:0:100";
2295 :     } elsif ($objective eq "cpd00166" && $additionalexchange !~ m/cpd01997/) {
2296 :     $additionalexchange .= ";cpd01997[c]:0:100;cpd03422[c]:0:100";
2297 :     }
2298 :    
2299 :     #Running MFAToolkit
2300 :     my $filename = $self->figmodel()->filename();
2301 :     my $command;
2302 :     if (defined($AllReversible) && $AllReversible == 1) {
2303 :     $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($filename,$self->id(),$media,["ProductionMFA"],{"use simple variable and constraint names"=>1,"Make all reactions reversible in MFA"=>1, "Recursive MILP solution limit" => $solutionnum,"Generate pathways to objective" => 1,"MFASolver" => "SCIP","objective" => $objectivestring,"exchange species" => $additionalexchange},"MinimalPathways-".$media."-".$self->id().$self->selected_version().".txt",undef,$self->selected_version());
2304 :     } else {
2305 :     $command = $self->figmodel()->GenerateMFAToolkitCommandLineCall($filename,$self->id(),$media,["ProductionMFA"],{"use simple variable and constraint names"=>1,"Make all reactions reversible in MFA"=>0, "Recursive MILP solution limit" => $solutionnum,"Generate pathways to objective" => 1,"MFASolver" => "SCIP","objective" => $objectivestring,"exchange species" => $additionalexchange},"MinimalPathways-".$media."-".$self->id().$self->selected_version().".txt",undef,$self->selected_version());
2306 :     }
2307 :     print $command."\n";
2308 :     system($command);
2309 :    
2310 :     #Loading problem report
2311 :     my $results = $self->figmodel()->LoadProblemReport($filename);
2312 :     #Clearing output
2313 :     $self->figmodel()->clearing_output($filename,"MinimalPathways-".$media."-".$self->id()."-".$objective.".txt");
2314 :     if (!defined($results)) {
2315 :     print STDERR $objective." pathway results not found!\n";
2316 :     return;
2317 :     }
2318 :    
2319 :     #Parsing output
2320 :     my @Array;
2321 :     my $row = $results->get_row(1);
2322 :     if (defined($row->{"Notes"}->[0])) {
2323 :     $_ = $row->{"Notes"}->[0];
2324 :     @Array = /\d+:([^\|]+)\|/g;
2325 :     }
2326 :    
2327 :     #Writing output to file
2328 :     $self->figmodel()->database()->print_array_to_file($self->directory()."MinimalPathways-".$media."-".$objective."-".$self->id()."-".$AllReversible."-".$self->selected_version().".txt",[join("|",@Array)]);
2329 :     }
2330 :    
2331 :     sub combine_minimal_pathways {
2332 :     my ($self) = @_;
2333 :    
2334 :     my $tbl;
2335 :     if (-e $self->directory()."MinimalPathwayTable-".$self->id().$self->selected_version().".tbl") {
2336 :     $tbl = FIGMODELTable::load_table($self->directory()."MinimalPathwayTable-".$self->id().$self->selected_version().".tbl",";","|",0,["Objective","Media","Reversible"]);
2337 :     } else {
2338 :     $tbl = FIGMODELTable->new(["Objective","Media","Reactions","Reversible","Shortest path","Number of essentials","Essentials","Length"],$self->directory()."MinimalPathwayTable-".$self->id().$self->selected_version().".tbl",["Objective","Media","Reversible"],";","|");
2339 :     }
2340 :     my @files = glob($self->directory()."MinimalPathways-*");
2341 :     for (my $i=0; $i < @files;$i++) {
2342 :     if ($files[$i] =~ m/MinimalPathways\-(\S+)\-(cpd\d\d\d\d\d)\-(\w+)\-(\d)\-/ || $files[$i] =~ m/MinimalPathways\-(\S+)\-(ENERGY)\-(\w+)\-(\d)\-/) {
2343 :     my $reactions = $self->figmodel()->database()->load_single_column_file($files[$i],"");
2344 :     if (defined($reactions) && @{$reactions} > 0 && length($reactions->[0]) > 0) {
2345 :     my $newrow = {"Objective"=>[$2],"Media"=>[$1],"Reversible"=>[$4]};
2346 :     my $row = $tbl->get_table_by_key($newrow->{"Objective"}->[0],"Objective")->get_table_by_key($newrow->{"Media"}->[0],"Media")->get_row_by_key($newrow->{"Reversible"}->[0],"Reversible");
2347 :     if (!defined($row)) {
2348 :     $row = $tbl->add_row($newrow);
2349 :     }
2350 :     $row->{Reactions} = $self->figmodel()->database()->load_single_column_file($files[$i],"");
2351 :     delete($row->{"Shortest path"});
2352 :     delete($row->{"Number of essentials"});
2353 :     delete($row->{"Essentials"});
2354 :     delete($row->{"Length"});
2355 :     for (my $j=0; $j < @{$row->{Reactions}}; $j++) {
2356 :     my @array = split(/,/,$row->{Reactions}->[$j]);
2357 :     $row->{"Length"}->[$j] = @array;
2358 :     if (!defined($row->{"Shortest path"}->[0]) || $row->{"Length"}->[$j] < $row->{"Shortest path"}->[0]) {
2359 :     $row->{"Shortest path"}->[0] = $row->{"Length"}->[$j];
2360 :     }
2361 :     $row->{"Number of essentials"}->[0] = 0;
2362 :     for (my $k=0; $k < @array;$k++) {
2363 :     if ($array[$k] =~ m/(rxn\d\d\d\d\d)/) {
2364 :     my $class = $self->get_reaction_class($1,1);
2365 :     my $temp = $row->{Media}->[0].":Essential";
2366 :     if ($class =~ m/$temp/) {
2367 :     $row->{"Number of essentials"}->[$j]++;
2368 :     if (!defined($row->{"Essentials"}->[$j]) && length($row->{"Essentials"}->[$j]) > 0) {
2369 :     $row->{"Essentials"}->[$j] = $array[$k];
2370 :     } else {
2371 :     $row->{"Essentials"}->[$j] .= ",".$array[$k];
2372 :     }
2373 :     }
2374 :     }
2375 :     }
2376 :     }
2377 :     }
2378 :     }
2379 :     }
2380 :     $tbl->save();
2381 :     }
2382 :    
2383 :     =head3 calculate_growth
2384 :     Definition:
2385 :     string::growth = FIGMODELmodel->calculate_growth(string:media);
2386 :     Description:
2387 :     Calculating growth in the input media
2388 :     =cut
2389 :     sub calculate_growth {
2390 :     my ($self,$Media) = @_;
2391 :     my $UniqueFilename = $self->figmodel()->filename();
2392 :     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()));
2393 :     my $ProblemReport = $self->figmodel()->LoadProblemReport($UniqueFilename);
2394 :     my $Result;
2395 :     if (defined($ProblemReport)) {
2396 :     my $Row = $ProblemReport->get_row(0);
2397 :     if (defined($Row) && defined($Row->{"Objective"}->[0])) {
2398 :     if ($Row->{"Objective"}->[0] < 0.00000001) {
2399 :     $Result = "NOGROWTH:".$Row->{"Individual metabolites with zero production"}->[0];
2400 :     } else {
2401 :     $Result = $Row->{"Objective"}->[0];
2402 :     }
2403 :     }
2404 :     }
2405 :     return $Result;
2406 :     }
2407 :    
2408 :     =head3 classify_model_reactions
2409 :     Definition:
2410 :     (FIGMODELTable:Reaction classes,FIGMODELTable:Compound classes) = FIGMODELmodel->classify_model_reactions(string:media);
2411 :     Description:
2412 :     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.
2413 :     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".
2414 :     Possible values for "Class" include:
2415 :     1.) Positive: these reactions are essential in the forward direction.
2416 :     2.) Negative: these reactions are essential in the reverse direction.
2417 :     3.) Positive variable: these reactions are nonessential, but they only ever proceed in the forward direction.
2418 :     4.) Negative variable: these reactions are nonessential, but they only ever proceed in the reverse direction.
2419 :     5.) Variable: these reactions are nonessential and proceed in the forward or reverse direction.
2420 :     6.) Blocked: these reactions never carry any flux at all in the media condition tested.
2421 :     7.) Dead: these reactions are disconnected from the network.
2422 :     =cut
2423 :     sub classify_model_reactions {
2424 :     my ($self,$Media,$SaveChanges) = @_;
2425 :    
2426 :     #Getting unique file for printing model output
2427 :     my $UniqueFilename = $self->figmodel()->filename();
2428 :     #Running the MFAToolkit
2429 :     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()));
2430 :     #Reading in the output bounds file
2431 :     my ($ReactionTB,$CompoundTB,$DeadCompounds,$DeadEndCompounds,$DeadReactions);
2432 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsReactionData0.txt") {
2433 :     $ReactionTB = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsReactionData0.txt",";","|",1,["DATABASE ID"]);
2434 :     }
2435 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsCompoundData0.txt") {
2436 :     $CompoundTB = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsCompoundData0.txt",";","|",1,["DATABASE ID"]);
2437 :     }
2438 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadReactions.txt") {
2439 :     $DeadReactions = $self->figmodel()->put_array_in_hash($self->figmodel()->database()->load_single_column_file($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadReactions.txt",""));
2440 :     }
2441 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadMetabolites.txt") {
2442 :     $DeadCompounds = $self->figmodel()->put_array_in_hash($self->figmodel()->database()->load_single_column_file($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadMetabolites.txt",""));
2443 :     }
2444 :     if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadEndMetabolites.txt") {
2445 :     $DeadEndCompounds = $self->figmodel()->put_array_in_hash($self->figmodel()->database()->load_single_column_file($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/DeadEndMetabolites.txt",""));
2446 :     }
2447 :     if (!defined($ReactionTB) && !defined($CompoundTB)) {
2448 :     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";
2449 :     return (undef,undef);
2450 :     }
2451 :    
2452 :     #Clearing output
2453 :     $self->figmodel()->clearing_output($UniqueFilename,"Classify-".$self->id().$self->selected_version()."-".$UniqueFilename.".log");
2454 :     #Creating the table objects that will hold the results of the reaction classification
2455 :     my $rxnclasstable = $self->reaction_class_table();
2456 :     my $cpdclasstable = $self->compound_class_table();
2457 :     #Loading the compound table
2458 :     if (defined($CompoundTB)) {
2459 :     for (my $i=0; $i < $CompoundTB->size(); $i++) {
2460 :     my $Row = $CompoundTB->get_row($i);
2461 :     if (defined($Row->{"DATABASE ID"})) {
2462 :     #Getting the compound row
2463 :     my $CpdRow = $cpdclasstable->get_row_by_key($Row->{"DATABASE ID"}->[0].$Row->{COMPARTMENT}->[0],"COMPOUND",1);
2464 :     #Setting row values
2465 :     my $Max = 0;
2466 :     my $Min = 0;
2467 :     my $Class = "Unknown";
2468 :     if (defined($DeadCompounds) && defined($DeadCompounds->{$Row->{"DATABASE ID"}->[0]})) {
2469 :     $Class = "Dead";
2470 :     } elsif (defined($DeadEndCompounds) && defined($DeadEndCompounds->{$Row->{"DATABASE ID"}->[0]})) {
2471 :     $Class = "Deadend";
2472 :     } elsif (defined($Row->{"Min DRAIN_FLUX"}) && defined($Row->{"Max DRAIN_FLUX"}) && $Row->{"Min DRAIN_FLUX"}->[0] ne "1e+07") {
2473 :     $Max = $Row->{"Max DRAIN_FLUX"}->[0];
2474 :     $Min = $Row->{"Min DRAIN_FLUX"}->[0];
2475 :     if ($Row->{"Min DRAIN_FLUX"}->[0] > 0.00000001) {
2476 :     $Class = "Positive";
2477 :     } elsif ($Row->{"Max DRAIN_FLUX"}->[0] < -0.00000001) {
2478 :     $Class = "Negative";
2479 :     } elsif ($Row->{"Min DRAIN_FLUX"}->[0] < -0.00000001) {
2480 :     if ($Row->{"Max DRAIN_FLUX"}->[0] > 0.00000001) {
2481 :     $Class = "Variable";
2482 :     } else {
2483 :     $Max = 0;
2484 :     $Class = "Negative variable";
2485 :     }
2486 :     } elsif ($Row->{"Max DRAIN_FLUX"}->[0] > 0.00000001) {
2487 :     $Min = 0;
2488 :     $Class = "Positive variable";
2489 :     } else {
2490 :     $Min = 0;
2491 :     $Max = 0;
2492 :     $Class = "Blocked";
2493 :     }
2494 :     }
2495 :     my $index = 0;
2496 :     if (defined($CpdRow->{MEDIA})) {
2497 :     for (my $i=0; $i < @{$CpdRow->{MEDIA}};$i++) {
2498 :     $index++;
2499 :     if ($CpdRow->{MEDIA}->[$i] eq $Media) {
2500 :     $index = $i;
2501 :     last;
2502 :     }
2503 :     }
2504 :     }
2505 :     $CpdRow->{MIN}->[$index] = $Min;
2506 :     $CpdRow->{MAX}->[$index] = $Max;
2507 :     $CpdRow->{CLASS}->[$index] = $Class;
2508 :     $CpdRow->{MEDIA}->[$index] = $Media;
2509 :     }
2510 :     }
2511 :     if (!defined($SaveChanges) || $SaveChanges == 1) {
2512 :     $cpdclasstable->save();
2513 :     }
2514 :     }
2515 :     if (defined($ReactionTB)) {
2516 :     for (my $i=0; $i < $ReactionTB->size(); $i++) {
2517 :     my $Row = $ReactionTB->get_row($i);
2518 :     if (defined($Row->{"DATABASE ID"})) {
2519 :     #Getting the compound row
2520 :     my $Compartment = "c";
2521 :     if (defined($Row->{COMPARTMENT}->[0])) {
2522 :     $Compartment = $Row->{COMPARTMENT}->[0];
2523 :     }
2524 :     my $RxnRow = $rxnclasstable->get_row_by_key($Row->{"DATABASE ID"}->[0],"REACTION",1);
2525 :     my $Max = 0;
2526 :     my $Min = 0;
2527 :     my $Class = "Unknown";
2528 :     if (defined($DeadReactions) && defined($DeadReactions->{$Row->{"DATABASE ID"}->[0]})) {
2529 :     $Class = "Dead";
2530 :     } elsif (defined($Row->{"Min FLUX"}) && defined($Row->{"Max FLUX"})) {
2531 :     $Max = $Row->{"Max FLUX"}->[0];
2532 :     $Min = $Row->{"Min FLUX"}->[0];
2533 :     if ($Row->{"Min FLUX"}->[0] > 0.00000001) {
2534 :     $Class = "Positive";
2535 :     } elsif ($Row->{"Max FLUX"}->[0] < -0.00000001) {
2536 :     $Class = "Negative";
2537 :     } elsif ($Row->{"Min FLUX"}->[0] < -0.00000001) {
2538 :     if ($Row->{"Max FLUX"}->[0] > 0.00000001) {
2539 :     $Class = "Variable";
2540 :     } else {
2541 :     $Max = 0;
2542 :     $Class = "Negative variable";
2543 :     }
2544 :     } elsif ($Row->{"Max FLUX"}->[0] > 0.00000001) {
2545 :     $Min = 0;
2546 :     $Class = "Positive variable";
2547 :     } else {
2548 :     $Min = 0;
2549 :     $Max = 0;
2550 :     $Class = "Blocked";
2551 :     }
2552 :     }
2553 :     my $index = 0;
2554 :     if (defined($RxnRow->{MEDIA})) {
2555 :     for (my $i=0; $i < @{$RxnRow->{MEDIA}};$i++) {
2556 :     $index++;
2557 :     if ($RxnRow->{MEDIA}->[$i] eq $Media) {
2558 :     $index = $i;
2559 :     last;
2560 :     }
2561 :     }
2562 :     }
2563 :     $RxnRow->{MIN}->[$index] = $Min;
2564 :     $RxnRow->{MAX}->[$index] = $Max;
2565 :     $RxnRow->{CLASS}->[$index] = $Class;
2566 :     $RxnRow->{MEDIA}->[$index] = $Media;
2567 :     }
2568 :     }
2569 :     if (!defined($SaveChanges) || $SaveChanges == 1) {
2570 :     $rxnclasstable->save();
2571 :     }
2572 :     }
2573 :     return ($rxnclasstable,$cpdclasstable);
2574 :     }
2575 :    
2576 :     =head3 RunAllStudiesWithDataFast
2577 :     Definition:
2578 :     (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);
2579 :     Description:
2580 :     Simulates every experimental condition currently available for the model.
2581 :     =cut
2582 :    
2583 :     sub RunAllStudiesWithDataFast {
2584 :     my ($self,$Experiment,$PrintResults) = @_;
2585 :    
2586 :     #Printing lp and key file for model
2587 :     if (!-e $self->directory()."FBA-".$self->id().$self->selected_version().".lp") {
2588 :     $self->PrintModelLPFile();
2589 :     }
2590 :     my $UniqueFilename = $self->figmodel()->filename();
2591 :    
2592 :     #Determing the simulations that need to be run
2593 :     my $ExperimentalDataTable = $self->figmodel()->GetExperimentalDataTable($self->genome(),$Experiment);
2594 :     #Creating the table of jobs to submit
2595 :     my $JobArray = $self->GetSimulationJobTable($ExperimentalDataTable,$Experiment,$UniqueFilename);
2596 :     #Printing the job file
2597 :     if (!-d $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/") {
2598 :     system("mkdir ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/");
2599 :     }
2600 :     $JobArray->save();
2601 :    
2602 :     #Running simulations
2603 :     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");
2604 :     #Parsing the results
2605 :     my $Results = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt",";","\\|",0,undef);
2606 :     if (!defined($Results)) {
2607 :     $self->figmodel()->error_message("FIGMODELmodel:RunAllStudiesWithDataFast:Could not find simulation results: ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt");
2608 :     return undef;
2609 :     }
2610 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector,$SimulationResults) = $self->EvaluateSimulationResults($Results,$ExperimentalDataTable);
2611 :     #Printing results to file
2612 :     $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");
2613 :     $self->figmodel()->clearing_output($UniqueFilename);
2614 :    
2615 :     return ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector);
2616 :     }
2617 :    
2618 :     =head3 GetSimulationJobTable
2619 :     Definition:
2620 :     my $JobTable = $model->GetSimulationJobTable($Experiment,$PrintResults,$Version);
2621 :     Description:
2622 :     =cut
2623 :    
2624 :     sub GetSimulationJobTable {
2625 :     my ($self,$SimulationTable,$Experiment,$Folder) = @_;
2626 :    
2627 :     #Determing the simulations that need to be run
2628 :     if (!defined($SimulationTable)) {
2629 :     $SimulationTable = $self->figmodel()->GetExperimentalDataTable($self->genome(),$Experiment);
2630 :     if (!defined($SimulationTable)) {
2631 :     return undef;
2632 :     }
2633 :     }
2634 :    
2635 :     #Creating the job table
2636 :     my $JobTable = $self->figmodel()->CreateJobTable($Folder);
2637 :     for (my $i=0; $i < $SimulationTable->size(); $i++) {
2638 :     if ($SimulationTable->get_row($i)->{"Heading"}->[0] =~ m/Gene\sKO/) {
2639 :     my $Row = $JobTable->get_row_by_key("Gene KO","LABEL",1);
2640 :     $JobTable->add_data($Row,"MEDIA",$SimulationTable->get_row($i)->{"Media"}->[0],1);
2641 :     } elsif ($SimulationTable->get_row($i)->{"Heading"}->[0] =~ m/Media\sgrowth/) {
2642 :     my $Row = $JobTable->get_row_by_key("Growth phenotype","LABEL",1);
2643 :     $JobTable->add_data($Row,"MEDIA",$SimulationTable->get_row($i)->{"Media"}->[0],1);
2644 :     } elsif ($SimulationTable->get_row($i)->{"Heading"}->[0] =~ m/Interval\sKO/) {
2645 :     my $Row = $JobTable->get_row_by_key($SimulationTable->get_row($i)->{"Heading"}->[0],"LABEL",1);
2646 :     $JobTable->add_data($Row,"MEDIA",$SimulationTable->get_row($i)->{"Media"}->[0],1);
2647 :     $JobTable->add_data($Row,"GENE KO",$SimulationTable->get_row($i)->{"Experiment type"}->[0],1);
2648 :     }
2649 :     }
2650 :    
2651 :     #Filling in model specific elements of the job table
2652 :     for (my $i=0; $i < $JobTable->size(); $i++) {
2653 :     if ($JobTable->get_row($i)->{"LABEL"}->[0] =~ m/Gene\sKO/) {
2654 :     $JobTable->get_row($i)->{"RUNTYPE"}->[0] = "SINGLEKO";
2655 :     $JobTable->get_row($i)->{"SAVE NONESSENTIALS"}->[0] = 1;
2656 :     } else {
2657 :     $JobTable->get_row($i)->{"RUNTYPE"}->[0] = "GROWTH";
2658 :     $JobTable->get_row($i)->{"SAVE NONESSENTIALS"}->[0] = 0;
2659 :     }
2660 :     $JobTable->get_row($i)->{"LP FILE"}->[0] = $self->directory()."FBA-".$self->id().$self->selected_version();
2661 :     $JobTable->get_row($i)->{"MODEL"}->[0] = $self->directory().$self->id().$self->selected_version().".txt";
2662 :     $JobTable->get_row($i)->{"SAVE FLUXES"}->[0] = 0;
2663 :     }
2664 :    
2665 :     return $JobTable;
2666 :     }
2667 :    
2668 :     =head3 EvaluateSimulationResults
2669 :     Definition:
2670 :     (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);
2671 :     Description:
2672 :     Compares simulation results with experimental data to produce a table indicating where predictions are incorrect.
2673 :     =cut
2674 :    
2675 :     sub EvaluateSimulationResults {
2676 :     my ($self,$Results,$ExperimentalDataTable) = @_;
2677 :    
2678 :     #Comparing experimental results with simulation results
2679 :     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);
2680 :     my $FalsePostives = 0;
2681 :     my $FalseNegatives = 0;
2682 :     my $CorrectNegatives = 0;
2683 :     my $CorrectPositives = 0;
2684 :     my @Errorvector;
2685 :     my @HeadingVector;
2686 :     my $ReactionKOWithGeneHash;
2687 :     for (my $i=0; $i < $Results->size(); $i++) {
2688 :     if ($Results->get_row($i)->{"LABEL"}->[0] eq "Gene KO") {
2689 :     if (defined($Results->get_row($i)->{"REACTION KO WITH GENES"})) {
2690 :     for (my $j=0; $j < @{$Results->get_row($i)->{"REACTION KO WITH GENES"}}; $j++) {
2691 :     my @Temp = split(/:/,$Results->get_row($i)->{"REACTION KO WITH GENES"}->[$j]);
2692 :     if (defined($Temp[1]) && length($Temp[1]) > 0) {
2693 :     $ReactionKOWithGeneHash->{$Temp[0]} = $Temp[1];
2694 :     }
2695 :     }
2696 :     }
2697 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] == 0) {
2698 :     for (my $j=0; $j < @{$Results->get_row($i)->{"NONESSENTIALGENES"}}; $j++) {
2699 :     my $Row = $ExperimentalDataTable->get_row_by_key("Gene KO:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j],"Heading");
2700 :     if (defined($Row)) {
2701 :     my $KOReactions = "none";
2702 :     if (defined($ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]})) {
2703 :     $KOReactions = $ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]};
2704 :     }
2705 :     push(@HeadingVector,$Row->{"Heading"}->[0].":".$KOReactions);
2706 :     my $Status = "Unknown";
2707 :     if ($Row->{"Growth"}->[0] > 0) {
2708 :     $Status = "False negative";
2709 :     $FalseNegatives++;
2710 :     push(@Errorvector,3);
2711 :     } else {
2712 :     $Status = "False positive";
2713 :     $FalsePostives++;
2714 :     push(@Errorvector,2);
2715 :     }
2716 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Gene KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => [$KOReactions]});
2717 :     }
2718 :     }
2719 :     } else {
2720 :     for (my $j=0; $j < @{$Results->get_row($i)->{"ESSENTIALGENES"}}; $j++) {
2721 :     #print $j."\t".$Results->get_row($i)->{"ESSENTIALGENES"}->[$j]."\n";
2722 :     my $Row = $ExperimentalDataTable->get_row_by_key("Gene KO:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"ESSENTIALGENES"}->[$j],"Heading");
2723 :     if (defined($Row)) {
2724 :     my $KOReactions = "none";
2725 :     if (defined($ReactionKOWithGeneHash->{$Results->get_row($i)->{"ESSENTIALGENES"}->[$j]})) {
2726 :     $KOReactions = $ReactionKOWithGeneHash->{$Results->get_row($i)->{"ESSENTIALGENES"}->[$j]};
2727 :     }
2728 :     push(@HeadingVector,$Row->{"Heading"}->[0].":".$KOReactions);
2729 :     my $Status = "Unknown";
2730 :     if ($Row->{"Growth"}->[0] > 0) {
2731 :     $Status = "False negative";
2732 :     $FalseNegatives++;
2733 :     push(@Errorvector,3);
2734 :     } else {
2735 :     $Status = "Correct negative";
2736 :     $CorrectNegatives++;
2737 :     push(@Errorvector,1);
2738 :     }
2739 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Gene KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => [$KOReactions]});
2740 :     }
2741 :     }
2742 :     for (my $j=0; $j < @{$Results->get_row($i)->{"NONESSENTIALGENES"}}; $j++) {
2743 :     my $Row = $ExperimentalDataTable->get_row_by_key("Gene KO:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j],"Heading");
2744 :     if (defined($Row)) {
2745 :     my $KOReactions = "none";
2746 :     if (defined($ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]})) {
2747 :     $KOReactions = $ReactionKOWithGeneHash->{$Results->get_row($i)->{"NONESSENTIALGENES"}->[$j]};
2748 :     }
2749 :     push(@HeadingVector,$Row->{"Heading"}->[0].":".$KOReactions);
2750 :     my $Status = "Unknown";
2751 :     if ($Row->{"Growth"}->[0] > 0) {
2752 :     $Status = "Correct positive";
2753 :     $CorrectPositives++;
2754 :     push(@Errorvector,0);
2755 :     } else {
2756 :     $Status = "False positive";
2757 :     $FalsePostives++;
2758 :     push(@Errorvector,2);
2759 :     }
2760 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Gene KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => [$KOReactions]});
2761 :     }
2762 :     }
2763 :     }
2764 :     } elsif ($Results->get_row($i)->{"LABEL"}->[0] eq "Growth phenotype") {
2765 :     my $Row = $ExperimentalDataTable->get_row_by_key("Media growth:".$Results->get_row($i)->{"MEDIA"}->[0].":".$Results->get_row($i)->{"MEDIA"}->[0],"Heading");
2766 :     if (defined($Row)) {
2767 :     push(@HeadingVector,$Row->{"Heading"}->[0].":none");
2768 :     my $Status = "Unknown";
2769 :     if ($Row->{"Growth"}->[0] > 0) {
2770 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2771 :     $Status = "Correct positive";
2772 :     $CorrectPositives++;
2773 :     push(@Errorvector,0);
2774 :     } else {
2775 :     $Status = "False negative";
2776 :     $FalseNegatives++;
2777 :     push(@Errorvector,3);
2778 :     }
2779 :     } else {
2780 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2781 :     $Status = "False positive";
2782 :     $FalsePostives++;
2783 :     push(@Errorvector,2);
2784 :     } else {
2785 :     $Status = "Correct negative";
2786 :     $CorrectNegatives++;
2787 :     push(@Errorvector,1);
2788 :     }
2789 :     }
2790 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Media growth"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Media"}->[0]],"Reactions knocked out" => ["none"]});
2791 :     }
2792 :     } elsif ($Results->get_row($i)->{"LABEL"}->[0] =~ m/Interval\sKO/ && defined($Results->get_row($i)->{"KOGENES"}->[0])) {
2793 :     my $Row = $ExperimentalDataTable->get_row_by_key($Results->get_row($i)->{"LABEL"}->[0],"Heading");
2794 :     if (defined($Row)) {
2795 :     my $Status = "Unknown";
2796 :     if ($Row->{"Growth"}->[0] > 0) {
2797 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2798 :     $Status = "Correct positive";
2799 :     $CorrectPositives++;
2800 :     push(@Errorvector,0);
2801 :     } else {
2802 :     $Status = "False negative";
2803 :     $FalseNegatives++;
2804 :     push(@Errorvector,3);
2805 :     }
2806 :     } else {
2807 :     if ($Results->get_row($i)->{"OBJECTIVE"}->[0] > 0) {
2808 :     $Status = "False positive";
2809 :     $FalsePostives++;
2810 :     push(@Errorvector,2);
2811 :     } else {
2812 :     $Status = "Correct negative";
2813 :     $CorrectNegatives++;
2814 :     push(@Errorvector,1);
2815 :     }
2816 :     }
2817 :     $SimulationResults->add_row({"Run result" => [$Status],"Experiment type" => ["Interval KO"],"Media" => [$Row->{"Media"}->[0]],"Experiment ID" => [$Row->{"Experiment ID"}->[0]],"Reactions knocked out" => ["none"]});
2818 :     }
2819 :     }
2820 :     }
2821 :    
2822 :     return ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,join(";",@Errorvector),join(";",@HeadingVector),$SimulationResults);
2823 :     }
2824 :    
2825 :     =head3 InspectSolution
2826 :     Definition:
2827 :     $model->InspectSolution(string::gene knocked out,string::media condition,[string]::list of reactions);
2828 :     Description:
2829 :     =cut
2830 :    
2831 :     sub InspectSolution {
2832 :     my ($self,$GeneKO,$Media,$ReactionList) = @_;
2833 :    
2834 :     #Getting a directory for the results
2835 :     my $UniqueFilename = $self->figmodel()->filename();
2836 :     system("mkdir ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/");
2837 :     my $TempVersion = "V".$UniqueFilename;
2838 :    
2839 :     #Setting gene ko to none if no genes are to be knocked out
2840 :     if ($GeneKO !~ m/^peg\./) {
2841 :     $GeneKO = "none";
2842 :     }
2843 :    
2844 :     #Implementing the input solution in the test model
2845 :     my $ReactionArray;
2846 :     my $DirectionArray;
2847 :     my %SolutionHash;
2848 :     for (my $k=0; $k < @{$ReactionList}; $k++) {
2849 :     if ($ReactionList->[$k] =~ m/(.+)(rxn\d\d\d\d\d)/) {
2850 :     my $Reaction = $2;
2851 :     my $Sign = $1;
2852 :     if (defined($SolutionHash{$Reaction})) {
2853 :     $SolutionHash{$Reaction} = "<=>";
2854 :     } elsif ($Sign eq "-") {
2855 :     $SolutionHash{$Reaction} = "<=";
2856 :     } elsif ($Sign eq "+") {
2857 :     $SolutionHash{$Reaction} = "=>";
2858 :     } else {
2859 :     $SolutionHash{$Reaction} = $Sign;
2860 :     }
2861 :     }
2862 :     }
2863 :     my @TempList = keys(%SolutionHash);
2864 :     for (my $k=0; $k < @TempList; $k++) {
2865 :     push(@{$ReactionArray},$TempList[$k]);
2866 :     push(@{$DirectionArray},$SolutionHash{$TempList[$k]});
2867 :     }
2868 :    
2869 :     print "Integrating solution!\n";
2870 :     $self->figmodel()->IntegrateGrowMatchSolution($self->id().$self->selected_version(),$self->directory().$self->id().$TempVersion.".txt",$ReactionArray,$DirectionArray,"SolutionInspection",1,1);
2871 :    
2872 :     #Printing lp and key file for model
2873 :     $self->PrintModelLPFile();
2874 :    
2875 :     #Running FBA on the test model
2876 :     my $JobTable = $self->figmodel()->CreateJobTable($UniqueFilename);
2877 :     $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]});
2878 :     $JobTable->save();
2879 :    
2880 :     #Running simulations
2881 :     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");
2882 :    
2883 :     #Parsing the results
2884 :     my $Results = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt",";","\\|",0,undef);
2885 :     if (!defined($Results)) {
2886 :     $self->figmodel()->error_message("FIGMODELmodel:InspectSolution:Could not load problem report ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/Output.txt");
2887 :     return undef;
2888 :     }
2889 :    
2890 :     #Making sure that the model grew with all reactions present
2891 :     my $Found = 0;
2892 :     for (my $i=0; $i < $Results->size(); $i++) {
2893 :     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) {
2894 :     $Found = 1;
2895 :     }
2896 :     }
2897 :     if ($Found == 0) {
2898 :     print "Solution no longer valid\n";
2899 :     return undef;
2900 :     }
2901 :    
2902 :     #Making sure all of the reactions added are still necessary
2903 :     my $FinalReactionList;
2904 :     for (my $k=0; $k < $Results->size(); $k++) {
2905 :     if (defined($Results->get_row($k)->{"KOGENES"}->[0]) && $Results->get_row($k)->{"KOGENES"}->[0] eq $GeneKO) {
2906 :     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) {
2907 :     push(@{$FinalReactionList},$Results->get_row($k)->{"KOREACTIONS"}->[0]);
2908 :     }
2909 :     }
2910 :     }
2911 :    
2912 :     #Deleting extra files created
2913 :     unlink($self->directory()."FBA-".$self->id().$TempVersion.".lp");
2914 :     unlink($self->directory()."FBA-".$self->id().$TempVersion.".key");
2915 :     unlink($self->directory().$self->id().$TempVersion.".txt");
2916 :    
2917 :     #Deleting the test model and the MFA folder
2918 :     $self->figmodel()->clearing_output($UniqueFilename);
2919 :    
2920 :     return $FinalReactionList;
2921 :     }
2922 :    
2923 :     =head3 GapFillingAlgorithm
2924 :    
2925 :     Definition:
2926 :     FIGMODELmodel->GapFillingAlgorithm();
2927 :    
2928 :     Description:
2929 :     This is a wrapper for running the gap filling algorithm on any model in the database.
2930 :     The algorithm performs a gap filling for any false negative prediction of the avialable experimental data.
2931 :     This function is threaded to improve efficiency: one thread does nothing but using the MFAToolkit to fill gaps for every false negative prediction.
2932 :     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.
2933 :     This function prints two important output files in the Model directory:
2934 :     1.) GapFillingOutput.txt: this is a summary of the results of the gap filling analysis
2935 :     2.) GapFillingErrorMatrix.txt: this lists the correct and incorrect predictions for each gapfilling solution implemented in a test model.
2936 :     =cut
2937 :    
2938 :     sub GapFillingAlgorithm {
2939 :     my ($self) = @_;
2940 :    
2941 :     #First the input model version and model filename should be simulated and the false negatives identified
2942 :     my ($FalsePostives,$FalseNegatives,$CorrectNegatives,$CorrectPositives,$Errorvector,$HeadingVector) = $self->RunAllStudiesWithDataFast("All");
2943 :    
2944 :     #Getting the filename
2945 :     my $UniqueFilename = $self->figmodel()->filename();
2946 :    
2947 :     #Printing the original performance vector
2948 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-OPEM".".txt",[$HeadingVector,$Errorvector]);
2949 :    
2950 :     my $PreviousGapFilling;
2951 :     if (-e $self->directory().$self->id().$self->selected_version()."-GFS.txt") {
2952 :     #Backing up the old solution file
2953 :     system("cp ".$self->directory().$self->id().$self->selected_version()."-GFS.txt ".$self->directory().$self->id().$self->selected_version()."-OldGFS.txt");
2954 :     unlink($self->directory().$self->id().$self->selected_version()."-GFS.txt");
2955 :     }
2956 :     if (-e $self->directory().$self->id().$self->selected_version()."-OldGFS.txt") {
2957 :     #Reading in the solution file from the previous gap filling if it exists
2958 :     $PreviousGapFilling = $self->figmodel()->database()->load_table($self->directory().$self->id().$self->selected_version()."-OldGFS.txt",";",",",0,["Experiment"]);
2959 :     }
2960 :    
2961 :     #Now we use the simulation output to make the gap filling run data
2962 :     my @Errors = split(/;/,$Errorvector);
2963 :     my @Headings = split(/;/,$HeadingVector);
2964 :     my $GapFillingRunSpecs = "";
2965 :     my $Count = 0;
2966 :     my $RescuedPreviousResults;
2967 :     my $RunCount = 0;
2968 :     my $SolutionExistedCount = 0;
2969 :     my $AcceptedSolutions = 0;
2970 :     my $RejectedSolutions = 0;
2971 :     my $NoExistingSolutions = 0;
2972 :     for (my $i=0; $i < @Errors; $i++) {
2973 :     if ($Errors[$i] == 3) {
2974 :     my @HeadingDataArray = split(/:/,$Headings[$i]);
2975 :     if ($HeadingDataArray[2] !~ m/^peg\./ || $HeadingDataArray[3] ne "none") {
2976 :     my $SolutionFound = 0;
2977 :     if (defined($PreviousGapFilling) && defined($PreviousGapFilling->get_row_by_key($HeadingDataArray[2],"Experiment"))) {
2978 :     my @Rows = $PreviousGapFilling->get_rows_by_key($HeadingDataArray[2],"Experiment");
2979 :     for (my $j=0; $j < @Rows; $j++) {
2980 :     if ($HeadingDataArray[2] =~ m/^peg\./) {
2981 :     my $ReactionList = $self->InspectSolution($HeadingDataArray[2],$HeadingDataArray[1],$Rows[$j]->{"Solution reactions"});
2982 :     if (defined($ReactionList)) {
2983 :     print join(",",@{$Rows[$j]->{"Solution reactions"}})."\t".join(",",@{$ReactionList})."\n";
2984 :     $SolutionFound++;
2985 :     push(@{$RescuedPreviousResults},$Rows[$j]->{"Experiment"}->[0].";".$Rows[$j]->{"Solution index"}->[0].";".$Rows[$j]->{"Solution cost"}->[0].";".join(",",@{$ReactionList}));
2986 :     $AcceptedSolutions++;
2987 :     } else {
2988 :     $RejectedSolutions++;
2989 :     }
2990 :     } else {
2991 :     my $ReactionList = $self->InspectSolution($HeadingDataArray[2],$HeadingDataArray[1],$Rows[$j]->{"Solution reactions"});
2992 :     if (defined($ReactionList)) {
2993 :     print join(",",@{$Rows[$j]->{"Solution reactions"}})."\t".join(",",@{$ReactionList})."\n";
2994 :     $SolutionFound++;
2995 :     push(@{$RescuedPreviousResults},$Rows[$j]->{"Experiment"}->[0].";".$Rows[$j]->{"Solution index"}->[0].";".$Rows[$j]->{"Solution cost"}->[0].";".join(",",@{$ReactionList}));
2996 :     $AcceptedSolutions++;
2997 :     } else {
2998 :     $RejectedSolutions++;
2999 :     }
3000 :     }
3001 :     }
3002 :     } else {
3003 :     $NoExistingSolutions++;
3004 :     }
3005 :     if ($SolutionFound == 0) {
3006 :     $RunCount++;
3007 :     if (length($GapFillingRunSpecs) > 0) {
3008 :     $GapFillingRunSpecs .= ";";
3009 :     }
3010 :     $GapFillingRunSpecs .= $HeadingDataArray[2].":".$HeadingDataArray[1].":".$HeadingDataArray[3];
3011 :     } else {
3012 :     $SolutionExistedCount++;
3013 :     }
3014 :     }
3015 :     $Count++;
3016 :     }
3017 :     }
3018 :    
3019 :     #Updating the growmatch progress table
3020 :     my $Row = $self->figmodel()->database()->get_row_by_key("GROWMATCH TABLE",$self->genome(),"ORGANISM",1);
3021 :     $Row->{"INITIAL FP"}->[0] = $FalsePostives;
3022 :     $Row->{"INITIAL FN"}->[0] = $FalseNegatives;
3023 :     $Row->{"GF TIMING"}->[0] = time()."-";
3024 :     $Row->{"FN WITH SOL"}->[0] = $FalseNegatives-$NoExistingSolutions;
3025 :     $Row->{"FN WITH ACCEPTED SOL"}->[0] = $SolutionExistedCount;
3026 :     $Row->{"TOTAL ACCEPTED GF SOL"}->[0] = $AcceptedSolutions;
3027 :     $Row->{"TOTAL REJECTED GF SOL"}->[0] = $RejectedSolutions;
3028 :     $Row->{"FN WITH NO SOL"}->[0] = $NoExistingSolutions+$RejectedSolutions;
3029 :     $self->figmodel()->database()->update_row("GROWMATCH TABLE",$Row,"ORGANISM");
3030 :    
3031 :     #Running the gap filling once to correct all false negative errors
3032 :     my $SolutionsFound = 0;
3033 :     my $GapFillingArray;
3034 :     push(@{$GapFillingArray},split(/;/,$GapFillingRunSpecs));
3035 :     my $GapFillingResults = $self->datagapfill($GapFillingArray,"GFS");
3036 :     if (defined($GapFillingResults)) {
3037 :     $SolutionsFound = 1;
3038 :     }
3039 :    
3040 :     if (defined($RescuedPreviousResults) && @{$RescuedPreviousResults} > 0) {
3041 :     #Printing previous solutions to GFS file
3042 :     $self->figmodel()->database()->print_array_to_file($self->directory().$self->id().$self->selected_version()."-GFS.txt",$RescuedPreviousResults,1);
3043 :     $SolutionsFound = 1;
3044 :     }
3045 :    
3046 :     #Recording the finishing of the gapfilling
3047 :     $Row = $self->figmodel()->database()->get_row_by_key("GROWMATCH TABLE",$self->genome(),"ORGANISM",1);
3048 :     $Row->{"GF TIMING"}->[0] .= time();
3049 :     $self->figmodel()->database()->update_row("GROWMATCH TABLE",$Row,"ORGANISM");
3050 :    
3051 :     if ($SolutionsFound == 1) {
3052 :     #Scheduling solution testing
3053 :     $self->figmodel()->add_job_to_queue("testsolutions?".$self->id().$self->selected_version()."?-1?GF","QSUB","fast","master","BACK");
3054 :     } else {
3055 :     $self->figmodel()->error_message("No false negative predictions found. Data gap filling not necessary!");
3056 :     }
3057 :    
3058 :     return $self->success();
3059 :     }
3060 :    
3061 :     =head3 SolutionReconciliation
3062 :     Definition:
3063 :     FIGMODELmodel->SolutionReconciliation();
3064 :     Description:
3065 :     This is a wrapper for running the solution reconciliation algorithm on any model in the database.
3066 :     The algorithm performs a reconciliation of any gap filling solutions to identify the combination of solutions that results in the optimal model.
3067 :     This function prints out one output file in the Model directory: ReconciliationOutput.txt: this is a summary of the results of the reconciliation analysis
3068 :     =cut
3069 :    
3070 :     sub SolutionReconciliation {
3071 :     my ($self,$GapFill,$Stage) = @_;
3072 :    
3073 :     #Setting the output filenames
3074 :     my $OutputFilename;
3075 :     my $OutputFilenameTwo;
3076 :     if ($GapFill == 1) {
3077 :     $OutputFilename = $self->directory().$self->id().$self->selected_version()."-GFReconciliation.txt";
3078 :     $OutputFilenameTwo = $self->directory().$self->id().$self->selected_version()."-GFSRS.txt";
3079 :     } else {
3080 :     $OutputFilename = $self->directory().$self->id().$self->selected_version()."-GGReconciliation.txt";
3081 :     $OutputFilenameTwo = $self->directory().$self->id().$self->selected_version()."-GGSRS.txt";
3082 :     }
3083 :    
3084 :     #In stage one, we run the reconciliation and create a test file to check combined solution performance
3085 :     if (!defined($Stage) || $Stage == 1) {
3086 :     my $GrowMatchTable = $self->figmodel()->database()->LockDBTable("GROWMATCH TABLE");
3087 :     my $Row = $GrowMatchTable->get_row_by_key($self->genome(),"ORGANISM",1);
3088 :     $Row->{"GF RECONCILATION TIMING"}->[0] = time()."-";
3089 :     $GrowMatchTable->save();
3090 :     $self->figmodel()->database()->UnlockDBTable("GROWMATCH TABLE");
3091 :    
3092 :     #Getting a unique filename
3093 :     my $UniqueFilename = $self->figmodel()->filename();
3094 :    
3095 :     #Copying over the necessary files
3096 :     if ($GapFill == 1) {
3097 :     if (!-e $self->directory().$self->id().$self->selected_version()."-GFEM.txt") {
3098 :     print STDERR "FIGMODEL:SolutionReconciliation:".$self->directory().$self->id().$self->selected_version()."-GFEM.txt file not found. Could not reconcile!";
3099 :     return 0;
3100 :     }
3101 :     if (!-e $self->directory().$self->id().$self->selected_version()."-OPEM.txt") {
3102 :     print STDERR "FIGMODEL:SolutionReconciliation:".$self->directory().$self->id().$self->selected_version()."-OPEM.txt file not found. Could not reconcile!";
3103 :     return 0;
3104 :     }
3105 :     system("cp ".$self->directory().$self->id().$self->selected_version()."-GFEM.txt ".$self->figmodel()->config("MFAToolkit input files")->[0].$UniqueFilename."-GFEM.txt");
3106 :     system("cp ".$self->directory().$self->id().$self->selected_version()."-OPEM.txt ".$self->figmodel()->config("MFAToolkit input files")->[0].$UniqueFilename."-OPEM.txt");
3107 :     #Backing up and deleting the existing reconciliation file
3108 :     if (-e $OutputFilename) {
3109 :     system("cp ".$OutputFilename." ".$self->directory().$self->id().$self->selected_version()."-OldGFReconciliation.txt");
3110 :     unlink($OutputFilename);
3111 :     }
3112 :     } else {
3113 :     if (!-e $self->directory().$self->id().$self->selected_version()."-GGEM.txt") {
3114 :     print STDERR "FIGMODEL:SolutionReconciliation:".$self->directory().$self->id().$self->selected_version()."-GGEM.txt file not found. Could not reconcile!";
3115 :     return 0;
3116 :     }
3117 :     if (!-e $self->directory().$self->id().$self->selected_version()."-GGOPEM.txt") {
3118 :     print STDERR "FIGMODEL:SolutionReconciliation:".$self->directory().$self->id().$self->selected_version()."-GGOPEM.txt file not found. Could not reconcile!";
3119 :     return 0;
3120 :     }
3121 :     system("cp ".$self->directory().$self->id().$self->selected_version()."-GGEM.txt ".$self->figmodel()->config("MFAToolkit input files")->[0].$UniqueFilename."-GGEM.txt");
3122 :     system("cp ".$self->directory().$self->id().$self->selected_version()."-GGOPEM.txt ".$self->figmodel()->config("MFAToolkit input files")->[0].$UniqueFilename."-OPEM.txt");
3123 :     #Backing up and deleting the existing reconciliation file
3124 :     if (-e $OutputFilename) {
3125 :     system("cp ".$OutputFilename." ".$self->directory().$self->id().$self->selected_version()."-OldGGReconciliation.txt");
3126 :     unlink($OutputFilename);
3127 :     }
3128 :     }
3129 :    
3130 :     #Running the reconciliation
3131 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),"NONE",["SolutionReconciliation"],{"Solution data for model optimization" => $UniqueFilename},"Reconciliation".$UniqueFilename.".log",undef,$self->selected_version()));
3132 :     $GrowMatchTable = $self->figmodel()->database()->LockDBTable("GROWMATCH TABLE");
3133 :     $Row = $GrowMatchTable->get_row_by_key($self->genome(),"ORGANISM",1);
3134 :     $Row->{"GF RECONCILATION TIMING"}->[0] .= time();
3135 :     $GrowMatchTable->save();
3136 :     $self->figmodel()->database()->UnlockDBTable("GROWMATCH TABLE");
3137 :    
3138 :     #Loading the problem report from the reconciliation run
3139 :     my $ReconciliatonOutput = $self->figmodel()->LoadProblemReport($UniqueFilename);
3140 :     print $UniqueFilename."\n";
3141 :     #Clearing output files
3142 :     $self->figmodel()->clearing_output($UniqueFilename,"Reconciliation".$UniqueFilename.".log");
3143 :     $ReconciliatonOutput->save("/home/chenry/Test.txt");
3144 :    
3145 :     #Checking the a problem report was found and was loaded
3146 :     if (!defined($ReconciliatonOutput) || $ReconciliatonOutput->size() < 1 || !defined($ReconciliatonOutput->get_row(0)->{"Notes"}->[0])) {
3147 :     print STDERR "FIGMODEL:SolutionReconciliation: MFAToolkit output from SolutionReconciliation of ".$self->id()." not found!\n\n";
3148 :     return 0;
3149 :     }
3150 :    
3151 :     #Processing the solutions
3152 :     my $SolutionCount = 0;
3153 :     my $ReactionSetHash;
3154 :     my $SingleReactionHash;
3155 :     my $ReactionDataHash;
3156 :     for (my $n=0; $n < $ReconciliatonOutput->size(); $n++) {
3157 :     if (defined($ReconciliatonOutput->get_row($n)->{"Notes"}->[0]) && $ReconciliatonOutput->get_row($n)->{"Notes"}->[0] =~ m/^Recursive\sMILP\s([^;]+)/) {
3158 :     #Breaking up the solution into reaction sets
3159 :     my @ReactionSets = split(/\|/,$1);
3160 :     #Creating reaction lists for each set
3161 :     my $SolutionHash;
3162 :     for (my $i=0; $i < @ReactionSets; $i++) {
3163 :     if (length($ReactionSets[$i]) > 0) {
3164 :     my @Alternatives = split(/:/,$ReactionSets[$i]);
3165 :     for (my $j=1; $j < @Alternatives; $j++) {
3166 :     if (length($Alternatives[$j]) > 0) {
3167 :     push(@{$SolutionHash->{$Alternatives[$j]}},$Alternatives[0]);
3168 :     }
3169 :     }
3170 :     if (@Alternatives == 1) {
3171 :     $SingleReactionHash->{$Alternatives[0]}->{$SolutionCount} = 1;
3172 :     if (!defined($SingleReactionHash->{$Alternatives[0]}->{"COUNT"})) {
3173 :     $SingleReactionHash->{$Alternatives[0]}->{"COUNT"} = 0;
3174 :     }
3175 :     $SingleReactionHash->{$Alternatives[0]}->{"COUNT"}++;
3176 :     }
3177 :     }
3178 :     }
3179 :     #Identifying reactions sets and storing the sets in the reactions set hash
3180 :     foreach my $Solution (keys(%{$SolutionHash})) {
3181 :     my $SetKey = join(",",sort(@{$SolutionHash->{$Solution}}));
3182 :     if (!defined($ReactionSetHash->{$SetKey}->{$SetKey}->{$SolutionCount})) {
3183 :     $ReactionSetHash->{$SetKey}->{$SetKey}->{$SolutionCount} = 1;
3184 :     if (!defined($ReactionSetHash->{$SetKey}->{$SetKey}->{"COUNT"})) {
3185 :     $ReactionSetHash->{$SetKey}->{$SetKey}->{"COUNT"} = 0;
3186 :     }
3187 :     $ReactionSetHash->{$SetKey}->{$SetKey}->{"COUNT"}++;
3188 :     }
3189 :     $ReactionSetHash->{$SetKey}->{$Solution}->{$SolutionCount} = 1;
3190 :     if (!defined($ReactionSetHash->{$SetKey}->{$Solution}->{"COUNT"})) {
3191 :     $ReactionSetHash->{$SetKey}->{$Solution}->{"COUNT"} = 0;
3192 :     }
3193 :     $ReactionSetHash->{$SetKey}->{$Solution}->{"COUNT"}++;
3194 :     }
3195 :     $SolutionCount++;
3196 :     }
3197 :     }
3198 :    
3199 :     #Handling the scenario where no solutions were found
3200 :     if ($SolutionCount == 0) {
3201 :     print STDERR "FIGMODEL:SolutionReconciliation: Reconciliation unsuccessful. No solution found.\n\n";
3202 :     return 0;
3203 :     }
3204 :    
3205 :     #Printing results without solution performance figures. Also printing solution test file
3206 :     open (RECONCILIATION, ">$OutputFilename");
3207 :     #Printing the file heading
3208 :     print RECONCILIATION "DATABASE;DEFINITION;REVERSIBLITY;DELTAG;DIRECTION;NUMBER OF SOLUTIONS";
3209 :     for (my $i=0; $i < $SolutionCount; $i++) {
3210 :     print RECONCILIATION ";Solution ".$i;
3211 :     }
3212 :     print RECONCILIATION "\n";
3213 :     #Printing the singlet reactions first
3214 :     my $Solutions;
3215 :     print RECONCILIATION "SINGLET REACTIONS\n";
3216 :     my @SingletReactions = keys(%{$SingleReactionHash});
3217 :     for (my $j=0; $j < $SolutionCount; $j++) {
3218 :     $Solutions->[$j]->{"BASE"} = $j;
3219 :     }
3220 :     for (my $i=0; $i < @SingletReactions; $i++) {
3221 :     my $ReactionData;
3222 :     if (defined($ReactionDataHash->{$SingletReactions[$i]})) {
3223 :     $ReactionData = $ReactionDataHash->{$SingletReactions[$i]};
3224 :     } else {
3225 :     my $Direction = substr($SingletReactions[$i],0,1);
3226 :     if ($Direction eq "+") {
3227 :     $Direction = "=>";
3228 :     } else {
3229 :     $Direction = "<=";
3230 :     }
3231 :     my $Reaction = substr($SingletReactions[$i],1);
3232 :     $ReactionData = FIGMODELObject->load($self->figmodel()->config("reaction directory")->[0].$Reaction,"\t");
3233 :     $ReactionData->{"DIRECTIONS"}->[0] = $Direction;
3234 :     $ReactionData->{"REACTIONS"}->[0] = $Reaction;
3235 :     if (!defined($ReactionData->{"DEFINITION"}->[0])) {
3236 :     $ReactionData->{"DEFINITION"}->[0] = "UNKNOWN";
3237 :     }
3238 :     if (!defined($ReactionData->{"THERMODYNAMIC REVERSIBILITY"}->[0])) {
3239 :     $ReactionData->{"THERMODYNAMIC REVERSIBILITY"}->[0] = "UNKNOWN";
3240 :     }
3241 :     if (!defined($ReactionData->{"DELTAG"}->[0])) {
3242 :     $ReactionData->{"DELTAG"}->[0] = "UNKNOWN";
3243 :     }
3244 :     $ReactionDataHash->{$SingletReactions[$i]} = $ReactionData;
3245 :     }
3246 :     print RECONCILIATION $ReactionData->{"REACTIONS"}->[0].";".$ReactionData->{"DEFINITION"}->[0].";".$ReactionData->{"THERMODYNAMIC REVERSIBILITY"}->[0].";".$ReactionData->{"DELTAG"}->[0].";".$ReactionData->{"DIRECTIONS"}->[0].";".$SingleReactionHash->{$SingletReactions[$i]}->{"COUNT"};
3247 :     for (my $j=0; $j < $SolutionCount; $j++) {
3248 :     print RECONCILIATION ";";
3249 :     if (defined($SingleReactionHash->{$SingletReactions[$i]}->{$j})) {
3250 :     $Solutions->[$j]->{$SingletReactions[$i]} = 1;
3251 :     $Solutions->[$j]->{"BASE"} = $j;
3252 :     print RECONCILIATION "|".$j."|";
3253 :     }
3254 :     }
3255 :     print RECONCILIATION "\n";
3256 :     }
3257 :     #Printing the reaction sets with alternatives
3258 :     print RECONCILIATION "Reaction sets with alternatives\n";
3259 :     my @ReactionSets = keys(%{$ReactionSetHash});
3260 :     foreach my $ReactionSet (@ReactionSets) {
3261 :     my $NewSolutions;
3262 :     my $BaseReactions;
3263 :     my $AltList = [$ReactionSet];
3264 :     push(@{$AltList},keys(%{$ReactionSetHash->{$ReactionSet}}));
3265 :     for (my $j=0; $j < @{$AltList}; $j++) {
3266 :     my $CurrentNewSolutions;
3267 :     my $Index;
3268 :     if ($j == 0) {
3269 :     print RECONCILIATION "NEW SET\n";
3270 :     } elsif ($AltList->[$j] ne $ReactionSet) {
3271 :     print RECONCILIATION "ALTERNATIVE SET\n";
3272 :     #For each base solution in which this set is represented, we copy the base solution to the new solution
3273 :     my $NewSolutionCount = 0;
3274 :     for (my $k=0; $k < $SolutionCount; $k++) {
3275 :     if (defined($ReactionSetHash->{$ReactionSet}->{$AltList->[$j]}->{$k})) {
3276 :     if (defined($Solutions)) {
3277 :     $Index->{$k} = @{$Solutions} + $NewSolutionCount;
3278 :     } else {
3279 :     $Index->{$k} = $NewSolutionCount;
3280 :     }
3281 :     if (defined($NewSolutions) && @{$NewSolutions} > 0) {
3282 :     $Index->{$k} += @{$NewSolutions};
3283 :     }
3284 :     $CurrentNewSolutions->[$NewSolutionCount] = {};
3285 :     foreach my $Reaction (keys(%{$Solutions->[$k]})) {
3286 :     $CurrentNewSolutions->[$NewSolutionCount]->{$Reaction} = $Solutions->[$k]->{$Reaction};
3287 :     }
3288 :     $NewSolutionCount++;
3289 :     }
3290 :     }
3291 :     }
3292 :     if ($j == 0 || $AltList->[$j] ne $ReactionSet) {
3293 :     my @SingletReactions = split(/,/,$AltList->[$j]);
3294 :     for (my $i=0; $i < @SingletReactions; $i++) {
3295 :     #Adding base reactions to base solutions and set reactions the new solutions
3296 :     if ($j == 0) {
3297 :     push(@{$BaseReactions},$SingletReactions[$i]);
3298 :     } else {
3299 :     for (my $k=0; $k < @{$CurrentNewSolutions}; $k++) {
3300 :     $CurrentNewSolutions->[$k]->{$SingletReactions[$i]} = 1;
3301 :     }
3302 :     }
3303 :     #Getting reaction data and printing reaction in output file
3304 :     my $ReactionData;
3305 :     if (defined($ReactionDataHash->{$SingletReactions[$i]})) {
3306 :     $ReactionData = $ReactionDataHash->{$SingletReactions[$i]};
3307 :     } else {
3308 :     my $Direction = substr($SingletReactions[$i],0,1);
3309 :     if ($Direction eq "+") {
3310 :     $Direction = "=>";
3311 :     } else {
3312 :     $Direction = "<=";
3313 :     }
3314 :     my $Reaction = substr($SingletReactions[$i],1);
3315 :     $ReactionData = FIGMODELObject->load($self->figmodel()->config("reaction directory")->[0].$Reaction,"\t");
3316 :     $ReactionData->{"DIRECTIONS"}->[0] = $Direction;
3317 :     $ReactionData->{"REACTIONS"}->[0] = $Reaction;
3318 :     if (!defined($ReactionData->{"DEFINITION"}->[0])) {
3319 :     $ReactionData->{"DEFINITION"}->[0] = "UNKNOWN";
3320 :     }
3321 :     if (!defined($ReactionData->{"THERMODYNAMIC REVERSIBILITY"}->[0])) {
3322 :     $ReactionData->{"THERMODYNAMIC REVERSIBILITY"}->[0] = "UNKNOWN";
3323 :     }
3324 :     if (!defined($ReactionData->{"DELTAG"}->[0])) {
3325 :     $ReactionData->{"DELTAG"}->[0] = "UNKNOWN";
3326 :     }
3327 :     $ReactionDataHash->{$SingletReactions[$i]} = $ReactionData;
3328 :     }
3329 :     print RECONCILIATION $ReactionData->{"REACTIONS"}->[0].";".$ReactionData->{"DEFINITION"}->[0].";".$ReactionData->{"THERMODYNAMIC REVERSIBILITY"}->[0].";".$ReactionData->{"DELTAG"}->[0].";".$ReactionData->{"DIRECTIONS"}->[0].";".$ReactionSetHash->{$ReactionSet}->{$AltList->[$j]}->{"COUNT"};
3330 :     for (my $k=0; $k < $SolutionCount; $k++) {
3331 :     print RECONCILIATION ";";
3332 :     if (defined($ReactionSetHash->{$ReactionSet}->{$AltList->[$j]}->{$k})) {
3333 :     if ($j == 0) {
3334 :     print RECONCILIATION "|".$k."|";
3335 :     } else {
3336 :     print RECONCILIATION "|".$Index->{$k}."|";
3337 :     }
3338 :     }
3339 :     }
3340 :     print RECONCILIATION "\n";
3341 :     }
3342 :     #Adding the current new solutions to the new solutions array
3343 :     if (defined($CurrentNewSolutions) && @{$CurrentNewSolutions} > 0) {
3344 :     push(@{$NewSolutions},@{$CurrentNewSolutions});
3345 :     }
3346 :     }
3347 :     }
3348 :     #Adding the base reactions to all existing solutions
3349 :     for (my $j=0; $j < @{$Solutions}; $j++) {
3350 :     if (defined($ReactionSetHash->{$ReactionSet}->{$ReactionSet}->{$Solutions->[$j]->{"BASE"}})) {
3351 :     foreach my $SingleReaction (@{$BaseReactions}) {
3352 :     $Solutions->[$j]->{$SingleReaction} = 1;
3353 :     }
3354 :     }
3355 :     }
3356 :     #Adding the new solutions to the set of existing solutions
3357 :     push(@{$Solutions},@{$NewSolutions});
3358 :     }
3359 :     close(RECONCILIATION);
3360 :     #Now printing a file that defines all of the solutions in a format the testsolutions function understands
3361 :     open (RECONCILIATION, ">$OutputFilenameTwo");
3362 :     print RECONCILIATION "Experiment;Solution index;Solution cost;Solution reactions\n";
3363 :     for (my $i=0; $i < @{$Solutions}; $i++) {
3364 :     delete($Solutions->[$i]->{"BASE"});
3365 :     print RECONCILIATION "SR".$i.";".$i.";10;".join(",",keys(%{$Solutions->[$i]}))."\n";
3366 :     }
3367 :     close(RECONCILIATION);
3368 :    
3369 :     $GrowMatchTable = $self->figmodel()->database()->LockDBTable("GROWMATCH TABLE");
3370 :     $Row = $GrowMatchTable->get_row_by_key($self->genome(),"ORGANISM",1);
3371 :     $Row->{"GF RECON TESTING TIMING"}->[0] = time()."-";
3372 :     $Row->{"GF RECON SOLUTIONS"}->[0] = @{$Solutions};
3373 :     $GrowMatchTable->save();
3374 :     $self->figmodel()->database()->UnlockDBTable("GROWMATCH TABLE");
3375 :    
3376 :     #Scheduling the solution testing
3377 :     if ($GapFill == 1) {
3378 :     system($self->figmodel()->config("scheduler executable")->[0]." \"add:testsolutions?".$self->id().$self->selected_version()."?-1?GFSR:BACK:fast:QSUB\"");
3379 :     } else {
3380 :     system($self->figmodel()->config("scheduler executable")->[0]." \"add:testsolutions?".$self->id().$self->selected_version()."?-1?GGSR:BACK:fast:QSUB\"");
3381 :     }
3382 :     } else {
3383 :     #Reading in the solution testing results
3384 :     my $Data;
3385 :     if ($GapFill == 1) {
3386 :     $Data = $self->figmodel()->database()->load_single_column_file($self->directory().$self->id().$self->selected_version()."-GFSREM.txt","");
3387 :     } else {
3388 :     $Data = $self->figmodel()->database()->load_single_column_file($self->directory().$self->id().$self->selected_version()."-GGSREM.txt","");
3389 :     }
3390 :    
3391 :     #Reading in the preliminate reconciliation report
3392 :     my $OutputData = $self->figmodel()->database()->load_single_column_file($OutputFilename,"");
3393 :     #Replacing the file tags with actual performance data
3394 :     my $Count = 0;
3395 :     for (my $i=0; $i < @{$Data}; $i++) {
3396 :     if ($Data->[$i] =~ m/^SR(\d+);.+;(\d+\/\d+);/) {
3397 :     my $Index = $1;
3398 :     my $Performance = $Index."/".$2;
3399 :     for (my $j=0; $j < @{$OutputData}; $j++) {
3400 :     $OutputData->[$j] =~ s/\|$Index\|/$Performance/g;
3401 :     }
3402 :     }
3403 :     }
3404 :     $self->figmodel()->database()->print_array_to_file($OutputFilename,$OutputData);
3405 :    
3406 :     my $GrowMatchTable = $self->figmodel()->database()->LockDBTable("GROWMATCH TABLE");
3407 :     my $Row = $GrowMatchTable->get_row_by_key($self->genome(),"ORGANISM",1);
3408 :     $Row->{"GF RECON TESTING TIMING"}->[0] .= time();
3409 :     $GrowMatchTable->save();
3410 :     $self->figmodel()->database()->UnlockDBTable("GROWMATCH TABLE");
3411 :     }
3412 :    
3413 :     return 1;
3414 :     }
3415 :    
3416 :     =head3 BuildSpecificBiomassReaction
3417 :     Definition:
3418 :     FIGMODELmodel->BuildSpecificBiomassReaction();
3419 :     Description:
3420 :     =cut
3421 :     sub BuildSpecificBiomassReaction {
3422 :     my ($self) = @_;
3423 :    
3424 :     my $biomassrxn;
3425 :     my $OrganismID = $self->genome();
3426 :     #Checking for a biomass override
3427 :     if (defined($self->config("biomass reaction override")->{$OrganismID})) {
3428 :     my $biomassID = $self->config("biomass reaction override")->{$OrganismID};
3429 :     my $tbl = $self->database()->get_table("BIOMASS",1);
3430 :     $biomassrxn = $tbl->get_row_by_key($biomassID,"DATABASE");
3431 :     print "Overriding biomass template and selecting ".$biomassID." for ".$OrganismID.".\n";
3432 :     } else {#Creating biomass reaction from the template
3433 :     #Getting the genome stats
3434 :     my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
3435 :     my $Class = $genomestats->{CLASS}->[0];
3436 :     my $Name = $genomestats->{NAME}->[0];
3437 :    
3438 :     #Checking for phoenix variants
3439 :     my $PhoenixVariantTable = $self->figmodel()->database()->GetDBTable("Phoenix variants table");
3440 :     my $Phoenix = 0;
3441 :     my @Rows = $PhoenixVariantTable->get_rows_by_key($OrganismID,"GENOME");
3442 :     my $VariantHash;
3443 :     for (my $i=0; $i < @Rows; $i++) {
3444 :     $Phoenix = 1;
3445 :     if (defined($Rows[$i]->{"SUBSYSTEM"}) && defined($Rows[$i]->{"VARIANT"})) {
3446 :     $VariantHash->{$Rows[$i]->{"SUBSYSTEM"}->[0]} = $Rows[$i]->{"VARIANT"}->[0];
3447 :     }
3448 :     }
3449 :    
3450 :     #Collecting genome data
3451 :     my $RoleHash;
3452 :     my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
3453 :     for (my $i=0; $i < $FeatureTable->size(); $i++) {
3454 :     for (my $j=0; $j < @{$FeatureTable->get_row($i)->{"ROLES"}}; $j++) {
3455 :     $RoleHash->{$FeatureTable->get_row($i)->{"ROLES"}->[$j]} = 1;
3456 :     my $Subsystems = $self->figmodel()->subsystems_of_role($FeatureTable->get_row($i)->{"ROLES"}->[$j]);
3457 :     if (defined($Subsystems)) {
3458 :     for (my $k=0; $k < @{$Subsystems}; $k++) {
3459 :     if ($Phoenix == 0) {
3460 :     $VariantHash->{$Subsystems->[$k]} = 1;
3461 :     }
3462 :     }
3463 :     }
3464 :     }
3465 :     }
3466 :    
3467 :     #Scanning through the template item by item and determinine which biomass components should be added
3468 :     my $ComponentTypes;
3469 :     my $EquationData;
3470 :     my $BiomassReactionTemplateTable = $self->figmodel()->database()->GetDBTable("BIOMASS TEMPLATE");
3471 :     for (my $i=0; $i < $BiomassReactionTemplateTable->size(); $i++) {
3472 :     my $Row = $BiomassReactionTemplateTable->get_row($i);
3473 :     if (defined($Row->{"INCLUSION CRITERIA"}->[0]) && $Row->{"INCLUSION CRITERIA"}->[0] eq "UNIVERSAL") {
3474 :     push(@{$EquationData},$Row);
3475 :     $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;
3476 :     } else {
3477 :     my $Criteria = $Row->{"INCLUSION CRITERIA"}->[0];
3478 :     my $End = 0;
3479 :     while ($End == 0) {
3480 :     if ($Criteria =~ m/^(.+)(AND)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(AND)\{([^{^}]+)\}$/ || $Criteria =~ m/^(.+)(OR)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(OR)\{([^{^}]+)\}$/) {
3481 :     print $Criteria."\n";
3482 :     my $Start = "";
3483 :     my $End = "";
3484 :     my $Condition = $1;
3485 :     my $Data = $2;
3486 :     if ($1 ne "AND" && $1 ne "OR") {
3487 :     $Start = $1;
3488 :     $End = $4;
3489 :     $Condition = $2;
3490 :     $Data = $3;
3491 :     }
3492 :     my $Result = "YES";
3493 :     if ($Condition eq "OR") {
3494 :     $Result = "NO";
3495 :     }
3496 :     my @Array = split(/\|/,$Data);
3497 :     for (my $j=0; $j < @Array; $j++) {
3498 :     if ($Array[$j] eq "YES" && $Condition eq "OR") {
3499 :     $Result = "YES";
3500 :     last;
3501 :     } elsif ($Array[$j] eq "NO" && $Condition eq "AND") {
3502 :     $Result = "NO";
3503 :     last;
3504 :     } elsif ($Array[$j] =~ m/^COMPOUND:(.+)/) {
3505 :     my $Match = 0;
3506 :     for (my $k=0; $k < @{$EquationData}; $k++) {
3507 :     if ($EquationData->[$k]->{"ID"}->[0] eq $1) {
3508 :     $Match = 1;
3509 :     last;
3510 :     }
3511 :     }
3512 :     if ($Match == 1 && $Condition eq "OR") {
3513 :     $Result = "YES";
3514 :     last;
3515 :     } elsif ($Match != 1 && $Condition eq "AND") {
3516 :     $Result = "NO";
3517 :     last;
3518 :     }
3519 :     } elsif ($Array[$j] =~ m/^!COMPOUND:(.+)/) {
3520 :     my $Match = 0;
3521 :     for (my $k=0; $k < @{$EquationData}; $k++) {
3522 :     if ($EquationData->[$k]->{"ID"}->[0] eq $1) {
3523 :     $Match = 1;
3524 :     last;
3525 :     }
3526 :     }
3527 :     if ($Match != 1 && $Condition eq "OR") {
3528 :     $Result = "YES";
3529 :     last;
3530 :     } elsif ($Match == 1 && $Condition eq "AND") {
3531 :     $Result = "NO";
3532 :     last;
3533 :     }
3534 :     } elsif ($Array[$j] =~ m/^NAME:(.+)/) {
3535 :     my $Comparison = $1;
3536 :     if ((!defined($Comparison) || !defined($Name) || $Name =~ m/$Comparison/) && $Condition eq "OR") {
3537 :     $Result = "YES";
3538 :     last;
3539 :     } elsif (defined($Comparison) && defined($Name) && $Name !~ m/$Comparison/ && $Condition eq "AND") {
3540 :     $Result = "NO";
3541 :     last;
3542 :     }
3543 :     } elsif ($Array[$j] =~ m/^!NAME:(.+)/) {
3544 :     my $Comparison = $1;
3545 :     if ((!defined($Comparison) || !defined($Name) || $Name !~ m/$Comparison/) && $Condition eq "OR") {
3546 :     $Result = "YES";
3547 :     last;
3548 :     } elsif (defined($Comparison) && defined($Name) && $Name =~ m/$Comparison/ && $Condition eq "AND") {
3549 :     $Result = "NO";
3550 :     last;
3551 :     }
3552 :     } elsif ($Array[$j] =~ m/^SUBSYSTEM:(.+)/) {
3553 :     my @SubsystemArray = split(/`/,$1);
3554 :     if (@SubsystemArray == 1) {
3555 :     if (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} ne -1 && $Condition eq "OR") {
3556 :     $Result = "YES";
3557 :     last;
3558 :     } elsif ((!defined($VariantHash->{$SubsystemArray[0]}) || $VariantHash->{$SubsystemArray[0]} eq -1) && $Condition eq "AND") {
3559 :     $Result = "NO";
3560 :     last;
3561 :     }
3562 :     } else {
3563 :     my $Match = 0;
3564 :     for (my $k=1; $k < @SubsystemArray; $k++) {
3565 :     if (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} eq $SubsystemArray[$k]) {
3566 :     $Match = 1;
3567 :     last;
3568 :     }
3569 :     }
3570 :     if ($Match == 1 && $Condition eq "OR") {
3571 :     $Result = "YES";
3572 :     last;
3573 :     } elsif ($Match != 1 && $Condition eq "AND") {
3574 :     $Result = "NO";
3575 :     last;
3576 :     }
3577 :     }
3578 :     } elsif ($Array[$j] =~ m/^!SUBSYSTEM:(.+)/) {
3579 :     my @SubsystemArray = split(/`/,$1);
3580 :     if (@SubsystemArray == 1) {
3581 :     if ((!defined($VariantHash->{$SubsystemArray[0]}) || $VariantHash->{$SubsystemArray[0]} eq -1) && $Condition eq "OR") {
3582 :     $Result = "YES";
3583 :     last;
3584 :     } elsif (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} ne -1 && $Condition eq "AND") {
3585 :     $Result = "NO";
3586 :     last;
3587 :     }
3588 :     } else {
3589 :     my $Match = 0;
3590 :     for (my $k=1; $k < @SubsystemArray; $k++) {
3591 :     if (defined($VariantHash->{$SubsystemArray[0]}) && $VariantHash->{$SubsystemArray[0]} eq $SubsystemArray[$k]) {
3592 :     $Match = 1;
3593 :     last;
3594 :     }
3595 :     }
3596 :     if ($Match != 1 && $Condition eq "OR") {
3597 :     $Result = "YES";
3598 :     last;
3599 :     } elsif ($Match == 1 && $Condition eq "AND") {
3600 :     $Result = "NO";
3601 :     last;
3602 :     }
3603 :     }
3604 :     } elsif ($Array[$j] =~ m/^ROLE:(.+)/) {
3605 :     if (defined($RoleHash->{$1}) && $Condition eq "OR") {
3606 :     $Result = "YES";
3607 :     last;
3608 :     } elsif (!defined($RoleHash->{$1}) && $Condition eq "AND") {
3609 :     $Result = "NO";
3610 :     last;
3611 :     }
3612 :     } elsif ($Array[$j] =~ m/^!ROLE:(.+)/) {
3613 :     if (!defined($RoleHash->{$1}) && $Condition eq "OR") {
3614 :     $Result = "YES";
3615 :     last;
3616 :     } elsif (defined($RoleHash->{$1}) && $Condition eq "AND") {
3617 :     $Result = "NO";
3618 :     last;
3619 :     }
3620 :     } elsif ($Array[$j] =~ m/^CLASS:(.+)/) {
3621 :     if ($Class eq $1 && $Condition eq "OR") {
3622 :     $Result = "YES";
3623 :     last;
3624 :     } elsif ($Class ne $1 && $Condition eq "AND") {
3625 :     $Result = "NO";
3626 :     last;
3627 :     }
3628 :     } elsif ($Array[$j] =~ m/^!CLASS:(.+)/) {
3629 :     if ($Class ne $1 && $Condition eq "OR") {
3630 :     $Result = "YES";
3631 :     last;
3632 :     } elsif ($Class eq $1 && $Condition eq "AND") {
3633 :     $Result = "NO";
3634 :     last;
3635 :     }
3636 :     }
3637 :     }
3638 :     $Criteria = $Start.$Result.$End;
3639 :     print $Criteria."\n";
3640 :     } else {
3641 :     $End = 1;
3642 :     last;
3643 :     }
3644 :     }
3645 :     if ($Criteria eq "YES") {
3646 :     push(@{$EquationData},$Row);
3647 :     $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;
3648 :     }
3649 :     }
3650 :     }
3651 :    
3652 :     #Building biomass equation
3653 :     my %Reactants;
3654 :     my %Products;
3655 :     foreach my $EquationRow (@{$EquationData}) {
3656 :     #First determine what the coefficient should be
3657 :     my $Coefficient;
3658 :     if (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^[0-9\.]+$/) {
3659 :     $Coefficient = $EquationRow->{"COEFFICIENT"}->[0];
3660 :     } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/cpd\d\d\d\d\d/) {
3661 :     $Coefficient = 0;
3662 :     my @CompoundList = split(/,/,$EquationRow->{"COEFFICIENT"}->[0]);
3663 :     foreach my $Compound (@CompoundList) {
3664 :     if (defined($Reactants{$Compound})) {
3665 :     $Coefficient += $Reactants{$Compound};
3666 :     }
3667 :     }
3668 :     } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^([0-9\.]+)\/(.+)$/) {
3669 :     my @Keys = keys(%{$ComponentTypes->{$2}});
3670 :     my $MW = 1;
3671 :     my $CompoundData = $self->figmodel()->LoadObject($EquationRow->{"ID"}->[0]);
3672 :     if (defined($CompoundData->{"MASS"})) {
3673 :     $MW = $CompoundData->{"MASS"}->[0];
3674 :     }
3675 :     $Coefficient = $1/@Keys/$MW;
3676 :     }
3677 :     if (defined($EquationRow->{"REACTANT"}) && $EquationRow->{"REACTANT"}->[0] eq "YES") {
3678 :     if (defined($Reactants{$EquationRow->{"ID"}->[0]})) {
3679 :     $Reactants{$EquationRow->{"ID"}->[0]} += $Coefficient;
3680 :     } elsif (defined($Products{$EquationRow->{"ID"}->[0]}) && $Products{$EquationRow->{"ID"}->[0]} > $Coefficient) {
3681 :     $Products{$EquationRow->{"ID"}->[0]} -= $Coefficient;
3682 :     } elsif (defined($Products{$EquationRow->{"ID"}->[0]}) && $Products{$EquationRow->{"ID"}->[0]} < $Coefficient) {
3683 :     $Reactants{$EquationRow->{"ID"}->[0]} = $Coefficient - $Products{$EquationRow->{"ID"}->[0]};
3684 :     delete $Products{$EquationRow->{"ID"}->[0]};
3685 :     } else {
3686 :     $Reactants{$EquationRow->{"ID"}->[0]} = $Coefficient;
3687 :     }
3688 :     } else {
3689 :     if (defined($Products{$EquationRow->{"ID"}->[0]})) {
3690 :     $Products{$EquationRow->{"ID"}->[0]} += $Coefficient;
3691 :     } elsif (defined($Reactants{$EquationRow->{"ID"}->[0]}) && $Reactants{$EquationRow->{"ID"}->[0]} > $Coefficient) {
3692 :     $Reactants{$EquationRow->{"ID"}->[0]} -= $Coefficient;
3693 :     } elsif (defined($Reactants{$EquationRow->{"ID"}->[0]}) && $Reactants{$EquationRow->{"ID"}->[0]} < $Coefficient) {
3694 :     $Products{$EquationRow->{"ID"}->[0]} = $Coefficient - $Reactants{$EquationRow->{"ID"}->[0]};
3695 :     delete $Reactants{$EquationRow->{"ID"}->[0]};
3696 :     } else {
3697 :     $Products{$EquationRow->{"ID"}->[0]} = $Coefficient;
3698 :     }
3699 :     }
3700 :     }
3701 :     my $Equation = "";
3702 :     my @ReactantList = sort(keys(%Reactants));
3703 :     for (my $i=0; $i < @ReactantList; $i++) {
3704 :     if (length($Equation) > 0) {
3705 :     $Equation .= " + ";
3706 :     }
3707 :     $Equation .= $self->figmodel()->format_coefficient($Reactants{$ReactantList[$i]})." ".$ReactantList[$i];
3708 :     }
3709 :     $Equation .= " => ";
3710 :     my $First = 1;
3711 :     @ReactantList = sort(keys(%Products));
3712 :     for (my $i=0; $i < @ReactantList; $i++) {
3713 :     if ($First == 0) {
3714 :     $Equation .= " + ";
3715 :     }
3716 :     $First = 0;
3717 :     $Equation .= $self->figmodel()->format_coefficient($Products{$ReactantList[$i]})." ".$ReactantList[$i];
3718 :     }
3719 :    
3720 :     #Adding the biomass equation to the biomass table
3721 :     my $NewRow = $self->figmodel()->add_biomass_reaction($Equation,$self->id(),"Template:".$self->genome());
3722 :     $biomassrxn = $NewRow;
3723 :     }
3724 :     return $biomassrxn;
3725 :     }
3726 :    
3727 :     =head3 PrintSBMLFile
3728 :     Definition:
3729 :     FIGMODELmodel->PrintSBMLFile();
3730 :     Description:
3731 :     Printing file with model data in SBML format
3732 :     =cut
3733 :     sub PrintSBMLFile {
3734 :     my($self) = @_;
3735 :    
3736 :     #Opening the SBML file for printing
3737 :     my $Filename = $self->config("SBML files")->[0].$self->id().".xml";
3738 :     if ($self->owner() ne "master") {
3739 :     if (!-d $self->config("SBML files")->[0].$self->owner()."/") {
3740 :     system("mkdir ".$self->config("SBML files")->[0].$self->owner()."/");
3741 :     }
3742 :     $Filename = $self->config("SBML files")->[0].$self->owner()."/".$self->id().".xml";
3743 :     }
3744 :     if (!open (SBMLOUTPUT, ">$Filename")) {
3745 :     return;
3746 :     }
3747 :    
3748 :     #Loading and parsing the model data
3749 :     my $ModelTable = $self->reaction_table();
3750 :     if (!defined($ModelTable) || !defined($ModelTable->{"array"})) {
3751 :     print "Failed to load ".$self->id()."\n";
3752 :     return;
3753 :     }
3754 :    
3755 :     #Adding intracellular metabolites that also need exchange fluxes to the exchange hash
3756 :     my $ExchangeHash = {"cpd11416" => "c"};
3757 :    
3758 :     my %CompartmentsPresent;
3759 :     $CompartmentsPresent{"c"} = 1;
3760 :     my %CompoundList;
3761 :     my @ReactionList;
3762 :     my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");
3763 :     for (my $i=0; $i < $ModelTable->size(); $i++) {
3764 :     my $Reaction = $ModelTable->get_row($i)->{"LOAD"}->[0];
3765 :     if (defined($ReactionTable->get_row_by_key($Reaction,"DATABASE")) && defined($ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0])) {
3766 :     push(@ReactionList,$Reaction);
3767 :     $_ = $ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0];
3768 :     my @MatchArray = /(cpd\d\d\d\d\d)/g;
3769 :     for (my $j=0; $j < @MatchArray; $j++) {
3770 :     $CompoundList{$MatchArray[$j]}->{"c"} = 1;
3771 :     }
3772 :     $_ = $ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0];
3773 :     @MatchArray = /(cpd\d\d\d\d\d\[\D\])/g;
3774 :     for (my $j=0; $j < @MatchArray; $j++) {
3775 :     if ($MatchArray[$j] =~ m/(cpd\d\d\d\d\d)\[(\D)\]/) {
3776 :     $CompartmentsPresent{lc($2)} = 1;
3777 :     $CompoundList{$1}->{lc($2)} = 1;
3778 :     }
3779 :     }
3780 :     }
3781 :     }
3782 :    
3783 :     #Printing header to SBML file
3784 :     my $ModelName = $self->id();
3785 :     $ModelName =~ s/\./_/;
3786 :     print SBMLOUTPUT '<?xml version="1.0" encoding="UTF-8"?>'."\n";
3787 :     print SBMLOUTPUT '<sbml xmlns="http://www.sbml.org/sbml/level2" level="2" version="1" xmlns:html="http://www.w3.org/1999/xhtml">' . "\n";
3788 :     if (defined($self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Organism name"}->[0])) {
3789 :     print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Organism name"}->[0].' SEED model">'."\n";
3790 :     } else {
3791 :     print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->id().' SEED model">'."\n";
3792 :     }
3793 :    
3794 :     #Printing the unit data
3795 :     print SBMLOUTPUT "<listOfUnitDefinitions>\n";
3796 :     print SBMLOUTPUT "\t<unitDefinition id=\"mmol_per_gDW_per_hr\">\n";
3797 :     print SBMLOUTPUT "\t\t<listOfUnits>\n";
3798 :     print SBMLOUTPUT "\t\t\t<unit kind=\"mole\" scale=\"-3\"/>\n";
3799 :     print SBMLOUTPUT "\t\t\t<unit kind=\"gram\" exponent=\"-1\"/>\n";
3800 :     print SBMLOUTPUT "\t\t\t<unit kind=\"second\" multiplier=\".00027777\" exponent=\"-1\"/>\n";
3801 :     print SBMLOUTPUT "\t\t</listOfUnits>\n";
3802 :     print SBMLOUTPUT "\t</unitDefinition>\n";
3803 :     print SBMLOUTPUT "</listOfUnitDefinitions>\n";
3804 :    
3805 :     #Printing compartments for SBML file
3806 :     print SBMLOUTPUT '<listOfCompartments>'."\n";
3807 :     foreach my $Compartment (keys(%CompartmentsPresent)) {
3808 :     if (defined($self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0])) {
3809 :     my @OutsideList = split(/\//,$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Outside"}->[0]);
3810 :     my $Printed = 0;
3811 :     foreach my $Outside (@OutsideList) {
3812 :     if (defined($CompartmentsPresent{$Outside}) && defined($self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Outside}->[0]->{"Name"}->[0])) {
3813 :     print SBMLOUTPUT '<compartment id="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'" outside="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Outside}->[0]->{"Name"}->[0].'"/>'."\n";
3814 :     $Printed = 1;
3815 :     last;
3816 :     }
3817 :     }
3818 :     if ($Printed eq 0) {
3819 :     print SBMLOUTPUT '<compartment id="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'"/>'."\n";
3820 :     }
3821 :     }
3822 :     }
3823 :     print SBMLOUTPUT '</listOfCompartments>'."\n";
3824 :    
3825 :     #Printing the list of metabolites involved in the model
3826 :     print SBMLOUTPUT '<listOfSpecies>'."\n";
3827 :     foreach my $Compound (keys(%CompoundList)) {
3828 :     my $Name = $Compound;
3829 :     my $Formula = "";
3830 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0])) {
3831 :     $Formula = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0];
3832 :     }
3833 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0])) {
3834 :     $Name = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0];
3835 :     $Name =~ s/\s/_/;
3836 :     $Name .= "_".$Formula;
3837 :     }
3838 :     $Name =~ s/[<>;&\*]//;
3839 :     my $Charge = 0;
3840 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0])) {
3841 :     $Charge = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0];
3842 :     }
3843 :     foreach my $Compartment (keys(%{$CompoundList{$Compound}})) {
3844 :     if ($Compartment eq "e") {
3845 :     $ExchangeHash->{$Compound} = "e";
3846 :     }
3847 :     print SBMLOUTPUT '<species id="'.$Compound.'_'.$Compartment.'" name="'.$Name.'" compartment="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'" charge="'.$Charge.'" boundaryCondition="false"/>'."\n";
3848 :     }
3849 :     }
3850 :     #Printing the boundary species
3851 :     foreach my $Compound (keys(%{$ExchangeHash})) {
3852 :     my $Name = $Compound;
3853 :     my $Formula = "";
3854 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0])) {
3855 :     $Formula = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0];
3856 :     }
3857 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0])) {
3858 :     $Name = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0];
3859 :     $Name =~ s/\s/_/;
3860 :     $Name .= "_".$Formula;
3861 :     }
3862 :     $Name =~ s/[<>;&\*]//;
3863 :     my $Charge = 0;
3864 :     if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0])) {
3865 :     $Charge = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0];
3866 :     }
3867 :     print SBMLOUTPUT '<species id="'.$Compound.'_b" name="'.$Name.'" compartment="Extracellular" charge="'.$Charge.'" boundaryCondition="true"/>'."\n";
3868 :     }
3869 :     print SBMLOUTPUT '</listOfSpecies>'."\n";
3870 :    
3871 :     #Printing the list of reactions involved in the model
3872 :     my $ObjectiveCoef;
3873 :     print SBMLOUTPUT '<listOfReactions>'."\n";
3874 :     foreach my $Reaction (@ReactionList) {
3875 :     $ObjectiveCoef = "0.0";
3876 :     if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"EQUATION"}->[0])) {
3877 :     if ($Reaction =~ m/^bio/) {
3878 :     $ObjectiveCoef = "1.0";
3879 :     }
3880 :     my $LowerBound = -10000;
3881 :     my $UpperBound = 10000;
3882 :     my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateData($Reaction);
3883 :     my $Name = $Reaction;
3884 :     if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"NAME"}->[0])) {
3885 :     $Name = $self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"NAME"}->[0];
3886 :     $Name =~ s/[<>;&]//g;
3887 :     }
3888 :     my $Reversibility = "true";
3889 :     if (defined($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0])) {
3890 :     if ($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0] ne "<=>") {
3891 :     $LowerBound = 0;
3892 :     $Reversibility = "false";
3893 :     }
3894 :     if ($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0] eq "<=") {
3895 :     my $Temp = $Products;
3896 :     $Products = $Reactants;
3897 :     $Reactants = $Temp;
3898 :     }
3899 :     }
3900 :     print SBMLOUTPUT '<reaction id="'.$Reaction.'" name="'.$Name.'" reversible="'.$Reversibility.'">'."\n";
3901 :     print SBMLOUTPUT "<notes>\n";
3902 :     my $ECData = "";
3903 :     if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"ENZYME"}->[0])) {
3904 :     $ECData = $self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"ENZYME"}->[0];
3905 :     }
3906 :     my $SubsystemData = "";
3907 :     if (defined($ModelTable->{$Reaction}->[0]->{"SUBSYSTEM"}->[0])) {
3908 :     $SubsystemData = $ModelTable->{$Reaction}->[0]->{"SUBSYSTEM"}->[0];
3909 :     }
3910 :     my $GeneAssociation = "";
3911 :     my $ProteinAssociation = "";
3912 :     if (defined($ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0])) {
3913 :     if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} == 1 && $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0] !~ m/\+/) {
3914 :     $GeneAssociation = $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0];
3915 :     } else {
3916 :     if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} > 1) {
3917 :     $GeneAssociation = "( ";
3918 :     }
3919 :     for (my $i=0; $i < @{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}}; $i++) {
3920 :     if ($i > 0) {
3921 :     $GeneAssociation .= " ) or ( ";
3922 :     }
3923 :     $GeneAssociation .= $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[$i];
3924 :     }
3925 :     if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} > 1) {
3926 :     $GeneAssociation .= " )";
3927 :     }
3928 :     }
3929 :     $GeneAssociation =~ s/\+/ and /g;
3930 :     if ($GeneAssociation =~ m/\sor\s/ || $GeneAssociation =~ m/\sand\s/) {
3931 :     $GeneAssociation = "( ".$GeneAssociation." )";
3932 :     }
3933 :     $ProteinAssociation = $GeneAssociation;
3934 :     if (defined($self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Genome ID"}->[0])) {
3935 :     $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]);
3936 :     }
3937 :     }
3938 :     print SBMLOUTPUT "<html:p>GENE_ASSOCIATION:".$GeneAssociation."</html:p>\n";
3939 :     print SBMLOUTPUT "<html:p>PROTEIN_ASSOCIATION:".$ProteinAssociation."</html:p>\n";
3940 :     print SBMLOUTPUT "<html:p>SUBSYSTEM: ".$SubsystemData."</html:p>\n";
3941 :     print SBMLOUTPUT "<html:p>PROTEIN_CLASS: ".$ECData."</html:p>\n";
3942 :     print SBMLOUTPUT "</notes>\n";
3943 :     print SBMLOUTPUT "<listOfReactants>\n";
3944 :     foreach my $Reactant (@{$Reactants}) {
3945 :     print SBMLOUTPUT '<speciesReference species="'.$Reactant->{"DATABASE"}->[0]."_".$Reactant->{"COMPARTMENT"}->[0].'" stoichiometry="'.$Reactant->{"COEFFICIENT"}->[0].'"/>'."\n";
3946 :     }
3947 :     print SBMLOUTPUT "</listOfReactants>\n";
3948 :     print SBMLOUTPUT "<listOfProducts>\n";
3949 :     foreach my $Product (@{$Products}) {
3950 :     print SBMLOUTPUT '<speciesReference species="'.$Product->{"DATABASE"}->[0]."_".$Product->{"COMPARTMENT"}->[0].'" stoichiometry="'.$Product->{"COEFFICIENT"}->[0].'"/>'."\n";
3951 :     }
3952 :     print SBMLOUTPUT "</listOfProducts>\n";
3953 :     print SBMLOUTPUT "<kineticLaw>\n";
3954 :     print SBMLOUTPUT "\t<math xmlns=\"http://www.w3.org/1998/Math/MathML\">\n";
3955 :     print SBMLOUTPUT "\t\t\t<ci> FLUX_VALUE </ci>\n";
3956 :     print SBMLOUTPUT "\t</math>\n";
3957 :     print SBMLOUTPUT "\t<listOfParameters>\n";
3958 :     print SBMLOUTPUT "\t\t<parameter id=\"LOWER_BOUND\" value=\"".$LowerBound."\" units=\"mmol_per_gDW_per_hr\"/>\n";
3959 :     print SBMLOUTPUT "\t\t<parameter id=\"UPPER_BOUND\" value=\"".$UpperBound."\" units=\"mmol_per_gDW_per_hr\"/>\n";
3960 :     print SBMLOUTPUT "\t\t<parameter id=\"OBJECTIVE_COEFFICIENT\" value=\"".$ObjectiveCoef."\"/>\n";
3961 :     print SBMLOUTPUT "\t\t<parameter id=\"FLUX_VALUE\" value=\"0.0\" units=\"mmol_per_gDW_per_hr\"/>\n";
3962 :     print SBMLOUTPUT "\t</listOfParameters>\n";
3963 :     print SBMLOUTPUT "</kineticLaw>\n";
3964 :     print SBMLOUTPUT '</reaction>'."\n";
3965 :     }
3966 :     }
3967 :    
3968 :     my @ExchangeList = keys(%{$ExchangeHash});
3969 :     foreach my $ExCompound (@ExchangeList) {
3970 :     my $ExCompoundName = $ExCompound;
3971 :     my $Row = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->get_row_by_key($ExCompound,"DATABASE");
3972 :     if (defined($Row) && defined($Row->{"NAME"}->[0])) {
3973 :     $ExCompoundName = $Row->{"NAME"}->[0];
3974 :     $ExCompoundName =~ s/[<>;&]//g;
3975 :     }
3976 :     $ObjectiveCoef = "0.0";
3977 :     print SBMLOUTPUT '<reaction id="EX_'.$ExCompound.'_'.$ExchangeHash->{$ExCompound}.'" name="EX_'.$ExCompoundName.'_'.$ExchangeHash->{$ExCompound}.'" reversible="true">'."\n";
3978 :     print SBMLOUTPUT "\t".'<notes>'."\n";
3979 :     print SBMLOUTPUT "\t\t".'<html:p>GENE_ASSOCIATION: </html:p>'."\n";
3980 :     print SBMLOUTPUT "\t\t".'<html:p>PROTEIN_ASSOCIATION: </html:p>'."\n";
3981 :     print SBMLOUTPUT "\t\t".'<html:p>SUBSYSTEM: S_</html:p>'."\n";
3982 :     print SBMLOUTPUT "\t\t".'<html:p>PROTEIN_CLASS: </html:p>'."\n";
3983 :     print SBMLOUTPUT "\t".'</notes>'."\n";
3984 :     print SBMLOUTPUT "\t".'<listOfReactants>'."\n";
3985 :     print SBMLOUTPUT "\t\t".'<speciesReference species="'.$ExCompound.'_'.$ExchangeHash->{$ExCompound}.'" stoichiometry="1.000000"/>'."\n";
3986 :     print SBMLOUTPUT "\t".'</listOfReactants>'."\n";
3987 :     print SBMLOUTPUT "\t".'<listOfProducts>'."\n";
3988 :     print SBMLOUTPUT "\t\t".'<speciesReference species="'.$ExCompound.'_b" stoichiometry="1.000000"/>'."\n";
3989 :     print SBMLOUTPUT "\t".'</listOfProducts>'."\n";
3990 :     print SBMLOUTPUT "\t".'<kineticLaw>'."\n";
3991 :     print SBMLOUTPUT "\t\t".'<math xmlns="http://www.w3.org/1998/Math/MathML">'."\n";
3992 :     print SBMLOUTPUT "\t\t\t\t".'<ci> FLUX_VALUE </ci>'."\n";
3993 :     print SBMLOUTPUT "\t\t".'</math>'."\n";
3994 :     print SBMLOUTPUT "\t\t".'<listOfParameters>'."\n";
3995 :     print SBMLOUTPUT "\t\t\t".'<parameter id="LOWER_BOUND" value="-10000.000000" units="mmol_per_gDW_per_hr"/>'."\n";
3996 :     print SBMLOUTPUT "\t\t\t".'<parameter id="UPPER_BOUND" value="10000.000000" units="mmol_per_gDW_per_hr"/>'."\n";
3997 :     print SBMLOUTPUT "\t\t\t".'<parameter id="OBJECTIVE_COEFFICIENT" value="'.$ObjectiveCoef.'"/>'."\n";
3998 :     print SBMLOUTPUT "\t\t\t".'<parameter id="FLUX_VALUE" value="0.000000" units="mmol_per_gDW_per_hr"/>'."\n";
3999 :     print SBMLOUTPUT "\t\t".'</listOfParameters>'."\n";
4000 :     print SBMLOUTPUT "\t".'</kineticLaw>'."\n";
4001 :     print SBMLOUTPUT '</reaction>'."\n";
4002 :     }
4003 :    
4004 :     #Closing out the file
4005 :     print SBMLOUTPUT '</listOfReactions>'."\n";
4006 :     print SBMLOUTPUT '</model>'."\n";
4007 :     print SBMLOUTPUT "</sbml>\n";
4008 :     close(SBMLOUTPUT);
4009 :     }
4010 :    
4011 :     =head3 PrintModelSimpleReactionTable
4012 :     Definition:
4013 :     success()/fail() FIGMODELmodel->PrintModelSimpleReactionTable();
4014 :     Description:
4015 :     Prints the table of model data
4016 :     =cut
4017 :     sub PrintModelSimpleReactionTable {
4018 :     my ($self) = @_;
4019 :    
4020 :     my $rxntbl = $self->reaction_table();
4021 :     my $tbl = $self->create_table_prototype("ModelSimpleReactionTable");
4022 :     for (my $i=0; $i < $rxntbl->size(); $i++) {
4023 :     my $row = $rxntbl->get_row($i);
4024 :     $row->{DATABASE} = $row->{LOAD};
4025 :     $tbl->add_row($row);
4026 :     }
4027 :     $tbl->save();
4028 :     system("cp ".$tbl->filename()." ".$self->figmodel()->config("Model table download directory")->[0].$self->figmodel()->config("ModelSimpleReactionTable")->{filename_prefix}->[0]."-".$self->id().".txt");
4029 :     }
4030 :    
4031 :     =head3 PrintModelLPFile
4032 :     Definition:
4033 :     success()/fail() FIGMODELmodel->PrintModelLPFile();
4034 :     Description:
4035 :     Prints the lp file needed to run the model using the mpifba program
4036 :     =cut
4037 :     sub PrintModelLPFile {
4038 :     my ($self) = @_;
4039 :     #Printing lp and key file for model
4040 :     my $UniqueFilename = $self->figmodel()->filename();
4041 :     #Printing the standard FBA file
4042 :     system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),"NoBounds",["ProdFullFBALP"],undef,$self->id().$self->selected_version()."-LPPrint.log",undef,$self->selected_version()));
4043 :     system("cp ".$self->config("MFAToolkit output directory")->[0].$UniqueFilename."/CurrentProblem.lp ".$self->directory()."FBA-".$self->id().$self->selected_version().".lp");
4044 :     my $KeyTable = FIGMODELTable::load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/VariableKey.txt",";","|",0,undef);
4045 :     if (!defined($KeyTable)) {
4046 :     print STDERR "FIGMODEL:RunAllStudiesWithDataFast: ".$self->id()." LP file could not be printed.\n";
4047 :     return 0;
4048 :     }
4049 :     $KeyTable->headings(["Variable type","Variable ID"]);
4050 :     $KeyTable->save($self->directory()."FBA-".$self->id().$self->selected_version().".key");
4051 :     unlink($self->config("database message file directory")->[0].$self->id().$self->selected_version()."-LPPrint.log");
4052 :     $self->figmodel()->clearing_output($UniqueFilename,"FBA-".$self->id().$self->selected_version().".lp");
4053 :     }
4054 :    
4055 :     =head3 patch_model
4056 :     Definition:
4057 :     FIGMODELTable:patch results = FIGMODELmodel->patch_model(FIGMODELTable:patch table);
4058 :     Description:
4059 :     =cut
4060 :     sub patch_model {
4061 :     my ($self,$tbl) = @_;
4062 :    
4063 :     #Instantiating table
4064 :     my $results = FIGMODELTable->new(["Reactions","New genes","Old genes","Genes","Roles","Status"],$self->directory()."PatchResults-".$self->id().$self->selected_version().".tbl",["Reaction"],"\t",";",undef);
4065 :     #Getting genome annotations
4066 :     my $features = $self->figmodel()->database()->get_genome_feature_table($self->genome());
4067 :     #Gettubg reaction table
4068 :     my $reactions = $self->reaction_table();
4069 :     #Checking for patched roles
4070 :     for (my $i=0; $i < $tbl->size(); $i++) {
4071 :     my $row = $tbl->get_row($i);
4072 :     my @genes = $features->get_rows_by_key($row->{ROLE}->[0],"ROLES");
4073 :     if (@genes > 0) {
4074 :     for (my $j=0; $j < @{$row->{REACTIONS}};$j++) {
4075 :     my $resultrxn = $results->get_row_by_key($row->{REACTIONS}->[$j],"Reactions");
4076 :     if (!defined($resultrxn)) {
4077 :     $resultrxn = $results->add_row({"Reactions"=>[$row->{REACTIONS}->[$j]],"Roles"=>[$row->{ROLE}->[0]]});
4078 :     }
4079 :     my $rxnrow = $reactions->get_row_by_key($row->{REACTIONS}->[$j],"LOAD");
4080 :     if (defined($rxnrow) && !defined($resultrxn->{"Old genes"})) {
4081 :     $resultrxn->{"Old genes"} = $rxnrow->{"ASSOCIATED PEG"};
4082 :     if ($resultrxn->{"Old genes"}->[0] !~ m/GAP|BOF|UNIVERSAL|SPONTANEOUS/) {
4083 :     push(@{$resultrxn->{"Genes"}},@{$resultrxn->{"Old genes"}});
4084 :     }
4085 :     }
4086 :     delete $resultrxn->{"Current gene set"};
4087 :     if (defined($resultrxn->{"Genes"})) {
4088 :     push(@{$resultrxn->{"Current gene set"}},@{$resultrxn->{"Genes"}});
4089 :     }
4090 :     for (my $k=0; $k < @genes; $k++) {
4091 :     if ($genes[$k]->{ID}->[0] =~ m/(peg\.\d+)/) {
4092 :     my $gene = $1;
4093 :     my $addgene = 1;
4094 :     if (defined($resultrxn->{"Old genes"})) {
4095 :     for (my $m=0; $m < @{$resultrxn->{"Old genes"}}; $m++) {
4096 :     if ($resultrxn->{"Old genes"}->[$m] =~ m/$gene/) {
4097 :     $addgene = 0;
4098 :     }
4099 :     }
4100 :     }
4101 :     if ($addgene == 1) {
4102 :     push(@{$resultrxn->{"New genes"}},$gene);
4103 :     if ($row->{COMPLEX}->[0] ne "0" && defined($resultrxn->{"Current gene set"})) {
4104 :     my $added = 0;
4105 :     for (my $m=0; $m < @{$resultrxn->{"Current gene set"}}; $m++) {
4106 :     if ($row->{COMPLEX}->[0] eq "1") {
4107 :     $resultrxn->{"Current gene set"}->[$m] = $resultrxn->{"Current gene set"}->[$m]."+".$gene;
4108 :     $added = 1;
4109 :     } else {
4110 :     my @geneset = split(/\+/,$resultrxn->{"Current gene set"}->[$m]);
4111 :     for (my $n=0; $n < @geneset;$n++) {
4112 :     if ($self->figmodel()->colocalized_genes($geneset[$n],$gene,$self->genome()) == 1) {
4113 :     $resultrxn->{"Current gene set"}->[$m] = $resultrxn->{"Current gene set"}->[$m]."+".$gene;
4114 :     $added = 1;
4115 :     last;
4116 :     }
4117 :     }
4118 :     }
4119 :     }
4120 :     if ($added == 0) {
4121 :     push(@{$resultrxn->{"Current gene set"}},$gene);
4122 :     }
4123 :     } else {
4124 :     push(@{$resultrxn->{"Current gene set"}},$gene);
4125 :     }
4126 :     }
4127 :     }
4128 :     }
4129 :     delete $resultrxn->{"Genes"};
4130 :     push(@{$resultrxn->{"Genes"}},@{$resultrxn->{"Current gene set"}});
4131 :     }
4132 :     }
4133 :     }
4134 :    
4135 :     #Ensuring that the old model is preserved
4136 :     $self->ArchiveModel();
4137 :     #Modifing the reaction list
4138 :     for (my $i=0; $i < $results->size();$i++) {
4139 :     my $row = $results->get_row($i);
4140 :     my $rxnrow = $reactions->get_row_by_key($row->{"Reactions"}->[0],"LOAD");
4141 :     if (defined($rxnrow)) {
4142 :     $rxnrow->{"ASSOCIATED PEG"} = $row->{"Genes"};
4143 :     } else {
4144 :     $reactions->add_row({LOAD=>[$row->{"Reactions"}->[0]],DIRECTIONALITY=>[$self->figmodel()->reversibility_of_reaction($row->{"Reactions"}->[0])],COMPARTMENT=>["c"],"ASSOCIATED PEG"=>$row->{"Genes"},SUBSYSTEM=>["NONE"],CONFIDENCE=>[2],REFERENCE=>["NONE"],NOTES=>["PATCH"]});
4145 :     }
4146 :     }
4147 :     $reactions->save();
4148 :     $results->save();
4149 :     $self->update_model_stats();
4150 :     $self->PrintModelLPFile();
4151 :     $self->run_default_model_predictions();
4152 :     #Returning results
4153 :     return $results;
4154 :     }
4155 :    
4156 :     =head3 translate_genes
4157 :     Definition:
4158 :     FIGMODELmodel->translate_genes();
4159 :     Description:
4160 :     =cut
4161 :     sub translate_genes {
4162 :     my ($self) = @_;
4163 :    
4164 :     #Loading gene translations
4165 :     if (!defined($self->{_gene_aliases})) {
4166 :     #Loading gene aliases from feature table
4167 :     my $tbl = $self->figmodel()->GetGenomeFeatureTable($self->genome());
4168 :     if (defined($tbl)) {
4169 :     for (my $i=0; $i < $tbl->size(); $i++) {
4170 :     my $row = $tbl->get_row($i);
4171 :     if ($row->{ID}->[0] =~ m/(peg\.\d+)/) {
4172 :     my $geneID = $1;
4173 :     for (my $j=0; $j < @{$row->{ALIASES}}; $j++) {
4174 :     $self->{_gene_aliases}->{$row->{ALIASES}->[$j]} = $geneID;
4175 :     }
4176 :     }
4177 :     }
4178 :     }
4179 :     #Loading additional gene aliases from the database
4180 :     if (-e $self->figmodel()->config("Translation directory")->[0]."AdditionalAliases/".$self->genome().".txt") {
4181 :     my $AdditionalAliases = $self->figmodel()->database()->load_multiple_column_file($self->figmodel()->config("Translation directory")->[0]."AdditionalAliases/".$self->genome().".txt","\t");
4182 :     for (my $i=0; $i < @{$AdditionalAliases}; $i++) {
4183 :     $self->{_gene_aliases}->{$AdditionalAliases->[$i]->[1]} = $AdditionalAliases->[$i]->[0];
4184 :     }
4185 :     }
4186 :     }
4187 :    
4188 :     #Cycling through reactions and translating genes
4189 :     for (my $i=0; $i < $self->reaction_table()->size(); $i++) {
4190 :     my $row = $self->reaction_table()->get_row($i);
4191 :     if (defined($row->{"ASSOCIATED PEG"})) {
4192 :     for (my $j=0; $j < @{$row->{"ASSOCIATED PEG"}}; $j++) {
4193 :     my $Original = $row->{"ASSOCIATED PEG"}->[$j];
4194 :     $Original =~ s/\sand\s/:/g;
4195 :     $Original =~ s/\sor\s/;/g;
4196 :     my @GeneNames = split(/[,\+\s\(\):;]/,$Original);
4197 :     foreach my $Gene (@GeneNames) {
4198 :     if (length($Gene) > 0 && defined($self->{_gene_aliases}->{$Gene})) {
4199 :     my $Replace = $self->{_gene_aliases}->{$Gene};
4200 :     $Original =~ s/([^\w])$Gene([^\w])/$1$Replace$2/g;
4201 :     $Original =~ s/^$Gene([^\w])/$Replace$1/g;
4202 :     $Original =~ s/([^\w])$Gene$/$1$Replace/g;
4203 :     $Original =~ s/^$Gene$/$Replace/g;
4204 :     }
4205 :     }
4206 :     $Original =~ s/:/ and /g;
4207 :     $Original =~ s/;/ or /g;
4208 :     $row->{"ASSOCIATED PEG"}->[$j] = $Original;
4209 :     }
4210 :     }
4211 :     }
4212 :    
4213 :     #Archiving model and saving reaction table
4214 :     $self->ArchiveModel();
4215 :     $self->reaction_table()->save();
4216 :     }
4217 :    
4218 :     =head3 feature_web_data
4219 :     Definition:
4220 :     string:web output for feature/model connection = FIGMODELmodel->feature_web_data(FIGMODELfeature:feature);
4221 :     Description:
4222 :     =cut
4223 :     sub feature_web_data {
4224 :     my ($self,$feature) = @_;
4225 :    
4226 :     #First checking if the feature is in the model
4227 :     my $featureGenome = $feature->{GENOME}->[0];
4228 :     if ($feature->{GENOME}->[0] =~ m/\>(\d+\.\d+)\</) {
4229 :     $featureGenome = $1;
4230 :     }
4231 :     if ($featureGenome ne $self->genome()) {
4232 :     return "Not in model";
4233 :     }
4234 :     my $data = $self->get_feature_data($feature->{ID}->[0]);
4235 :     if (!defined($data)) {
4236 :     return "Not in model";
4237 :     }
4238 :    
4239 :     my $output;
4240 :     #Printing predictions
4241 :     if (defined($data->{$self->id()."PREDICTIONS"})) {
4242 :     #Parsing essentiality data
4243 :     my $essentialityData;
4244 :     if (defined($feature->{ESSENTIALITY})) {
4245 :     for (my $i=0; $i < @{$feature->{ESSENTIALITY}};$i++) {
4246 :     my @temp = split(/:/,$feature->{ESSENTIALITY}->[$i]);
4247 :     $essentialityData->{$temp[0]} = $temp[1];
4248 :     }
4249 :     }
4250 :     my $predictionHash;
4251 :     for (my $i=0; $i < @{$data->{$self->id()."PREDICTIONS"}};$i++) {
4252 :     my @temp = split(/:/,$data->{$self->id()."PREDICTIONS"}->[$i]);
4253 :     if (defined($essentialityData->{$temp[0]})) {
4254 :     if ($temp[1] eq "essential" && $essentialityData->{$temp[0]} eq "essential") {
4255 :     push(@{$predictionHash->{"Correct negative"}},$temp[0]);
4256 :     } elsif ($temp[1] eq "nonessential" && $essentialityData->{$temp[0]} eq "essential") {
4257 :     push(@{$predictionHash->{"False positive"}},$temp[0]);
4258 :     } elsif ($temp[1] eq "essential" && $essentialityData->{$temp[0]} eq "nonessential") {
4259 :     push(@{$predictionHash->{"False negative"}},$temp[0]);
4260 :     } elsif ($temp[1] eq "nonessential" && $essentialityData->{$temp[0]} eq "nonessential") {
4261 :     push(@{$predictionHash->{"Correct positive"}},$temp[0]);
4262 :     }
4263 :     } else {
4264 :     push(@{$predictionHash->{$temp[1]}},$temp[0]);
4265 :     }
4266 :     }
4267 :     foreach my $key (keys(%{$predictionHash})) {
4268 :     my $string = $key;
4269 :     $string = ucfirst($string);
4270 :     push(@{$output},'<span title="'.join(", ",@{$predictionHash->{$key}}).'">'.$string.'</span>');
4271 :     }
4272 :     }
4273 :    
4274 :     #Printing reactions
4275 :     if (defined($data->{$self->id()."REACTIONS"})) {
4276 :     for (my $i=0; $i < @{$data->{$self->id()."REACTIONS"}};$i++) {
4277 :     my $rxnData = $self->get_reaction_data($data->{$self->id()."REACTIONS"}->[$i]);
4278 :     my $reactionString = $self->figmodel()->web()->create_reaction_link($data->{$self->id()."REACTIONS"}->[$i],join(" or ",@{$rxnData->{"ASSOCIATED PEG"}}),$self->id());
4279 :     if (defined($rxnData->{PREDICTIONS})) {
4280 :     my $predictionHash;
4281 :     for (my $i=0; $i < @{$rxnData->{PREDICTIONS}};$i++) {
4282 :     my @temp = split(/:/,$rxnData->{PREDICTIONS}->[$i]);
4283 :     push(@{$predictionHash->{$temp[1]}},$temp[0]);
4284 :     }
4285 :     $reactionString .= "(";
4286 :     foreach my $key (keys(%{$predictionHash})) {
4287 :     if ($key eq "Essential =>") {
4288 :     $reactionString .= '<span title="Essential in '.join(",",@{$predictionHash->{$key}}).'">E=></span>,';
4289 :     } elsif ($key eq "Essential <=") {
4290 :     $reactionString .= '<span title="Essential in '.join(",",@{$predictionHash->{$key}}).'">E<=</span>,';
4291 :     } elsif ($key eq "Active =>") {
4292 :     $reactionString .= '<span title="Active in '.join(",",@{$predictionHash->{$key}}).'">A=></span>,';
4293 :     } elsif ($key eq "Active <=") {
4294 :     $reactionString .= '<span title="Active in '.join(",",@{$predictionHash->{$key}}).'">A<=</span>,';
4295 :     } elsif ($key eq "Active <=>") {
4296 :     $reactionString .= '<span title="Active in '.join(",",@{$predictionHash->{$key}}).'">A</span>,';
4297 :     } elsif ($key eq "Inactive") {
4298 :     $reactionString .= '<span title="Inactive in '.join(",",@{$predictionHash->{$key}}).'">I</span>,';
4299 :     } elsif ($key eq "Dead") {
4300 :     $reactionString .= '<span title="Dead">D</span>,';
4301 :     }
4302 :     }
4303 :     $reactionString =~ s/,$/)/;
4304 :     }
4305 :     push(@{$output},$reactionString);
4306 :     }
4307 :     }
4308 :    
4309 :     #Returning output
4310 :     return join("<br>",@{$output});
4311 :     }
4312 :    
4313 :     =head3 remove_obsolete_reactions
4314 :     Definition:
4315 :     void FIGMODELmodel->remove_obsolete_reactions();
4316 :     Description:
4317 :     =cut
4318 :     sub remove_obsolete_reactions {
4319 :     my ($self) = @_;
4320 :    
4321 :     (my $dummy,my $translation) = $self->figmodel()->put_two_column_array_in_hash($self->figmodel()->database()->load_multiple_column_file($self->figmodel()->config("Translation directory")->[0]."ObsoleteRxnIDs.txt","\t"));
4322 :     my $rxnTbl = $self->reaction_table();
4323 :     if (defined($rxnTbl)) {
4324 :     for (my $i=0; $i < $rxnTbl->size(); $i++) {
4325 :     my $row = $rxnTbl->get_row($i);
4326 :     if (defined($translation->{$row->{LOAD}->[0]}) || defined($translation->{$row->{LOAD}->[0]."r"})) {
4327 :     my $direction = $row->{DIRECTION}->[0];
4328 :     my $newRxn;
4329 :     if (defined($translation->{$row->{LOAD}->[0]."r"})) {
4330 :     $newRxn = $translation->{$row->{LOAD}->[0]."r"};
4331 :     if ($direction eq "<=") {
4332 :     $direction = "=>";
4333 :     } elsif ($direction eq "=>") {
4334 :     $direction = "<=";
4335 :     }
4336 :     } else {
4337 :     $newRxn = $translation->{$row->{LOAD}->[0]};
4338 :     }
4339 :     #Checking if the new reaction is already in the model
4340 :     my $newRow = $rxnTbl->get_row_by_key($newRxn,"LOAD");
4341 :     if (defined($newRow)) {
4342 :     #Handling direction
4343 :     if ($newRow->{DIRECTION}->[0] ne $direction) {
4344 :     $newRow->{DIRECTION}->[0] = "<=>";
4345 :     }
4346 :     push(@{$row->{"ASSOCIATED PEG"}},@{$rxnTbl->get_row($i)->{"ASSOCIATED PEG"}});
4347 :     } else {
4348 :     $rxnTbl->get_row($i)->{LOAD}->[0] = $newRxn;
4349 :     $rxnTbl->get_row($i)->{DIRECTION}->[0] = $direction;
4350 :     }
4351 :     }
4352 :     }
4353 :     $rxnTbl->save();
4354 :     }
4355 :     }
4356 :    
4357 :     1;

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3