From be602badd6e6c53e7b08d27147577a3a971f09d8 Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Fri, 1 Jul 2022 19:23:53 -0400 Subject: [PATCH 1/6] polishing MSBuilder --- modelseedpy/core/msbuilder.py | 218 +++++++++++++++------------------- 1 file changed, 94 insertions(+), 124 deletions(-) diff --git a/modelseedpy/core/msbuilder.py b/modelseedpy/core/msbuilder.py index 54b854b1..a08caf2a 100755 --- a/modelseedpy/core/msbuilder.py +++ b/modelseedpy/core/msbuilder.py @@ -1,17 +1,17 @@ import logging -import itertools +import itertools # !!! this import is not used +logger = logging.getLogger(__name__) + import cobra -from modelseedpy.core.rast_client import RastClient +from modelseedpy.core.rast_client import RastClient # !!! this import is not used from modelseedpy.core.msgenome import normalize_role -from modelseedpy.core.msmodel import get_gpr_string, get_reaction_constraints_from_direction -from cobra.core import Gene, Metabolite, Model, Reaction -from modelseedpy.core import FBAHelper +from modelseedpy.core.msmodel import get_gpr_string, get_reaction_constraints_from_direction # !!! none of these imports are used +from cobra.core import Gene, Metabolite, Model, Reaction # !!! Gene is not used +from modelseedpy.core.fbahelper import FBAHelper from modelseedpy.fbapkg.mspackagemanager import MSPackageManager SBO_ANNOTATION = "sbo" -logger = logging.getLogger(__name__) - ### temp stuff ### core_biomass = { 'cpd00032_c0': -1.7867, @@ -251,18 +251,21 @@ def build_biomass(rxn_id, cobra_model, template, biomass_compounds, index='0'): def _aaaa(genome, ontology_term): - search_name_to_genes = {} - search_name_to_orginal = {} - for f in genome.features: - if ontology_term in f.ontology_terms: - functions = f.ontology_terms[ontology_term] - for function in functions: + """ + A search function + + """ + + search_name_to_genes, search_name_to_orginal = {}, {} + for feature in genome.features: + if ontology_term in feature.ontology_terms: + for function in feature.ontology_terms[ontology_term]: f_norm = normalize_role(function) if f_norm not in search_name_to_genes: search_name_to_genes[f_norm] = set() search_name_to_orginal[f_norm] = set() search_name_to_orginal[f_norm].add(function) - search_name_to_genes[f_norm].add(f.id) + search_name_to_genes[f_norm].add(feature.id) return search_name_to_genes, search_name_to_orginal @@ -294,16 +297,27 @@ def build_gpr2(cpx_sets): return ' or '.join(list_of_ors) return list_of_ors[0] -def build_gpr(cpx_gene_role): +def _reaction_sinks(self,model): + reactions_sinks = [] + for cpd_id in ['cpd02701_c0', 'cpd11416_c0', 'cpd15302_c0']: + if cpd_id in model.metabolites: + met = model.metabolites.get_by_id(cpd_id) + rxn_exchange = Reaction('SK_'+met.id, 'Sink for '+met.name, 'exchanges', 0, 1000) + rxn_exchange.add_metabolites({met: -1}) + rxn_exchange.annotation[SBO_ANNOTATION] = "SBO:0000627" + reactions_sinks.append(rxn_exchange) + return reactions_sinks + +def build_gpr(cpx_gene_role): #!!! Unused and redundant function """ example input: - {'sdh': [{'b0721': 'sdhC', 'b0722': 'sdhD', 'b0723': 'sdhA', 'b0724': 'sdhB'}]} + {'sdh': [{'b0721': 'sdhC', 'b0722': 'sdhD', 'b0723': 'sdhA', 'b0724': 'sdhB'}]} - (b0721 and b0722 and b0724 and b0723) + (b0721 and b0722 and b0724 and b0723) - {'cpx1': [{'g1': 'role1', 'g3': 'role2'}, {'g2': 'role1', 'g3': 'role2'}]} + {'cpx1': [{'g1': 'role1', 'g3': 'role2'}, {'g2': 'role1', 'g3': 'role2'}]} - (g1 and g3) or (g2 and g3) + (g1 and g3) or (g2 and g3) :param cpx_gene_role: :return: @@ -320,18 +334,15 @@ def build_gpr(cpx_gene_role): class MSBuilder: - def __init__(self, genome, template=None): """ for future methods with better customization """ - self.genome = genome - self.template = template + self.genome = genome; self.template = template self.search_name_to_genes, self.search_name_to_original = _aaaa(genome, 'RAST') def _get_template_reaction_complexes(self, template_reaction): """ - :param template_reaction: :return: """ @@ -339,12 +350,10 @@ def _get_template_reaction_complexes(self, template_reaction): for cpx in template_reaction.get_complexes(): template_reaction_complexes[cpx.id] = {} for role, (triggering, optional) in cpx.roles.items(): - sn = normalize_role(role.name) + rn_norm = normalize_role(role.name) template_reaction_complexes[cpx.id][role.id] = [ - sn, - triggering, - optional, - set() if sn not in self.search_name_to_genes else set(self.search_name_to_genes[sn]) + rn_norm, triggering, optional, + set() if rn_norm not in self.search_name_to_genes else set(self.search_name_to_genes[rn_norm]) ] return template_reaction_complexes @@ -356,11 +365,11 @@ def _build_reaction_complex_gpr_sets2(match_complex, allow_incomplete_complexes= roles = set() role_genes = {} for role_id in match_complex[cpx_id]: - t = match_complex[cpx_id][role_id] - complete &= len(t[3]) > 0 or not t[1] or t[2] - if len(t[3]) > 0: + complx = match_complex[cpx_id][role_id] + complete &= len(complx[3]) > 0 or not complx[1] or complx[2] + if len(complx[3]) > 0: roles.add(role_id) - role_genes[role_id] = t[3] + role_genes[role_id] = complx[3] #print(cpx_id, complete, roles) if len(roles) > 0 and (allow_incomplete_complexes or complete): complexes[cpx_id] = {} @@ -369,7 +378,7 @@ def _build_reaction_complex_gpr_sets2(match_complex, allow_incomplete_complexes= #print(role_id, role_genes[role_id]) return complexes - @staticmethod + @staticmethod #!!! Unused and redundant function def _build_reaction_complex_gpr_sets(match_complex, allow_incomplete_complexes=True): complexes = {} for cpx_id in match_complex: @@ -377,12 +386,12 @@ def _build_reaction_complex_gpr_sets(match_complex, allow_incomplete_complexes=T roles = set() role_genes = {} for role_id in match_complex[cpx_id]: - t = match_complex[cpx_id][role_id] - complete &= len(t[3]) > 0 or not t[1] or t[2] # true if has genes or is not triggering or is optional + complx = match_complex[cpx_id][role_id] + complete &= len(complx[3]) > 0 or not complx[1] or complx[2] # true if has genes or is not triggering or is optional # print(t[3]) - if len(t[3]) > 0: + if len(complx[3]) > 0: roles.add(role_id) - role_genes[role_id] = t[3] + role_genes[role_id] = complx[3] # print(t) # it is never complete if has no genes, only needed if assuming a complex can have all # roles be either non triggering or optional @@ -404,38 +413,31 @@ def get_gpr_from_template_reaction(self, template_reaction, allow_incomplete_com template_reaction_complexes = self._get_template_reaction_complexes(template_reaction) if len(template_reaction_complexes) == 0: return None - # self.map_gene(template_reaction_complexes) - # print(template_reaction_complexes) gpr_set = self._build_reaction_complex_gpr_sets2(template_reaction_complexes, allow_incomplete_complexes) return gpr_set @staticmethod - def _build_reaction(reaction_id, gpr_set, template, index='0', sbo=None): + def _build_reaction(reaction_id, template, gpr_set = [], index='0', sbo=None): template_reaction = template.reactions.get_by_id(reaction_id) - reaction_compartment = template_reaction.compartment metabolites = {} - for cpd, value in template_reaction.metabolites.items(): + for cpd, stoich in template_reaction.metabolites.items(): compartment = f"{cpd.compartment}{index}" - name = f"{cpd.name}_{compartment}" - cpd = Metabolite(cpd.id + str(index), cpd.formula, name, cpd.charge, compartment) - metabolites[cpd] = value + met = Metabolite(cpd.id+index, cpd.formula, f"{cpd.name}_{compartment}", cpd.charge, compartment) + metabolites[met] = stoich - reaction = Reaction( - "{}{}".format(template_reaction.id, index), - "{}_{}{}".format(template_reaction.name, reaction_compartment, index), - '', + reaction = Reaction(template_reaction.id+index, + f'{template_reaction.name}_{reaction_compartment}{index}', '', template_reaction.lower_bound, template_reaction.upper_bound ) - gpr_str = build_gpr2(gpr_set) if gpr_set else '' + gpr_str = build_gpr2(gpr_set) reaction.add_metabolites(metabolites) - if gpr_str and len(gpr_str) > 0: - reaction.gene_reaction_rule = gpr_str # get_gpr_string(gpr_ll) - reaction.annotation["seed.reaction"] = template_reaction.reference_id + if len(gpr_str) > 0: + reaction.gene_reaction_rule = gpr_str # get_gpr_string(gpr_ll) if sbo: reaction.annotation[SBO_ANNOTATION] = sbo return reaction @@ -449,29 +451,28 @@ def build_exchanges(model, extra_cell='e0'): :return: """ reactions_exchanges = [] - for m in model.metabolites: - if m.compartment == extra_cell: - rxn_exchange_id = 'EX_' + m.id + for met in model.metabolites: + if met.compartment == extra_cell: + rxn_exchange_id = 'EX_' + met.id if rxn_exchange_id not in model.reactions: - rxn_exchange = Reaction(rxn_exchange_id, 'Exchange for ' + m.name, 'exchanges', -1000, 1000) - rxn_exchange.add_metabolites({m: -1}) + rxn_exchange = Reaction(rxn_exchange_id, 'Exchange for '+met.name, 'exchanges', -1000, 1000) + rxn_exchange.add_metabolites({met: -1}) rxn_exchange.annotation[SBO_ANNOTATION] = "SBO:0000627" reactions_exchanges.append(rxn_exchange) model.add_reactions(reactions_exchanges) - return reactions_exchanges @staticmethod def build_biomasses(model, template, index): - res = [] + biomass_reactions = [] if template.name.startswith('CoreModel'): - res.append(build_biomass('bio1', model, template, core_biomass, index)) - res.append(build_biomass('bio2', model, template, core_atp, index)) + biomass_reactions.append(build_biomass('bio1', model, template, core_biomass, index)) + biomass_reactions.append(build_biomass('bio2', model, template, core_atp, index)) if template.name.startswith('GramNeg'): - res.append(build_biomass('bio1', model, template, gramneg, index)) + biomass_reactions.append(build_biomass('bio1', model, template, gramneg, index)) if template.name.startswith('GramPos'): - res.append(build_biomass('bio1', model, template, grampos, index)) - return res + biomass_reactions.append(build_biomass('bio1', model, template, grampos, index)) + return biomass_reactions def auto_select_template(self): """ @@ -480,9 +481,7 @@ def auto_select_template(self): """ from modelseedpy.helpers import get_template, get_classifier from modelseedpy.core.mstemplate import MSTemplateBuilder - genome_classifier = get_classifier('knn_ACNP_RAST_filter') - genome_class = genome_classifier.classify(self.genome) - + genome_class = get_classifier('knn_ACNP_RAST_filter').classify(self.genome) template_genome_scale_map = { 'A': 'template_gram_neg', 'C': 'template_gram_neg', @@ -499,8 +498,7 @@ def auto_select_template(self): if genome_class in template_genome_scale_map and genome_class in template_core_map: self.template = MSTemplateBuilder.from_dict(get_template(template_genome_scale_map[genome_class])).build() elif self.template is None: - raise Exception(f'unable to select template for {genome_class}') - + raise Exception(f'A template is not defined for the genome class {genome_class}.') return genome_class def build_metabolic_reactions(self, index='0', allow_incomplete_complexes=True): @@ -511,58 +509,41 @@ def build_metabolic_reactions(self, index='0', allow_incomplete_complexes=True): metabolic_reactions[template_reaction.id] = gpr_set logger.debug("[%s] gpr set: %s", template_reaction.id, gpr_set) - reactions = list(map( - lambda x: self._build_reaction(x[0], x[1], self.template, index, "SBO:0000176"), - metabolic_reactions.items())) - - return reactions + return [self._build_reaction(key, self.template, val, index, "SBO:0000176") + for key, val in metabolic_reactions.items()] - def build_non_metabolite_reactions(self, cobra_model, index='0', allow_all_non_grp_reactions=False): - reactions_no_gpr = [] - reactions_in_model = set(map(lambda x: x.id, cobra_model.reactions)) - metabolites_in_model = set(map(lambda x: x.id, cobra_model.metabolites)) + def build_non_metabolite_reactions(self, cobra_model, index='0', allow_all_non_grp_reactions=False): #!!! allow_all_non_grp_reactions in not used + reactions_sans_gpr = [] + model_rxns = set(x.id for x in cobra_model.reactions) + model_mets = set(x.id for x in cobra_model.metabolites) for rxn in self.template.reactions: - if rxn.type == 'universal' or rxn.type == 'spontaneous': - reaction = self._build_reaction(rxn.id, {}, self.template, index, "SBO:0000176") - reaction_metabolite_ids = set(map(lambda x: x.id, set(reaction.metabolites))) - if (len(metabolites_in_model & reaction_metabolite_ids) > 0 or allow_all_non_grp_reactions) and \ - reaction.id not in reactions_in_model: - reaction.annotation["seed.reaction"] = rxn.id - reactions_no_gpr.append(reaction) + if rxn.type in ['universal', 'spontaneous']: + reaction = self._build_reaction(rxn.id, self.template, [], index, "SBO:0000176") + if len(model_mets.intersection(x.id for x in set(reaction.metabolites))) > 0: + if reaction.id not in model_rxns: + reaction.annotation["seed.reaction"] = rxn.id + reactions_sans_gpr.append(reaction) - return reactions_no_gpr + return reactions_sans_gpr def build(self, model_id, index='0', allow_all_non_grp_reactions=False, annotate_with_rast=True): - if annotate_with_rast: - rast = RastClient() - res = rast.annotate_genome(self.genome) self.search_name_to_genes, self.search_name_to_original = _aaaa(self.genome, 'RAST') - # rxn_roles = aux_template(self.template) # needs to be fixed to actually reflect template GPR rules - if self.template is None: + if not self.template: self.auto_select_template() + # construct the model cobra_model = Model(model_id) cobra_model.add_reactions(self.build_metabolic_reactions(index=index)) - cobra_model.add_reactions(self.build_non_metabolite_reactions(cobra_model, index, allow_all_non_grp_reactions)) + if allow_all_non_grp_reactions: + cobra_model.add_reactions(self.build_non_metabolite_reactions(cobra_model, index)) self.build_exchanges(cobra_model) - if self.template.name.startswith('CoreModel') or \ - self.template.name.startswith('GramNeg') or self.template.name.startswith('GramPos'): + if any([self.template.name.startswith(x) for x in ('GramPos', 'CoreModel', 'GramNeg')]): cobra_model.add_reactions(self.build_biomasses(cobra_model, self.template, index)) cobra_model.objective = 'bio1' - - reactions_sinks = [] - for cpd_id in ['cpd02701_c0', 'cpd11416_c0', 'cpd15302_c0']: - if cpd_id in cobra_model.metabolites: - m = cobra_model.metabolites.get_by_id(cpd_id) - rxn_exchange = Reaction('SK_' + m.id, 'Sink for ' + m.name, 'exchanges', 0, 1000) - rxn_exchange.add_metabolites({m: -1}) - rxn_exchange.annotation[SBO_ANNOTATION] = "SBO:0000627" - reactions_sinks.append(rxn_exchange) - cobra_model.add_reactions(reactions_sinks) - + cobra_model.add_reactions(_reaction_sinks(cobra_model)) return cobra_model @staticmethod @@ -574,10 +555,10 @@ def build_full_template_model(template, model_id=None, index='0'): :param index: index for the metabolites :return: """ - model = Model(model_id if model_id else template.id) + model = Model(model_id or template.id) all_reactions = [] for rxn in template.reactions: - reaction = MSBuilder._build_reaction(rxn.id, {}, template, index, "SBO:0000176") + reaction = MSBuilder._build_reaction(rxn.id, template, [], index, "SBO:0000176") reaction.annotation["seed.reaction"] = rxn.id all_reactions.append(reaction) model.add_reactions(all_reactions) @@ -588,24 +569,15 @@ def build_full_template_model(template, model_id=None, index='0'): bio_rxn2 = build_biomass('bio2', model, template, core_atp, index) model.add_reactions([bio_rxn1, bio_rxn2]) model.objective = 'bio1' - if template.name.startswith('GramNeg'): + elif template.name.startswith('GramNeg'): bio_rxn1 = build_biomass('bio1', model, template, gramneg, index) model.add_reactions([bio_rxn1]) model.objective = 'bio1' - if template.name.startswith('GramPos'): + elif template.name.startswith('GramPos'): bio_rxn1 = build_biomass('bio1', model, template, grampos, index) model.add_reactions([bio_rxn1]) model.objective = 'bio1' - - reactions_sinks = [] - for cpd_id in ['cpd02701_c0', 'cpd11416_c0', 'cpd15302_c0']: - if cpd_id in model.metabolites: - m = model.metabolites.get_by_id(cpd_id) - rxn_exchange = Reaction('SK_' + m.id, 'Sink for ' + m.name, 'exchanges', 0, 1000) - rxn_exchange.add_metabolites({m: -1}) - rxn_exchange.annotation[SBO_ANNOTATION] = "SBO:0000627" - reactions_sinks.append(rxn_exchange) - model.add_reactions(reactions_sinks) + model.add_reactions(_reaction_sinks(model)) return model @staticmethod @@ -613,7 +585,6 @@ def build_metabolic_model(model_id, genome, gapfill_media=None, template=None, i allow_all_non_grp_reactions=False, annotate_with_rast=True, gapfill_model=True): builder = MSBuilder(genome, template) model = builder.build(model_id, index, allow_all_non_grp_reactions, annotate_with_rast) - # Gapfilling model if gapfill_model: model = MSBuilder.gapfill_model(model, 'bio1', builder.template, gapfill_media) return model @@ -621,7 +592,7 @@ def build_metabolic_model(model_id, genome, gapfill_media=None, template=None, i @staticmethod def gapfill_model(original_mdl, target_reaction, template, media): FBAHelper.set_objective_from_target_reaction(original_mdl, target_reaction) - model = cobra.io.json.from_json(cobra.io.json.to_json(original_mdl)) + model = cobra.io.json.from_json(cobra.io.json.to_json(original_mdl)) #!!! what is the benefit of this I/O processing? pkgmgr = MSPackageManager.get_pkg_mgr(model) pkgmgr.getpkg("GapfillingPkg").build_package({ "default_gapfill_templates": [template], @@ -632,7 +603,6 @@ def gapfill_model(original_mdl, target_reaction, template, media): pkgmgr.getpkg("KBaseMediaPkg").build_package(media) #with open('Gapfilling.lp', 'w') as out: # out.write(str(model.solver)) - sol = model.optimize() gfresults = pkgmgr.getpkg("GapfillingPkg").compute_gapfilled_solution() for rxnid in gfresults["reversed"]: rxn = original_mdl.reactions.get_by_id(rxnid) @@ -650,4 +620,4 @@ def gapfill_model(original_mdl, target_reaction, template, media): else: rxn.upper_bound = 0 rxn.lower_bound = -100 - return original_mdl \ No newline at end of file + return original_mdl From 29d40559e5929af0e2efff15c5882712fffe820f Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Fri, 1 Jul 2022 19:43:44 -0400 Subject: [PATCH 2/6] polishing MSBuilder --- modelseedpy/core/msbuilder.py | 36 +++++++++++++++-------------------- 1 file changed, 15 insertions(+), 21 deletions(-) diff --git a/modelseedpy/core/msbuilder.py b/modelseedpy/core/msbuilder.py index a08caf2a..1dd65f5f 100755 --- a/modelseedpy/core/msbuilder.py +++ b/modelseedpy/core/msbuilder.py @@ -249,13 +249,8 @@ def build_biomass(rxn_id, cobra_model, template, biomass_compounds, index='0'): bio_rxn.annotation[SBO_ANNOTATION] = "SBO:0000629" return bio_rxn - +# A search function def _aaaa(genome, ontology_term): - """ - A search function - - """ - search_name_to_genes, search_name_to_orginal = {}, {} for feature in genome.features: if ontology_term in feature.ontology_terms: @@ -284,10 +279,9 @@ def aux_template(template): def build_gpr2(cpx_sets): list_of_ors = [] - for cpx in cpx_sets: + for cpx, role_ids in cpx_sets.items(): list_of_ands = [] - for role_id in cpx_sets[cpx]: - gene_ors = cpx_sets[cpx][role_id] + for role_id, gene_ors in role_ids.items(): if len(gene_ors) > 1: list_of_ands.append('(' + ' or '.join(gene_ors) + ')') else: @@ -311,13 +305,13 @@ def _reaction_sinks(self,model): def build_gpr(cpx_gene_role): #!!! Unused and redundant function """ example input: - {'sdh': [{'b0721': 'sdhC', 'b0722': 'sdhD', 'b0723': 'sdhA', 'b0724': 'sdhB'}]} + {'sdh': [{'b0721': 'sdhC', 'b0722': 'sdhD', 'b0723': 'sdhA', 'b0724': 'sdhB'}]} - (b0721 and b0722 and b0724 and b0723) + (b0721 and b0722 and b0724 and b0723) - {'cpx1': [{'g1': 'role1', 'g3': 'role2'}, {'g2': 'role1', 'g3': 'role2'}]} + {'cpx1': [{'g1': 'role1', 'g3': 'role2'}, {'g2': 'role1', 'g3': 'role2'}]} - (g1 and g3) or (g2 and g3) + (g1 and g3) or (g2 and g3) :param cpx_gene_role: :return: @@ -334,6 +328,7 @@ def build_gpr(cpx_gene_role): #!!! Unused and redundant function class MSBuilder: + def __init__(self, genome, template=None): """ for future methods with better customization @@ -418,7 +413,7 @@ def get_gpr_from_template_reaction(self, template_reaction, allow_incomplete_com return gpr_set @staticmethod - def _build_reaction(reaction_id, template, gpr_set = [], index='0', sbo=None): + def _build_reaction(reaction_id, gpr_set, template, index='0', sbo=None): template_reaction = template.reactions.get_by_id(reaction_id) reaction_compartment = template_reaction.compartment metabolites = {} @@ -430,13 +425,12 @@ def _build_reaction(reaction_id, template, gpr_set = [], index='0', sbo=None): reaction = Reaction(template_reaction.id+index, f'{template_reaction.name}_{reaction_compartment}{index}', '', - template_reaction.lower_bound, template_reaction.upper_bound - ) + template_reaction.lower_bound, template_reaction.upper_bound) - gpr_str = build_gpr2(gpr_set) + gpr_str = build_gpr2(gpr_set) if gpr_set else '' reaction.add_metabolites(metabolites) reaction.annotation["seed.reaction"] = template_reaction.reference_id - if len(gpr_str) > 0: + if gpr_str and len(gpr_str) > 0: reaction.gene_reaction_rule = gpr_str # get_gpr_string(gpr_ll) if sbo: reaction.annotation[SBO_ANNOTATION] = sbo @@ -509,7 +503,7 @@ def build_metabolic_reactions(self, index='0', allow_incomplete_complexes=True): metabolic_reactions[template_reaction.id] = gpr_set logger.debug("[%s] gpr set: %s", template_reaction.id, gpr_set) - return [self._build_reaction(key, self.template, val, index, "SBO:0000176") + return [self._build_reaction(key, val, self.template, index, "SBO:0000176") for key, val in metabolic_reactions.items()] def build_non_metabolite_reactions(self, cobra_model, index='0', allow_all_non_grp_reactions=False): #!!! allow_all_non_grp_reactions in not used @@ -518,7 +512,7 @@ def build_non_metabolite_reactions(self, cobra_model, index='0', allow_all_non_g model_mets = set(x.id for x in cobra_model.metabolites) for rxn in self.template.reactions: if rxn.type in ['universal', 'spontaneous']: - reaction = self._build_reaction(rxn.id, self.template, [], index, "SBO:0000176") + reaction = self._build_reaction(rxn.id, {}, self.template, index, "SBO:0000176") if len(model_mets.intersection(x.id for x in set(reaction.metabolites))) > 0: if reaction.id not in model_rxns: reaction.annotation["seed.reaction"] = rxn.id @@ -558,7 +552,7 @@ def build_full_template_model(template, model_id=None, index='0'): model = Model(model_id or template.id) all_reactions = [] for rxn in template.reactions: - reaction = MSBuilder._build_reaction(rxn.id, template, [], index, "SBO:0000176") + reaction = MSBuilder._build_reaction(rxn.id, {}, template, index, "SBO:0000176") reaction.annotation["seed.reaction"] = rxn.id all_reactions.append(reaction) model.add_reactions(all_reactions) From 5679652db0c9269db6d7608e2e4b3c62337027a3 Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Fri, 1 Jul 2022 20:39:25 -0400 Subject: [PATCH 3/6] polishing MSBuilder --- modelseedpy/core/msbuilder.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modelseedpy/core/msbuilder.py b/modelseedpy/core/msbuilder.py index 1dd65f5f..f7e2e2b9 100755 --- a/modelseedpy/core/msbuilder.py +++ b/modelseedpy/core/msbuilder.py @@ -328,7 +328,7 @@ def build_gpr(cpx_gene_role): #!!! Unused and redundant function class MSBuilder: - + def __init__(self, genome, template=None): """ for future methods with better customization From 2364cb5f892c79e44bfb164b6f75b3c5dfd4098c Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Fri, 1 Jul 2022 20:49:11 -0400 Subject: [PATCH 4/6] polishing MSBuilder --- modelseedpy/core/msbuilder.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/modelseedpy/core/msbuilder.py b/modelseedpy/core/msbuilder.py index f7e2e2b9..8d2b9534 100755 --- a/modelseedpy/core/msbuilder.py +++ b/modelseedpy/core/msbuilder.py @@ -373,8 +373,8 @@ def _build_reaction_complex_gpr_sets2(match_complex, allow_incomplete_complexes= #print(role_id, role_genes[role_id]) return complexes - @staticmethod #!!! Unused and redundant function - def _build_reaction_complex_gpr_sets(match_complex, allow_incomplete_complexes=True): + @staticmethod + def _build_reaction_complex_gpr_sets(match_complex, allow_incomplete_complexes=True): #!!! Unused and redundant function complexes = {} for cpx_id in match_complex: complete = True @@ -522,6 +522,8 @@ def build_non_metabolite_reactions(self, cobra_model, index='0', allow_all_non_g def build(self, model_id, index='0', allow_all_non_grp_reactions=False, annotate_with_rast=True): if annotate_with_rast: + rast = RastClient() + res = rast.annotate_genome(self.genome) # !!! res is never used self.search_name_to_genes, self.search_name_to_original = _aaaa(self.genome, 'RAST') # rxn_roles = aux_template(self.template) # needs to be fixed to actually reflect template GPR rules if not self.template: @@ -530,8 +532,7 @@ def build(self, model_id, index='0', allow_all_non_grp_reactions=False, annotate # construct the model cobra_model = Model(model_id) cobra_model.add_reactions(self.build_metabolic_reactions(index=index)) - if allow_all_non_grp_reactions: - cobra_model.add_reactions(self.build_non_metabolite_reactions(cobra_model, index)) + cobra_model.add_reactions(self.build_non_metabolite_reactions(cobra_model, index, allow_all_non_grp_reactions)) self.build_exchanges(cobra_model) if any([self.template.name.startswith(x) for x in ('GramPos', 'CoreModel', 'GramNeg')]): @@ -597,6 +598,7 @@ def gapfill_model(original_mdl, target_reaction, template, media): pkgmgr.getpkg("KBaseMediaPkg").build_package(media) #with open('Gapfilling.lp', 'w') as out: # out.write(str(model.solver)) + sol = model.optimize() # !!! sol is never used gfresults = pkgmgr.getpkg("GapfillingPkg").compute_gapfilled_solution() for rxnid in gfresults["reversed"]: rxn = original_mdl.reactions.get_by_id(rxnid) From ee641ebde44bee38854ee2eb0e60bdb98acbd01a Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Fri, 22 Jul 2022 22:00:30 -0400 Subject: [PATCH 5/6] msgrowthphenotypes polishing --- modelseedpy/core/msbuilder.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/modelseedpy/core/msbuilder.py b/modelseedpy/core/msbuilder.py index 8d2b9534..2e3fc659 100755 --- a/modelseedpy/core/msbuilder.py +++ b/modelseedpy/core/msbuilder.py @@ -1,6 +1,5 @@ import logging import itertools # !!! this import is not used -logger = logging.getLogger(__name__) import cobra from modelseedpy.core.rast_client import RastClient # !!! this import is not used @@ -12,6 +11,8 @@ SBO_ANNOTATION = "sbo" +logger = logging.getLogger(__name__) + ### temp stuff ### core_biomass = { 'cpd00032_c0': -1.7867, @@ -526,7 +527,7 @@ def build(self, model_id, index='0', allow_all_non_grp_reactions=False, annotate res = rast.annotate_genome(self.genome) # !!! res is never used self.search_name_to_genes, self.search_name_to_original = _aaaa(self.genome, 'RAST') # rxn_roles = aux_template(self.template) # needs to be fixed to actually reflect template GPR rules - if not self.template: + if self.template is None: self.auto_select_template() # construct the model From bd8de86cc0729c0be930d2fb247e3fad108174f4 Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Thu, 15 Sep 2022 15:27:05 -0500 Subject: [PATCH 6/6] final polish --- modelseedpy/core/msbuilder.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/modelseedpy/core/msbuilder.py b/modelseedpy/core/msbuilder.py index 2e3fc659..36f3f50a 100755 --- a/modelseedpy/core/msbuilder.py +++ b/modelseedpy/core/msbuilder.py @@ -334,7 +334,8 @@ def __init__(self, genome, template=None): """ for future methods with better customization """ - self.genome = genome; self.template = template + self.genome = genome + self.template = template self.search_name_to_genes, self.search_name_to_original = _aaaa(genome, 'RAST') def _get_template_reaction_complexes(self, template_reaction):