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

Diff of /FigKernelPackages/FIGMODELmodel.pm

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.18, Mon May 31 18:05:08 2010 UTC revision 1.27, Fri Aug 20 14:35:38 2010 UTC
# Line 13  Line 13 
13          This is the constructor for the FIGMODELmodel object.          This is the constructor for the FIGMODELmodel object.
14  =cut  =cut
15  sub new {  sub new {
16          my ($class,$figmodel,$id) = @_;          my ($class,$figmodel,$id,$metagenome) = @_;
   
17          #Error checking first          #Error checking first
18          if (!defined($figmodel)) {          if (!defined($figmodel)) {
19                  print STDERR "FIGMODELmodel->new(undef,".$id."):figmodel must be defined to create a model object!\n";                  print STDERR "FIGMODELmodel->new(undef,".$id."):figmodel must be defined to create a model object!\n";
# Line 26  Line 25 
25                  $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,undef):id must be defined to create a model object");                  $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,undef):id must be defined to create a model object");
26                  return undef;                  return undef;
27          }          }
   
28          #Checking that the id exists          #Checking that the id exists
29            if (!defined($metagenome) || $metagenome != 1) {
30          my $modelHandle = $self->figmodel()->database()->get_object_manager("model");          my $modelHandle = $self->figmodel()->database()->get_object_manager("model");
31          if (!defined($modelHandle)) {                  if (defined($modelHandle)) {
                 $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not load MODELS table. Check database!");  
                 return undef;  
         }  
   
32          #If the id is a number, we get the model row by index          #If the id is a number, we get the model row by index
33          if ($id =~ m/^\d+$/) {          if ($id =~ m/^\d+$/) {
34                  my $objects = $modelHandle->get_objects();                  my $objects = $modelHandle->get_objects();
35                                    if (defined($objects->[$id])) {
36                  $self->{_data} = $objects->[$id];                  $self->{_data} = $objects->[$id];
37                  $self->figmodel()->{_models}->{$id} = $self;                  $self->figmodel()->{_models}->{$id} = $self;
38                                    }
39          } else {          } else {
40                  my $objects = $modelHandle->get_objects({id => $id});                  my $objects = $modelHandle->get_objects({id => $id});
41                  if (!defined($objects->[0]) && $id =~ m/(.+)(V[^V]+)/) {                                  if (defined($objects->[0])) {
42                                            $self->{_data} = $objects->[0];
43                                    } elsif (!defined($objects->[0]) && $id =~ m/(.+)(V[^V]+)/) {
44                          $objects = $modelHandle->get_objects({id => $1});                          $objects = $modelHandle->get_objects({id => $1});
45                          if (!defined($objects->[0])) {                                          if (defined($objects->[0])) {
                                 $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find model ".$id." in database!");  
                                 return undef;  
                         }  
46                          $self->{_selectedversion} = $2;                          $self->{_selectedversion} = $2;
47                  } elsif (!defined($objects->[0])) {                                                  $self->{_data} = $objects->[0];
48                          $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find model ".$id." in database!");                                          }
                         return undef;  
49                  }                  }
50                            }
51                    }
52                    if (defined($self->{_data})) {
53                            $self->{_modeltype} = "genome";
54                    }
55            }
56            if (!defined($self->{_data})) {
57                    my $modelHandle = $self->figmodel()->database()->get_object_manager("mgmodel");
58                    if (defined($modelHandle)) {
59                            #If the id is a number, we get the model row by index
60                            if ($id =~ m/^\d+$/) {
61                                    my $objects = $modelHandle->get_objects();
62                                    if (defined($objects->[$id])) {
63                                            $self->{_data} = $objects->[$id];
64                                            $self->figmodel()->{_models}->{"MG".$id} = $self;
65                                    }
66                            } else {
67                                    my $objects = $modelHandle->get_objects({id => $id});
68                                    if (defined($objects->[0])) {
69                                            $self->{_data} = $objects->[0];
70                                    } elsif (!defined($objects->[0]) && $id =~ m/(.+)(V[^V]+)/) {
71                                            $objects = $modelHandle->get_objects({id => $1});
72                                            if (defined($objects->[0])) {
73                                                    $self->{_selectedversion} = $2;
74                  $self->{_data} = $objects->[0];                  $self->{_data} = $objects->[0];
75          }          }
76                                    }
77                            }
78                    }
79                    if (defined($self->{_data})) {
80                            $self->{_modeltype} = "metagenome";
81                    }
82            }
83          if (!defined($self->{_data})) {          if (!defined($self->{_data})) {
84                  $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find specified id in database!");                  $self->figmodel()->error_message("FIGMODELmodel->new(figmodel,".$id."):could not find model ".$id." in database!");
85                  return undef;                  return undef;
86          }          }
87          $self->figmodel()->{_models}->{$self->id()} = $self;          $self->figmodel()->{_models}->{$self->id()} = $self;
   
88          return $self;          return $self;
89  }  }
90    
# Line 75  Line 100 
100          return $self->figmodel()->config($key);          return $self->figmodel()->config($key);
101  }  }
102    
103    =head3 error_message
104    Definition:
105            string:message text = FIGMODELmodel->error_message(string::message);
106    Description:
107    =cut
108    sub error_message {
109            my ($self,$message) = @_;
110            return $self->figmodel()->error_message("FIGMODELmodel:".$self->id().":".$message);
111    }
112    
113  =head3 fail  =head3 fail
114  Definition:  Definition:
115          -1 = FIGMODELmodel->fail();          -1 = FIGMODELmodel->fail();
# Line 166  Line 201 
201  =cut  =cut
202  sub mgdata {  sub mgdata {
203          my ($self) = @_;          my ($self) = @_;
204          if (!defined($self->{_mgdata}) && $self->source() =~ /^MGRAST/) {          if (!defined($self->{_mgdata})) {
205                  require MGRAST;                  my $mgrastdbhandle = $self->figmodel()->database()->get_object_manager("mgjob");
206                  $self->{_mgdata} = $self->figmodel()->mgrast()->Job->get_objects( { 'genome_id' => $self->genome() } )                  my $objects = $mgrastdbhandle->get_objects( { 'genome_id' => $self->genome() } );
207                    $self->{_mgdata} = $objects->[0];
208          }          }
209          return $self->{_mgdata};          return $self->{_mgdata};
210  }  }
# Line 184  Line 220 
220          return $self->{_data}->id();          return $self->{_data}->id();
221  }  }
222    
223    =head3 get_model_type
224    Definition:
225            string = FIGMODELmodel->get_model_type();
226    Description:
227            Returns the type of the model
228    =cut
229    sub get_model_type {
230            my ($self) = @_;
231            return $self->{_modeltype};
232    }
233    
234  =head3 owner  =head3 owner
235  Definition:  Definition:
236          string = FIGMODELmodel->owner();          string = FIGMODELmodel->owner();
# Line 195  Line 242 
242          return $self->{_data}->owner();          return $self->{_data}->owner();
243  }  }
244    
245    =head3 users
246    Definition:
247            string = FIGMODELmodel->users();
248    Description:
249    =cut
250    sub users {
251            my ($self) = @_;
252            return $self->{_data}->users();
253    }
254    
255  =head3 name  =head3 name
256  Definition:  Definition:
257          string = FIGMODELmodel->name();          string = FIGMODELmodel->name();
# Line 248  Line 305 
305                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
306                          } elsif ($ClassRow->{CLASS}->[0] eq "Negative") {                          } elsif ($ClassRow->{CLASS}->[0] eq "Negative") {
307                                  $class = "Essential <=";                                  $class = "Essential <=";
308                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$max)." to ".sprintf("%.3g",$min)."]<br>" : $class.="<br>[Flux: ".$max." to ".$min."]<br>";
309                          } elsif ($ClassRow->{CLASS}->[0] eq "Positive variable") {                          } elsif ($ClassRow->{CLASS}->[0] eq "Positive variable") {
310                                  $class = "Active =>";                                  $class = "Active =>";
311                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
312                          } elsif ($ClassRow->{CLASS}->[0] eq "Negative variable") {                          } elsif ($ClassRow->{CLASS}->[0] eq "Negative variable") {
313                                  $class = "Active <=";                                  $class = "Active <=";
314                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$max)." to ".sprintf("%.3g",$min)."]<br>" : $class.="<br>[Flux: ".$max." to ".$min."]<br>";
315                          } elsif ($ClassRow->{CLASS}->[0] eq "Variable") {                          } elsif ($ClassRow->{CLASS}->[0] eq "Variable") {
316                                  $class = "Active <=>";                                  $class = "Active <=>";
317                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $class.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $class.="<br>[Flux: ".$min." to ".$max."]<br>";
# Line 295  Line 352 
352                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
353                          } elsif ($ClassRow->{CLASS}->[$i] eq "Negative") {                          } elsif ($ClassRow->{CLASS}->[$i] eq "Negative") {
354                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Essential <=";                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Essential <=";
355                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$max)." to ".sprintf("%.3g",$min)."]<br>" : $NewClass.="<br>[Flux: ".$max." to ".$min."]<br>";
356                          } elsif ($ClassRow->{CLASS}->[$i] eq "Positive variable") {                          } elsif ($ClassRow->{CLASS}->[$i] eq "Positive variable") {
357                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Active =>";                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Active =>";
358                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
359                          } elsif ($ClassRow->{CLASS}->[$i] eq "Negative variable") {                          } elsif ($ClassRow->{CLASS}->[$i] eq "Negative variable") {
360                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=";                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=";
361                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$max)." to ".sprintf("%.3g",$min)."]<br>" : $NewClass.="<br>[Flux: ".$max." to ".$min."]<br>";
362                          } elsif ($ClassRow->{CLASS}->[$i] eq "Variable") {                          } elsif ($ClassRow->{CLASS}->[$i] eq "Variable") {
363                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=>";                                  $NewClass = $ClassRow->{MEDIA}->[$i].":Active <=>";
364                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";                                  $brief_flux ? $NewClass.="<br>[Flux: ".sprintf("%.3g",$min)." to ".sprintf("%.3g",$max)."]<br>" : $NewClass.="<br>[Flux: ".$min." to ".$max."]<br>";
# Line 352  Line 409 
409  =cut  =cut
410  sub get_reaction_data {  sub get_reaction_data {
411          my ($self,$reaction) = @_;          my ($self,$reaction) = @_;
   
412          if (!defined($self->reaction_table())) {          if (!defined($self->reaction_table())) {
413                  return undef;                  return undef;
414          }          }
415          if ($reaction =~ m/^\d+$/) {          if ($reaction =~ m/^\d+$/) {
416                  $self->reaction_table()->get_row($reaction);                  return $self->reaction_table()->get_row($reaction);
417          }          }
418          return $self->reaction_table()->get_row_by_key($reaction,"LOAD");          return $self->reaction_table()->get_row_by_key($reaction,"LOAD");
419  }  }
# Line 387  Line 443 
443  =cut  =cut
444  sub load_model_table {  sub load_model_table {
445          my ($self,$name,$refresh) = @_;          my ($self,$name,$refresh) = @_;
446          if (!defined($refresh)) {          if (defined($refresh) && $refresh == 1) {
                 $refresh = 1;  
         }  
         if ($refresh == 1) {  
447                  delete $self->{"_".$name};                  delete $self->{"_".$name};
448          }          }
449          if (!defined($self->{"_".$name})) {          if (!defined($self->{"_".$name})) {
450                  my $tbldef = $self->figmodel()->config("$name");                  my $tbldef = $self->figmodel()->config($name);
451                  if (!defined($tbldef)) {                  if (!defined($tbldef)) {
452                          return undef;                          return undef;
453                  }                  }
454                  my $filename = $self->directory().$name."-".$self->id().$self->selected_version().".tbl";                  my $itemDelim = "|";
455                  $self->{"_".$name} = $self->figmodel()->database()->load_table($filename,"\t","|",$tbldef->{headingline}->[0],$tbldef->{hashcolumns});                  if (defined($tbldef->{itemdelimiter}->[0])) {
456                  if (!defined($self->{"_".$name})) {                          $itemDelim = $tbldef->{itemdelimiter}->[0];
457                            if ($itemDelim eq "SC") {
458                                    $itemDelim = ";";
459                            }
460                    }
461                    my $columnDelim = "\t";
462                    if (defined($tbldef->{columndelimiter}->[0])) {
463                            $columnDelim = $tbldef->{columndelimiter}->[0];
464                            if ($columnDelim eq "SC") {
465                                    $columnDelim = ";";
466                            }
467                    }
468                    my $suffix = ".tbl";
469                    if (defined($tbldef->{filename_suffix}->[0])) {
470                            $suffix = $tbldef->{filename_suffix}->[0];
471                    }
472                    my $filename = $self->directory().$name."-".$self->id().$self->selected_version().$suffix;
473                    if (defined($tbldef->{filename_prefix}->[0])) {
474                            if ($tbldef->{filename_prefix}->[0] eq "NONE") {
475                                    $filename = $self->directory().$self->id().$self->selected_version().$suffix;
476                            } else {
477                                    $filename = $self->directory().$tbldef->{filename_prefix}->[0]."-".$self->id().$self->selected_version().$suffix;
478                            }
479                    }
480                    if (-e $filename) {
481                            $self->{"_".$name} = $self->figmodel()->database()->load_table($filename,$columnDelim,$itemDelim,$tbldef->{headingline}->[0],$tbldef->{hashcolumns});
482                    } else {
483                          if (defined($tbldef->{prefix})) {                          if (defined($tbldef->{prefix})) {
484                                  $self->{"_".$name} = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},"\t","|",join(@{$tbldef->{prefix}},"\n"));                                  $self->{"_".$name} = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},$columnDelim,$itemDelim,join(@{$tbldef->{prefix}},"\n"));
485                          } else {                          } else {
486                                  $self->{"_".$name} = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},"\t","|");                                  $self->{"_".$name} = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},$columnDelim,$itemDelim);
487                          }                          }
488                  }                  }
489          }          }
# Line 421  Line 500 
500  =cut  =cut
501  sub create_table_prototype {  sub create_table_prototype {
502          my ($self,$TableName) = @_;          my ($self,$TableName) = @_;
   
503          #Checking if the table definition exists in the FIGMODELconfig file          #Checking if the table definition exists in the FIGMODELconfig file
504          if (!defined($self->figmodel()->config($TableName))) {          my $tbldef = $self->figmodel()->config($TableName);
505            if (!defined($tbldef)) {
506                  $self->figmodel()->error_message("FIGMODELdatabase:create_table_prototype:Definition not found for ".$TableName);                  $self->figmodel()->error_message("FIGMODELdatabase:create_table_prototype:Definition not found for ".$TableName);
507                  return undef;                  return undef;
508          }          }
509          #Checking that this is a database table          #Checking that this is a database table
510          if (!defined($self->config($TableName)->{tabletype}) || $self->config($TableName)->{tabletype}->[0] ne "ModelTable") {          if (!defined($tbldef->{tabletype}) || $tbldef->{tabletype}->[0] ne "ModelTable") {
511                  $self->figmodel()->error_message("FIGMODELdatabase:create_table_prototype:".$TableName." is not a model table!");                  $self->figmodel()->error_message("FIGMODELdatabase:create_table_prototype:".$TableName." is not a model table!");
512                  return undef;                  return undef;
513          }          }
514          if (!defined($self->config($TableName)->{delimiter})) {          #Setting default values for table parameters
515                  $self->config($TableName)->{delimiter}->[0] = ";";          my $prefix;
516            if (defined($tbldef->{prefix})) {
517                    $prefix = join("\n",@{$self->config($TableName)->{prefix}})."\n";
518          }          }
519          if (!defined($self->config($TableName)->{itemdelimiter})) {          my $itemDelim = "|";
520                  $self->config($TableName)->{itemdelimiter}->[0] = "|";          if (defined($tbldef->{itemdelimiter}->[0])) {
521                    $itemDelim = $tbldef->{itemdelimiter}->[0];
522                    if ($itemDelim eq "SC") {
523                            $itemDelim = ";";
524          }          }
         my $prefix;  
         if (defined($self->config($TableName)->{prefix})) {  
                 $prefix = join("\n",@{$self->config($TableName)->{prefix}});  
525          }          }
526          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);          my $columnDelim = "\t";
527            if (defined($tbldef->{columndelimiter}->[0])) {
528                    $columnDelim = $tbldef->{columndelimiter}->[0];
529                    if ($columnDelim eq "SC") {
530                            $columnDelim = ";";
531                    }
532            }
533            my $suffix = ".tbl";
534            if (defined($tbldef->{filename_suffix}->[0])) {
535                    $suffix = $tbldef->{filename_suffix}->[0];
536            }
537            my $filename = $self->directory().$TableName."-".$self->id().$self->selected_version().$suffix;
538            if (defined($tbldef->{filename_prefix}->[0])) {
539                    if ($tbldef->{filename_prefix}->[0] eq "NONE") {
540                            $filename = $self->directory().$self->id().$self->selected_version().$suffix;
541                    } else {
542                            $filename = $self->directory().$tbldef->{filename_prefix}->[0]."-".$self->id().$self->selected_version().$suffix;
543                    }
544            }
545            #Creating the table prototype
546            my $tbl = FIGMODELTable->new($tbldef->{columns},$filename,$tbldef->{hashcolumns},$columnDelim,$itemDelim,$prefix);
547          return $tbl;          return $tbl;
548  }  }
549    
# Line 467  Line 568 
568          Returns FIGMODELTable with the reaction list for the model          Returns FIGMODELTable with the reaction list for the model
569  =cut  =cut
570  sub reaction_table {  sub reaction_table {
571          my ($self) = @_;          my ($self,$clear) = @_;
572            if (defined($clear) && $clear == 1) {
573          if (!defined($self->{_reaction_data})) {                  delete $self->{_reaction_table};
574                  $self->{_reaction_data} = $self->figmodel()->database()->GetDBModel($self->id());          }
575            if (!defined($self->{_reaction_table})) {
576                    $self->{_reaction_table} = $self->load_model_table("ModelReactions",$clear);
577                  my $classTbl = $self->reaction_class_table();                  my $classTbl = $self->reaction_class_table();
578                    if (defined($classTbl)) {
579                  for (my $i=0; $i < $classTbl->size(); $i++) {                  for (my $i=0; $i < $classTbl->size(); $i++) {
580                          my $row = $classTbl->get_row($i);                          my $row = $classTbl->get_row($i);
581                          my $rxnRow = $self->{_reaction_data}->get_row_by_key($row->{"REACTION"}->[0],"LOAD");                                  if (defined($row->{REACTION})) {
582                                            my $rxnRow = $self->{_reaction_table}->get_row_by_key($row->{"REACTION"}->[0],"LOAD");
583                                            if (defined($row->{MEDIA})) {
584                          for (my $j=0; $j < @{$row->{MEDIA}};$j++) {                          for (my $j=0; $j < @{$row->{MEDIA}};$j++) {
585                                  my $class = "Active <=>";                                  my $class = "Active <=>";
586                                  if ($row->{CLASS}->[$j] eq "Positive") {                                  if ($row->{CLASS}->[$j] eq "Positive") {
# Line 496  Line 602 
602                          }                          }
603                  }                  }
604          }          }
605                            }
606                    }
607            }
608            return $self->{_reaction_table};
609    }
610    
611    =head3 essentials_table
612    Definition:
613            FIGMODELTable = FIGMODELmodel->essentials_table();
614    Description:
615            Returns FIGMODELTable with the essential genes for the model
616    =cut
617    sub essentials_table {
618            my ($self,$clear) = @_;
619            my $tbl = $self->load_model_table("ModelEssentialGenes",$clear);
620            return $tbl;
621    }
622    
623          return $self->{_reaction_data};  =head3 model_history
624    Definition:
625            FIGMODELTable = FIGMODELmodel->model_history();
626    Description:
627            Returns FIGMODELTable with the history of model changes
628    =cut
629    sub model_history {
630            my ($self,$clear) = @_;
631            return $self->load_model_table("ModelHistory",$clear);
632  }  }
633    
634  =head3 feature_table  =head3 feature_table
# Line 551  Line 682 
682                  }                  }
683              }              }
684              #Loading predictions              #Loading predictions
685              my @files = glob($self->directory()."EssentialGenes-".$self->id()."-*");                  my $esstbl = $self->essentials_table();
             foreach my $file (@files) {  
                 if ($file =~ m/\-([^\-]+)\.tbl/) {  
                         my $media = $1;  
                         my $list = $self->figmodel()->database()->load_single_column_file($file,"");  
                         my $hash;  
                         foreach my $gene (@{$list}) {  
                                 $hash->{$gene} = 1;  
                         }  
686                          for (my $i=0; $i < $self->{_feature_data}->size(); $i++) {                          for (my $i=0; $i < $self->{_feature_data}->size(); $i++) {
687                                  my $Row = $self->{_feature_data}->get_row($i);                                  my $Row = $self->{_feature_data}->get_row($i);
688                                  if ($Row->{ID}->[0] =~ m/(peg\.\d+)/) {                                  if ($Row->{ID}->[0] =~ m/(peg\.\d+)/) {
689                                          my $gene = $1;                                          my $gene = $1;
690                                          if (defined($hash->{$gene})) {                                  my @rows = $esstbl->get_rows_by_key($gene,"ESSENTIAL GENES");
691                                                  push(@{$Row->{$self->id()."PREDICTIONS"}},$media.":essential");                                  my $mediahash;
692                                    for (my $j=0; $j < $esstbl->size(); $j++) {
693                                            $mediahash->{$esstbl->get_row($j)->{MEDIA}->[0]} = 0;
694                                    }
695                                    for (my $j=0; $j < @rows; $j++) {
696                                            $mediahash->{$rows[$j]->{MEDIA}->[0]} = 1;
697                                    }
698                                    my @mediaList = keys(%{$mediahash});
699                                    for (my $j=0; $j < @mediaList; $j++) {
700                                            if ($mediahash->{$mediaList[$j]} == 1) {
701                                                    push(@{$Row->{$self->id()."PREDICTIONS"}},$mediaList[$j].":essential");
702                                          } else {                                          } else {
703                                                  push(@{$Row->{$self->id()."PREDICTIONS"}},$media.":nonessential");                                                  push(@{$Row->{$self->id()."PREDICTIONS"}},$mediaList[$j].":nonessential");
                                         }  
704                                  }                                  }
705                          }                          }
706                  }                  }
# Line 584  Line 716 
716          Returns FIGMODELTable with the reaction class data, and creates the table file  if it does not exist          Returns FIGMODELTable with the reaction class data, and creates the table file  if it does not exist
717  =cut  =cut
718  sub reaction_class_table {  sub reaction_class_table {
719          my ($self) = @_;          my ($self,$clear) = @_;
720            return $self->load_model_table("ModelReactionClasses",$clear);
         if (!defined($self->{_reaction_class_table})) {  
                 if (-e $self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl") {  
                         $self->{_reaction_class_table} = $self->figmodel()->database()->load_table($self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl",";","|",0,["REACTION","CLASS","MEDIA"]);  
                 } else {  
                         $self->{_reaction_class_table} = FIGMODELTable->new(["REACTION","MEDIA","CLASS","MIN","MAX"],$self->directory()."ReactionClassification-".$self->id().$self->selected_version().".tbl",["REACTION","CLASS","MEDIA"],";","|",undef);  
                 }  
         }  
   
         return $self->{_reaction_class_table};  
721  }  }
722    
723  =head3 compound_class_table  =head3 compound_class_table
# Line 604  Line 727 
727          Returns FIGMODELTable with the compound class data, and creates the table file  if it does not exist          Returns FIGMODELTable with the compound class data, and creates the table file  if it does not exist
728  =cut  =cut
729  sub compound_class_table {  sub compound_class_table {
730          my ($self) = @_;          my ($self,$clear) = @_;
731            return $self->load_model_table("ModelCompoundClasses",$clear);
         if (!defined($self->{_compound_class_table})) {  
                 if (-e $self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl") {  
                         $self->{_compound_class_table} = $self->figmodel()->database()->load_table($self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl",";","|",0,["COMPOUND","CLASS","MEDIA"]);  
                 } else {  
                         $self->{_compound_class_table} = FIGMODELTable->new(["COMPOUND","MEDIA","CLASS","MIN","MAX"],$self->directory()."CompoundClassification-".$self->id().$self->selected_version().".tbl",["COMPOUND","CLASS","MEDIA"],";","|",undef);  
                 }  
         }  
   
         return $self->{_compound_class_table};  
732  }  }
733    
734  =head3 get_essential_genes  =head3 get_essential_genes
# Line 625  Line 739 
739  =cut  =cut
740  sub get_essential_genes {  sub get_essential_genes {
741          my ($self,$media) = @_;          my ($self,$media) = @_;
742            my $tbl = $self->essentials_table();
743          if (!defined($media)) {          my $row = $tbl->get_row_by_key($media,"MEDIA");
744                  $media = "Complete";          if (defined($row)) {
745          }                  return $row->{"ESSENTIAL GENES"};
         if (!defined($self->{_essential_genes}->{$media})) {  
                 $self->{_essential_genes}->{$media} = undef;  
                 if (-e $self->directory()."EssentialGenes-".$self->id().$self->selected_version()."-".$media.".tbl") {  
                         $self->{_essential_genes}->{$media} = $self->figmodel()->database()->load_single_column_file($self->directory()."EssentialGenes-".$self->id().$self->selected_version()."-".$media.".tbl","");  
746                  }                  }
747          }          return undef;
   
         return $self->{_essential_genes}->{$media};  
748  }  }
749    
750  =head3 compound_table  =head3 compound_table
# Line 649  Line 757 
757          my ($self) = @_;          my ($self) = @_;
758    
759          if (!defined($self->{_compound_table})) {          if (!defined($self->{_compound_table})) {
760                  $self->{_compound_table} = $self->figmodel()->database()->GetDBModelCompounds($self->id());                  $self->{_compound_table} = $self->create_table_prototype("ModelCompounds");
761                    #Loading the reactions
762                    my $ReactionTable = $self->figmodel()->database()->get_table("REACTIONS");
763                    my $BiomassTable = $self->figmodel()->database()->get_table("BIOMASS");
764                    #Loading the model
765                    my $ModelTable = $self->reaction_table();
766                    #Checking that the tables were loaded
767                    if (!defined($ModelTable) || !defined($ReactionTable)) {
768                            return undef;
769                    }
770                    #Finding the biomass reaction
771                    for (my $i=0; $i < $ModelTable->size(); $i++) {
772                            my $ID = $ModelTable->get_row($i)->{"LOAD"}->[0];
773                            my $Row = $ReactionTable->get_row_by_key($ID,"DATABASE");
774                            my $IsBiomass = 0;
775                            if (!defined($Row)) {
776                                    $Row = $BiomassTable->get_row_by_key($ID,"DATABASE");
777                                    $IsBiomass = 1;
778                            }
779                            if (defined($Row->{"EQUATION"}->[0])) {
780                                    $_ = $Row->{"EQUATION"}->[0];
781                                    my @OriginalArray = /(cpd\d\d\d\d\d[\[\w]*)/g;
782                                    foreach my $Compound (@OriginalArray) {
783                                            my $ID = substr($Compound,0,8);
784                                            my $NewRow = $self->{_compound_table}->get_row_by_key($ID,"DATABASE",1);
785                                            if ($IsBiomass == 1) {
786                                                    $self->{_compound_table}->add_data($NewRow,"BIOMASS",$Row->{"DATABASE"}->[0],1);
787                                            }
788                                            if (length($Compound) > 8) {
789                                                    #print $Compound."\t".$Row->{"EQUATION"}->[0]."\t".$Row->{"DATABASE"}->[0]."\n";
790                                                    my $Compartment = substr($Compound,8,1);
791                                                    $self->{_compound_table}->add_data($NewRow,"COMPARTMENTS",$Compartment,1);
792                                                    $self->{_compound_table}->add_data($NewRow,"TRANSPORTERS",$Row->{"DATABASE"}->[0],1);
793                                            }
794                                    }
795                            }
796                    }
797          }          }
798    
799          return $self->{_compound_table};          return $self->{_compound_table};
# Line 769  Line 913 
913                  }                  }
914                  my $source = $self->source();                  my $source = $self->source();
915                  if ($source =~ /^MGRAST/) {                  if ($source =~ /^MGRAST/) {
916                          $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";                          $self->{_directory} = $self->figmodel()->config("mgrast model directory")->[0].$userdirectory.$self->genome()."/";
917                  } elsif ($source =~ /^RAST/) {                  } elsif ($source =~ /^RAST/) {
918                          $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";                          $self->{_directory} = $self->figmodel()->config("organism directory")->[0].$userdirectory.$self->genome()."/";
919                  } elsif ($source =~ /^SEED/) {                  } elsif ($source =~ /^SEED/) {
# Line 909  Line 1053 
1053          return $self->{_data}->cellwalltype();          return $self->{_data}->cellwalltype();
1054  }  }
1055    
1056    sub autocompleteMedia {
1057            my ($self,$newMedia) = @_;
1058            if (defined($newMedia)) {
1059                    return $self->{_data}->autoCompleteMedia($newMedia);
1060            }
1061            return $self->{_data}->autoCompleteMedia();
1062    }
1063    
1064    sub biomassReaction {
1065            my ($self,$newBiomass) = @_;
1066            if (!defined($newBiomass)) {
1067                    return $self->{_data}->biomassReaction();
1068            } else {
1069                    #Figuring out what the old biomass is
1070                    my $oldBiomass = $self->{_data}->biomassReaction();
1071                    $self->{_data}->biomassReaction($newBiomass);
1072                    #Changing the biomass reaction in the model file
1073                    my $rxnTbl = $self->reaction_table();
1074                    for (my $i=0; $i < $rxnTbl->size(); $i++) {
1075                            my $row = $rxnTbl->get_row($i);
1076                            if ($row->{LOAD}->[0] =~ m/^bio/) {
1077                                    $row->{LOAD}->[0] = $newBiomass;
1078                            }
1079                    }
1080                    $rxnTbl->save();
1081                    if ($newBiomass ne $oldBiomass) {
1082                            #Figuring out if the new biomass exists
1083                            my $handle = $self->figmodel()->database()->get_object_manager("bof");
1084                            my $objects = $handle->get_objects({id=>$newBiomass});
1085                            if (!defined($objects) || !defined($objects->[0])) {
1086                                    print STDERR "Could not find new biomass reaction ".$newBiomass."\n";
1087                                    return $oldBiomass;
1088                            }
1089                    }
1090                    return $newBiomass;
1091            }
1092    }
1093    
1094  =head3 taxonomy  =head3 taxonomy
1095  Definition:  Definition:
1096          string = FIGMODELmodel->taxonomy();          string = FIGMODELmodel->taxonomy();
# Line 992  Line 1174 
1174    
1175          #Assuming complete media if none is provided          #Assuming complete media if none is provided
1176          if (!defined($Media)) {          if (!defined($Media)) {
1177                  $Media = "Complete";                  $Media = $self->autocompleteMedia();
1178          }          }
1179    
1180          #Predicting essentiality          #Predicting essentiality
1181          my $result = $self->figmodel()->RunFBASimulation($self->id(),"SINGLEKO",undef,undef,[$self->id()],[$Media]);          my $result = $self->figmodel()->RunFBASimulation($self->id(),"SINGLEKO",undef,undef,[$self->id()],[$Media]);
1182          #Checking that the table is defined and the output file exists          #Checking that the table is defined and the output file exists
1183          if (defined($result) && defined($result->get_row(0)->{"ESSENTIALGENES"})) {          if (defined($result) && defined($result->get_row(0)->{"ESSENTIALGENES"})) {
1184                  $self->figmodel()->database()->print_array_to_file($self->directory()."EssentialGenes-".$self->id()."-".$Media.".tbl",[join("\n",@{$result->get_row(0)->{"ESSENTIALGENES"}})]);                  my $tbl = $self->essentials_table();
1185                    my $row = $tbl->get_row_by_key($Media,"MEDIA",1);
1186                    $row->{"ESSENTIAL GENES"} = $result->get_row(0)->{"ESSENTIALGENES"};
1187                    $tbl->save();
1188          } else {          } else {
1189                  $self->figmodel()->error_message("FIGMODELmodel:run_default_model_predictions:could not identify essential reactions for model ".$self->id().$self->selected_version().".");                  $self->figmodel()->error_message("FIGMODELmodel:run_default_model_predictions:could not identify essential reactions for model ".$self->id().$self->selected_version().".");
1190                  return $self->figmodel()->fail();                  return $self->figmodel()->fail();
# Line 1028  Line 1213 
1213          $self->{_data}->modificationDate(time());          $self->{_data}->modificationDate(time());
1214          my $version = $self->{_data}->autocompleteVersion();          my $version = $self->{_data}->autocompleteVersion();
1215          $self->{_data}->autocompleteVersion($version+1);          $self->{_data}->autocompleteVersion($version+1);
1216            $self->update_model_stats();
1217  }  }
1218    
1219  =head3 update_stats_for_build  =head3 update_stats_for_build
# Line 1054  Line 1240 
1240          #Getting reaction table          #Getting reaction table
1241          my $rxntbl = $self->reaction_table();          my $rxntbl = $self->reaction_table();
1242          if (!defined($rxntbl)) {          if (!defined($rxntbl)) {
1243                  $self->figmodel()->error_message("FIGMODELmodel:update_model_stats:Could not load reaction list for ".$self->id());                  die $self->error_message("update_model_stats:Could not load reaction list!");
                 return undef;  
1244          }          }
1245          my $cpdtbl = $self->compound_table();          my $cpdtbl = $self->compound_table();
1246    
# Line 1109  Line 1294 
1294          #Setting the reaction count          #Setting the reaction count
1295          $self->{_data}->reactions($rxntbl->size());          $self->{_data}->reactions($rxntbl->size());
1296          #Setting the metabolite count          #Setting the metabolite count
1297          $self->{_data}->compounds($rxntbl->size());          $self->{_data}->compounds($cpdtbl->size());
1298          #Setting the gene count          #Setting the gene count
1299          my $geneCount = @genes + @othergenes;          my $geneCount = @genes + @othergenes;
1300          $self->{_data}->associatedGenes($geneCount);          $self->{_data}->associatedGenes($geneCount);
# Line 1159  Line 1344 
1344          my $UniqueFilename = $self->figmodel()->filename();          my $UniqueFilename = $self->figmodel()->filename();
1345          my $StartTime = time();          my $StartTime = time();
1346    
1347          #Archiving the existing model          #Reading original reaction table
1348          $self->ArchiveModel();          my $OriginalRxn = $self->reaction_table();
1349            #Clearing the table
1350            $self->reaction_table(1);
1351    
1352          #Removing any gapfilling reactions that may be currently present in the model          #Removing any gapfilling reactions that may be currently present in the model
1353          if (!defined($donotclear) || $donotclear != 1) {          if (!defined($donotclear) || $donotclear != 1) {
1354                  my $ModelTable = $self->reaction_table();                  my $ModelTable = $self->reaction_table();
1355                  for (my $i=0; $i < $ModelTable->size(); $i++) {                  for (my $i=0; $i < $ModelTable->size(); $i++) {
1356                          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/) {                          if (!defined($ModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0]) || $ModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] eq "AUTOCOMPLETION" || $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/) {
1357                                  $ModelTable->delete_row($i);                                  $ModelTable->delete_row($i);
1358                                  $i--;                                  $i--;
1359                          }                          }
# Line 1175  Line 1362 
1362          }          }
1363    
1364          #Calling the MFAToolkit to run the gap filling optimization          #Calling the MFAToolkit to run the gap filling optimization
1365          my $MinimalMediaTable = $self->figmodel()->database()->GetDBTable("MINIMAL MEDIA TABLE");          my $Media = $self->autocompleteMedia();
1366          my $Media = "Complete";          if ($Media eq "Complete") {
1367          if (defined($MinimalMediaTable->get_row_by_key($self->genome(),"Organism"))) {                  system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),undef,["GapFilling"],{"Reactions to knockout" => $self->config("permanently knocked out reactions")->[0]},"GapFill".$self->id().".log",undef));
1368                  $Media = $MinimalMediaTable->get_row_by_key($self->genome(),"Organism")->{"Minimal media"}->[0];          } else {
1369                  #Loading media, changing bounds, saving media as a test media                  #Loading media, changing bounds, saving media as a test media
1370                  my $MediaTable = FIGMODELTable::load_table($self->config("Media directory")->[0].$Media.".txt",";","",0,["VarName"]);                  my $MediaTable = FIGMODELTable::load_table($self->config("Media directory")->[0].$Media.".txt",";","",0,["VarName"]);
1371                  for (my $i=0; $i < $MediaTable->size(); $i++) {                  for (my $i=0; $i < $MediaTable->size(); $i++) {
# Line 1190  Line 1377 
1377                          }                          }
1378                  }                  }
1379                  $MediaTable->save($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");                  $MediaTable->save($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");
1380                    #print $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)."\n";
1381                  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));                  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));
1382                  unlink($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");                  unlink($self->config("Media directory")->[0].$UniqueFilename."TestMedia.txt");
         } else {  
                 system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),undef,["GapFilling"],{"Reactions to knockout" => $self->config("permanently knocked out reactions")->[0]},"GapFill".$self->id().".log",undef));  
