概要

libSBMLというSystemBiologyのためのマークアップ言語を扱うライブラリのPython APIの例がPython使いには読みにくかったので、勝手に書き直す。

元ページ：https://synonym.caltech.edu/software/libsbml/5.18.0/docs/formatted/python-api/libsbml-python-creating-model.html
元ファイル：https://synonym.caltech.edu/software/libsbml/5.18.0/docs/formatted/python-api/create_simple_model_8py-example.html

環境

$ python --version
Python 3.10.0
$ pip freeze | grep libsbml
python-libsbml==5.19.5

内容

元ファイル

コード内容（createSimpleModel.py）

createSimpleModel.py(抜粋)

 import sys
 from libsbml import *
 
 
 def check(value, message):
   """If 'value' is None, prints an error message constructed using
   'message' and then exits with status code 1.  If 'value' is an integer,
   it assumes it is a libSBML return status code.  If the code value is
   LIBSBML_OPERATION_SUCCESS, returns without further action; if it is not,
   prints an error message constructed using 'message' along with text from
   libSBML explaining the meaning of the code, and exits with status code 1.
   """
   if value is None:
     raise SystemExit('LibSBML returned a null value trying to ' + message + '.')
   elif type(value) is int:
     if value == LIBSBML_OPERATION_SUCCESS:
       return
     else:
       err_msg = 'Error encountered trying to ' + message + '.' \
                 + 'LibSBML returned error code ' + str(value) + ': "' \
                 + OperationReturnValue_toString(value).strip() + '"'
       raise SystemExit(err_msg)
   else:
     return
 
 
 def create_model():
   """Returns a simple but complete SBML Level 3 model for illustration."""
 
   # Create an empty SBMLDocument object.  It's a good idea to check for
   # possible errors.  Even when the parameter values are hardwired like
   # this, it is still possible for a failure to occur (e.g., if the
   # operating system runs out of memory).
 
   try:
     document = SBMLDocument(3, 1)
   except ValueError:
     raise SystemExit('Could not create SBMLDocument object')
 
   # Create the basic Model object inside the SBMLDocument object.  To
   # produce a model with complete units for the reaction rates, we need
   # to set the 'timeUnits' and 'extentUnits' attributes on Model.  We
   # set 'substanceUnits' too, for good measure, though it's not strictly
   # necessary here because we also set the units for individual species
   # in their definitions.
 
   model = document.createModel()
   check(model,                              'create model')
   check(model.setTimeUnits("second"),       'set model-wide time units')
   check(model.setExtentUnits("mole"),       'set model units of extent')
   check(model.setSubstanceUnits('mole'),    'set model substance units')
 
   # Create a unit definition we will need later.  Note that SBML Unit
   # objects must have all four attributes 'kind', 'exponent', 'scale'
   # and 'multiplier' defined.
 
   per_second = model.createUnitDefinition()
   check(per_second,                         'create unit definition')
   check(per_second.setId('per_second'),     'set unit definition id')
   unit = per_second.createUnit()
   check(unit,                               'create unit on per_second')
   check(unit.setKind(UNIT_KIND_SECOND),     'set unit kind')
   check(unit.setExponent(-1),               'set unit exponent')
   check(unit.setScale(0),                   'set unit scale')
   check(unit.setMultiplier(1),              'set unit multiplier')
 
   # Create a compartment inside this model, and set the required
   # attributes for an SBML compartment in SBML Level 3.
 
   c1 = model.createCompartment()
   check(c1,                                 'create compartment')
   check(c1.setId('c1'),                     'set compartment id')
   check(c1.setConstant(True),               'set compartment "constant"')
   check(c1.setSize(1),                      'set compartment "size"')
   check(c1.setSpatialDimensions(3),         'set compartment dimensions')
   check(c1.setUnits('litre'),               'set compartment size units')
 
   # Create two species inside this model, set the required attributes
   # for each species in SBML Level 3 (which are the 'id', 'compartment',
   # 'constant', 'hasOnlySubstanceUnits', and 'boundaryCondition'
   # attributes), and initialize the amount of the species along with the
   # units of the amount.
 
   s1 = model.createSpecies()
   check(s1,                                 'create species s1')
   check(s1.setId('s1'),                     'set species s1 id')
   check(s1.setCompartment('c1'),            'set species s1 compartment')
   check(s1.setConstant(False),              'set "constant" attribute on s1')
   check(s1.setInitialAmount(5),             'set initial amount for s1')
   check(s1.setSubstanceUnits('mole'),       'set substance units for s1')
   check(s1.setBoundaryCondition(False),     'set "boundaryCondition" on s1')
   check(s1.setHasOnlySubstanceUnits(False), 'set "hasOnlySubstanceUnits" on s1')
 
   s2 = model.createSpecies()
   check(s2,                                 'create species s2')
   check(s2.setId('s2'),                     'set species s2 id')
   check(s2.setCompartment('c1'),            'set species s2 compartment')
   check(s2.setConstant(False),              'set "constant" attribute on s2')
   check(s2.setInitialAmount(0),             'set initial amount for s2')
   check(s2.setSubstanceUnits('mole'),       'set substance units for s2')
   check(s2.setBoundaryCondition(False),     'set "boundaryCondition" on s2')
   check(s2.setHasOnlySubstanceUnits(False), 'set "hasOnlySubstanceUnits" on s2')
 
   # Create a parameter object inside this model, set the required
   # attributes 'id' and 'constant' for a parameter in SBML Level 3, and
   # initialize the parameter with a value along with its units.
 
   k = model.createParameter()
   check(k,                                  'create parameter k')
   check(k.setId('k'),                       'set parameter k id')
   check(k.setConstant(True),                'set parameter k "constant"')
   check(k.setValue(1),                      'set parameter k value')
   check(k.setUnits('per_second'),           'set parameter k units')
 
   # Create a reaction inside this model, set the reactants and products,
   # and set the reaction rate expression (the SBML "kinetic law").  We
   # set the minimum required attributes for all of these objects.  The
   # units of the reaction rate are determined from the 'timeUnits' and
   # 'extentUnits' attributes on the Model object.
 
   r1 = model.createReaction()
   check(r1,                                 'create reaction')
   check(r1.setId('r1'),                     'set reaction id')
   check(r1.setReversible(False),            'set reaction reversibility flag')
   check(r1.setFast(False),                  'set reaction "fast" attribute')
 
   species_ref1 = r1.createReactant()
   check(species_ref1,                       'create reactant')
   check(species_ref1.setSpecies('s1'),      'assign reactant species')
   check(species_ref1.setConstant(True),     'set "constant" on species ref 1')
 
   species_ref2 = r1.createProduct()
   check(species_ref2,                       'create product')
   check(species_ref2.setSpecies('s2'),      'assign product species')
   check(species_ref2.setConstant(True),     'set "constant" on species ref 2')
 
   math_ast = parseL3Formula('k * s1 * c1')
   check(math_ast,                           'create AST for rate expression')
 
   kinetic_law = r1.createKineticLaw()
   check(kinetic_law,                        'create kinetic law')
   check(kinetic_law.setMath(math_ast),      'set math on kinetic law')
 
   # And we're done creating the basic model.
   # Now return a text string containing the model in XML format.
 
   return writeSBMLToString(document)
 
 
 if __name__ == '__main__':
   print(create_model())

