Skip to content
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
192 changes: 65 additions & 127 deletions modelseedpy/core/mseditorapi.py
Original file line number Diff line number Diff line change
@@ -1,9 +1,10 @@
import logging
logger = logging.getLogger(__name__)
from cobra.core import Reaction
import cobra
import re

logger = logging.getLogger(__name__)

class MSEditorAPI:

@staticmethod
Expand All @@ -22,7 +23,7 @@ def remove_reactions(model, rxn_id_list = []):
# ASSUMPTIONS:
# an arrow will exist in the program, either =>, <=>, or <=
@staticmethod
def edit_reaction(model, rxn_id, direction=None, gpr=None, genome=None):
def edit_reaction(model, rxn_id, direction=None, gpr=None, genome=None): # !!! genome is never used
# Direction: =>, <=, or <=>
if direction is not None:
lower_bound = model.reactions.get_by_id(rxn_id).lower_bound
Expand Down Expand Up @@ -51,97 +52,80 @@ def edit_reaction(model, rxn_id, direction=None, gpr=None, genome=None):
if gpr is not None:
model.reactions.get_by_id(rxn_id).gene_reaction_rule = gpr
except:
raise Exception('invalid gpr statement, check parentheses') # not working, unsure exactly why
raise Exception(f'The gpr {gpr} is invalid. Perhaps check parentheses.') # not working, unsure exactly why

@staticmethod
def edit_biomass_compound(model,biomass_id,cpd_id,new_coef,rescale = 1):
def edit_biomass_compound(model, biomass_id, cpd_id, new_coef):
if cpd_id not in model.metabolites:
raise Exception('Metabolite', cpd_id, ' is not in the model.')
if biomass_id in model.reactions:
if cpd_id in model.metabolites:
model.reactions.get_by_id(biomass_id).add_metabolites({model.metabolites.get_by_id(cpd_id): new_coef},
combine=False)
else:
raise Exception('Metabolite', cpd_id, ' is not in the model.')
model.reactions.get_by_id(biomass_id).add_metabolites(
{model.metabolites.get_by_id(cpd_id): new_coef}, combine=False)
else: # if there is no biomass reaction
biomass_rxn = Reaction(biomass_id)
model.add_reaction(biomass_rxn)
if cpd_id in model.metabolites:
biomass_rxn.add_metabolites({model.metabolites.get_by_id(cpd_id): new_coef})
else:
raise Exception('Metabolite ', cpd_id, ' is not in the model.')
biomass_rxn.add_metabolites({model.metabolites.get_by_id(cpd_id): new_coef})

@staticmethod
def compute_molecular_weight(model, metabolite_id):
if metabolite_id not in model.metabolites:
raise Exception('Error, metabolite', metabolite_id, 'not found in the model')
raise ValueError(f'The metabolite {metabolite_id} is not defined in the model')
return model.metabolites.get_by_id(metabolite_id).formula_weight

@staticmethod

def add_custom_reaction(model,rxn_id,MSEquation,gpr = None,genome = None):
def add_custom_reaction(model, rxn_id, MSEquation, gpr=None, genome=None): # !!! genome is never used
new_rxn = Reaction(id=rxn_id)
# going on the assumption that all metabolites are present in the model
metabolites = {}
for key in MSEquation.equation:
met_id = key[0] + '_' + key[1]
for key, stoich in MSEquation.equation.items():
met_id = key[0]+'_'+key[1]
if met_id in model.metabolites:
metabolites[met_id] = MSEquation.equation[key]
metabolites[met_id] = stoich
else:
raise Exception("Error,", met_id, "not in model metabolites list")
raise ValueError(f"The {met_id} is not in model.")
model.add_reaction(new_rxn)
# new_rxn.gene_reaction_rule = gpr
new_rxn.add_metabolites(metabolites)
if gpr:
new_rxn.gene_reaction_rule = gpr

