Skip to content

Instantly share code, notes, and snippets.

@lukauskas
Last active June 7, 2021 04:12
Show Gist options
  • Save lukauskas/d1e30bdccc5b801d341d to your computer and use it in GitHub Desktop.
Save lukauskas/d1e30bdccc5b801d341d to your computer and use it in GitHub Desktop.
Parse SBML stoichiometry matrix
from __future__ import print_function
import libsbml
import argparse
def _parser():
parser = argparse.ArgumentParser(description="Parse stoichiometry matrix of SBML file")
parser.add_argument('file', metavar="filename", type=argparse.FileType('r'),
help="Filename of SBML file to parse")
return parser
def _main():
parser = _parser()
args = parser.parse_args()
file_ = args.file
species, reactions, stoichiometry_matrix = parse_file(file_)
_print_stoichiometry_matrix(species, reactions, stoichiometry_matrix)
def parse_file(open_file_):
reader = libsbml.SBMLReader()
document = reader.readSBML(open_file_.name)
sbml_model = document.getModel()
species = [s.getName() for s in sbml_model.getListOfSpecies()]
reactions = [r.getId() for r in sbml_model.getListOfReactions()]
stoichiometry_matrix = {}
for reaction_index, reaction in enumerate(sbml_model.getListOfReactions()):
reactants = {r.getSpecies(): r.getStoichiometry() for r in reaction.getListOfReactants()}
products = {p.getSpecies(): p.getStoichiometry() for p in reaction.getListOfProducts()}
for species_index, species_node in enumerate(sbml_model.getListOfSpecies()):
species_id = species_node.getId()
net_stoichiometry = int(reactants.get(species_id, 0)) - int(products.get(species_id, 0))
stoichiometry_matrix[species_index, reaction_index] = net_stoichiometry
return species, reactions, stoichiometry_matrix
def _print_stoichiometry_matrix(species, reactions, stoichiometry_matrix):
print('\t'.join(['---'] + reactions))
for species_ix, species_label in enumerate(species):
to_print = [species_label]
for reaction_ix in range(len(reactions)):
stoichiometry = stoichiometry_matrix[species_ix, reaction_ix]
to_print.append(stoichiometry)
print('\t'.join(map(str, to_print)))
if __name__ == '__main__':
_main()
@lukauskas
Copy link
Author

@qacwnfq It's been a while since I worked on this so not sure how to fix it.
Would you suggest making it into a float instead?

@qacwnfq
Copy link

qacwnfq commented Sep 28, 2019

@lukauskas
Yes, that could work :)
The main reason behind my comment was just documentation for any one else using this gist, because it worked very well on some sbml files and not on others for me.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment