[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.23, Mon Jul 26 05:53:06 2010 UTC revision 1.28, Wed Aug 25 20:39:16 2010 UTC
# Line 100  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 232  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 285  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 332  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 1035  Line 1055 
1055    
1056  sub autocompleteMedia {  sub autocompleteMedia {
1057          my ($self,$newMedia) = @_;          my ($self,$newMedia) = @_;
1058            if (defined($newMedia)) {
1059          return $self->{_data}->autoCompleteMedia($newMedia);          return $self->{_data}->autoCompleteMedia($newMedia);
1060  }  }
1061            return $self->{_data}->autoCompleteMedia();
1062    }
1063    
1064  sub biomassReaction {  sub biomassReaction {
1065          my ($self,$newBiomass) = @_;          my ($self,$newBiomass) = @_;
# Line 1045  Line 1068 
1068          } else {          } else {
1069                  #Figuring out what the old biomass is                  #Figuring out what the old biomass is
1070                  my $oldBiomass = $self->{_data}->biomassReaction();                  my $oldBiomass = $self->{_data}->biomassReaction();
1071                  if ($newBiomass ne $oldBiomass) {                  $self->{_data}->biomassReaction($newBiomass);
                         #Figuring out if the new biomass exists  
                         my $handle = $self->database()->get_object_manager("bof");  
                         my $objects = $handle->get_objects({"DATABASE" => [$newBiomass]});  
                         if (!defined($objects) || !defined($objects->[0])) {  
                                 print STDERR "Could not find new biomass reaction ".$newBiomass."\n";  
                                 return $oldBiomass;  
                         }  
1072                          #Changing the biomass reaction in the model file                          #Changing the biomass reaction in the model file
1073                          my $rxnTbl = $self->reaction_table();                          my $rxnTbl = $self->reaction_table();
1074                          for (my $i=0; $i < $rxnTbl->size(); $i++) {                          for (my $i=0; $i < $rxnTbl->size(); $i++) {
# Line 1062  Line 1078 
1078                                  }                                  }
1079                          }                          }
1080                          $rxnTbl->save();                          $rxnTbl->save();
1081                          #Changing the biomass reaction in the database                  if ($newBiomass ne $oldBiomass) {
1082                          return  $self->{_data}->biomassReaction($newBiomass);                          #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    
# Line 1204  Line 1227 
1227          $self->{_data}->modificationDate(time());          $self->{_data}->modificationDate(time());
1228          my $version = $self->{_data}->version();          my $version = $self->{_data}->version();
1229          $self->{_data}->version($version+1);          $self->{_data}->version($version+1);
1230            $self->update_model_stats();
1231  }  }
1232    
1233  =head3 update_model_stats  =head3 update_model_stats
# Line 1217  Line 1241 
1241          #Getting reaction table          #Getting reaction table
1242          my $rxntbl = $self->reaction_table();          my $rxntbl = $self->reaction_table();
1243          if (!defined($rxntbl)) {          if (!defined($rxntbl)) {
1244                  $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;  
1245          }          }
1246          my $cpdtbl = $self->compound_table();          my $cpdtbl = $self->compound_table();
1247    
# Line 2036  Line 2059 
2059                  }                  }
2060          }          }
2061    
2062          #Checking if a biomass reaction already exists          #Creating biomass reaction for model
2063          my $BiomassReactionRow = $self->BuildSpecificBiomassReaction();          my $biomassID = $self->BuildSpecificBiomassReaction();
2064          if (!defined($BiomassReactionRow)) {          if ($biomassID !~ m/bio\d\d\d\d\d/) {
2065                  $self->set_status(-2,"Preliminary reconstruction failed: could not generate biomass reaction");                  $self->set_status(-2,"Preliminary reconstruction failed: could not generate biomass reaction");
2066                  $self->figmodel()->error_message("FIGMODEL:CreateModelReactionList: Could not generate biomass function for ".$self->id()."!");                  die $self->error_message("CreateSingleGenomeReactionList: Could not generate biomass reaction!");
2067                  return $self->fail();          }
2068            #Getting the biomass reaction PPO object
2069            my $bioRxn = $self->figmodel()->database()->get_object("bof",{id=>$biomassID});
2070            if (!defined($bioRxn)) {
2071                    die $self->error_message("CreateSingleGenomeReactionList: Could not find biomass reaction ".$biomassID."!");
2072            }
2073            #Getting the list of essential reactions for biomass reaction
2074            my $ReactionList;
2075            my $essentialReactions = $bioRxn->essentialRxn();
2076            if (defined($essentialReactions) && $essentialReactions =~ m/rxn\d\d\d\d\d/) {
2077                    push(@{$ReactionList},split(/\|/,$essentialReactions));
2078                    if ($essentialReactions !~ m/$biomassID/) {
2079                            push(@{$ReactionList},$biomassID);
2080                    }
2081          }          }
         my $ReactionList = $BiomassReactionRow->{"ESSENTIAL REACTIONS"};  
         push(@{$ReactionList},$BiomassReactionRow->{DATABASE}->[0]);  