# adjust the bounds based on the arrow direction -1000, 1000, 0
if MSEquation.direction == 'left':
new_rxn.lower_bound = 0
new_rxn.upper_bound = 1000
if MSEquation.direction == '<':
new_rxn.lower_bound = -1000
new_rxn.upper_bound = 0
if MSEquation.direction == 'reversable':
elif MSEquation.direction == '=':
new_rxn.lower_bound = -1000
new_rxn.upper_bound = 1000

@staticmethod
def add_ms_reaction(model, rxn_id, modelseed, compartment_equivalents = {'0':'c0', '1':'e0'}, direction = 'forward'):#Andrew
def add_ms_reaction(model, rxn_id, modelseed, compartment_equivalents = {'0':'c0', '1':'e0'}, direction = '>'):#Andrew
''' Add a reaction with ModelSEED parameters to the FBA simulation
"model" (COBRA object): The metabolic model that is defined by COBRApy
"rxn_id" (Python object, string): The ModelSEED reaction id that will be added to the model
"Compartment_equivalents" (Python object, dictionary): The compartment codes that are used in formating the passed reactions
"direction" (Python object, string): The direction of the defined reaction
"modelseed" (ModelSEED object): The ModelSEED database that describes the reaction and metabolites of the argument
'''
modelseed_reaction = modelseed.get_seed_reaction(rxn_id)
reaction_stoich = modelseed_reaction.cstoichiometry
cobra_reaction = cobra.Reaction(rxn_id)
modelseed_reaction = modelseed.get_seed_reaction(rxn_id)
cobra_reaction.name = modelseed_reaction.data['name']
reaction_stoich = modelseed_reaction.cstoichiometry

metabolites_to_add = {}
for metabolite, stoich in reaction_stoich.items():
id = metabolite[0]
compound = modelseed.get_seed_compound(id).data
compartment_number = metabolite[1]
compartment_string = compartment_equivalents[compartment_number]

metabolites_to_add[cobra.Metabolite(id, name = compound['name'], compartment = compartment_string)] = stoich
metabolites_to_add[cobra.Metabolite(
metabolite[0], name = modelseed.get_seed_compound(metabolite[0]).data['name'],
compartment = compartment_equivalents[metabolite[1]])
] = stoich

cobra_reaction.add_metabolites(metabolites_to_add)
cobra_reaction.reaction

if direction == 'reversible':
cobra_reaction.lower_bound = 0
cobra_reaction.upper_bound = 1000
if direction == '=':
cobra_reaction.lower_bound = -1000
elif direction == 'backward':
elif direction == '<':
cobra_reaction.lower_bound = -1000
cobra_reaction.upper_bound = 0
elif direction == 'forward':
pass
else:
directions = ['forward', 'backward', 'reversible']
print('ERROR: The \'direction\' argument is not among the accepted {}, {}, and {} values.'.format(directions[0], directions[1], directions[2]))
direction = input('What is the direction of the reaction?')
while direction not in directions:
print('ERROR: The \'direction\' argument is not among the accepted {}, {}, and {} values.'.format(directions[0], directions[1], directions[2]))
direction = input('What is the direction of the reaction?')
MSEditorAPI.add_ms_reaction(model, rxn_id, modelseed, direction = direction)


model.add_reactions([cobra_reaction])

@staticmethod
Expand All @@ -150,113 +134,67 @@ def copy_model_reactions(model,source_model,rxn_id_list = []):
if rxnid in source_model.reactions:
model.add_reactions([source_model.reactions.get_by_id(rxnid)])
else:
raise Exception('The reaction', rxnid, 'in the reaction list is not in the model.')
raise ValueError(f'The {rxnid} reaction ID is not in the source model, and thus cannot be added to the model.')

@staticmethod
def copy_all_model_reactions(model,source_model): #new method that copies all reactions, may not be necessary
for rxn in source_model.reactions:
model.add_reactions([source_model.reactions.get_by_id(rxn.id)])

model.add_reactions([source_model.reactions.get_by_id(rxn.id) for rxn in source_model.reactions if rxn not in model.reactions])

class MSEquation:

def __init__(self, stoichiometry, d):
def __init__(self, stoichiometry, direction = None):
self.equation = stoichiometry
self.direction = d
self.direction = direction