1383          }          }
1384    
1385          #Parse the solutions file for the model and get the reaction list from it          #Parse the solutions file for the model and get the reaction list from it
1386          my $SolutionData = $self->figmodel()->LoadProblemReport($UniqueFilename);          my $SolutionData = $self->figmodel()->LoadProblemReport($UniqueFilename);
1387    
1388          #Clearing the mfatoolkit output and log file          #Clearing the mfatoolkit output and log file
1389          $self->figmodel()->clearing_output($UniqueFilename,"GapFill".$self->id().".log");          $self->figmodel()->clearing_output($UniqueFilename,"GapFill".$self->id().".log");
1390          if (!defined($SolutionData) || $SolutionData->size() == 0) {          if (!defined($SolutionData) || $SolutionData->size() == 0) {
# Line 1205  Line 1392 
1392                  $self->figmodel()->error_message("FIGMODEL:GapFillModel: Gap filling report not found!");                  $self->figmodel()->error_message("FIGMODEL:GapFillModel: Gap filling report not found!");
1393                  return $self->fail();                  return $self->fail();
1394          }          }
         $SolutionData->save("/home/chenry/Solution.txt");  
1395    
1396          #Looking for the last printed recursive MILP solution          #Looking for the last printed recursive MILP solution
1397          for (my $i=($SolutionData->size()-1); $i >=0; $i--) {          for (my $i=($SolutionData->size()-1); $i >=0; $i--) {
# Line 1233  Line 1419 
1419                                                  push(@{$ReactionList},$ID);                                                  push(@{$ReactionList},$ID);
1420                                          }                                          }
1421                                  }                                  }
1422                                  $self->figmodel()->IntegrateGrowMatchSolution($self->id(),undef,$ReactionList,$DirectionList,"GAP FILLING",0,1);                                  $self->figmodel()->IntegrateGrowMatchSolution($self->id(),undef,$ReactionList,$DirectionList,"AUTOCOMPLETION",0,1);
1423                          }                          }
1424                          last;                          last;
1425                  }                  }
1426          }          }
1427    
1428          #Updating model stats with gap filling results          #Updating model stats with gap filling results
1429          my $ElapsedTime = time() - $StartTime;          my $ElapsedTime = time() - $StartTime;
1430          $self->figmodel()->database()->ClearDBModel($self->id(),1);          $self->reaction_table(1);
1431            $self->calculate_model_changes($OriginalRxn,"AUTOCOMPLETION");
1432    
1433          #Determining why each gap filling reaction was added          #Determining why each gap filling reaction was added
1434          $self->figmodel()->IdentifyDependancyOfGapFillingReactions($self->id(),$Media);          $self->figmodel()->IdentifyDependancyOfGapFillingReactions($self->id(),$Media);
1435          if (!defined($donotclear) || $donotclear != 1) {          if (!defined($donotclear) || $donotclear != 1) {
                 $self->figmodel()->AddBiologTransporters($self->id());  
1436                  if ($self->id() !~ m/MGRast/) {                  if ($self->id() !~ m/MGRast/) {
1437                          $self->update_stats_for_gap_filling($ElapsedTime);                          $self->update_stats_for_gap_filling($ElapsedTime);
1438                  }                  }
# Line 1260  Line 1448 
1448          return $self->success();          return $self->success();
1449  }  }
1450    
1451    =head3 calculate_model_changes
1452    Definition:
1453            FIGMODELmodel->calculate_model_changes(FIGMODELTable:original reaction table,string:modification cause);
1454    Description:
1455    
1456    =cut
1457    
1458    sub calculate_model_changes {
1459            my ($self,$originalReactions,$cause,$tbl,$version) = @_;
1460            my $modTime = time();
1461            if (!defined($version)) {
1462                    $version = $self->selected_version();
1463            }
1464            my $user = $self->figmodel()->user();
1465            #Getting the history table
1466            my $histTbl = $self->model_history();
1467            #Checking for differences
1468            if (!defined($tbl)) {
1469                    $tbl = $self->reaction_table();
1470            }
1471            for (my $i=0; $i < $tbl->size(); $i++) {
1472                    my $row = $tbl->get_row($i);
1473                    my $orgRow = $originalReactions->get_row_by_key($row->{LOAD}->[0],"LOAD");
1474                    if (!defined($orgRow)) {
1475                            $histTbl->add_row({Reaction => [$row->{LOAD}->[0]], DirectionChange => $row->{DIRECTIONALITY}, GeneChange => $row->{"ASSOCIATED PEG"}, Action => ["ADDED"], ModificationTime => [$modTime], ModifcationCause => [$cause], User => [$user], Version => [$version]});
1476                    } else {
1477                            my $geneChanges;
1478                            my $directionChange;
1479                            if ($orgRow->{"DIRECTIONALITY"}->[0] ne $row->{"DIRECTIONALITY"}->[0]) {
1480                                    $directionChange = $orgRow->{"DIRECTIONALITY"}->[0]." to ".$row->{"DIRECTIONALITY"}->[0];
1481                            }
1482                            for (my $j=0; $j < @{$row->{"ASSOCIATED PEG"}}; $j++) {
1483                                    my $match = 0;
1484                                    if (defined($orgRow->{"ASSOCIATED PEG"})) {
1485                                            for (my $k=0; $k < @{$orgRow->{"ASSOCIATED PEG"}}; $k++) {
1486                                                    if ($row->{"ASSOCIATED PEG"}->[$j] eq $orgRow->{"ASSOCIATED PEG"}->[$k]) {
1487                                                            $match = 1;
1488                                                    }
1489                                            }
1490                                    }
1491                                    if ($match == 0) {
1492                                            push(@{$geneChanges},"Added ".$row->{"ASSOCIATED PEG"}->[$j]);
1493                                    }
1494                            }
1495                            if (defined($orgRow->{"ASSOCIATED PEG"})) {
1496                                    for (my $k=0; $k < @{$orgRow->{"ASSOCIATED PEG"}}; $k++) {
1497                                            my $match = 0;
1498                                            if (defined($row->{"ASSOCIATED PEG"})) {
1499                                                    for (my $j=0; $j < @{$row->{"ASSOCIATED PEG"}}; $j++) {
1500                                                            if ($row->{"ASSOCIATED PEG"}->[$j] eq $orgRow->{"ASSOCIATED PEG"}->[$k]) {
1501                                                                    $match = 1;
1502                                                            }
1503                                                    }
1504                                            }
1505                                            if ($match == 0) {
1506                                                    push(@{$geneChanges},"Removed ".$orgRow->{"ASSOCIATED PEG"}->[$k]);
1507                                            }
1508                                    }
1509                            }
1510                            if ((defined($directionChange) && length($directionChange) > 0) || defined($geneChanges) && @{$geneChanges} > 0) {
1511                                    $histTbl->add_row({Reaction => [$row->{LOAD}->[0]], DirectionChange => [$directionChange], GeneChange => $geneChanges, Action => ["CHANGE"], ModificationTime => [$modTime], ModifcationCause => [$cause], User => [$user], Version => [$version]});
1512                            }
1513                    }
1514            }
1515            #Looking for removed reactions
1516            for (my $i=0; $i < $originalReactions->size(); $i++) {
1517                    my $row = $originalReactions->get_row($i);
1518                    my $orgRow = $tbl->get_row_by_key($row->{LOAD}->[0],"LOAD");
1519                    if (!defined($orgRow)) {
1520                            $histTbl->add_row({Reaction => [$row->{LOAD}->[0]], DirectionChange => $row->{DIRECTIONALITY}, GeneChange => $row->{"ASSOCIATED PEG"}, Action => ["REMOVED"], ModificationTime => [$modTime], ModifcationCause => [$cause], User => [$user], Version => [$version]});
1521                    }
1522            }
1523            $histTbl->save();
1524    }
1525    
1526  =head3 GapGenModel  =head3 GapGenModel
1527  Definition:  Definition:
1528          FIGMODELmodel->GapGenModel();          FIGMODELmodel->GapGenModel();
# Line 1497  Line 1760 
1760  =cut  =cut
1761  sub status {  sub status {
1762          my ($self) = @_;          my ($self) = @_;
1763          return $self->{_data}->{status}->[0];          return $self->{_data}->status();
1764  }  }
1765    
1766  =head3 message  =head3 message
# Line 1508  Line 1771 
1771  =cut  =cut
1772  sub message {  sub message {
1773          my ($self) = @_;          my ($self) = @_;
1774          return $self->{_data}->{message}->[0];          return $self->{_data}->message();
1775  }  }
1776    
1777  =head3 set_status  =head3 set_status
# Line 1523  Line 1786 
1786  =cut  =cut
1787  sub set_status {  sub set_status {
1788          my ($self,$NewStatus,$Message) = @_;          my ($self,$NewStatus,$Message) = @_;
1789            $self->{_data}->status($NewStatus);
1790          #Getting the model row from the MODELS table          $self->{_data}->message($Message);
         $self->{_data}->{status}->[0] = $NewStatus;  
         $self->{_data}->{message}->[0] = $Message;  
         $self->figmodel()->database()->update_row("MODELS",$self->{_data},"id");  
1791          return $self->config("SUCCESS")->[0];          return $self->config("SUCCESS")->[0];
1792  }  }
1793    
# Line 1566  Line 1826 
1826          }          }
1827          #Setting up needed variables          #Setting up needed variables
1828          my $OriginalModelTable = undef;          my $OriginalModelTable = undef;
1829          #Locking model status table          if ($self->status() == 0) {
         my $ModelTable = $self->figmodel()->database()->LockDBTable("MODELS");  
         my $Row = $ModelTable->get_row_by_key($self->id(),"id");  
         if (!defined($Row) || !defined($Row->{status}->[0])) {  
                 $self->figmodel()->database()->UnlockDBTable("MODELS");  
                 $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList: ".$self->id()." has no row in models table!");  
                 return $self->fail();  
         } elsif ($Row->{status}->[0] == 0) {  
1830                  $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList:Model is already being built. Canceling current build.");                  $self->figmodel()->error_message("FIGMODEL:CreateSingleGenomeReactionList:Model is already being built. Canceling current build.");
                 $self->figmodel()->database()->UnlockDBTable("MODELS");  
1831                  return $self->fail();                  return $self->fail();
1832          }elsif ($Row->{status}->[0] == 1) {          }elsif ($self->status() == 1) {
1833                  $OriginalModelTable = $self->reaction_table();                  $OriginalModelTable = $self->reaction_table();
1834                  $self->ArchiveModel();                  $self->set_status(0,"Rebuilding preliminary reconstruction");
                 $Row->{status}->[0] = 0;  
                 $Row->{message}->[0] = "Rebuilding preliminary reconstruction";  
1835          } else {          } else {
1836                  $Row->{status}->[0] = 0;                  $self->set_status(0,"Preliminary reconstruction");
                 $Row->{message}->[0] = "Preliminary reconstruction";  
1837          }          }
         #Updating the status table  
         $self->figmodel()->database()->save_table($ModelTable);  
         $self->figmodel()->database()->UnlockDBTable("MODELS");  