2082    
2083          #Adding biomass reactions to the model table          #Adding biomass reactions to the model table
2084          foreach my $BOFReaction (@{$ReactionList}) {          foreach my $BOFReaction (@{$ReactionList}) {
# Line 2149  Line 2183 
2183          #Updating the model stats table          #Updating the model stats table
2184          $self->update_stats_for_build();          $self->update_stats_for_build();
2185          $self->PrintSBMLFile();          $self->PrintSBMLFile();
2186            $self->PrintModelLPFile();
2187            $self->PrintModelSimpleReactionTable();
2188          if (defined($OriginalModelTable)) {          if (defined($OriginalModelTable)) {
2189                  $self->calculate_model_changes($OriginalModelTable,"REBUILD");                  $self->calculate_model_changes($OriginalModelTable,"REBUILD");
2190          }          }
# Line 3751  Line 3787 
3787          return 1;          return 1;
3788  }  }
3789    
3790  =head3 BuildSpecificBiomassReaction  =head3 DetermineCofactorLipidCellWallComponents
3791  Definition:  Definition:
3792          FIGMODELmodel->BuildSpecificBiomassReaction();          {cofactor=>{string:compound id=>float:coefficient},lipid=>...cellWall=>} = FIGMODELmodel->DetermineCofactorLipidCellWallComponents();
3793  Description:  Description:
3794  =cut  =cut
3795  sub BuildSpecificBiomassReaction {  sub DetermineCofactorLipidCellWallComponents {
3796          my ($self) = @_;          my ($self) = @_;
3797            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  
3798                  my $genomestats = $self->figmodel()->get_genome_stats($self->genome());                  my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
3799                  my $Class = $genomestats->{CLASS}->[0];          my $Class = $self->{_data}->cellwalltype();
3800                  my $Name = $genomestats->{NAME}->[0];          my $Name = $self->name();
3801            my $translation = {COFACTOR=>"cofactor",LIPIDS=>"lipid","CELL WALL"=>"cellWall"};
3802                  #Checking for phoenix variants                  #Checking for phoenix variants
3803                  my $PhoenixVariantTable = $self->figmodel()->database()->GetDBTable("Phoenix variants table");                  my $PhoenixVariantTable = $self->figmodel()->database()->GetDBTable("Phoenix variants table");
3804                  my $Phoenix = 0;                  my $Phoenix = 0;
3805                  my @Rows = $PhoenixVariantTable->get_rows_by_key($OrganismID,"GENOME");          my @Rows = $PhoenixVariantTable->get_rows_by_key($self->genome(),"GENOME");
3806                  my $VariantHash;                  my $VariantHash;
3807                  for (my $i=0; $i < @Rows; $i++) {                  for (my $i=0; $i < @Rows; $i++) {
3808                          $Phoenix = 1;                          $Phoenix = 1;
# Line 3784  Line 3810 
3810                                  $VariantHash->{$Rows[$i]->{"SUBSYSTEM"}->[0]} = $Rows[$i]->{"VARIANT"}->[0];                                  $VariantHash->{$Rows[$i]->{"SUBSYSTEM"}->[0]} = $Rows[$i]->{"VARIANT"}->[0];
3811                          }                          }
3812                  }                  }
   
3813                  #Collecting genome data                  #Collecting genome data
3814                  my $RoleHash;                  my $RoleHash;
3815                  my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());                  my $FeatureTable = $self->figmodel()->GetGenomeFeatureTable($self->genome());
# Line 3801  Line 3826 
3826                                  }                                  }
3827                          }                          }
3828                  }                  }
   
3829                  #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
3830                  my $ComponentTypes;          my $includedHash;
3831                  my $EquationData;          my $BiomassReactionTemplateTable = $self->figmodel()->database()->get_table("BIOMASSTEMPLATE");
                 my $BiomassReactionTemplateTable = $self->figmodel()->database()->GetDBTable("BIOMASS TEMPLATE");  