盆栽ed

bonsaied.py

from libsbml import SBMLDocument, writeSBMLToString, UNIT_KIND_SECOND, parseL3Formula


def create_sbml_model_string():
    document = SBMLDocument(3, 1)  # create an empty SBMLDocument object

    model = document.createModel()  # create model
    model.setTimeUnits("second")
    model.setExtentUnits("mole")
    model.setSubstanceUnits('mole')

    per_second = model.createUnitDefinition()  # create unit definition
    per_second.setId('per_second')

    unit = per_second.createUnit()  # create unit on per_second
    unit.setKind(UNIT_KIND_SECOND)
    unit.setExponent(-1)
    unit.setScale(0)
    unit.setMultiplier(1)

    c1 = model.createCompartment()  # create compartment
    c1.setId('c1')
    c1.setConstant(True)
    c1.setSize(1)
    c1.setSpatialDimensions(3)
    c1.setUnits('litre')


    s1 = model.createSpecies()  # create species s1
    s1.setId('s1')
    s1.setCompartment('c1')
    s1.setConstant(False)
    s1.setInitialAmount(5)
    s1.setSubstanceUnits('mole')
    s1.setBoundaryCondition(False)
    s1.setHasOnlySubstanceUnits(False)

    s2 = model.createSpecies()  # create species s2
    s2.setId('s2')
    s2.setCompartment('c1')
    s2.setConstant(False)
    s2.setInitialAmount(0)
    s2.setSubstanceUnits('mole')
    s2.setBoundaryCondition(False)
    s2.setHasOnlySubstanceUnits(False)

    k = model.createParameter()  # create parameter k
    k.setId('k')
    k.setConstant(True)
    k.setValue(1)
    k.setUnits('per_second')


    r1 = model.createReaction()  # create reaction
    r1.setId('r1')
    r1.setReversible(False)
    r1.setFast(False)

    species_ref1 = r1.createReactant()  # create reactant of r1
    species_ref1.setSpecies('s1')
    species_ref1.setConstant(True)

    species_ref2 = r1.createProduct()  # create product of r1
    species_ref2.setSpecies('s2')
    species_ref2.setConstant(True)

    math_ast = parseL3Formula('k * s1 * c1')  # create AST for rate expression

    kinetic_law = r1.createKineticLaw()  # create kinetic law
    kinetic_law.setMath(math_ast)


    return writeSBMLToString(document)