1838          #Sorting GenomeData by gene location on the chromosome          #Sorting GenomeData by gene location on the chromosome
1839            my $ftrTbl = $self->figmodel()->database()->get_table("ROLERXNMAPPING");
1840          $FeatureTable->sort_rows("MIN LOCATION");          $FeatureTable->sort_rows("MIN LOCATION");
1841          my ($ComplexHash,$SuggestedMappings,$UnrecognizedReactions,$ReactionHash);          my ($ComplexHash,$SuggestedMappings,$UnrecognizedReactions,$ReactionHash);
1842          my %LastGenePosition;          my %LastGenePosition;
# Line 1661  Line 1908 
1908                          }                          }
1909                  }                  }
1910          }          }
   
1911          #Creating nonadjacent complexes          #Creating nonadjacent complexes
1912          my @ReactionList = keys(%{$ReactionHash});          my @ReactionList = keys(%{$ReactionHash});
1913          foreach my $Reaction (@ReactionList) {          foreach my $Reaction (@ReactionList) {
# Line 1767  Line 2013 
2013                  $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"]});                  $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"]});
2014          }          }
2015    
2016            #Getting feature rows for features that are lumped
2017            my @rows = $ftrTbl->get_rows_by_key("2","MASTER");
2018            for (my $i=0; $i < @rows; $i++) {
2019                    my $rxn = $rows[$i]->{REACTION}->[0];
2020                    my $role = $rows[$i]->{ROLE}->[0];
2021                    my @orgrows = $FeatureTable->get_rows_by_key($role,"ROLES");
2022                    my $genes;
2023                    for (my $j=0; $j < @orgrows; $j++) {
2024                            if ($orgrows[$j]->{ID}->[0] =~ m/(peg\.\d+)/) {
2025                                    push(@{$genes},$1);
2026                            }
2027                    }
2028                    if (defined($genes) && @{$genes} > 0) {
2029                            my $row = $NewModelTable->get_row_by_key($rxn,"LOAD");
2030                            my $newGeneAssociations;
2031                            for (my $k=0; $k < @{$genes}; $k++) {
2032                                    for (my $j=0; $j < @{$row->{"ASSOCIATED PEG"}}; $j++) {
2033                                            my $peg = $genes->[$k];
2034                                            if ($row->{"ASSOCIATED PEG"}->[$j] !~ m/$peg/) {
2035                                                    push(@{$newGeneAssociations},$row->{"ASSOCIATED PEG"}->[$j]."+".$peg);
2036                                            }
2037                                    }
2038                            }
2039                            $row->{"ASSOCIATED PEG"} = $newGeneAssociations;
2040                    }
2041            }
2042    
2043          #Adding the spontaneous and universal reactions          #Adding the spontaneous and universal reactions
2044          foreach my $ReactionID (@{$self->config("spontaneous reactions")}) {          foreach my $ReactionID (@{$self->config("spontaneous reactions")}) {
2045                  #Getting the thermodynamic reversibility from the database                  #Getting the thermodynamic reversibility from the database
# Line 1785  Line 2058 
2058                  }                  }
2059          }          }
2060    
2061          #Checking if a biomass reaction already exists          #Creating biomass reaction for model
2062          my $BiomassReactionRow = $self->BuildSpecificBiomassReaction();          my $biomassID = $self->BuildSpecificBiomassReaction();
2063          if (!defined($BiomassReactionRow)) {          if ($biomassID !~ m/bio\d\d\d\d\d/) {
2064                  $self->set_status(-2,"Preliminary reconstruction failed: could not generate biomass reaction");                  $self->set_status(-2,"Preliminary reconstruction failed: could not generate biomass reaction");
2065                  $self->figmodel()->error_message("FIGMODEL:CreateModelReactionList: Could not generate biomass function for ".$self->id()."!");                  die $self->error_message("CreateSingleGenomeReactionList: Could not generate biomass reaction!");
2066                  return $self->fail();          }
2067            #Getting the biomass reaction PPO object
2068            my $bioRxn = $self->figmodel()->database()->get_object("bof",{id=>$biomassID});
2069            if (!defined($bioRxn)) {
2070                    die $self->error_message("CreateSingleGenomeReactionList: Could not find biomass reaction ".$biomassID."!");
2071            }
2072            #Getting the list of essential reactions for biomass reaction
2073            my $ReactionList;
2074            my $essentialReactions = $bioRxn->essentialRxn();
2075            if (defined($essentialReactions) && $essentialReactions =~ m/rxn\d\d\d\d\d/) {
2076                    push(@{$ReactionList},split(/\|/,$essentialReactions));
2077                    if ($essentialReactions !~ m/$biomassID/) {
2078                            push(@{$ReactionList},$biomassID);
2079                    }
2080          }          }
         my $ReactionList = $BiomassReactionRow->{"ESSENTIAL REACTIONS"};  
         push(@{$ReactionList},$BiomassReactionRow->{DATABASE}->[0]);  