@staticmethod
def build_from_palsson_string(equation_string, default_group='c'): # add default group

def clean_ends(lst):
"""
FIXME: use .strip
:param lst:
:return:
"""
for i in range(len(lst)):
# remove whitespace from the front
while lst[i][0] == ' ':
lst[i] = lst[i][1:]
# remove whitespace from the back
while lst[i][-1] == ' ':
lst[i] = lst[i][:-1]
return lst
return [i.strip() for i in lst]

def get_coef_and_group(lst, return_dict, side):
# for side variable, -1 is left side, 1 is right side, for coeficients
for reactant in lst:
for reagent in lst:
coeficient = side
if reactant.find('(') != -1 or reactant.find(
')') != -1: # if one is present, check to make sure both there
if equation_string.find(')') == -1:
raise Exception("Error, ')' character missing in string", reactant)
if equation_string.find('(') == -1:
raise Exception("Error, '(' character missing in string", reactant)
identifier = default_group
if '(' in reagent and ')' in reagent:
number = ''
position = 1
while reactant[position] != ')':
number += reactant[position]
while reagent[position] != ')':
number += reagent[position]
position += 1
coeficient = side * float(number)
reactant = reactant[position + 1:]

identifier = default_group
if reactant.find('[') != -1 or reactant.find(
']') != -1: # if one is present, check to make sure both there
# check to see both are present
if equation_string.find(']') == -1:
raise Exception("Error, ']' character missing in string", reactant)
if equation_string.find('[') == -1:
raise Exception("Error, '[' character missing in string", reactant)
reagent = reagent[position+1: ]
elif '[' in reagent and ']' in reagent:
s = ''
position = -2
while reactant[position] != '[':
s = reactant[position] + s
while reagent[position] != '[':
s = reagent[position] + s
position -= 1
identifier = s
reactant = reactant[:position]

return_dict[(reactant.strip(), identifier)] = coeficient
reagent = reagent[:position]
elif any([x in reagent for x in ['(', ')', '[', ']']]):
raise ValueError("A closing or opening parentheses or bracket is missing in the reaction string", reagent)
return_dict[(reagent.strip(), identifier)] = coeficient
return return_dict

# check for the '=' character, throw exception otherwise
if equation_string.find('=') == -1:
raise Exception("Error, '=' character missing, unable to split string", equation_string)

# check direction
reversible = False
right = False
left = False
ret_str = ''
reversible = equation_string.find('<=>') != -1
if reversible:
ret_str = '='
else: # not two ways, so check right
right = equation_string.find('=>') != -1
if right:
ret_str = '>'
else: # if not right, check left
left = equation_string.find('<=') != -1
if left:
ret_str = '<'
else: # if not left, error
ret_str = "?"
if '=' not in equation_string:
raise ValueError(f"Error: The '=' character is missing; hence, the reaction string {equation_string} cannot be split.")
if '<=>' in equation_string:
direction = '='
elif '=>' in equation_string:
direction = '>'
elif '<=' in equation_string:
direction = '<'
else:
direction = '?'

# cpd00001 + cpd00002[e] => (2)cpd00003 + cpd00004
# get substrings for either side of the euqation
# get substrings for either side of the equation
reactants_substring_list = equation_string[0:equation_string.find('=') - 1].split('+')
products_substring_list = equation_string[equation_string.find('=') + 2:len(equation_string)].split('+')

# clean up our substrings:
reactants_substring_list = list(map(lambda x: x.strip(), reactants_substring_list))
products_substring_list = list(map(lambda x: x.strip(), products_substring_list))

variables = {}
# add reactants to the dictionary

get_coef_and_group(reactants_substring_list, variables, -1)
get_coef_and_group(products_substring_list, variables, 1)

ret_mse = MSEquation(variables, ret_str)

return ret_mse
reactant_dict = get_coef_and_group([x.strip() for x in reactants_substring_list], {}, -1)
products_dict = get_coef_and_group([x.strip() for x in products_substring_list], reactant_dict, 1)
return MSEquation(reactant_dict.update(products_dict), direction)