Create & Annotate An SBML Model¶
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 | # in this file, we show a minimal example of constructing
# an SBML model containing composite annotations
# we want to describe "cytosolic glucose in a pancreatic beta cell"
import antimony
# load a simple SBML model with MIRIAM-style annotations
# for 'glucose' and the 'cytosol'
sb='''
// Created by libAntimony v2.11.0
model *beta_cell_model()
// Compartments and Species:
compartment cytosol;
species glucose in cytosol;
// Reactions:
glucose_import: => glucose; glucose_import_glucose_import_rate;
// Species initializations:
glucose = 0;
glucose has molar_per_volume;
// Compartment initializations:
cytosol = 1;
// Variable initializations:
glucose_import_glucose_import_rate = 1;
// Other declarations:
const cytosol;
// Unit definitions:
unit molar = 1e1 mole / 1e1 litre;
unit molar_per_volume = 1e1 mole / (1e1 litre * litre);
// CV terms:
cytosol identity "https://identifiers.org/GO:0005829"
glucose identity "http://identifiers.org/CHEBI:17234"
end
'''
antimony.loadAntimonyString(sb)
sbml = antimony.getSBMLString('beta_cell_model')
# print the raw SBML content
print('**********************************************************')
print('SBML output:')
print('**********************************************************')
print(sbml)
# import the model into libSemSim
import semsim
from semsim import OPB, CL, PhysicalProperty, bqb
# importing the model into libSemSim will automatically
# add metaids for any annotatable SBML elements that lack them
importer = semsim.SBMLImporter(sbml)
model = importer.getSBMLModel()
# print the extracted annotations as RDF
print('**********************************************************')
print('RDF serialization of annotations:')
print('**********************************************************\n')
# the RDF always needs the location of the original SBML file
# (which is actually a URI)
# for COMBINE archives, this will be a relative path inside the archive
print(model.getRDF('./my-sbml-file.xml', 'turtle'))
# this isn't bad - the RDF already contains a composite annotation
# telling is that the species "glucose" is inside the "cytosol"
# we currently have "cytosolic glucose"
# however - what we want is "cytosolic glucose in a pancreatic beta cell"
# to describe this, we need to annotate the compartment (not the species!)
# and describe the cell type it is part of
# the compartment currently has a *singular annotation*, which can't
# represent membership in a "pancreatic beta cell", so we change this to
# a composite annotation
c = model.getComponentForId('cytosol')
c.setAnnotation(c.getAnnotation().makeComposite(PhysicalProperty(OPB.get(523))))
# now add a descriptor term specifying that the cytoplasm
#is part of a pancreatic beta cell (CL:0000169)
c.getCompositeAnnotation().addTerm(
bqb.isPartOf, # the relation (the cytoplasm *is part of* ...)
CL.get(169) # the resource (pancreatic beta cell)
)
print('**********************************************************')
print('RDF serialization with cell type annotation:')
print('**********************************************************\n')
# the RDF always needs the location of the original SBML file
# (which is actually a URI)
# for COMBINE archives, this will be a relative path inside the archive
print(model.getRDF('./my-sbml-file.xml', 'turtle'))
print('**********************************************************')
print('Humanized version of the annotation information')
print('(ontology terms replaced by human-readable names):')
print('**********************************************************\n')
print(model.humanize())
|
Expected Output¶
**********************************************************
SBML output:
**********************************************************
<?xml version="1.0" encoding="UTF-8"?>
<!-- Created by libAntimony version v2.11.0 with libSBML version 5.17.1. -->
<sbml xmlns="http://www.sbml.org/sbml/level3/version1/core" level="3" version="1">
<model metaid="beta_cell_model" id="beta_cell_model">
<listOfUnitDefinitions>
<unitDefinition id="molar_per_volume">
<listOfUnits>
<unit kind="mole" exponent="1" scale="1" multiplier="1"/>
<unit kind="litre" exponent="-1" scale="1" multiplier="1"/>
<unit kind="litre" exponent="-1" scale="0" multiplier="1"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="molar">
<listOfUnits>
<unit kind="mole" exponent="1" scale="1" multiplier="1"/>
<unit kind="litre" exponent="-1" scale="1" multiplier="1"/>
</listOfUnits>
</unitDefinition>
</listOfUnitDefinitions>
<listOfCompartments>
<compartment metaid="beta_cell_model.cytosol" id="cytosol" spatialDimensions="3" size="1" constant="true">
<annotation>
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#" xmlns:dcterms="http://purl.org/dc/terms/" xmlns:vCard="http://www.w3.org/2001/vcard-rdf/3.0#" xmlns:vCard4="http://www.w3.org/2006/vcard/ns#" xmlns:bqbiol="http://biomodels.net/biology-qualifiers/" xmlns:bqmodel="http://biomodels.net/model-qualifiers/">
<rdf:Description rdf:about="#beta_cell_model.cytosol">
<bqbiol:is>
<rdf:Bag>
<rdf:li rdf:resource="https://identifiers.org/GO:0005829"/>
</rdf:Bag>
</bqbiol:is>
</rdf:Description>
</rdf:RDF>
</annotation>
</compartment>
</listOfCompartments>
<listOfSpecies>
<species metaid="beta_cell_model.glucose" id="glucose" compartment="cytosol" initialConcentration="0" substanceUnits="molar" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false">
<annotation>
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#" xmlns:dcterms="http://purl.org/dc/terms/" xmlns:vCard="http://www.w3.org/2001/vcard-rdf/3.0#" xmlns:vCard4="http://www.w3.org/2006/vcard/ns#" xmlns:bqbiol="http://biomodels.net/biology-qualifiers/" xmlns:bqmodel="http://biomodels.net/model-qualifiers/">
<rdf:Description rdf:about="#beta_cell_model.glucose">
<bqbiol:is>
<rdf:Bag>
<rdf:li rdf:resource="http://identifiers.org/CHEBI:17234"/>
</rdf:Bag>
</bqbiol:is>
</rdf:Description>
</rdf:RDF>
</annotation>
</species>
</listOfSpecies>
<listOfParameters>
<parameter id="glucose_import_glucose_import_rate" value="1" constant="true"/>
</listOfParameters>
<listOfReactions>
<reaction id="glucose_import" reversible="false" fast="false">
<listOfProducts>
<speciesReference species="glucose" stoichiometry="1" constant="true"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> glucose_import_glucose_import_rate </ci>
</math>
</kineticLaw>
</reaction>
</listOfReactions>
</model>
</sbml>
**********************************************************
RDF serialization of annotations:
**********************************************************
@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqb: <http://biomodels.net/biology-qualifiers/> .
<#beta_cell_model.glucose_entity>
bqb:is <http://identifiers.org/CHEBI:17234> ;
bqb:occursIn <./my-sbml-file.xml#beta_cell_model.cytosol> .
<./my-sbml-file.xml#beta_cell_model.cytosol>
bqb:is <https://identifiers.org/GO:0005829> .
<./my-sbml-file.xml#beta_cell_model.glucose>
bqb:is <http://identifiers.org/opb/OPB_00425> ;
bqb:isPropertyOf <#beta_cell_model.glucose_entity> .
**********************************************************
RDF serialization with cell type annotation:
**********************************************************
@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqb: <http://biomodels.net/biology-qualifiers/> .
<#beta_cell_model.cytosol_entity>
bqb:is <https://identifiers.org/GO:0005829> ;
bqb:isPartOf <#beta_cell_model.cytosol_entity_term1> .
<#beta_cell_model.cytosol_entity_term1>
bqb:is <https://identifiers.org/CL:0000169> .
<#beta_cell_model.glucose_entity>
bqb:is <http://identifiers.org/CHEBI:17234> ;
bqb:occursIn <./my-sbml-file.xml#beta_cell_model.cytosol> .
<./my-sbml-file.xml#beta_cell_model.cytosol>
bqb:is <http://identifiers.org/opb/OPB_00523> ;
bqb:isPropertyOf <#beta_cell_model.cytosol_entity> .
<./my-sbml-file.xml#beta_cell_model.glucose>
bqb:is <http://identifiers.org/opb/OPB_00425> ;
bqb:isPropertyOf <#beta_cell_model.glucose_entity> .
**********************************************************
Humanized version of the annotation information
(ontology terms replaced by human-readable names):
**********************************************************
OPB:spatial volume -> (isPropertyOf) -> #beta_cell_model.cytosol -> (is) -> GO:cytoplasm -> (isPartOf) -> CL:type B pancreatic cell
OPB:chemical concentration -> (isPropertyOf) -> #beta_cell_model.glucose -> (is) -> CHEBI:glucose -> (occursIn) -> #beta_cell_model.cytosol