2081    
2082          #Adding biomass reactions to the model table          #Adding biomass reactions to the model table
2083          foreach my $BOFReaction (@{$ReactionList}) {          foreach my $BOFReaction (@{$ReactionList}) {
# Line 1823  Line 2107 
2107                  }                  }
2108          }          }
2109    
2110            #If an original model exists, we copy the gap filling solution from that model
2111            if (defined($OriginalModelTable)) {
2112                    for (my $i=0; $i < $OriginalModelTable->size(); $i++) {
2113                            if ($OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] =~ m/GAP/ && $OriginalModelTable->get_row($i)->{"ASSOCIATED PEG"}->[0] ne "INITIAL GAP FILLING") {
2114                                    my $Row = $NewModelTable->get_row_by_key($OriginalModelTable->get_row($i)->{"LOAD"}->[0],"LOAD");
2115                                    if (!defined($Row)) {
2116                                            $NewModelTable->add_row($OriginalModelTable->get_row($i));
2117                                    }
2118                            }
2119                    }
2120            }
2121    
2122          #Now we compare the model to the previous model to determine if any differences exist that aren't gap filling reactions          #Now we compare the model to the previous model to determine if any differences exist that aren't gap filling reactions
2123          if (defined($OriginalModelTable)) {          if (defined($OriginalModelTable)) {
2124                  my $PerfectMatch = 1;                  my $PerfectMatch = 1;
# Line 1882  Line 2178 
2178          #Saving the new model to file          #Saving the new model to file
2179          $self->set_status(1,"Preliminary reconstruction complete");          $self->set_status(1,"Preliminary reconstruction complete");
2180          $self->figmodel()->database()->save_table($NewModelTable);          $self->figmodel()->database()->save_table($NewModelTable);
2181          $self->{_reaction_data} = $NewModelTable;          $self->reaction_table(1);
         #Clearing the previous model from the cache  
         $self->figmodel()->database()->ClearDBModel($self->id(),1);  
2182          #Updating the model stats table          #Updating the model stats table
2183          $self->update_stats_for_build();          $self->update_stats_for_build();
2184          $self->PrintSBMLFile();          $self->PrintSBMLFile();
2185            if (defined($OriginalModelTable)) {
2186                    $self->calculate_model_changes($OriginalModelTable,"REBUILD");
2187            }
2188    
2189          #Adding model to gapfilling queue          #Adding model to gapfilling queue
2190          if (defined($RunGapFilling) && $RunGapFilling == 1) {          if (defined($RunGapFilling) && $RunGapFilling == 1) {
# Line 1906  Line 2203 
2203    
2204  sub CreateMetaGenomeReactionList {  sub CreateMetaGenomeReactionList {
2205          my ($self) = @_;          my ($self) = @_;
   
2206          #Checking if the metagenome file exists          #Checking if the metagenome file exists
2207          if (!-e $self->config("raw MGRAST directory")->[0].$self->genome().".summary") {          if (!-e $self->config("raw MGRAST directory")->[0].$self->genome().".summary") {
2208                  $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());                  $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());
2209                    return $self->fail();
2210          }          }
2211          #Loading metagenome data          #Loading metagenome data
2212          my $MGRASTData = $self->figmodel()->database()->load_multiple_column_file($self->config("raw MGRAST directory")->[0].$self->genome().".summary","\t");          my $MGRASTData = $self->figmodel()->database()->load_multiple_column_file($self->config("raw MGRAST directory")->[0].$self->genome().".summary","\t");
2213          if (!defined($MGRASTData)) {          if (!defined($MGRASTData)) {
2214                  $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());                  $self->error_message("FIGMODEL:CreateMetaGenomeReactionList: could not find raw data file for metagenome ".$self->genome());
2215                    return $self->fail();
2216          }          }
   
2217          #Setting up needed variables          #Setting up needed variables
2218          my $OriginalModelTable = undef;          my $OriginalModelTable = undef;
   
2219          #Checking status          #Checking status
2220          if ($self->status() < 0) {          if ($self->status() < 0) {
2221                  $self->set_status(0,"Preliminary reconstruction");                  $self->set_status(0,"Preliminary reconstruction");
# Line 1931  Line 2227 
2227                  $self->ArchiveModel();                  $self->ArchiveModel();
2228                  $self->set_status(0,"Rebuilding preliminary reconstruction");                  $self->set_status(0,"Rebuilding preliminary reconstruction");
2229          }          }
2230            #Creating a hash of escores and pegs associated with each role
2231            my $rolePegHash;
2232            my $roleEscores;
2233            for (my $i=0; $i < @{$MGRASTData};$i++) {
2234                    #MD5,PEG,number of sims,role,sim e-scores,max escore,min escore,ave escore,stdev escore,ave exponent,stddev exponent
2235                    $rolePegHash->{$MGRASTData->[$i]->[3]}->{substr($MGRASTData->[$i]->[1],4)} = 1;
2236                    push(@{$roleEscores->{$MGRASTData->[$i]->[3]}},split(/;/,$MGRASTData->[$i]->[4]));
2237            }
2238          #Getting the reaction table          #Getting the reaction table
2239          my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");          my $ReactionTable = $self->figmodel()->database()->get_table("REACTIONS");
2240          #Creating model table          #Creating model table
2241          my $ModelTable = FIGMODELTable->new(["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG","SUBSYSTEM","CONFIDENCE","REFERENCE","NOTES"],$self->directory().$self->id().".txt",["LOAD"],";","|","REACTIONS\n");          my $ModelTable = $self->create_table_prototype("ModelReactions");
2242          for (my $i=0; $i < @{$MGRASTData};$i++) {          print $ModelTable->filename();
2243                  #MD5,PEG,number of sims,role,sim e-scores          my @roles = keys(%{$rolePegHash});
2244                  my $Role = $MGRASTData->[$i]->[3];          for (my $i=0; $i < @roles; $i++) {
2245                  my $MD5 = $MGRASTData->[$i]->[0];                  my $min = -1;
2246                  my $peg = $MGRASTData->[$i]->[1];                  my $max = -1;
2247                  my $sims = $MGRASTData->[$i]->[4];                  my $count = @{$roleEscores->{$roles[$i]}};
2248                  $sims =~ s/;/,/g;                  my $ave = 0;
2249                    my $stdev = 0;
2250                    my $aveexp = 0;
2251                    my $stdevexp = 0;
2252                    for (my $j=0; $j < @{$roleEscores->{$roles[$i]}}; $j++) {
2253                            if ($roleEscores->{$roles[$i]} < $min || $min == -1) {
2254                                    $min = $roleEscores->{$roles[$i]};
2255                            }
2256                            if ($roleEscores->{$roles[$i]} > $max || $max == -1) {
2257                                    $max = $roleEscores->{$roles[$i]};
2258                            }
2259                            $ave += $roleEscores->{$roles[$i]}->[$j];
2260                            if ($roleEscores->{$roles[$i]}->[$j] =~ m/e(-\d+$)/) {
2261                                    $aveexp += $1;
2262                            }
2263                    }
2264                    $ave = $ave/$count;
2265                    $aveexp = $aveexp/$count;
2266                    for (my $j=0; $j < @{$roleEscores->{$roles[$i]}}; $j++) {
2267                            $stdev += ($roleEscores->{$roles[$i]}->[$j]-$ave)*($roleEscores->{$roles[$i]}->[$j]-$ave);
2268                            if ($roleEscores->{$roles[$i]}->[$j] =~ m/e(-\d+$)/) {
2269                                    $stdevexp += ($1-$aveexp)*($1-$aveexp);
2270                            }
2271                    }
2272                    $stdev = sqrt($stdev/$count);
2273                    $stdevexp = sqrt($stdevexp/$count);
2274                  #Checking for subsystems                  #Checking for subsystems
2275                  my $GeneSubsystems = $self->figmodel()->subsystems_of_role($Role);                  my $GeneSubsystems = $self->figmodel()->subsystems_of_role($roles[$i]);
2276                  #Checking if there are reactions associated with this role                  #Checking if there are reactions associated with this role
2277                  my $ReactionHashArrayRef = $self->figmodel()->reactions_of_role($Role);                  my $ReactionHashArrayRef = $self->figmodel()->reactions_of_role($roles[$i]);
2278                  if ($ReactionHashArrayRef != 0) {                  if ($ReactionHashArrayRef != 0) {
2279                          foreach my $Mapping (@{$ReactionHashArrayRef}) {                          foreach my $Mapping (@{$ReactionHashArrayRef}) {
2280                                  if (defined($Mapping->{"REACTION"}) && defined($Mapping->{"MASTER"}) && defined($Mapping->{"SUBSYSTEM"}) && defined($Mapping->{"SOURCE"})) {                                  if (defined($Mapping->{"REACTION"}) && defined($Mapping->{"MASTER"}) && defined($Mapping->{"SUBSYSTEM"}) && defined($Mapping->{"SOURCE"})) {
# Line 1959  Line 2286 
2286                                                                  $ReactionRow = {"LOAD" => [$Mapping->{"REACTION"}->[0]],"DIRECTIONALITY" => [$self->figmodel()->reversibility_of_reaction($Mapping->{"REACTION"}->[0])],"COMPARTMENT" => ["c"]};                                                                  $ReactionRow = {"LOAD" => [$Mapping->{"REACTION"}->[0]],"DIRECTIONALITY" => [$self->figmodel()->reversibility_of_reaction($Mapping->{"REACTION"}->[0])],"COMPARTMENT" => ["c"]};
2287                                                                  $ModelTable->add_row($ReactionRow);                                                                  $ModelTable->add_row($ReactionRow);
2288                                                          }                                                          }
2289                                                          push(@{$ReactionRow->{"ASSOCIATED PEG"}},substr($peg,4));                                                          my %pegHash = %{$rolePegHash->{$roles[$i]}};
2290                                                          push(@{$ReactionRow->{"REFERENCE"}},$MD5.":".$Role);                                                          if (defined($ReactionRow->{"ASSOCIATED PEG"})) {
2291                                                          push(@{$ReactionRow->{"CONFIDENCE"}},$sims);                                                                  for (my $j=0; $j < @{$ReactionRow->{"ASSOCIATED PEG"}}; $j++) {
2292                                                                            $pegHash{$ReactionRow->{"ASSOCIATED PEG"}->[$j]} = 1;
2293                                                                    }
2294                                                            }
2295                                                            delete $ReactionRow->{"ASSOCIATED PEG"};
2296                                                            push(@{$ReactionRow->{"ASSOCIATED PEG"}},keys(%pegHash));
2297                                                            push(@{$ReactionRow->{"REFERENCE"}},$count.":".$ave.":".$stdev.":".$aveexp.":".$stdevexp.":".$min.":".$max);
2298                                                          if (defined($GeneSubsystems)) {                                                          if (defined($GeneSubsystems)) {
2299                                                                  push(@{$ReactionRow->{"SUBSYSTEM"}},@{$GeneSubsystems});                                                                  push(@{$ReactionRow->{"SUBSYSTEM"}},@{$GeneSubsystems});
2300                                                          }                                                          }
# Line 2387  Line 2720 
2720          Calculating growth in the input media          Calculating growth in the input media
2721  =cut  =cut
2722  sub calculate_growth {  sub calculate_growth {
2723          my ($self,$Media) = @_;          my ($self,$Media,$outputDirectory,$InParameters,$saveLPFile) = @_;
2724            #Setting the Media
2725            if (!defined($Media) || length($Media) == 0) {
2726                    $Media = $self->autocompleteMedia();
2727            }
2728            #Setting parameters for the run
2729            my $DefaultParameters = $self->figmodel()->defaultParameters();
2730            if (defined($InParameters)) {
2731                    my @parameters = keys(%{$InParameters});
2732                    for (my $i=0; $i < @parameters; $i++) {
2733                            $DefaultParameters->{$parameters[$i]} = $InParameters->{$parameters[$i]};
2734                    }
2735            }
2736            $DefaultParameters->{"optimize metabolite production if objective is zero"} = 1;
2737            #Setting filenames
2738          my $UniqueFilename = $self->figmodel()->filename();          my $UniqueFilename = $self->figmodel()->filename();
2739          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()));          if (!defined($outputDirectory)) {
2740                    $outputDirectory = $self->config("database message file directory")->[0];
2741            }
2742            my $fluxFilename = $outputDirectory."Fluxes-".$self->id()."-".$Media.".txt";
2743            my $cpdFluxFilename = $outputDirectory."CompoundFluxes-".$self->id()."-".$Media.".txt";
2744            #Running FBA
2745            #print $self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$Media,["ProductionMFA"],$DefaultParameters,$self->id()."-".$Media."-GrowthTest.txt",undef,$self->selected_version())."\n";
2746            system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$Media,["ProductionMFA"],$DefaultParameters,$self->id()."-".$Media."-GrowthTest.txt",undef,$self->selected_version()));
2747            #Saving LP file if requested
2748            if (defined($saveLPFile) && $saveLPFile == 1 && -e $self->figmodel()->{"MFAToolkit output directory"}->[0].$UniqueFilename."/CurrentProblem.lp") {
2749                    system("cp ".$self->figmodel()->config("MFAToolkit output directory")->[0].$UniqueFilename."/CurrentProblem.lp ".$self->directory().$self->id().".lp");
2750            }
2751          my $ProblemReport = $self->figmodel()->LoadProblemReport($UniqueFilename);          my $ProblemReport = $self->figmodel()->LoadProblemReport($UniqueFilename);
2752          my $Result;          my $Result;
2753          if (defined($ProblemReport)) {          if (defined($ProblemReport)) {
2754                  my $Row = $ProblemReport->get_row(0);                  my $Row = $ProblemReport->get_row(0);
2755                  if (defined($Row) && defined($Row->{"Objective"}->[0])) {                  if (defined($Row) && defined($Row->{"Objective"}->[0])) {
2756                          if ($Row->{"Objective"}->[0] < 0.00000001) {                          if ($Row->{"Objective"}->[0] < 0.00000001 || $Row->{"Objective"}->[0] == 1e7) {
2757                                  $Result = "NOGROWTH:".$Row->{"Individual metabolites with zero production"}->[0];                                  $Result = "NOGROWTH";
2758                                    if (defined($Row->{"Individual metabolites with zero production"}->[0]) && $Row->{"Individual metabolites with zero production"}->[0] =~ m/cpd\d\d\d\d\d/) {
2759                                            $Result .= ":".$Row->{"Individual metabolites with zero production"}->[0];
2760                                    }
2761                          } else {                          } else {
2762                                    if (-e $self->figmodel()->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/SolutionReactionData0.txt") {
2763                                            system("cp ".$self->figmodel()->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/SolutionReactionData0.txt ".$fluxFilename);
2764                                            system("cp ".$self->figmodel()->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/SolutionCompoundData0.txt ".$cpdFluxFilename);
2765                                    }
2766                                  $Result = $Row->{"Objective"}->[0];                                  $Result = $Row->{"Objective"}->[0];
2767                          }                          }
2768                  }                  }
2769          }          }
2770            #Deleting files if necessary
2771            if ($self->figmodel()->config("preserve all log files")->[0] ne "yes") {
2772                    $self->figmodel()->cleardirectory($UniqueFilename);
2773                    unlink($self->figmodel()->config("database message file directory")->[0].$self->id()."-".$Media."-GrowthTest.txt");
2774            }
2775            #Returning result
2776          return $Result;          return $Result;
2777  }  }
2778    
# Line 3413  Line 3784 
3784          return 1;          return 1;
3785  }  }
3786    
3787  =head3 BuildSpecificBiomassReaction  =head3 DetermineCofactorLipidCellWallComponents
3788  Definition:  Definition:
3789          FIGMODELmodel->BuildSpecificBiomassReaction();          {cofactor=>{string:compound id=>float:coefficient},lipid=>...cellWall=>} = FIGMODELmodel->DetermineCofactorLipidCellWallComponents();
3790  Description:  Description:
3791  =cut  =cut
3792  sub BuildSpecificBiomassReaction {  sub DetermineCofactorLipidCellWallComponents {
3793          my ($self) = @_;          my ($self) = @_;
3794            my $templateResults;
         my $biomassrxn;  
         my $OrganismID = $self->genome();  
         #Checking for a biomass override  
         if (defined($self->config("biomass reaction override")->{$OrganismID})) {  
                 my $biomassID = $self->config("biomass reaction override")->{$OrganismID};  
                 my $tbl = $self->database()->get_table("BIOMASS",1);  
                 $biomassrxn = $tbl->get_row_by_key($biomassID,"DATABASE");  
                 print "Overriding biomass template and selecting ".$biomassID." for ".$OrganismID.".\n";  
         } else {#Creating biomass reaction from the template  
                 #Getting the genome stats  
3795                  my $genomestats = $self->figmodel()->get_genome_stats($self->genome());                  my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
3796                  my $Class = $genomestats->{CLASS}->[0];          my $Class = $self->{_data}->cellwalltype();
3797                  my $Name = $genomestats->{NAME}->[0];          my $Name = $self->name();
3798            my $translation = {COFACTOR=>"cofactor",LIPIDS=>"lipid","CELL WALL"=>"cellWall"};
3799                  #Checking for phoenix variants                  #Checking for phoenix variants
3800                  my $PhoenixVariantTable = $self->figmodel()->database()->GetDBTable("Phoenix variants table");                  my $PhoenixVariantTable = $self->figmodel()->database()->GetDBTable("Phoenix variants table");
3801                  my $Phoenix = 0;                  my $Phoenix = 0;
3802                  my @Rows = $PhoenixVariantTable->get_rows_by_key($OrganismID,"GENOME");          my @Rows = $PhoenixVariantTable->get_rows_by_key($self->genome(),"GENOME");
3803                  my $VariantHash;                  my $VariantHash;
3804                  for (my $i=0; $i < @Rows; $i++) {                  for (my $i=0; $i < @Rows; $i++) {
3805                          $Phoenix = 1;                          $Phoenix = 1;
# Line 3446  Line 3807 
3807                                  $VariantHash->{$Rows[$i]->{"SUBSYSTEM"}->[0]} = $Rows[$i]->{"VARIANT"}->[0];                                  $VariantHash->{$Rows[$i]->{"SUBSYSTEM"}->[0]} = $Rows[$i]->{"VARIANT"}->[0];
3808                          }                          }
3809                  }                  }
   
3810                  #Collecting genome data                  #Collecting genome data
3811                  my $RoleHash;                  my $RoleHash;
3812                  my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());                  my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
# Line 3463  Line 3823 
3823                                  }                                  }
3824                          }                          }
3825                  }                  }
   
3826                  #Scanning through the template item by item and determinine which biomass components should be added                  #Scanning through the template item by item and determinine which biomass components should be added
3827                  my $ComponentTypes;          my $includedHash;
3828                  my $EquationData;          my $BiomassReactionTemplateTable = $self->figmodel()->database()->get_table("BIOMASSTEMPLATE");
                 my $BiomassReactionTemplateTable = $self->figmodel()->database()->GetDBTable("BIOMASS TEMPLATE");  
3829                  for (my $i=0; $i < $BiomassReactionTemplateTable->size(); $i++) {                  for (my $i=0; $i < $BiomassReactionTemplateTable->size(); $i++) {
3830                          my $Row = $BiomassReactionTemplateTable->get_row($i);                          my $Row = $BiomassReactionTemplateTable->get_row($i);
3831                    if (defined($translation->{$Row->{CLASS}->[0]})) {
3832                            my $coef = -1;
3833                            if ($Row->{"REACTANT"}->[0] eq "NO") {
3834                                    $coef = 1;
3835                                    if ($Row->{"COEFFICIENT"}->[0] =~ m/cpd/) {
3836                                            $coef = $Row->{"COEFFICIENT"}->[0];
3837                                    }
3838                            }
3839                          if (defined($Row->{"INCLUSION CRITERIA"}->[0]) && $Row->{"INCLUSION CRITERIA"}->[0] eq "UNIVERSAL") {                          if (defined($Row->{"INCLUSION CRITERIA"}->[0]) && $Row->{"INCLUSION CRITERIA"}->[0] eq "UNIVERSAL") {
3840                                  push(@{$EquationData},$Row);                                  $includedHash->{$Row->{"ID"}->[0]} = 1;
3841                                  $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;                                  $templateResults->{$translation->{$Row->{CLASS}->[0]}}->{$Row->{"ID"}->[0]} = $coef;
3842                          } else {                          } elsif (defined($Row->{"INCLUSION CRITERIA"}->[0])) {
3843                                  my $Criteria = $Row->{"INCLUSION CRITERIA"}->[0];                                  my $Criteria = $Row->{"INCLUSION CRITERIA"}->[0];
3844                                  my $End = 0;                                  my $End = 0;
3845                                  while ($End == 0) {                                  while ($End == 0) {
3846                                          if ($Criteria =~ m/^(.+)(AND)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(AND)\{([^{^}]+)\}$/ || $Criteria =~ m/^(.+)(OR)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(OR)\{([^{^}]+)\}$/) {                                          if ($Criteria =~ m/^(.+)(AND)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(AND)\{([^{^}]+)\}$/ || $Criteria =~ m/^(.+)(OR)\{([^{^}]+)\}(.+)$/ || $Criteria =~ m/^(OR)\{([^{^}]+)\}$/) {
3847                                                  print $Criteria."\n";                                                  print $Criteria." : ";
3848                                                  my $Start = "";                                                  my $Start = "";
3849                                                  my $End = "";                                                  my $End = "";
3850                                                  my $Condition = $1;                                                  my $Condition = $1;
# Line 3502  Line 3868 
3868                                                                  $Result = "NO";                                                                  $Result = "NO";
3869                                                                  last;                                                                  last;
3870                                                          } elsif ($Array[$j] =~ m/^COMPOUND:(.+)/) {                                                          } elsif ($Array[$j] =~ m/^COMPOUND:(.+)/) {
3871                                                                  my $Match = 0;                                                                  if (defined($includedHash->{$1}) && $Condition eq "OR") {
                                                                 for (my $k=0; $k < @{$EquationData}; $k++) {  
                                                                         if ($EquationData->[$k]->{"ID"}->[0] eq $1) {  
                                                                                 $Match = 1;  
                                                                                 last;  
                                                                         }  
                                                                 }  
                                                                 if ($Match == 1 && $Condition eq "OR") {  
3872                                                                          $Result = "YES";                                                                          $Result = "YES";
3873                                                                          last;                                                                          last;
3874                                                                  } elsif ($Match != 1 && $Condition eq "AND") {                                                                  } elsif (!defined($includedHash->{$1}) && $Condition eq "AND") {
3875                                                                          $Result = "NO";                                                                          $Result = "NO";
3876                                                                          last;                                                                          last;
3877                                                                  }                                                                  }
3878                                                          } elsif ($Array[$j] =~ m/^!COMPOUND:(.+)/) {                                                          } elsif ($Array[$j] =~ m/^!COMPOUND:(.+)/) {
3879                                                                  my $Match = 0;                                                                  if (!defined($includedHash->{$1}) && $Condition eq "OR") {
                                                                 for (my $k=0; $k < @{$EquationData}; $k++) {  
                                                                         if ($EquationData->[$k]->{"ID"}->[0] eq $1) {  
                                                                                 $Match = 1;  
                                                                                 last;  
                                                                         }  
                                                                 }  
                                                                 if ($Match != 1 && $Condition eq "OR") {  
3880                                                                          $Result = "YES";                                                                          $Result = "YES";
3881                                                                          last;                                                                          last;
3882                                                                  } elsif ($Match == 1 && $Condition eq "AND") {                                                                  } elsif (defined($includedHash->{$1}) && $Condition eq "AND") {
3883                                                                          $Result = "NO";                                                                          $Result = "NO";
3884                                                                          last;                                                                          last;
3885                                                                  }                                                                  }
# Line 3643  Line 3995 
3995                                          }                                          }
3996                                  }                                  }
3997                                  if ($Criteria eq "YES") {                                  if ($Criteria eq "YES") {
3998                                          push(@{$EquationData},$Row);                                          $templateResults->{$translation->{$Row->{CLASS}->[0]}}->{$Row->{"ID"}->[0]} = $coef;
3999                                          $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;                                          $includedHash->{$Row->{"ID"}->[0]} = 1;
4000                                  }                                  }
4001                          }                          }
4002                  }                  }
4003            }
4004                  #Building biomass equation          my $types = ["cofactor","lipid","cellWall"];
4005                  my %Reactants;          my $cpdMgr = $self->figmodel()->database()->get_object_manager("compound");
4006                  my %Products;          for (my $i=0; $i < @{$types}; $i++) {
4007                  foreach my $EquationRow (@{$EquationData}) {                  my @list =      keys(%{$templateResults->{$types->[$i]}});
4008                          #First determine what the coefficient should be                  my $entries = 0;
4009                          my $Coefficient;                  for (my $j=0; $j < @list; $j++) {
4010                          if (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^[0-9\.]+$/) {                          if ($templateResults->{$types->[$i]}->{$list[$j]} eq "-1") {
4011                                  $Coefficient = $EquationRow->{"COEFFICIENT"}->[0];                                  my $objs = $cpdMgr->get_objects({id=>$list[$j]});
4012                          } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/cpd\d\d\d\d\d/) {                                  if (!defined($objs->[0]) || $objs->[0]->mass() == 0) {
4013                                  $Coefficient = 0;                                          $templateResults->{$types->[$i]}->{$list[$j]} = -1e-5;
                                 my @CompoundList = split(/,/,$EquationRow->{"COEFFICIENT"}->[0]);  
                                 foreach my $Compound (@CompoundList) {  
                                         if (defined($Reactants{$Compound})) {  
                                                 $Coefficient += $Reactants{$Compound};  
                                         }  
                                 }  
                         } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^([0-9\.]+)\/(.+)$/) {  
                                 my @Keys = keys(%{$ComponentTypes->{$2}});  
                                 my $MW = 1;  
                                 my $CompoundData = $self->figmodel()->LoadObject($EquationRow->{"ID"}->[0]);  
                                 if (defined($CompoundData->{"MASS"})) {  
                                         $MW = $CompoundData->{"MASS"}->[0];  
                                 }  
                                 $Coefficient = $1/@Keys/$MW;  
                         }  
                         if (defined($EquationRow->{"REACTANT"}) && $EquationRow->{"REACTANT"}->[0] eq "YES") {  
                                 if (defined($Reactants{$EquationRow->{"ID"}->[0]})) {  
                                         $Reactants{$EquationRow->{"ID"}->[0]} += $Coefficient;  
                                 } elsif (defined($Products{$EquationRow->{"ID"}->[0]}) && $Products{$EquationRow->{"ID"}->[0]} > $Coefficient) {  
                                         $Products{$EquationRow->{"ID"}->[0]} -= $Coefficient;  
                                 } elsif (defined($Products{$EquationRow->{"ID"}->[0]}) && $Products{$EquationRow->{"ID"}->[0]} < $Coefficient) {  
                                         $Reactants{$EquationRow->{"ID"}->[0]} = $Coefficient - $Products{$EquationRow->{"ID"}->[0]};  
                                         delete $Products{$EquationRow->{"ID"}->[0]};  
                                 } else {  
                                         $Reactants{$EquationRow->{"ID"}->[0]} = $Coefficient;  
                                 }  
                         } else {  
                                 if (defined($Products{$EquationRow->{"ID"}->[0]})) {  
                                         $Products{$EquationRow->{"ID"}->[0]} += $Coefficient;  
                                 } elsif (defined($Reactants{$EquationRow->{"ID"}->[0]}) && $Reactants{$EquationRow->{"ID"}->[0]} > $Coefficient) {  
                                         $Reactants{$EquationRow->{"ID"}->[0]} -= $Coefficient;  
                                 } elsif (defined($Reactants{$EquationRow->{"ID"}->[0]}) && $Reactants{$EquationRow->{"ID"}->[0]} < $Coefficient) {  
                                         $Products{$EquationRow->{"ID"}->[0]} = $Coefficient - $Reactants{$EquationRow->{"ID"}->[0]};  
                                         delete $Reactants{$EquationRow->{"ID"}->[0]};  
4014                                  } else {                                  } else {
4015                                          $Products{$EquationRow->{"ID"}->[0]} = $Coefficient;                                          $entries++;
4016                                    }
4017                            }
4018                    }
4019                    for (my $j=0; $j < @list; $j++) {
4020                            if ($templateResults->{$types->[$i]}->{$list[$j]} eq "-1") {
4021                                    $templateResults->{$types->[$i]}->{$list[$j]} = -1/$entries;
4022                            } elsif ($templateResults->{$types->[$i]}->{$list[$j]} =~ m/cpd/) {
4023                                    my $netCoef = 0;
4024                                    my @allcpd = split(/,/,$templateResults->{$types->[$i]}->{$list[$j]});
4025                                    for (my $k=0; $k < @allcpd; $k++) {
4026                                            if (defined($templateResults->{$types->[$i]}->{$allcpd[$k]}) && $templateResults->{$types->[$i]}->{$allcpd[$k]} ne "-1e-5") {
4027                                                    $netCoef += (1/$entries);
4028                                            } elsif (defined($templateResults->{$types->[$i]}->{$allcpd[$k]}) && $templateResults->{$types->[$i]}->{$allcpd[$k]} eq "-1e-5") {
4029                                                    $netCoef += 1e-5;
4030                                            }
4031                                    }
4032                                    $templateResults->{$types->[$i]}->{$list[$j]} = $netCoef;
4033                            }
4034                    }
4035            }
4036            return $templateResults;
4037    }
4038    
4039    =head3 BuildSpecificBiomassReaction
4040    Definition:
4041            FIGMODELmodel->BuildSpecificBiomassReaction();
4042    Description:
4043    =cut
4044    sub BuildSpecificBiomassReaction {
4045            my ($self) = @_;
4046            #Getting the database handle for biomass reactions
4047            my $bioMgr = $self->figmodel()->database()->get_object_manager("bof");
4048            #Checking if the current biomass reaction appears in more than on model, if not, this biomass reaction is conserved for this model
4049            my $biomassID = $self->biomassReaction();
4050            if ($biomassID =~ m/bio\d\d\d\d\d/) {
4051                    my $mdlMgr = $self->figmodel()->database()->get_object_manager("model");
4052                    my $mdlObs = $mdlMgr->get_objects({biomassReaction=>$biomassID});
4053                    if (defined($mdlObs->[1])) {
4054                            $biomassID = "NONE";
4055                    }
4056            }
4057            #If the biomass ID is "NONE", then we create a new biomass reaction for the model
4058            my $bioObj;
4059            my $originalPackages = "";
4060            my $originalEssReactions = "";
4061            if ($biomassID eq "NONE") {
4062                    #Getting the current largest ID
4063                    $biomassID = $self->figmodel()->database()->check_out_new_id("bof");
4064                    $bioObj = $bioMgr->create({
4065                            id=>$biomassID,owner=>$self->owner(),users=>$self->users(),name=>"Biomass",equation=>"NONE",protein=>"0",
4066                            energy=>"0",DNA=>"0",RNA=>"0",lipid=>"0",cellWall=>"0",cofactor=>"0",
4067                            modificationDate=>time(),creationDate=>time(),
4068                            cofactorPackage=>"NONE",lipidPackage=>"NONE",cellWallPackage=>"NONE",
4069                            DNACoef=>"NONE",RNACoef=>"NONE",proteinCoef=>"NONE",lipidCoef=>"NONE",
4070                            cellWallCoef=>"NONE",cofactorCoef=>"NONE",essentialRxn=>"NONE"});
4071                    if (!defined($bioObj)) {
4072                            die $self->error_message("BuildSpecificBiomassReaction():Could not create new biomass reaction ".$biomassID."!");
4073                    }
4074            } else {
4075                    #Getting the biomass DB handler from the database
4076                    my $objs = $bioMgr->get_objects({id=>$biomassID});
4077                    if (!defined($objs->[0])) {
4078                            die $self->error_message("BuildSpecificBiomassReaction():Could not find biomass reaction ".$biomassID." in database!");
4079                    }
4080                    $bioObj = $objs->[0];
4081                    $bioObj->owner($self->owner());
4082                    $bioObj->users($self->users());
4083                    if (defined($bioObj->essentialRxn())) {
4084                            $originalEssReactions = $bioObj->essentialRxn();
4085                            $originalPackages = $bioObj->cofactorPackage().$bioObj->lipidPackage().$bioObj->cellWallPackage();
4086                                  }                                  }
4087                          }                          }
4088            #Getting genome stats
4089            my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
4090            my $Class = $self->{_data}->cellwalltype();
4091            #Setting global coefficients based on cell wall type
4092            my $biomassCompounds;
4093            my $compounds;
4094            if ($Class eq "Gram positive") {
4095                    $compounds->{RNA} = {cpd00002=>-0.262,cpd00012=>1,cpd00038=>-0.323,cpd00052=>-0.199,cpd00062=>-0.215};
4096                    $compounds->{protein} = {cpd00001=>1,cpd00023=>-0.0637,cpd00033=>-0.0999,cpd00035=>-0.0653,cpd00039=>-0.0790,cpd00041=>-0.0362,cpd00051=>-0.0472,cpd00053=>-0.0637,cpd00054=>-0.0529,cpd00060=>-0.0277,cpd00065=>-0.0133,cpd00066=>-0.0430,cpd00069=>-0.0271,cpd00084=>-0.0139,cpd00107=>-0.0848,cpd00119=>-0.0200,cpd00129=>-0.0393,cpd00132=>-0.0362,cpd00156=>-0.0751,cpd00161=>-0.0456,cpd00322=>-0.0660};
4097                    $bioObj->protein("0.5284");
4098                    $bioObj->DNA("0.026");
4099                    $bioObj->RNA("0.0655");
4100                    $bioObj->lipid("0.075");
4101                    $bioObj->cellWall("0.25");
4102                    $bioObj->cofactor("0.10");
4103            } else {
4104                    $compounds->{RNA} = {cpd00002=>-0.262,cpd00012=>1,cpd00038=>-0.322,cpd00052=>-0.2,cpd00062=>-0.216};
4105                    $compounds->{protein} = {cpd00001=>1,cpd00023=>-0.0492,cpd00033=>-0.1145,cpd00035=>-0.0961,cpd00039=>-0.0641,cpd00041=>-0.0451,cpd00051=>-0.0554,cpd00053=>-0.0492,cpd00054=>-0.0403,cpd00060=>-0.0287,cpd00065=>-0.0106,cpd00066=>-0.0347,cpd00069=>-0.0258,cpd00084=>-0.0171,cpd00107=>-0.0843,cpd00119=>-0.0178,cpd00129=>-0.0414,cpd00132=>-0.0451,cpd00156=>-0.0791,cpd00161=>-0.0474,cpd00322=>-0.0543};
4106                    $bioObj->protein("0.563");
4107                    $bioObj->DNA("0.031");
4108                    $bioObj->RNA("0.21");
4109                    $bioObj->lipid("0.093");
4110                    $bioObj->cellWall("0.177");
4111                    $bioObj->cofactor("0.039");
4112            }
4113            #Setting energy coefficient for all reactions
4114            $bioObj->energy("40");
4115            $compounds->{energy} = {cpd00002=>-1,cpd00001=>-1,cpd00008=>1,cpd00009=>1,cpd00067=>1};
4116            #Setting DNA coefficients based on GC content
4117            my $gc = $self->figmodel()->get_genome_gc_content($self->genome());
4118            $compounds->{DNA} = {cpd00012=>1,cpd00115=>0.5*(1-$gc),cpd00241=>0.5*$gc,cpd00356=>0.5*$gc,cpd00357=>0.5*(1-$gc)};
4119            #Setting Lipid,cell wall,and cofactor coefficients based on biomass template
4120            my $templateResults = $self->DetermineCofactorLipidCellWallComponents();
4121            $compounds->{cofactor} = $templateResults->{cofactor};
4122            $compounds->{lipid} = $templateResults->{lipid};
4123            $compounds->{cellWall} = $templateResults->{cellWall};
4124            #Getting package number for cofactor, lipid, and cell wall
4125            my $packages;
4126            my $cpdgrpMgr = $self->figmodel()->database()->get_object_manager("cpdgrp");
4127            my $packageTypes = ["Cofactor","Lipid","CellWall"];
4128            my $translation = {"Cofactor"=>"cofactor","Lipid"=>"lipid","CellWall"=>"cellWall"};
4129            for (my $i=0; $i < @{$packageTypes}; $i++) {
4130                    my @cpdList = keys(%{$compounds->{$translation->{$packageTypes->[$i]}}});
4131                    my $function = $translation->{$packageTypes->[$i]}."Package";
4132                    if (@cpdList == 0) {
4133                            $bioObj->$function("NONE");
4134                    } else {
4135                            my $cpdgrpObs = $cpdgrpMgr->get_objects({type=>$packageTypes->[$i]."Package"});
4136                            for (my $j=0; $j < @{$cpdgrpObs}; $j++) {
4137                                    $packages->{$packageTypes->[$i]}->{$cpdgrpObs->[$j]->grouping()}->{$cpdgrpObs->[$j]->COMPOUND()} = 1;
4138                            }
4139                            my @packageList = keys(%{$packages->{$packageTypes->[$i]}});
4140                            my $packageHash;
4141                            for (my $j=0; $j < @packageList; $j++) {
4142                                    $packageHash->{join("|",sort(keys(%{$packages->{$packageTypes->[$i]}->{$packageList[$j]}})))} = $packageList[$j];
4143                            }
4144                            if (defined($packageHash->{join("|",sort(keys(%{$compounds->{$translation->{$packageTypes->[$i]}}})))})) {
4145                                    $bioObj->$function($packageHash->{join("|",sort(keys(%{$compounds->{$translation->{$packageTypes->[$i]}}})))});
4146                            } else {
4147                                    my $newPackageID = $self->figmodel()->database()->check_out_new_id($packageTypes->[$i]."Pkg");
4148                                    $bioObj->$function($newPackageID);
4149                                    my @cpdList = keys(%{$compounds->{$translation->{$packageTypes->[$i]}}});
4150                                    for (my $j=0; $j < @cpdList; $j++) {
4151                                            $cpdgrpMgr = $self->figmodel()->database()->get_object_manager("cpdgrp");
4152                                            $cpdgrpMgr->create({COMPOUND=>$cpdList[$j],grouping=>$newPackageID,type=>$packageTypes->[$i]."Package"});
4153                                    }
4154                            }
4155                    }
4156            }
4157            #Filling in coefficient terms in database and calculating global reaction coefficients based on classification abundancies
4158            my $equationCompounds;
4159            my $types = ["RNA","DNA","protein","lipid","cellWall","cofactor","energy"];
4160            my $cpdMgr = $self->figmodel()->database()->get_object_manager("compound");
4161            for (my $i=0; $i < @{$types}; $i++) {
4162                    my $coefString = "";
4163                    my @compounds = sort(keys(%{$compounds->{$types->[$i]}}));
4164                    #Building coefficient strings and determining net mass for component types
4165                    my $netMass = 0;
4166                    for (my $j=0; $j < @compounds; $j++) {
4167                            my $objs = $cpdMgr->get_objects({id=>$compounds[$j]});
4168                            my $mass = 0;
4169                            if (defined($objs->[0]) && $objs->[0]->mass() != 0) {
4170                                    $mass = $objs->[0]->mass();
4171                                    $netMass += -$compounds->{$types->[$i]}->{$compounds[$j]}*$objs->[0]->mass();
4172                            }
4173                            if (!defined($equationCompounds->{$compounds[$j]})) {
4174                                    $equationCompounds->{$compounds[$j]}->{"coef"} = 0;
4175                                    $equationCompounds->{$compounds[$j]}->{"type"} = $types->[$i];
4176                                    $equationCompounds->{$compounds[$j]}->{"mass"} = $mass;
4177                            }
4178                            $coefString .= $compounds->{$types->[$i]}->{$compounds[$j]}."|";
4179                    }
4180                    $netMass = 0.001*$netMass;
4181                    #Calculating coefficients for all component compounds
4182                    for (my $j=0; $j < @compounds; $j++) {
4183                            #Normalizing certain type coefficients by mass
4184                            my $function = $types->[$i];
4185                            my $fraction = $bioObj->$function();
4186                            if ($types->[$i] ne "energy") {
4187                                    $fraction = $fraction/$netMass;
4188                            }
4189                            if ($compounds->{$types->[$i]}->{$compounds[$j]} eq 1e-5) {
4190                                    $fraction = 1;
4191                            }
4192                            $equationCompounds->{$compounds[$j]}->{"coef"} += $fraction*$compounds->{$types->[$i]}->{$compounds[$j]};
4193                    }
4194                    chop($coefString);
4195                    if (length($coefString) == 0) {
4196                            $coefString = "NONE";
4197                    }
4198                    my $function = $types->[$i]."Coef";
4199                    if ($types->[$i] ne "energy") {
4200                            $bioObj->$function($coefString);
4201                    }
4202            }
4203            #Adding biomass to compound list
4204            $equationCompounds->{cpd11416}->{coef} = 1;
4205            $equationCompounds->{cpd11416}->{type} = "macromolecule";
4206            #Building equation from hash and populating compound biomass table
4207            my @compoundList = keys(%{$equationCompounds});
4208            my ($reactants,$products);
4209            #Deleting existing Biomass Compound info
4210            my $cpdbofMgr = $self->figmodel()->database()->get_object_manager("cpdbof");
4211            my $matchingObjs = $cpdbofMgr->get_objects({BIOMASS=>$biomassID});
4212            for (my $i=0; $i < @{$matchingObjs}; $i++) {
4213                    $matchingObjs->[$i]->delete();
4214            }
4215            my $typeCategories = {"macromolecule"=>"M","RNA"=>"R","DNA"=>"D","protein"=>"P","lipid"=>"L","cellWall"=>"W","cofactor"=>"C","energy"=>"E"};
4216            my $productmass = 0;
4217            my $reactantmass = 0;
4218            my $totalmass = 0;
4219            foreach my $compound (@compoundList) {
4220                    if (defined($equationCompounds->{$compound}->{coef}) && defined($equationCompounds->{$compound}->{mass})) {
4221                            $totalmass += $equationCompounds->{$compound}->{coef}*0.001*$equationCompounds->{$compound}->{mass};
4222                    }
4223                    if ($equationCompounds->{$compound}->{coef} < 0) {
4224                            if (defined($equationCompounds->{$compound}->{coef}) && defined($equationCompounds->{$compound}->{mass})) {
4225                                    $reactantmass += $equationCompounds->{$compound}->{coef}*0.001*$equationCompounds->{$compound}->{mass};
4226                            }
4227                            $reactants->{$compound} = $self->figmodel()->format_coefficient(-1*$equationCompounds->{$compound}->{coef});
4228                    } else {
4229                            if (defined($equationCompounds->{$compound}->{coef}) && defined($equationCompounds->{$compound}->{mass})) {
4230                                    $productmass += $equationCompounds->{$compound}->{coef}*0.001*$equationCompounds->{$compound}->{mass};
4231                            }
4232                            $products->{$compound} = $self->figmodel()->format_coefficient($equationCompounds->{$compound}->{coef});
4233                    }
4234                    #Adding biomass reaction compounds to the biomass compound table
4235                    $cpdbofMgr = $self->figmodel()->database()->get_object_manager("cpdbof");
4236                    $cpdbofMgr->create({COMPOUND=>$compound,BIOMASS=>$biomassID,coefficient=>$equationCompounds->{$compound}->{coef},compartment=>"c",category=>$typeCategories->{$equationCompounds->{$compound}->{type}}});
4237                  }                  }
4238            print "Total mass = ".$totalmass.", Reactant mass = ".$reactantmass.", Product mass = ".$productmass."\n";
4239                  my $Equation = "";                  my $Equation = "";
4240                  my @ReactantList = sort(keys(%Reactants));          my @ReactantList = sort(keys(%{$reactants}));
4241                  for (my $i=0; $i < @ReactantList; $i++) {                  for (my $i=0; $i < @ReactantList; $i++) {
4242                          if (length($Equation) > 0) {                          if (length($Equation) > 0) {
4243                                  $Equation .= " + ";                                  $Equation .= " + ";
4244                          }                          }
4245                          $Equation .= $self->figmodel()->format_coefficient($Reactants{$ReactantList[$i]})." ".$ReactantList[$i];                  $Equation .= "(".$reactants->{$ReactantList[$i]}.") ".$ReactantList[$i];
4246                  }                  }
4247                  $Equation .= " => ";                  $Equation .= " => ";
4248                  my $First = 1;                  my $First = 1;
4249                  @ReactantList = sort(keys(%Products));          @ReactantList = sort(keys(%{$products}));
4250                  for (my $i=0; $i < @ReactantList; $i++) {                  for (my $i=0; $i < @ReactantList; $i++) {
4251                          if ($First == 0) {                          if ($First == 0) {
4252                                  $Equation .= " + ";                                  $Equation .= " + ";
4253                          }                          }
4254                          $First = 0;                          $First = 0;
4255                          $Equation .= $self->figmodel()->format_coefficient($Products{$ReactantList[$i]})." ".$ReactantList[$i];                  $Equation .= "(".$products->{$ReactantList[$i]}.") ".$ReactantList[$i];
4256            }
4257            $bioObj->equation($Equation);
4258            #Setting the biomass reaction of this model
4259            $self->biomassReaction($biomassID);
4260            $self->figmodel()->print_biomass_reaction_file($biomassID);
4261            #Checking if the biomass reaction remained unchanged
4262            if ($originalPackages ne "" && $originalPackages eq $bioObj->cofactorPackage().$bioObj->lipidPackage().$bioObj->cellWallPackage()) {
4263                    print "UNCHANGED!\n";
4264                    $bioObj->essentialRxn($originalEssReactions);
4265            } else {
4266                    #Copying essential reaction lists if the packages in this biomasses reaction exactly match those in another biomass reaction
4267                    my $matches = $bioMgr->get_objects({cofactorPackage=>$bioObj->cofactorPackage(),lipidPackage=>$bioObj->lipidPackage(),cellWallPackage=>$bioObj->cellWallPackage()});
4268                    my $matchFound = 0;
4269                    for (my $i=0; $i < @{$matches}; $i++) {
4270                            if ($matches->[$i]->id() ne $biomassID && defined($matches->[$i]->essentialRxn()) && length($matches->[$i]->essentialRxn())) {
4271                                    $bioObj->essentialRxn($matches->[$i]->essentialRxn());
4272                                    print "MATCH!\n";
4273                                    $matchFound = 1;
4274                                    last;
4275                  }                  }
   
                 #Adding the biomass equation to the biomass table  
                 my $NewRow = $self->figmodel()->add_biomass_reaction($Equation,$self->id(),"Template:".$self->genome());  
                 $biomassrxn = $NewRow;  
4276          }          }
4277          return $biomassrxn;                  #Otherwise, we calculate essential reactions
4278                    if ($matchFound == 0) {
4279                            print "NOMATCH!\n";
4280                            $self->figmodel()->determine_biomass_essential_reactions($biomassID);
4281                    }
4282            }
4283            return $biomassID;
4284  }  }
4285    
4286  =head3 PrintSBMLFile  =head3 PrintSBMLFile
# Line 3732  Line 4291 
4291  =cut  =cut
4292  sub PrintSBMLFile {  sub PrintSBMLFile {
4293          my($self) = @_;          my($self) = @_;
   
4294          #Opening the SBML file for printing          #Opening the SBML file for printing
4295          my $Filename = $self->config("SBML files")->[0].$self->id().".xml";          my $Filename = $self->directory().$self->id().".xml";
         if ($self->owner() ne "master") {  
                 if (!-d $self->config("SBML files")->[0].$self->owner()."/") {  
                         system("mkdir ".$self->config("SBML files")->[0].$self->owner()."/");  
                 }  
                 $Filename = $self->config("SBML files")->[0].$self->owner()."/".$self->id().".xml";  
         }  
4296          if (!open (SBMLOUTPUT, ">$Filename")) {          if (!open (SBMLOUTPUT, ">$Filename")) {
4297                  return;                  return;
4298          }          }
   
4299          #Loading and parsing the model data          #Loading and parsing the model data
4300          my $ModelTable = $self->reaction_table();          my $mdlTbl = $self->reaction_table();
4301          if (!defined($ModelTable) || !defined($ModelTable->{"array"})) {          if (!defined($mdlTbl) || !defined($mdlTbl->{"array"})) {
4302                  print "Failed to load ".$self->id()."\n";                  return $self->fail();
                 return;  
4303          }          }
4304            my $rxnMgr = $self->figmodel()->database()->get_object_manager("reaction");
4305            my $cmpTbl = $self->figmodel()->database()->get_table("COMPARTMENTS");
4306            my $cpdMgr = $self->figmodel()->database()->get_object_manager("compound");
4307            my $bioMgr = $self->figmodel()->database()->get_object_manager("bof");
4308          #Adding intracellular metabolites that also need exchange fluxes to the exchange hash          #Adding intracellular metabolites that also need exchange fluxes to the exchange hash
4309          my $ExchangeHash = {"cpd11416" => "c"};          my $ExchangeHash = {"cpd11416" => "c"};
   
4310          my %CompartmentsPresent;          my %CompartmentsPresent;
4311          $CompartmentsPresent{"c"} = 1;          $CompartmentsPresent{"c"} = 1;
4312          my %CompoundList;          my %CompoundList;
4313          my @ReactionList;          my @ReactionList;
4314          my $ReactionTable = $self->figmodel()->database()->GetDBTable("REACTIONS");          for (my $i=0; $i < $mdlTbl->size(); $i++) {
4315          for (my $i=0; $i < $ModelTable->size(); $i++) {                  my $rxnObj;
4316                  my $Reaction = $ModelTable->get_row($i)->{"LOAD"}->[0];                  if ($mdlTbl->get_row($i)->{"LOAD"}->[0] =~ m/rxn\d\d\d\d\d/) {
4317                  if (defined($ReactionTable->get_row_by_key($Reaction,"DATABASE")) && defined($ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0])) {                          $rxnObj = $rxnMgr->get_objects({id=>$mdlTbl->get_row($i)->{"LOAD"}->[0]})->[0];
4318                          push(@ReactionList,$Reaction);                  } elsif ($mdlTbl->get_row($i)->{"LOAD"}->[0] =~ m/bio\d\d\d\d\d/) {
4319                          $_ = $ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0];                          $rxnObj = $bioMgr->get_objects({id=>$mdlTbl->get_row($i)->{"LOAD"}->[0]})->[0];
4320                    }
4321                    if (!defined($rxnObj)) {
4322                            next;
4323                    }
4324                    push(@ReactionList,$rxnObj);
4325                    $_ = $rxnObj->equation();
4326                          my @MatchArray = /(cpd\d\d\d\d\d)/g;                          my @MatchArray = /(cpd\d\d\d\d\d)/g;
4327                          for (my $j=0; $j < @MatchArray; $j++) {                          for (my $j=0; $j < @MatchArray; $j++) {
4328                                  $CompoundList{$MatchArray[$j]}->{"c"} = 1;                                  $CompoundList{$MatchArray[$j]}->{"c"} = 1;
4329                          }                          }
4330                          $_ = $ReactionTable->get_row_by_key($Reaction,"DATABASE")->{"EQUATION"}->[0];                  $_ = $rxnObj->equation();
4331                          @MatchArray = /(cpd\d\d\d\d\d\[\D\])/g;                          @MatchArray = /(cpd\d\d\d\d\d\[\D\])/g;
4332                          for (my $j=0; $j < @MatchArray; $j++) {                          for (my $j=0; $j < @MatchArray; $j++) {
4333                                  if ($MatchArray[$j] =~ m/(cpd\d\d\d\d\d)\[(\D)\]/) {                                  if ($MatchArray[$j] =~ m/(cpd\d\d\d\d\d)\[(\D)\]/) {
# Line 3778  Line 4336 
4336                                  }                                  }
4337                          }                          }
4338                  }                  }
         }  