3832                  for (my $i=0; $i < $BiomassReactionTemplateTable->size(); $i++) {                  for (my $i=0; $i < $BiomassReactionTemplateTable->size(); $i++) {
3833                          my $Row = $BiomassReactionTemplateTable->get_row($i);                          my $Row = $BiomassReactionTemplateTable->get_row($i);
3834                    if (defined($translation->{$Row->{CLASS}->[0]})) {
3835                            my $coef = -1;
3836                            if ($Row->{"REACTANT"}->[0] eq "NO") {
3837                                    $coef = 1;
3838                                    if ($Row->{"COEFFICIENT"}->[0] =~ m/cpd/) {
3839                                            $coef = $Row->{"COEFFICIENT"}->[0];
3840                                    }
3841                            }
3842                          if (defined($Row->{"INCLUSION CRITERIA"}->[0]) && $Row->{"INCLUSION CRITERIA"}->[0] eq "UNIVERSAL") {                          if (defined($Row->{"INCLUSION CRITERIA"}->[0]) && $Row->{"INCLUSION CRITERIA"}->[0] eq "UNIVERSAL") {
3843                                  push(@{$EquationData},$Row);                                  $includedHash->{$Row->{"ID"}->[0]} = 1;
3844                                  $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;                                  $templateResults->{$translation->{$Row->{CLASS}->[0]}}->{$Row->{"ID"}->[0]} = $coef;
3845                          } else {                          } elsif (defined($Row->{"INCLUSION CRITERIA"}->[0])) {
3846                                  my $Criteria = $Row->{"INCLUSION CRITERIA"}->[0];                                  my $Criteria = $Row->{"INCLUSION CRITERIA"}->[0];
3847                                  my $End = 0;                                  my $End = 0;
3848                                  while ($End == 0) {                                  while ($End == 0) {
3849                                          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)\{([^{^}]+)\}$/) {
3850                                                  print $Criteria."\n";                                                  print $Criteria." : ";
3851                                                  my $Start = "";                                                  my $Start = "";
3852                                                  my $End = "";                                                  my $End = "";
3853                                                  my $Condition = $1;                                                  my $Condition = $1;
# Line 3840  Line 3871 
3871                                                                  $Result = "NO";                                                                  $Result = "NO";
3872                                                                  last;                                                                  last;
3873                                                          } elsif ($Array[$j] =~ m/^COMPOUND:(.+)/) {                                                          } elsif ($Array[$j] =~ m/^COMPOUND:(.+)/) {
3874                                                                  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") {  
3875                                                                          $Result = "YES";                                                                          $Result = "YES";
3876                                                                          last;                                                                          last;
3877                                                                  } elsif ($Match != 1 && $Condition eq "AND") {                                                                  } elsif (!defined($includedHash->{$1}) && $Condition eq "AND") {
3878                                                                          $Result = "NO";                                                                          $Result = "NO";
3879                                                                          last;                                                                          last;
3880                                                                  }                                                                  }
3881                                                          } elsif ($Array[$j] =~ m/^!COMPOUND:(.+)/) {                                                          } elsif ($Array[$j] =~ m/^!COMPOUND:(.+)/) {
3882                                                                  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") {  
3883                                                                          $Result = "YES";                                                                          $Result = "YES";
3884                                                                          last;                                                                          last;
3885                                                                  } elsif ($Match == 1 && $Condition eq "AND") {                                                                  } elsif (defined($includedHash->{$1}) && $Condition eq "AND") {
3886                                                                          $Result = "NO";                                                                          $Result = "NO";
3887                                                                          last;                                                                          last;
3888                                                                  }                                                                  }
# Line 3981  Line 3998 
3998                                          }                                          }
3999                                  }                                  }
4000                                  if ($Criteria eq "YES") {                                  if ($Criteria eq "YES") {
4001                                          push(@{$EquationData},$Row);                                          $templateResults->{$translation->{$Row->{CLASS}->[0]}}->{$Row->{"ID"}->[0]} = $coef;
4002                                          $ComponentTypes->{$Row->{"CLASS"}->[0]}->{$Row->{"ID"}->[0]} = 1;                                          $includedHash->{$Row->{"ID"}->[0]} = 1;
4003                                  }                                  }
4004                          }                          }
4005                  }                  }
4006            }
4007                  #Building biomass equation          my $types = ["cofactor","lipid","cellWall"];
4008                  my %Reactants;          my $cpdMgr = $self->figmodel()->database()->get_object_manager("compound");
4009                  my %Products;          for (my $i=0; $i < @{$types}; $i++) {
4010                  foreach my $EquationRow (@{$EquationData}) {                  my @list =      keys(%{$templateResults->{$types->[$i]}});
4011                          #First determine what the coefficient should be                  my $entries = 0;
4012                          my $Coefficient;                  for (my $j=0; $j < @list; $j++) {
4013                          if (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/^[0-9\.]+$/) {                          if ($templateResults->{$types->[$i]}->{$list[$j]} eq "-1") {
4014                                  $Coefficient = $EquationRow->{"COEFFICIENT"}->[0];                                  my $objs = $cpdMgr->get_objects({id=>$list[$j]});
4015                          } elsif (defined($EquationRow->{"COEFFICIENT"}) && $EquationRow->{"COEFFICIENT"}->[0] =~ m/cpd\d\d\d\d\d/) {                                  if (!defined($objs->[0]) || $objs->[0]->mass() == 0) {
4016                                  $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]};  
4017                                  } else {                                  } else {
4018                                          $Products{$EquationRow->{"ID"}->[0]} = $Coefficient;                                          $entries++;
4019                                    }
4020                            }
4021                    }
4022                    for (my $j=0; $j < @list; $j++) {
4023                            if ($templateResults->{$types->[$i]}->{$list[$j]} eq "-1") {
4024                                    $templateResults->{$types->[$i]}->{$list[$j]} = -1/$entries;
4025                            } elsif ($templateResults->{$types->[$i]}->{$list[$j]} =~ m/cpd/) {
4026                                    my $netCoef = 0;
4027                                    my @allcpd = split(/,/,$templateResults->{$types->[$i]}->{$list[$j]});
4028                                    for (my $k=0; $k < @allcpd; $k++) {
4029                                            if (defined($templateResults->{$types->[$i]}->{$allcpd[$k]}) && $templateResults->{$types->[$i]}->{$allcpd[$k]} ne "-1e-5") {
4030                                                    $netCoef += (1/$entries);
4031                                            } elsif (defined($templateResults->{$types->[$i]}->{$allcpd[$k]}) && $templateResults->{$types->[$i]}->{$allcpd[$k]} eq "-1e-5") {
4032                                                    $netCoef += 1e-5;
4033                                            }
4034                                    }
4035                                    $templateResults->{$types->[$i]}->{$list[$j]} = $netCoef;
4036                            }
4037                    }
4038            }
4039            return $templateResults;
4040    }
4041    
4042    =head3 BuildSpecificBiomassReaction
4043    Definition:
4044            FIGMODELmodel->BuildSpecificBiomassReaction();
4045    Description:
4046    =cut
4047    sub BuildSpecificBiomassReaction {
4048            my ($self) = @_;
4049            #Getting the database handle for biomass reactions
4050            my $bioMgr = $self->figmodel()->database()->get_object_manager("bof");
4051            #Checking if the current biomass reaction appears in more than on model, if not, this biomass reaction is conserved for this model
4052            my $biomassID = $self->biomassReaction();
4053            if ($biomassID =~ m/bio\d\d\d\d\d/) {
4054                    my $mdlMgr = $self->figmodel()->database()->get_object_manager("model");
4055                    my $mdlObs = $mdlMgr->get_objects({biomassReaction=>$biomassID});
4056                    if (defined($mdlObs->[1])) {
4057                            $biomassID = "NONE";
4058                    }
4059            }
4060            #If the biomass ID is "NONE", then we create a new biomass reaction for the model
4061            my $bioObj;
4062            my $originalPackages = "";
4063            my $originalEssReactions = "";
4064            if ($biomassID !~ m/bio\d\d\d\d\d/) {
4065                    #Getting the current largest ID
4066                    $biomassID = $self->figmodel()->database()->check_out_new_id("bof");
4067                    $bioObj = $bioMgr->create({
4068                            id=>$biomassID,owner=>$self->owner(),users=>$self->users(),name=>"Biomass",equation=>"NONE",protein=>"0",
4069                            energy=>"0",DNA=>"0",RNA=>"0",lipid=>"0",cellWall=>"0",cofactor=>"0",
4070                            modificationDate=>time(),creationDate=>time(),
4071                            cofactorPackage=>"NONE",lipidPackage=>"NONE",cellWallPackage=>"NONE",
4072                            DNACoef=>"NONE",RNACoef=>"NONE",proteinCoef=>"NONE",lipidCoef=>"NONE",
4073                            cellWallCoef=>"NONE",cofactorCoef=>"NONE",essentialRxn=>"NONE"});
4074                    if (!defined($bioObj)) {
4075                            die $self->error_message("BuildSpecificBiomassReaction():Could not create new biomass reaction ".$biomassID."!");
4076                    }
4077            } else {
4078                    #Getting the biomass DB handler from the database
4079                    my $objs = $bioMgr->get_objects({id=>$biomassID});
4080                    if (!defined($objs->[0])) {
4081                            die $self->error_message("BuildSpecificBiomassReaction():Could not find biomass reaction ".$biomassID." in database!");
4082                    }
4083                    $bioObj = $objs->[0];
4084                    $bioObj->owner($self->owner());
4085                    $bioObj->users($self->users());
4086                    if (defined($bioObj->essentialRxn())) {
4087                            $originalEssReactions = $bioObj->essentialRxn();
4088                            $originalPackages = $bioObj->cofactorPackage().$bioObj->lipidPackage().$bioObj->cellWallPackage();
4089                                  }                                  }
4090                          }                          }
4091            #Getting genome stats
4092            my $genomestats = $self->figmodel()->get_genome_stats($self->genome());
4093            my $Class = $self->{_data}->cellwalltype();
4094            #Setting global coefficients based on cell wall type
4095            my $biomassCompounds;
4096            my $compounds;
4097            if ($Class eq "Gram positive") {
4098                    $compounds->{RNA} = {cpd00002=>-0.262,cpd00012=>1,cpd00038=>-0.323,cpd00052=>-0.199,cpd00062=>-0.215};
4099                    $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};
4100                    $bioObj->protein("0.5284");
4101                    $bioObj->DNA("0.026");
4102                    $bioObj->RNA("0.0655");
4103                    $bioObj->lipid("0.075");
4104                    $bioObj->cellWall("0.25");
4105                    $bioObj->cofactor("0.10");
4106            } else {
4107                    $compounds->{RNA} = {cpd00002=>-0.262,cpd00012=>1,cpd00038=>-0.322,cpd00052=>-0.2,cpd00062=>-0.216};
4108                    $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};
4109                    $bioObj->protein("0.563");
4110                    $bioObj->DNA("0.031");
4111                    $bioObj->RNA("0.21");
4112                    $bioObj->lipid("0.093");
4113                    $bioObj->cellWall("0.177");
4114                    $bioObj->cofactor("0.039");
4115            }
4116            #Setting energy coefficient for all reactions
4117            $bioObj->energy("40");
4118            $compounds->{energy} = {cpd00002=>-1,cpd00001=>-1,cpd00008=>1,cpd00009=>1,cpd00067=>1};
4119            #Setting DNA coefficients based on GC content
4120            my $gc = $self->figmodel()->get_genome_gc_content($self->genome());
4121            $compounds->{DNA} = {cpd00012=>1,cpd00115=>0.5*(1-$gc),cpd00241=>0.5*$gc,cpd00356=>0.5*$gc,cpd00357=>0.5*(1-$gc)};
4122            #Setting Lipid,cell wall,and cofactor coefficients based on biomass template
4123            my $templateResults = $self->DetermineCofactorLipidCellWallComponents();
4124            $compounds->{cofactor} = $templateResults->{cofactor};
4125            $compounds->{lipid} = $templateResults->{lipid};
4126            $compounds->{cellWall} = $templateResults->{cellWall};
4127            #Getting package number for cofactor, lipid, and cell wall
4128            my $packages;
4129            my $cpdgrpMgr = $self->figmodel()->database()->get_object_manager("cpdgrp");
4130            my $packageTypes = ["Cofactor","Lipid","CellWall"];
4131            my $translation = {"Cofactor"=>"cofactor","Lipid"=>"lipid","CellWall"=>"cellWall"};
4132            for (my $i=0; $i < @{$packageTypes}; $i++) {
4133                    my @cpdList = keys(%{$compounds->{$translation->{$packageTypes->[$i]}}});
4134                    my $function = $translation->{$packageTypes->[$i]}."Package";
4135                    if (@cpdList == 0) {
4136                            $bioObj->$function("NONE");
4137                    } else {
4138                            my $cpdgrpObs = $cpdgrpMgr->get_objects({type=>$packageTypes->[$i]."Package"});
4139                            for (my $j=0; $j < @{$cpdgrpObs}; $j++) {
4140                                    $packages->{$packageTypes->[$i]}->{$cpdgrpObs->[$j]->grouping()}->{$cpdgrpObs->[$j]->COMPOUND()} = 1;
4141                            }
4142                            my @packageList = keys(%{$packages->{$packageTypes->[$i]}});
4143                            my $packageHash;
4144                            for (my $j=0; $j < @packageList; $j++) {
4145                                    $packageHash->{join("|",sort(keys(%{$packages->{$packageTypes->[$i]}->{$packageList[$j]}})))} = $packageList[$j];
4146                            }
4147                            if (defined($packageHash->{join("|",sort(keys(%{$compounds->{$translation->{$packageTypes->[$i]}}})))})) {
4148                                    $bioObj->$function($packageHash->{join("|",sort(keys(%{$compounds->{$translation->{$packageTypes->[$i]}}})))});
4149                            } else {
4150                                    my $newPackageID = $self->figmodel()->database()->check_out_new_id($packageTypes->[$i]."Pkg");
4151                                    $bioObj->$function($newPackageID);
4152                                    my @cpdList = keys(%{$compounds->{$translation->{$packageTypes->[$i]}}});
4153                                    for (my $j=0; $j < @cpdList; $j++) {
4154                                            $cpdgrpMgr = $self->figmodel()->database()->get_object_manager("cpdgrp");
4155                                            $cpdgrpMgr->create({COMPOUND=>$cpdList[$j],grouping=>$newPackageID,type=>$packageTypes->[$i]."Package"});
4156                                    }
4157                            }
4158                    }
4159            }
4160            #Filling in coefficient terms in database and calculating global reaction coefficients based on classification abundancies
4161            my $equationCompounds;
4162            my $types = ["RNA","DNA","protein","lipid","cellWall","cofactor","energy"];
4163            my $cpdMgr = $self->figmodel()->database()->get_object_manager("compound");
4164            for (my $i=0; $i < @{$types}; $i++) {
4165                    my $coefString = "";
4166                    my @compounds = sort(keys(%{$compounds->{$types->[$i]}}));
4167                    #Building coefficient strings and determining net mass for component types
4168                    my $netMass = 0;
4169                    for (my $j=0; $j < @compounds; $j++) {
4170                            my $objs = $cpdMgr->get_objects({id=>$compounds[$j]});
4171                            my $mass = 0;
4172                            if (defined($objs->[0]) && $objs->[0]->mass() != 0) {
4173                                    $mass = $objs->[0]->mass();
4174                                    $netMass += -$compounds->{$types->[$i]}->{$compounds[$j]}*$objs->[0]->mass();
4175                            }
4176                            if (!defined($equationCompounds->{$compounds[$j]})) {
4177                                    $equationCompounds->{$compounds[$j]}->{"coef"} = 0;
4178                                    $equationCompounds->{$compounds[$j]}->{"type"} = $types->[$i];
4179                                    $equationCompounds->{$compounds[$j]}->{"mass"} = $mass;
4180                            }
4181                            $coefString .= $compounds->{$types->[$i]}->{$compounds[$j]}."|";
4182                    }
4183                    $netMass = 0.001*$netMass;
4184                    #Calculating coefficients for all component compounds
4185                    for (my $j=0; $j < @compounds; $j++) {
4186                            #Normalizing certain type coefficients by mass
4187                            my $function = $types->[$i];
4188                            my $fraction = $bioObj->$function();
4189                            if ($types->[$i] ne "energy") {
4190                                    $fraction = $fraction/$netMass;
4191                            }
4192                            if ($compounds->{$types->[$i]}->{$compounds[$j]} eq 1e-5) {
4193                                    $fraction = 1;
4194                            }
4195                            $equationCompounds->{$compounds[$j]}->{"coef"} += $fraction*$compounds->{$types->[$i]}->{$compounds[$j]};
4196                    }
4197                    chop($coefString);
4198                    if (length($coefString) == 0) {
4199                            $coefString = "NONE";
4200                    }
4201                    my $function = $types->[$i]."Coef";
4202                    if ($types->[$i] ne "energy") {
4203                            $bioObj->$function($coefString);
4204                    }
4205            }
4206            #Adding biomass to compound list
4207            $equationCompounds->{cpd11416}->{coef} = 1;
4208            $equationCompounds->{cpd11416}->{type} = "macromolecule";
4209            #Building equation from hash and populating compound biomass table
4210            my @compoundList = keys(%{$equationCompounds});
4211            my ($reactants,$products);
4212            #Deleting existing Biomass Compound info
4213            my $cpdbofMgr = $self->figmodel()->database()->get_object_manager("cpdbof");
4214            my $matchingObjs = $cpdbofMgr->get_objects({BIOMASS=>$biomassID});
4215            for (my $i=0; $i < @{$matchingObjs}; $i++) {
4216                    $matchingObjs->[$i]->delete();
4217            }
4218            my $typeCategories = {"macromolecule"=>"M","RNA"=>"R","DNA"=>"D","protein"=>"P","lipid"=>"L","cellWall"=>"W","cofactor"=>"C","energy"=>"E"};
4219            my $productmass = 0;
4220            my $reactantmass = 0;
4221            my $totalmass = 0;
4222            foreach my $compound (@compoundList) {
4223                    if (defined($equationCompounds->{$compound}->{coef}) && defined($equationCompounds->{$compound}->{mass})) {
4224                            $totalmass += $equationCompounds->{$compound}->{coef}*0.001*$equationCompounds->{$compound}->{mass};
4225                    }
4226                    if ($equationCompounds->{$compound}->{coef} < 0) {
4227                            if (defined($equationCompounds->{$compound}->{coef}) && defined($equationCompounds->{$compound}->{mass})) {
4228                                    $reactantmass += $equationCompounds->{$compound}->{coef}*0.001*$equationCompounds->{$compound}->{mass};
4229                            }
4230                            $reactants->{$compound} = $self->figmodel()->format_coefficient(-1*$equationCompounds->{$compound}->{coef});
4231                    } else {
4232                            if (defined($equationCompounds->{$compound}->{coef}) && defined($equationCompounds->{$compound}->{mass})) {
4233                                    $productmass += $equationCompounds->{$compound}->{coef}*0.001*$equationCompounds->{$compound}->{mass};
4234                            }
4235                            $products->{$compound} = $self->figmodel()->format_coefficient($equationCompounds->{$compound}->{coef});
4236                    }
4237                    #Adding biomass reaction compounds to the biomass compound table
4238                    $cpdbofMgr = $self->figmodel()->database()->get_object_manager("cpdbof");
4239                    $cpdbofMgr->create({COMPOUND=>$compound,BIOMASS=>$biomassID,coefficient=>$equationCompounds->{$compound}->{coef},compartment=>"c",category=>$typeCategories->{$equationCompounds->{$compound}->{type}}});
4240                  }                  }
4241            print "Total mass = ".$totalmass.", Reactant mass = ".$reactantmass.", Product mass = ".$productmass."\n";
4242                  my $Equation = "";                  my $Equation = "";
4243                  my @ReactantList = sort(keys(%Reactants));          my @ReactantList = sort(keys(%{$reactants}));
4244                  for (my $i=0; $i < @ReactantList; $i++) {                  for (my $i=0; $i < @ReactantList; $i++) {
4245                          if (length($Equation) > 0) {                          if (length($Equation) > 0) {
4246                                  $Equation .= " + ";                                  $Equation .= " + ";
4247                          }                          }
4248                          $Equation .= $self->figmodel()->format_coefficient($Reactants{$ReactantList[$i]})." ".$ReactantList[$i];                  $Equation .= "(".$reactants->{$ReactantList[$i]}.") ".$ReactantList[$i];
4249                  }                  }
4250                  $Equation .= " => ";                  $Equation .= " => ";
4251                  my $First = 1;                  my $First = 1;
4252                  @ReactantList = sort(keys(%Products));          @ReactantList = sort(keys(%{$products}));
4253                  for (my $i=0; $i < @ReactantList; $i++) {                  for (my $i=0; $i < @ReactantList; $i++) {
4254                          if ($First == 0) {                          if ($First == 0) {
4255                                  $Equation .= " + ";                                  $Equation .= " + ";
4256                          }                          }
4257                          $First = 0;                          $First = 0;
4258                          $Equation .= $self->figmodel()->format_coefficient($Products{$ReactantList[$i]})." ".$ReactantList[$i];                  $Equation .= "(".$products->{$ReactantList[$i]}.") ".$ReactantList[$i];
4259            }
4260            $bioObj->equation($Equation);
4261            #Setting the biomass reaction of this model
4262            $self->biomassReaction($biomassID);
4263            $self->figmodel()->print_biomass_reaction_file($biomassID);
4264            #Checking if the biomass reaction remained unchanged
4265            if ($originalPackages ne "" && $originalPackages eq $bioObj->cofactorPackage().$bioObj->lipidPackage().$bioObj->cellWallPackage()) {
4266                    print "UNCHANGED!\n";
4267                    $bioObj->essentialRxn($originalEssReactions);
4268            } else {
4269                    #Copying essential reaction lists if the packages in this biomasses reaction exactly match those in another biomass reaction
4270                    my $matches = $bioMgr->get_objects({cofactorPackage=>$bioObj->cofactorPackage(),lipidPackage=>$bioObj->lipidPackage(),cellWallPackage=>$bioObj->cellWallPackage()});
4271                    my $matchFound = 0;
4272                    for (my $i=0; $i < @{$matches}; $i++) {
4273                            if ($matches->[$i]->id() ne $biomassID && defined($matches->[$i]->essentialRxn()) && length($matches->[$i]->essentialRxn())) {
4274                                    $bioObj->essentialRxn($matches->[$i]->essentialRxn());
4275                                    print "MATCH!\n";
4276                                    $matchFound = 1;
4277                                    last;
4278                            }
4279                    }
4280                    #Otherwise, we calculate essential reactions
4281                    if ($matchFound == 0) {
4282                            print "NOMATCH!\n";
4283                            $self->figmodel()->determine_biomass_essential_reactions($biomassID);
4284                  }                  }
   
                 #Adding the biomass equation to the biomass table  
                 my $NewRow = $self->figmodel()->add_biomass_reaction($Equation,$self->id(),"Template:".$self->genome());  
                 $biomassrxn = $NewRow;  
4285          }          }
4286          return $biomassrxn;          return $biomassID;
4287  }  }
4288    
4289  =head3 PrintSBMLFile  =head3 PrintSBMLFile
# Line 4070  Line 4294 
4294  =cut  =cut
4295  sub PrintSBMLFile {  sub PrintSBMLFile {
4296          my($self) = @_;          my($self) = @_;
   
4297          #Opening the SBML file for printing          #Opening the SBML file for printing
4298          my $Filename = $self->directory().$self->id().".xml";          my $Filename = $self->directory().$self->id().".xml";
4299          if (!open (SBMLOUTPUT, ">$Filename")) {          if (!open (SBMLOUTPUT, ">$Filename")) {
4300                  return;                  return;
4301          }          }
   
4302          #Loading and parsing the model data          #Loading and parsing the model data
4303          my $mdlTbl = $self->reaction_table();          my $mdlTbl = $self->reaction_table();
4304          if (!defined($mdlTbl) || !defined($mdlTbl->{"array"})) {          if (!defined($mdlTbl) || !defined($mdlTbl->{"array"})) {
4305                  return $self->fail();                  return $self->fail();
4306          }          }
4307          my $rxnTbl = $self->figmodel()->database()->get_table("REACTIONS");          my $rxnMgr = $self->figmodel()->database()->get_object_manager("reaction");
         my $bioTbl = $self->figmodel()->database()->get_table("BIOMASS");  
         my $cpdTbl = $self->figmodel()->database()->get_table("COMPOUNDS");  
4308          my $cmpTbl = $self->figmodel()->database()->get_table("COMPARTMENTS");          my $cmpTbl = $self->figmodel()->database()->get_table("COMPARTMENTS");
4309            my $cpdMgr = $self->figmodel()->database()->get_object_manager("compound");
4310            my $bioMgr = $self->figmodel()->database()->get_object_manager("bof");
4311          #Adding intracellular metabolites that also need exchange fluxes to the exchange hash          #Adding intracellular metabolites that also need exchange fluxes to the exchange hash
4312          my $ExchangeHash = {"cpd11416" => "c"};          my $ExchangeHash = {"cpd11416" => "c"};
4313          my %CompartmentsPresent;          my %CompartmentsPresent;
# Line 4094  Line 4315 
4315          my %CompoundList;          my %CompoundList;
4316          my @ReactionList;          my @ReactionList;
4317          for (my $i=0; $i < $mdlTbl->size(); $i++) {          for (my $i=0; $i < $mdlTbl->size(); $i++) {
4318                  my $Reaction = $mdlTbl->get_row($i)->{"LOAD"}->[0];                  my $rxnObj;
4319                  my $row = $rxnTbl->get_row_by_key($Reaction,"DATABASE");                  if ($mdlTbl->get_row($i)->{"LOAD"}->[0] =~ m/rxn\d\d\d\d\d/) {
4320                  if (!defined($row)) {                          $rxnObj = $rxnMgr->get_objects({id=>$mdlTbl->get_row($i)->{"LOAD"}->[0]})->[0];
4321                          $row = $bioTbl->get_row_by_key($Reaction,"DATABASE");                  } elsif ($mdlTbl->get_row($i)->{"LOAD"}->[0] =~ m/bio\d\d\d\d\d/) {
4322                          if (!defined($row)) {                          $rxnObj = $bioMgr->get_objects({id=>$mdlTbl->get_row($i)->{"LOAD"}->[0]})->[0];
                                 next;  
                         }  
4323                  }                  }
4324                  if (!defined($row->{"EQUATION"}->[0])) {                  if (!defined($rxnObj)) {
4325                          next;                          next;
4326                  }                  }
4327                  push(@ReactionList,$Reaction);                  push(@ReactionList,$rxnObj);
4328                  $_ = $row->{"EQUATION"}->[0];                  $_ = $rxnObj->equation();
4329                  my @MatchArray = /(cpd\d\d\d\d\d)/g;                  my @MatchArray = /(cpd\d\d\d\d\d)/g;
4330                  for (my $j=0; $j < @MatchArray; $j++) {                  for (my $j=0; $j < @MatchArray; $j++) {
4331                          $CompoundList{$MatchArray[$j]}->{"c"} = 1;                          $CompoundList{$MatchArray[$j]}->{"c"} = 1;
4332                  }                  }
4333                  $_ = $row->{"EQUATION"}->[0];                  $_ = $rxnObj->equation();
4334                  @MatchArray = /(cpd\d\d\d\d\d\[\D\])/g;                  @MatchArray = /(cpd\d\d\d\d\d\[\D\])/g;
4335                  for (my $j=0; $j < @MatchArray; $j++) {                  for (my $j=0; $j < @MatchArray; $j++) {
4336                          if ($MatchArray[$j] =~ m/(cpd\d\d\d\d\d)\[(\D)\]/) {                          if ($MatchArray[$j] =~ m/(cpd\d\d\d\d\d)\[(\D)\]/) {
# Line 4168  Line 4387 
4387          #Printing the list of metabolites involved in the model          #Printing the list of metabolites involved in the model
4388          print SBMLOUTPUT '<listOfSpecies>'."\n";          print SBMLOUTPUT '<listOfSpecies>'."\n";
4389          foreach my $Compound (keys(%CompoundList)) {          foreach my $Compound (keys(%CompoundList)) {
4390                  my $row = $cpdTbl->get_row_by_key($Compound,"DATABASE");                  my $cpdObj = $cpdMgr->get_objects({id=>$Compound})->[0];
4391                  if (!defined($row)) {                  if (!defined($cpdObj)) {
4392                          next;                          next;
4393                  }                  }
                 my $Name = $Compound;  
4394                  my $Formula = "";                  my $Formula = "";
4395                  if (defined($row->{"FORMULA"}->[0])) {                  if (defined($cpdObj->formula())) {
4396                          $Formula = $row->{"FORMULA"}->[0];                          $Formula = $cpdObj->formula();
4397                  }                  }
4398                  if (defined($row->{"NAME"}->[0])) {                  my $Name = $cpdObj->name();
                         $Name = $row->{"NAME"}->[0];  
4399                          $Name =~ s/\s/_/;                          $Name =~ s/\s/_/;
4400                          $Name .= "_".$Formula;                          $Name .= "_".$Formula;
                 }  
4401                  $Name =~ s/[<>;&\*]//;                  $Name =~ s/[<>;&\*]//;
4402                  my $Charge = 0;                  my $Charge = 0;
4403                  if (defined($row->{"CHARGE"}->[0])) {                  if (defined($cpdObj->charge())) {
4404                          $Charge = $row->{"CHARGE"}->[0];                          $Charge = $cpdObj->charge();
4405                  }                  }
4406                  foreach my $Compartment (keys(%{$CompoundList{$Compound}})) {                  foreach my $Compartment (keys(%{$CompoundList{$Compound}})) {
4407                          if ($Compartment eq "e") {                          if ($Compartment eq "e") {
# Line 4198  Line 4414 
4414    
4415          #Printing the boundary species          #Printing the boundary species
4416          foreach my $Compound (keys(%{$ExchangeHash})) {          foreach my $Compound (keys(%{$ExchangeHash})) {
4417                  my $Name = $Compound;                  my $cpdObj = $cpdMgr->get_objects({id=>$Compound})->[0];
4418                  my $Formula = "";                  if (!defined($cpdObj)) {
                 my $row = $cpdTbl->get_row_by_key($Compound,"DATABASE");  
                 if (!defined($row)) {  
4419                          next;                          next;
4420                  }                  }
4421                  if (defined($row->{"FORMULA"}->[0])) {                  my $Formula = "";
4422                          $Formula = $row->{"FORMULA"}->[0];                  if (defined($cpdObj->formula())) {
4423                            $Formula = $cpdObj->formula();
4424                  }                  }
4425                  if (defined($row->{"NAME"}->[0])) {                  my $Name = $cpdObj->name();
                         $Name = $row->{"NAME"}->[0];  
4426                          $Name =~ s/\s/_/;                          $Name =~ s/\s/_/;
4427                          $Name .= "_".$Formula;                          $Name .= "_".$Formula;
                 }  
4428                  $Name =~ s/[<>;&\*]//;                  $Name =~ s/[<>;&\*]//;
4429                  my $Charge = 0;                  my $Charge = 0;
4430                  if (defined($row->{"CHARGE"}->[0])) {                  if (defined($cpdObj->charge())) {
4431                          $Charge = $row->{"CHARGE"}->[0];                          $Charge = $cpdObj->charge();
4432                  }                  }
4433                  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";
4434          }          }
# Line 4225  Line 4438 
4438          my $ObjectiveCoef;          my $ObjectiveCoef;
4439          print SBMLOUTPUT '<listOfReactions>'."\n";          print SBMLOUTPUT '<listOfReactions>'."\n";
4440    
4441          foreach my $Reaction (@ReactionList) {          foreach my $rxnObj (@ReactionList) {
4442                  $ObjectiveCoef = "0.0";                  $ObjectiveCoef = "0.0";
4443                  my $mdlrow = $mdlTbl->get_row_by_key($Reaction,"LOAD");                  my $mdlrow = $mdlTbl->get_row_by_key($rxnObj->id(),"LOAD");
4444                  my $Reaction = $mdlrow->{"LOAD"}->[0];                  if ($rxnObj->id() =~ m/^bio/) {
                 my $row = $rxnTbl->get_row_by_key($Reaction,"DATABASE");  
                 if (!defined($row)) {  
                         $row = $bioTbl->get_row_by_key($Reaction,"DATABASE");  
                         if (!defined($row)) {  
                                 next;  
                         }  
                 }  
                 if (!defined($row->{"EQUATION"}->[0])) {  
                         next;  
                 }  
                 if ($Reaction =~ m/^bio/) {  
4445                          $ObjectiveCoef = "1.0";                          $ObjectiveCoef = "1.0";
4446                  }                  }
4447                  my $LowerBound = -10000;                  my $LowerBound = -10000;
4448                  my $UpperBound = 10000;                  my $UpperBound = 10000;
4449                  my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateDataFromEquation($row->{"EQUATION"}->[0]);                  my ($Reactants,$Products) = $self->figmodel()->GetReactionSubstrateDataFromEquation($rxnObj->equation());
4450                  my $Name = $Reaction;                  my $Name = $rxnObj->name();
                 if (defined($row->{"NAME"}->[0])) {  
                         $Name = $row->{"NAME"}->[0];  
4451                          $Name =~ s/[<>;&]//g;                          $Name =~ s/[<>;&]//g;
                 }  
4452                  my $Reversibility = "true";                  my $Reversibility = "true";
4453                  if (defined($mdlrow->{"DIRECTIONALITY"}->[0])) {                  if (defined($mdlrow->{"DIRECTIONALITY"}->[0])) {
4454                          if ($mdlrow->{"DIRECTIONALITY"}->[0] ne "<=>") {                          if ($mdlrow->{"DIRECTIONALITY"}->[0] ne "<=>") {
# Line 4262  Line 4461 
4461                                  $Reactants = $Temp;                                  $Reactants = $Temp;
4462                          }                          }
4463                  }                  }
4464                  print SBMLOUTPUT '<reaction id="'.$Reaction.'" name="'.$Name.'" reversible="'.$Reversibility.'">'."\n";                  print SBMLOUTPUT '<reaction id="'.$rxnObj->id().'" name="'.$Name.'" reversible="'.$Reversibility.'">'."\n";
4465                  print SBMLOUTPUT "<notes>\n";                  print SBMLOUTPUT "<notes>\n";
4466                  my $ECData = "";                  my $ECData = "";
4467                  if (defined($row->{"ENZYME"}->[0])) {                  if ($rxnObj->id() !~ m/^bio/) {
4468                          $ECData = $row->{"ENZYME"}->[0];                          if (defined($rxnObj->enzyme())) {
4469                                    my @ecList = split(/\|/,$rxnObj->enzyme());
4470                                    if (defined($ecList[1])) {
4471                                            $ECData = $ecList[1];
4472                                    }
4473                            }
4474                  }                  }
4475                  my $SubsystemData = "";                  my $SubsystemData = "";
4476                  if (defined($mdlrow->{"SUBSYSTEM"}->[0])) {                  if (defined($mdlrow->{"SUBSYSTEM"}->[0])) {
# Line 4331  Line 4535 
4535    
4536          my @ExchangeList = keys(%{$ExchangeHash});          my @ExchangeList = keys(%{$ExchangeHash});
4537          foreach my $ExCompound (@ExchangeList) {          foreach my $ExCompound (@ExchangeList) {
4538                  my $ExCompoundName = $ExCompound;                  my $cpdObj = $cpdMgr->get_objects({id=>$ExCompound})->[0];
4539                  my $Row = $cpdTbl->get_row_by_key($ExCompound,"DATABASE");                  if (!defined($cpdObj)) {
4540                  if (defined($Row) && defined($Row->{"NAME"}->[0])) {                          next;
                         $ExCompoundName = $Row->{"NAME"}->[0];  
                         $ExCompoundName =~ s/[<>;&]//g;  
4541                  }                  }
4542                    my $ExCompoundName = $cpdObj->name();
4543                    $ExCompoundName =~ s/[<>;&]//g;
4544                  $ObjectiveCoef = "0.0";                  $ObjectiveCoef = "0.0";
4545                  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";
4546                  print SBMLOUTPUT "\t".'<notes>'."\n";                  print SBMLOUTPUT "\t".'<notes>'."\n";
# Line 4417  Line 4621 
4621    
4622  =head3 patch_model  =head3 patch_model
4623  Definition:  Definition:
4624          FIGMODELTable:patch results = FIGMODELmodel->patch_model(FIGMODELTable:patch table);          FIGMODELmodel->patch_model([] -or- {} of patch arguments);
4625  Description:  Description:
4626  =cut  =cut
4627  sub patch_model {  sub patch_model {
4628          my ($self,$tbl) = @_;          my ($self,$arguments) = @_;
4629            if (($self->owner() eq "master" || $self->owner() eq "chenry") && $self->id() =~ m/Seed/) {
4630          #Instantiating table                  print "Model qualifies.";
4631          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();
4632          #Getting genome annotations                  my $reactionTable = $self->reaction_table();
4633          my $features = $self->figmodel()->database()->get_genome_feature_table($self->genome());                  my $row = $reactionTable->get_row_by_key("rxn05294","LOAD");
4634          #Gettubg reaction table                  if (defined($row)) {
4635          my $reactions = $self->reaction_table();                          $reactionTable->delete_row($reactionTable->row_index($row));
         #Checking for patched roles  
         for (my $i=0; $i < $tbl->size(); $i++) {  
                 my $row = $tbl->get_row($i);  
                 my @genes = $features->get_rows_by_key($row->{ROLE}->[0],"ROLES");  
                 if (@genes > 0) {  
                         for (my $j=0; $j < @{$row->{REACTIONS}};$j++) {  
                                 my $resultrxn = $results->get_row_by_key($row->{REACTIONS}->[$j],"Reactions");  
                                 if (!defined($resultrxn)) {  
                                         $resultrxn = $results->add_row({"Reactions"=>[$row->{REACTIONS}->[$j]],"Roles"=>[$row->{ROLE}->[0]]});  
                                 }  
                                 my $rxnrow = $reactions->get_row_by_key($row->{REACTIONS}->[$j],"LOAD");  
                                 if (defined($rxnrow) && !defined($resultrxn->{"Old genes"})) {  
                                         $resultrxn->{"Old genes"} = $rxnrow->{"ASSOCIATED PEG"};  
                                         if ($resultrxn->{"Old genes"}->[0] !~ m/GAP|BOF|UNIVERSAL|SPONTANEOUS/) {  
                                                 push(@{$resultrxn->{"Genes"}},@{$resultrxn->{"Old genes"}});  
                                         }  
                                 }  
                                 delete $resultrxn->{"Current gene set"};  
                                 if (defined($resultrxn->{"Genes"})) {  
                                         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"}});  
4636                          }                          }
4637                    $row = $reactionTable->get_row_by_key("rxn05295","LOAD");
4638                    if (defined($row)) {
4639                            $reactionTable->delete_row($reactionTable->row_index($row));
4640                  }                  }
4641                    $row = $reactionTable->get_row_by_key("rxn05296","LOAD");
4642                    if (defined($row)) {
4643                            $reactionTable->delete_row($reactionTable->row_index($row));
4644          }          }
4645                    $reactionTable->save();
4646          #Ensuring that the old model is preserved                  my $growth = $self->calculate_growth();
4647          $self->ArchiveModel();                  if ($growth =~ m/^NOGROWTH/ ||  $growth < 1e-5) {
4648          #Modifing the reaction list                          $self->error_message("patch_model:model did not grow, gapfilling!");
4649          for (my $i=0; $i < $results->size();$i++) {                          $self->GapFillModel(1);
4650                  my $row = $results->get_row($i);                          $growth = $self->calculate_growth();
4651                  my $rxnrow = $reactions->get_row_by_key($row->{"Reactions"}->[0],"LOAD");                          if ($growth =~ m/^NOGROWTH/ ||  $growth < 1e-5) {
4652                  if (defined($rxnrow)) {                                  $self->error_message("patch_model:model still did not grow!");
4653                          $rxnrow->{"ASSOCIATED PEG"} = $row->{"Genes"};                                  print "Patch failed.\n";
4654                  } else {                                  return;
                         $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"]});  
4655                  }                  }
4656          }          }
4657          $reactions->save();                  $self->PrintSBMLFile();
         $results->save();  
         $self->update_model_stats();  
4658          $self->PrintModelLPFile();          $self->PrintModelLPFile();
4659                    $self->PrintModelSimpleReactionTable();
4660          $self->run_default_model_predictions();          $self->run_default_model_predictions();
4661          #Returning results                  $self->update_model_stats();
4662          return $results;                  print "Patch complete.\n";
4663            }
4664  }  }
4665    
4666  =head3 translate_genes  =head3 translate_genes

Legend:
Removed from v.1.23  
changed lines
  Added in v.1.28

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3