if __name__ == '__main__':
    print(create_sbml_model_string())

checkは見にくいので削除。Errorをキャッチしてreraiseするのは、この際インタープリターに任せてしまう。(手堅く書いた方が良いよとのコメントもあるのだけれど)
インデントは半角スペース4つに。
関数名は、結局文字列を返すので、create_sbml_model_stringに🤔
sysはもうimportしない。
import内容も明示する。
列数overはひとまず気にしない。。

setの乱用が気になるが、他のAPIを知らないのでこのぐらいが関の山な気がする。

出力

<?xml version="1.0" encoding="UTF-8"?>
<sbml xmlns="http://www.sbml.org/sbml/level3/version1/core" level="3" version="1">
  <model substanceUnits="mole" timeUnits="second" extentUnits="mole">
    <listOfUnitDefinitions>
      <unitDefinition id="per_second">
        <listOfUnits>
          <unit kind="second" exponent="-1" scale="0" multiplier="1"/>
        </listOfUnits>
      </unitDefinition>
    </listOfUnitDefinitions>
    <listOfCompartments>
      <compartment id="c1" spatialDimensions="3" size="1" units="litre" constant="true"/>
    </listOfCompartments>
    <listOfSpecies>
      <species id="s1" compartment="c1" initialAmount="5" substanceUnits="mole" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="s2" compartment="c1" initialAmount="0" substanceUnits="mole" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
    </listOfSpecies>
    <listOfParameters>
      <parameter id="k" value="1" units="per_second" constant="true"/>
    </listOfParameters>
    <listOfReactions>
      <reaction id="r1" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="s1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="s2" constant="true"/>
        </listOfProducts>
        <kineticLaw>
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k </ci>
              <ci> s1 </ci>
              <ci> c1 </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
    </listOfReactions>
  </model>
</sbml>

感想

理解してから見てみると、元のコードもそこまで複雑ではないから、
やはり(この投稿は)下手の横好き、下手の盆栽だったのだなあ、と少し肩を落とす。