4339    
4340          #Printing header to SBML file          #Printing header to SBML file
4341          my $ModelName = $self->id();          my $ModelName = $self->id().$self->selected_version();
4342          $ModelName =~ s/\./_/;          $ModelName =~ s/\./_/;
4343          print SBMLOUTPUT '<?xml version="1.0" encoding="UTF-8"?>'."\n";          print SBMLOUTPUT '<?xml version="1.0" encoding="UTF-8"?>'."\n";
4344      print SBMLOUTPUT '<sbml xmlns="http://www.sbml.org/sbml/level2" level="2" version="1" xmlns:html="http://www.w3.org/1999/xhtml">' . "\n";      print SBMLOUTPUT '<sbml xmlns="http://www.sbml.org/sbml/level2" level="2" version="1" xmlns:html="http://www.w3.org/1999/xhtml">' . "\n";
4345      if (defined($self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Organism name"}->[0])) {          if (defined($self->name())) {
4346          print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Organism name"}->[0].' SEED model">'."\n";                  print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->name().' SEED model">'."\n";
4347      } else {      } else {
4348          print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->id().' SEED model">'."\n";                  print SBMLOUTPUT '<model id="'.$ModelName.'" name="'.$self->id().$self->selected_version().' SEED model">'."\n";
4349      }      }
4350    
4351          #Printing the unit data          #Printing the unit data
# Line 3805  Line 4362 
4362      #Printing compartments for SBML file      #Printing compartments for SBML file
4363      print SBMLOUTPUT '<listOfCompartments>'."\n";      print SBMLOUTPUT '<listOfCompartments>'."\n";
4364      foreach my $Compartment (keys(%CompartmentsPresent)) {      foreach my $Compartment (keys(%CompartmentsPresent)) {
4365          if (defined($self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0])) {                  my $row = $cmpTbl->get_row_by_key($Compartment,"Abbreviation");
4366                  my @OutsideList = split(/\//,$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Outside"}->[0]);                  if (!defined($row) && !defined($row->{"Name"}->[0])) {
4367                            next;
4368                    }
4369                    my @OutsideList = split(/\//,$row->{"Outside"}->[0]);
4370                  my $Printed = 0;                  my $Printed = 0;
4371                  foreach my $Outside (@OutsideList) {                  foreach my $Outside (@OutsideList) {
4372                                  if (defined($CompartmentsPresent{$Outside}) && defined($self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Outside}->[0]->{"Name"}->[0])) {                          if (defined($CompartmentsPresent{$Outside}) && defined($row->{"Name"}->[0])) {
4373                                  print SBMLOUTPUT '<compartment id="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'" outside="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Outside}->[0]->{"Name"}->[0].'"/>'."\n";                                  print SBMLOUTPUT '<compartment id="'.$row->{"Name"}->[0].'" outside="'.$row->{"Name"}->[0].'"/>'."\n";
4374                                  $Printed = 1;                                  $Printed = 1;
4375                                  last;                                  last;
4376                          }                          }
4377                  }                  }
4378                  if ($Printed eq 0) {                  if ($Printed eq 0) {
4379                          print SBMLOUTPUT '<compartment id="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'"/>'."\n";                          print SBMLOUTPUT '<compartment id="'.$row->{"Name"}->[0].'"/>'."\n";
                         }  
4380          }          }
4381      }      }
4382      print SBMLOUTPUT '</listOfCompartments>'."\n";      print SBMLOUTPUT '</listOfCompartments>'."\n";
# Line 3825  Line 4384 
4384      #Printing the list of metabolites involved in the model      #Printing the list of metabolites involved in the model
4385          print SBMLOUTPUT '<listOfSpecies>'."\n";          print SBMLOUTPUT '<listOfSpecies>'."\n";
4386      foreach my $Compound (keys(%CompoundList)) {      foreach my $Compound (keys(%CompoundList)) {
4387          my $Name = $Compound;                  my $cpdObj = $cpdMgr->get_objects({id=>$Compound})->[0];
4388                    if (!defined($cpdObj)) {
4389                            next;
4390                    }
4391                  my $Formula = "";                  my $Formula = "";
4392                  if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0])) {                  if (defined($cpdObj->formula())) {
4393                          $Formula = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0];                          $Formula = $cpdObj->formula();
4394                  }                  }
4395                  if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0])) {                  my $Name = $cpdObj->name();
                         $Name = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0];  
4396                          $Name =~ s/\s/_/;                          $Name =~ s/\s/_/;
4397                          $Name .= "_".$Formula;                          $Name .= "_".$Formula;
                 }  
4398                  $Name =~ s/[<>;&\*]//;                  $Name =~ s/[<>;&\*]//;
4399                  my $Charge = 0;                  my $Charge = 0;
4400                  if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0])) {                  if (defined($cpdObj->charge())) {
4401                          $Charge = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0];                          $Charge = $cpdObj->charge();
4402                  }                  }
4403                  foreach my $Compartment (keys(%{$CompoundList{$Compound}})) {                  foreach my $Compartment (keys(%{$CompoundList{$Compound}})) {
4404                  if ($Compartment eq "e") {                  if ($Compartment eq "e") {
4405                          $ExchangeHash->{$Compound} = "e";                          $ExchangeHash->{$Compound} = "e";
4406                  }                  }
4407                  print SBMLOUTPUT '<species id="'.$Compound.'_'.$Compartment.'" name="'.$Name.'" compartment="'.$self->figmodel()->database()->GetDBTable("COMPARTMENTS")->{$Compartment}->[0]->{"Name"}->[0].'" charge="'.$Charge.'" boundaryCondition="false"/>'."\n";                          my $cmprow = $cmpTbl->get_row_by_key($Compartment,"Abbreviation");
4408                            print SBMLOUTPUT '<species id="'.$Compound.'_'.$Compartment.'" name="'.$Name.'" compartment="'.$cmprow->{"Name"}->[0].'" charge="'.$Charge.'" boundaryCondition="false"/>'."\n";
4409          }          }
4410      }      }
4411    
4412          #Printing the boundary species          #Printing the boundary species
4413          foreach my $Compound (keys(%{$ExchangeHash})) {          foreach my $Compound (keys(%{$ExchangeHash})) {
4414                  my $Name = $Compound;                  my $cpdObj = $cpdMgr->get_objects({id=>$Compound})->[0];
4415                    if (!defined($cpdObj)) {
4416                            next;
4417                    }
4418                  my $Formula = "";                  my $Formula = "";
4419                  if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0])) {                  if (defined($cpdObj->formula())) {
4420                          $Formula = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"FORMULA"}->[0];                          $Formula = $cpdObj->formula();
4421                  }                  }
4422                  if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0])) {                  my $Name = $cpdObj->name();
                         $Name = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"NAME"}->[0];  
4423                          $Name =~ s/\s/_/;                          $Name =~ s/\s/_/;
4424                          $Name .= "_".$Formula;                          $Name .= "_".$Formula;
                 }  
4425                  $Name =~ s/[<>;&\*]//;                  $Name =~ s/[<>;&\*]//;
4426                  my $Charge = 0;                  my $Charge = 0;
4427                  if (defined($self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0])) {                  if (defined($cpdObj->charge())) {
4428                          $Charge = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->{$Compound}->[0]->{"CHARGE"}->[0];                          $Charge = $cpdObj->charge();
4429                  }                  }
4430                  print SBMLOUTPUT '<species id="'.$Compound.'_b" name="'.$Name.'" compartment="Extracellular" charge="'.$Charge.'" boundaryCondition="true"/>'."\n";                  print SBMLOUTPUT '<species id="'.$Compound.'_b" name="'.$Name.'" compartment="Extracellular" charge="'.$Charge.'" boundaryCondition="true"/>'."\n";
4431          }          }
# Line 3871  Line 4434 
4434      #Printing the list of reactions involved in the model      #Printing the list of reactions involved in the model
4435          my $ObjectiveCoef;          my $ObjectiveCoef;
4436      print SBMLOUTPUT '<listOfReactions>'."\n";      print SBMLOUTPUT '<listOfReactions>'."\n";
4437          foreach my $Reaction (@ReactionList) {  
4438            foreach my $rxnObj (@ReactionList) {
4439          $ObjectiveCoef = "0.0";          $ObjectiveCoef = "0.0";
4440                  if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"EQUATION"}->[0])) {                  my $mdlrow = $mdlTbl->get_row_by_key($rxnObj->id(),"LOAD");
4441                  if ($Reaction =~ m/^bio/) {                  if ($rxnObj->id() =~ m/^bio/) {
4442                                  $ObjectiveCoef = "1.0";                                  $ObjectiveCoef = "1.0";
4443                          }                          }
4444                          my $LowerBound = -10000;                          my $LowerBound = -10000;
4445                  my $UpperBound = 10000;                  my $UpperBound = 10000;
4446                          my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateData($Reaction);                  my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateDataFromEquation($rxnObj->equation());
4447                  my $Name = $Reaction;                  my $Name = $rxnObj->name();
                 if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"NAME"}->[0])) {  
                         $Name = $self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"NAME"}->[0];  
4448                                  $Name =~ s/[<>;&]//g;                                  $Name =~ s/[<>;&]//g;
                 }  
4449                  my $Reversibility = "true";                  my $Reversibility = "true";
4450                  if (defined($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0])) {                  if (defined($mdlrow->{"DIRECTIONALITY"}->[0])) {
4451                          if ($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0] ne "<=>") {                          if ($mdlrow->{"DIRECTIONALITY"}->[0] ne "<=>") {
4452                                  $LowerBound = 0;                                  $LowerBound = 0;
4453                                          $Reversibility = "false";                                          $Reversibility = "false";
4454                          }                          }
4455                          if ($ModelTable->{$Reaction}->[0]->{"DIRECTIONALITY"}->[0] eq "<=") {                          if ($mdlrow->{"DIRECTIONALITY"}->[0] eq "<=") {
4456                                  my $Temp = $Products;                                  my $Temp = $Products;
4457                                  $Products = $Reactants;                                  $Products = $Reactants;
4458                                  $Reactants = $Temp;                                  $Reactants = $Temp;
4459                          }                          }
4460                  }                  }
4461                          print SBMLOUTPUT '<reaction id="'.$Reaction.'" name="'.$Name.'" reversible="'.$Reversibility.'">'."\n";                  print SBMLOUTPUT '<reaction id="'.$rxnObj->id().'" name="'.$Name.'" reversible="'.$Reversibility.'">'."\n";
4462                          print SBMLOUTPUT "<notes>\n";                          print SBMLOUTPUT "<notes>\n";
4463                          my $ECData = "";                          my $ECData = "";
4464                          if (defined($self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"ENZYME"}->[0])) {                  if ($rxnObj->id() !~ m/^bio/) {
4465                                  $ECData = $self->figmodel()->database()->GetDBTable("REACTIONS")->{$Reaction}->[0]->{"ENZYME"}->[0];                          if (defined($rxnObj->enzyme())) {
4466                                    my @ecList = split(/\|/,$rxnObj->enzyme());
4467                                    if (defined($ecList[1])) {
4468                                            $ECData = $ecList[1];
4469                                    }
4470                            }
4471                          }                          }
4472                          my $SubsystemData = "";                          my $SubsystemData = "";
4473                          if (defined($ModelTable->{$Reaction}->[0]->{"SUBSYSTEM"}->[0])) {                  if (defined($mdlrow->{"SUBSYSTEM"}->[0])) {
4474                                  $SubsystemData = $ModelTable->{$Reaction}->[0]->{"SUBSYSTEM"}->[0];                          $SubsystemData = $mdlrow->{"SUBSYSTEM"}->[0];
4475                          }                          }
4476                          my $GeneAssociation = "";                          my $GeneAssociation = "";
4477                          my $ProteinAssociation = "";                          my $ProteinAssociation = "";
4478                          if (defined($ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0])) {                  if (defined($mdlrow->{"ASSOCIATED PEG"}->[0])) {
4479                                  if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} == 1 && $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0] !~ m/\+/) {                          if (@{$mdlrow->{"ASSOCIATED PEG"}} == 1 && $mdlrow->{"ASSOCIATED PEG"}->[0] !~ m/\+/) {
4480                                          $GeneAssociation = $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[0];                                  $GeneAssociation = $mdlrow->{"ASSOCIATED PEG"}->[0];
4481                                  } else {                                  } else {
4482                                          if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} > 1) {                                  if (@{$mdlrow->{"ASSOCIATED PEG"}} > 1) {
4483                                                  $GeneAssociation = "( ";                                                  $GeneAssociation = "( ";
4484                                          }                                          }
4485                                          for (my $i=0; $i < @{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}}; $i++) {                                  for (my $i=0; $i < @{$mdlrow->{"ASSOCIATED PEG"}}; $i++) {
4486                                                  if ($i > 0) {                                                  if ($i > 0) {
4487                                                          $GeneAssociation .= " )  or  ( ";                                                          $GeneAssociation .= " )  or  ( ";
4488                                                  }                                                  }
4489                                                  $GeneAssociation .= $ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}->[$i];                                          $GeneAssociation .= $mdlrow->{"ASSOCIATED PEG"}->[$i];
4490                                          }                                          }
4491                                          if (@{$ModelTable->{$Reaction}->[0]->{"ASSOCIATED PEG"}} > 1) {                                  if (@{$mdlrow->{"ASSOCIATED PEG"}} > 1) {
4492                                                  $GeneAssociation .= " )";                                                  $GeneAssociation .= " )";
4493                                          }                                          }
4494                                  }                                  }
# Line 3931  Line 4497 
4497                                          $GeneAssociation = "( ".$GeneAssociation." )";                                          $GeneAssociation = "( ".$GeneAssociation." )";
4498                                  }                                  }
4499                                  $ProteinAssociation = $GeneAssociation;                                  $ProteinAssociation = $GeneAssociation;
4500                                  if (defined($self->figmodel()->database()->GetDBTable("MODEL STATS")->{$self->id()}->[0]->{"Genome ID"}->[0])) {                          if (defined($self->genome())) {
4501                                          $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]);                                  $ProteinAssociation = $self->figmodel()->translate_gene_to_protein($mdlrow->{"ASSOCIATED PEG"}->[0],$self->genome());
4502                                  }                                  }
4503                          }                          }
4504                          print SBMLOUTPUT "<html:p>GENE_ASSOCIATION:".$GeneAssociation."</html:p>\n";                          print SBMLOUTPUT "<html:p>GENE_ASSOCIATION:".$GeneAssociation."</html:p>\n";
# Line 3963  Line 4529 
4529              print SBMLOUTPUT "</kineticLaw>\n";              print SBMLOUTPUT "</kineticLaw>\n";
4530                          print SBMLOUTPUT '</reaction>'."\n";                          print SBMLOUTPUT '</reaction>'."\n";
4531                  }                  }
         }  
4532    
4533          my @ExchangeList = keys(%{$ExchangeHash});          my @ExchangeList = keys(%{$ExchangeHash});
4534          foreach my $ExCompound (@ExchangeList) {          foreach my $ExCompound (@ExchangeList) {
4535                  my $ExCompoundName = $ExCompound;                  my $cpdObj = $cpdMgr->get_objects({id=>$ExCompound})->[0];
4536                  my $Row = $self->figmodel()->database()->GetDBTable("COMPOUNDS")->get_row_by_key($ExCompound,"DATABASE");                  if (!defined($cpdObj)) {
4537                  if (defined($Row) && defined($Row->{"NAME"}->[0])) {                          next;
                         $ExCompoundName = $Row->{"NAME"}->[0];  
                         $ExCompoundName =~ s/[<>;&]//g;  
4538                  }                  }
4539                    my $ExCompoundName = $cpdObj->name();
4540                    $ExCompoundName =~ s/[<>;&]//g;
4541                  $ObjectiveCoef = "0.0";                  $ObjectiveCoef = "0.0";
4542                  print SBMLOUTPUT '<reaction id="EX_'.$ExCompound.'_'.$ExchangeHash->{$ExCompound}.'" name="EX_'.$ExCompoundName.'_'.$ExchangeHash->{$ExCompound}.'" reversible="true">'."\n";                  print SBMLOUTPUT '<reaction id="EX_'.$ExCompound.'_'.$ExchangeHash->{$ExCompound}.'" name="EX_'.$ExCompoundName.'_'.$ExchangeHash->{$ExCompound}.'" reversible="true">'."\n";
4543                  print SBMLOUTPUT "\t".'<notes>'."\n";                  print SBMLOUTPUT "\t".'<notes>'."\n";
# Line 4025  Line 4590 
4590                  $tbl->add_row($row);                  $tbl->add_row($row);
4591          }          }
4592          $tbl->save();          $tbl->save();
         system("cp ".$tbl->filename()." ".$self->figmodel()->config("Model table download directory")->[0].$self->figmodel()->config("ModelSimpleReactionTable")->{filename_prefix}->[0]."-".$self->id().".txt");  
4593  }  }
4594    
4595  =head3 PrintModelLPFile  =head3 PrintModelLPFile
# Line 4054  Line 4618 
4618    
4619  =head3 patch_model  =head3 patch_model
4620  Definition:  Definition:
4621          FIGMODELTable:patch results = FIGMODELmodel->patch_model(FIGMODELTable:patch table);          FIGMODELmodel->patch_model([] -or- {} of patch arguments);
4622  Description:  Description:
4623  =cut  =cut
4624  sub patch_model {  sub patch_model {
4625          my ($self,$tbl) = @_;          my ($self,$arguments) = @_;
4626            if (($self->owner() eq "master" || $self->owner() eq "chenry") && $self->id() =~ m/Seed/) {
4627          #Instantiating table                  print "Model qualifies.";
4628          my $results = FIGMODELTable->new(["Reactions","New genes","Old genes","Genes","Roles","Status"],$self->directory()."PatchResults-".$self->id().$self->selected_version().".tbl",["Reaction"],"\t",";",undef);                  $self->BuildSpecificBiomassReaction();
4629          #Getting genome annotations                  my $reactionTable = $self->reaction_table();
4630          my $features = $self->figmodel()->database()->get_genome_feature_table($self->genome());                  my $row = $reactionTable->get_row_by_key("rxn05294","LOAD");
4631          #Gettubg reaction table                  if (defined($row)) {
4632          my $reactions = $self->reaction_table();                          $reactionTable->delete_row($reactionTable->row_index($row));
4633          #Checking for patched roles                  }
4634          for (my $i=0; $i < $tbl->size(); $i++) {                  $row = $reactionTable->get_row_by_key("rxn05295","LOAD");
4635                  my $row = $tbl->get_row($i);                  if (defined($row)) {
4636                  my @genes = $features->get_rows_by_key($row->{ROLE}->[0],"ROLES");                          $reactionTable->delete_row($reactionTable->row_index($row));
4637                  if (@genes > 0) {                  }
4638                          for (my $j=0; $j < @{$row->{REACTIONS}};$j++) {                  $row = $reactionTable->get_row_by_key("rxn05296","LOAD");
4639                                  my $resultrxn = $results->get_row_by_key($row->{REACTIONS}->[$j],"Reactions");                  if (defined($row)) {
4640                                  if (!defined($resultrxn)) {                          $reactionTable->delete_row($reactionTable->row_index($row));
4641                                          $resultrxn = $results->add_row({"Reactions"=>[$row->{REACTIONS}->[$j]],"Roles"=>[$row->{ROLE}->[0]]});                  }
4642                                  }                  $reactionTable->save();
4643                                  my $rxnrow = $reactions->get_row_by_key($row->{REACTIONS}->[$j],"LOAD");                  my $growth = $self->calculate_growth();
4644                                  if (defined($rxnrow) && !defined($resultrxn->{"Old genes"})) {                  if ($growth =~ m/^NOGROWTH/ ||  $growth < 1e-5) {
4645                                          $resultrxn->{"Old genes"} = $rxnrow->{"ASSOCIATED PEG"};                          $self->error_message("patch_model:model did not grow, gapfilling!");
4646                                          if ($resultrxn->{"Old genes"}->[0] !~ m/GAP|BOF|UNIVERSAL|SPONTANEOUS/) {                          $self->GapFillModel(1);
4647                                                  push(@{$resultrxn->{"Genes"}},@{$resultrxn->{"Old genes"}});                          $growth = $self->calculate_growth();
4648                                          }                          if ($growth =~ m/^NOGROWTH/ ||  $growth < 1e-5) {
4649                                  }                                  $self->error_message("patch_model:model still did not grow!");
4650                                  delete $resultrxn->{"Current gene set"};                                  print "Patch failed.\n";
4651                                  if (defined($resultrxn->{"Genes"})) {                                  return;
                                         push(@{$resultrxn->{"Current gene set"}},@{$resultrxn->{"Genes"}});  
                                 }  
                                 for (my $k=0; $k < @genes; $k++) {  
                                         if ($genes[$k]->{ID}->[0] =~ m/(peg\.\d+)/) {  
                                                 my $gene = $1;  
                                                 my $addgene = 1;  
                                                 if (defined($resultrxn->{"Old genes"})) {  
                                                         for (my $m=0; $m < @{$resultrxn->{"Old genes"}}; $m++) {  
                                                                 if ($resultrxn->{"Old genes"}->[$m] =~ m/$gene/) {  
                                                                         $addgene = 0;  
                                                                 }  
                                                         }  
                                                 }  
                                                 if ($addgene == 1) {  
                                                         push(@{$resultrxn->{"New genes"}},$gene);  
                                                         if ($row->{COMPLEX}->[0] ne "0" && defined($resultrxn->{"Current gene set"})) {  
                                                                 my $added = 0;  
                                                                 for (my $m=0; $m < @{$resultrxn->{"Current gene set"}}; $m++) {  
                                                                         if ($row->{COMPLEX}->[0] eq "1") {  
                                                                                 $resultrxn->{"Current gene set"}->[$m] = $resultrxn->{"Current gene set"}->[$m]."+".$gene;  
                                                                                 $added = 1;  
                                                                         } else {  
                                                                                 my @geneset = split(/\+/,$resultrxn->{"Current gene set"}->[$m]);  
                                                                                 for (my $n=0; $n < @geneset;$n++) {  
                                                                                         if ($self->figmodel()->colocalized_genes($geneset[$n],$gene,$self->genome()) == 1) {  
                                                                                                 $resultrxn->{"Current gene set"}->[$m] = $resultrxn->{"Current gene set"}->[$m]."+".$gene;  
                                                                                                 $added = 1;  
                                                                                                 last;  
                                                                                         }  
                                                                                 }  
                                                                         }  
                                                                 }  
                                                                 if ($added == 0) {  
                                                                         push(@{$resultrxn->{"Current gene set"}},$gene);  
                                                                 }  
                                                         } else {  
                                                                 push(@{$resultrxn->{"Current gene set"}},$gene);  
                                                         }  
                                                 }  
                                         }  
                                 }  
                                 delete $resultrxn->{"Genes"};  
                                 push(@{$resultrxn->{"Genes"}},@{$resultrxn->{"Current gene set"}});  
                         }  
                 }  
         }  
   
         #Ensuring that the old model is preserved  
         $self->ArchiveModel();  
         #Modifing the reaction list  
         for (my $i=0; $i < $results->size();$i++) {  
                 my $row = $results->get_row($i);  
                 my $rxnrow = $reactions->get_row_by_key($row->{"Reactions"}->[0],"LOAD");  
                 if (defined($rxnrow)) {  
                         $rxnrow->{"ASSOCIATED PEG"} = $row->{"Genes"};  
                 } else {  
                         $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"]});  
4652                  }                  }
4653          }          }
4654          $reactions->save();                  $self->PrintSBMLFile();
         $results->save();  
         $self->update_model_stats();  
4655          $self->PrintModelLPFile();          $self->PrintModelLPFile();
4656                    $self->PrintModelSimpleReactionTable();
4657          $self->run_default_model_predictions();          $self->run_default_model_predictions();
4658          #Returning results                  $self->update_model_stats();
4659          return $results;                  print "Patch complete.\n";
4660            }
4661  }  }
4662    
4663  =head3 translate_genes  =head3 translate_genes
# Line 4354  Line 4861 
4861          }          }
4862  }  }
4863    
4864    =pod
4865    
4866    =item * [string]:I<list of essential genes> = B<run_geneKO_slow> (string:I<media>,0/1:I<max growth>,0/1:I<save results>);
4867    
4868    =cut
4869    
4870    sub run_geneKO_slow {
4871            my ($self,$media,$maxGrowth,$save) = @_;
4872            my $output;
4873            my $UniqueFilename = $self->figmodel()->filename();
4874            if (defined($maxGrowth) && $maxGrowth == 1) {
4875                    system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$media,["ProductionMFA"],{"perform single KO experiments" => 1,"MFASolver" => "GLPK","Constrain objective to this fraction of the optimal value" => 0.999},"SlowGeneKO-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
4876            } else {
4877                    system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$media,["ProductionMFA"],{"perform single KO experiments" => 1,"MFASolver" => "GLPK","Constrain objective to this fraction of the optimal value" => 0.1},"SlowGeneKO-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
4878            }
4879            if (!-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."DeletionStudyResults.txt") {
4880                    print "Deletion study file not found!.\n";
4881                    return undef;
4882            }
4883            my $deltbl = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."DeletionStudyResults.txt",";","|",1,["Experiment"]);
4884            for (my $i=0; $i < $deltbl->size(); $i++) {
4885                    my $row = $deltbl->get_row($i);
4886                    if ($row->{"Insilico growth"}->[0] < 0.0000001) {
4887                            push(@{$output},$row->{Experiment}->[0]);
4888                    }
4889            }
4890            if (defined($output)) {
4891                    if (defined($save) && $save == 1) {
4892                            my $tbl = $self->essentials_table();
4893                            my $row = $tbl->get_row_by_key($media,"MEDIA",1);
4894                            $row->{"ESSENTIAL GENES"} = $output;
4895                            $tbl->save();
4896                    }
4897            }
4898            return $output;
4899    }
4900    
4901    =pod
4902    
4903    =item * [string]:I<list of minimal genes> = B<run_gene_minimization> (string:I<media>,0/1:I<max growth>,0/1:I<save results>);
4904    
4905    =cut
4906    
4907    sub run_gene_minimization {
4908            my ($self,$media,$maxGrowth,$save) = @_;
4909            my $output;
4910    
4911            #Running the MFAToolkit
4912            my $UniqueFilename = $self->figmodel()->filename();
4913            if (defined($maxGrowth) && $maxGrowth == 1) {
4914                    system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$media,["ProductionMFA"],{"optimize organism genes" => 1,"MFASolver" => "CPLEX","Constrain objective to this fraction of the optimal value" => 0.999},"MinimizeGenes-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
4915            } else {
4916                    system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$media,["ProductionMFA"],{"optimize organism genes" => 1,"MFASolver" => "CPLEX","Constrain objective to this fraction of the optimal value" => 0.1},"MinimizeGenes-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
4917            }
4918            my $tbl = $self->figmodel()->LoadProblemReport($UniqueFilename);
4919            if (!defined($tbl)) {
4920                    return undef;
4921            }
4922            for (my $i=0; $i < $tbl->size(); $i++) {
4923                    my $row = $tbl->get_row($i);
4924                    if ($row->{Notes}->[0] =~ m/Recursive\sMILP\sGENE_USE\soptimization/) {
4925                            my @array = split(/\|/,$row->{Notes}->[0]);
4926                            my $solution = $array[0];
4927                            $_ = $solution;
4928                            my @OriginalArray = /(peg\.\d+)/g;
4929                            push(@{$output},@OriginalArray);
4930                            last;
4931                    }
4932            }
4933    
4934            if (defined($output)) {
4935                    if (defined($save) && $save == 1) {
4936                            my $tbl = $self->load_model_table("MinimalGenes");
4937                            my $row = $tbl->get_table_by_key("MEDIA",$media)->get_row_by_key("MAXGROWTH",$maxGrowth);
4938                            if (defined($row)) {
4939                                    $row->{GENES} = $output;
4940                            } else {
4941                                    $tbl->add_row({GENES => $output,MEDIA => [$media],MAXGROWTH => [$maxGrowth]});
4942                            }
4943                            $tbl->save();
4944                    }
4945            }
4946            return $output;
4947    }
4948    
4949    =pod
4950    
4951    =item * [string]:I<list of inactive genes> = B<identify_inactive_genes> (string:I<media>,0/1:I<max growth>,0/1:I<save results>);
4952    
4953    =cut
4954    
4955    sub identify_inactive_genes {
4956            my ($self,$media,$maxGrowth,$save) = @_;
4957            my $output;
4958            #Running the MFAToolkit
4959            my $UniqueFilename = $self->figmodel()->filename();
4960            if (defined($maxGrowth) && $maxGrowth == 1) {
4961                    system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$media,["ProductionMFA"],{"find tight bounds" => 1,"MFASolver" => "GLPK","Constrain objective to this fraction of the optimal value" => 0.999},"Classify-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
4962            } else {
4963                    system($self->figmodel()->GenerateMFAToolkitCommandLineCall($UniqueFilename,$self->id(),$media,["ProductionMFA"],{"find tight bounds" => 1,"MFASolver" => "GLPK","Constrain objective to this fraction of the optimal value" => 0.1},"Classify-".$self->id().$self->selected_version()."-".$UniqueFilename.".log",undef,$self->selected_version()));
4964            }
4965            #Reading in the output bounds file
4966            my $ReactionTB;
4967            if (-e $self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsReactionData0.txt") {
4968                    $ReactionTB = $self->figmodel()->database()->load_table($self->config("MFAToolkit output directory")->[0].$UniqueFilename."/MFAOutput/TightBoundsReactionData0.txt",";","|",1,["DATABASE ID"]);
4969            }
4970            if (!defined($ReactionTB)) {
4971                    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";
4972                    return undef;
4973            }
4974            #Clearing output
4975            $self->figmodel()->clearing_output($UniqueFilename,"Classify-".$self->id().$self->selected_version()."-".$UniqueFilename.".log");
4976            my $geneHash;
4977            my $activeGeneHash;
4978            for (my $i=0; $i < $ReactionTB->size(); $i++) {
4979                    my $Row = $ReactionTB->get_row($i);
4980                    if (defined($Row->{"Min FLUX"}) && defined($Row->{"Max FLUX"}) && defined($Row->{"DATABASE ID"}) && $Row->{"DATABASE ID"}->[0] =~ m/rxn\d\d\d\d\d/) {
4981                            my $data = $self->get_reaction_data($Row->{"DATABASE ID"}->[0]);
4982                            if (defined($data->{"ASSOCIATED PEG"})) {
4983                                    my $active = 0;
4984                                    if ($Row->{"Min FLUX"}->[0] > 0.00000001 || $Row->{"Max FLUX"}->[0] < -0.00000001 || ($Row->{"Max FLUX"}->[0]-$Row->{"Min FLUX"}->[0]) > 0.00000001) {
4985                                            $active = 1;
4986                                    }
4987                                    for (my $j=0; $j < @{$data->{"ASSOCIATED PEG"}}; $j++) {
4988                                            $_ = $data->{"ASSOCIATED PEG"}->[$j];
4989                                            my @OriginalArray = /(peg\.\d+)/g;
4990                                            for (my $k=0; $k < @OriginalArray; $k++) {
4991                                                    if ($active == 1) {
4992                                                            $activeGeneHash->{$OriginalArray[$k]} = 1;
4993                                                    }
4994                                                    $geneHash->{$OriginalArray[$k]} = 1;
4995                                            }
4996                                    }
4997                            }
4998                    }
4999            }
5000            my @allGenes = keys(%{$geneHash});
5001            for (my $i=0; $i < @allGenes; $i++) {
5002                    if (!defined($activeGeneHash->{$allGenes[$i]})) {
5003                            push(@{$output},$allGenes[$i]);
5004                    }
5005            }
5006            if (defined($output)) {
5007                    if (defined($save) && $save == 1) {
5008                            my $tbl = $self->load_model_table("InactiveGenes");
5009                            my $row = $tbl->get_table_by_key("MEDIA",$media)->get_row_by_key("MAXGROWTH",$maxGrowth);
5010                            if (defined($row)) {
5011                                    $row->{GENES} = $output;
5012                            } else {
5013                                    $tbl->add_row({GENES => $output,MEDIA => [$media],MAXGROWTH => [$maxGrowth]});
5014                            }
5015                            $tbl->save();
5016                    }
5017            }
5018            return $output;
5019    }
5020    
5021    sub ConvertVersionsToHistoryFile {
5022            my ($self) = @_;
5023            my $vone = 0;
5024            my $vtwo = 0;
5025            my $continue = 1;
5026            my $lastTable;
5027            my $currentTable;
5028            my $cause;
5029            my $lastChanged = 0;
5030            my $noHitCount = 0;
5031            while ($continue == 1) {
5032                    $cause = "NONE";
5033                    $currentTable = undef;
5034                    if (-e $self->directory().$self->id()."V".($vone+1).".".$vtwo.".txt") {
5035                            $noHitCount = 0;
5036                            $lastChanged = 0;
5037                            $vone = $vone+1;
5038                            $currentTable = $self->figmodel()->database()->load_table($self->directory().$self->id()."V".$vone.".".$vtwo.".txt",";","|",1,["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG"]);
5039                            $cause = "RECONSTRUCTION";
5040                    } elsif (-e $self->directory().$self->id()."V".$vone.".".($vtwo+1).".txt") {
5041                            $noHitCount = 0;
5042                            $lastChanged = 0;
5043                            $vtwo = $vtwo+1;
5044                            $currentTable = $self->figmodel()->database()->load_table($self->directory().$self->id()."V".$vone.".".$vtwo.".txt",";","|",1,["LOAD","DIRECTIONALITY","COMPARTMENT","ASSOCIATED PEG"]);
5045                            $cause = "AUTOCOMPLETION";
5046                    } elsif ($lastChanged == 0) {
5047                            $lastChanged = 1;
5048                            $vone = $vone+1;
5049                            $cause = "RECONSTRUCTION";
5050                    } elsif ($lastChanged == 1) {
5051                            $lastChanged = 2;
5052                            $vone = $vone-1;
5053                            $vtwo = $vtwo+1;
5054                            $cause = "AUTOCOMPLETION";
5055                    } elsif ($lastChanged == 2) {
5056                            $lastChanged = 0;
5057                            $vone = $vone+1;
5058                            $cause = "RECONSTRUCTION";
5059                    }
5060                    if (defined($currentTable)) {
5061                            if (defined($lastTable)) {
5062                                    print $cause."\t".$self->directory().$self->id()."V".$vone.".".$vtwo.".txt\n";
5063                                    $self->calculate_model_changes($lastTable,$cause,$currentTable,"V".$vone.".".$vtwo);
5064                            }
5065                            $lastTable = $currentTable;
5066                    } else {
5067                            $noHitCount++;
5068                            if ($noHitCount >= 40) {
5069                                    last;
5070                            }
5071                    }
5072            }
5073    }
5074    
5075  1;  1;

Legend:
Removed from v.1.18  
changed lines
  Added in v.1.27

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3