Composite Annotations: SBML Examples

Annotations on a Physical Entity

In this case we build a new model using tellurium and antimony in Python. For simplicity, the string was copied from Python and stored as a string for the C++ and C examples.

In the example we create a PhysicalEntity type annotation of A, which has the metid SpeciesA. Our annotation states that we are simulating the chemical concetration (OPB:OPB_00340) of Smad3 (uniprot/P84022’) in the dermal fibroblasts `fma/FMA:63877 of the dermal papilla fma:FMA:70737 in the skin FMA_7163 of the surface or the arm FMA:24178.

Note

You can add an arbitrary number of anatomical locations, but the order you add them is meaningful. The first anatomical location you add will be at the smallest level of biological organization which increases when you add more locations.

Demonstrate the creation of a physical entity composite annotation in Python
from pyomexmeta import RDF, eUriType
sbml = """<?xml version="1.1" encoding="UTF-8"?>
<!-- Created by libAntimony version v2.12.0 with libSBML version 5.18.1. -->
<sbml xmlns="http://www.sbml.org/sbml/level3/version1/core" level="3" version="1">
  <model metaid="SimpleRegulation" id="SimpleRegulation">
    <listOfCompartments>
      <compartment id="cell" spatialDimensions="3" size="1" constant="true"/>
    </listOfCompartments>
    <listOfSpecies>
      <species id="A" metaid="SpeciesA" compartment="cell" initialConcentration="10" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="B" metaid="SpeciesB" compartment="cell" initialConcentration="0" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="C" metaid="SpeciesC" compartment="cell" initialConcentration="10" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="D" metaid="SpeciesD" compartment="cell" initialConcentration="0" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
    </listOfSpecies>
    <listOfParameters>
      <parameter id="S"  metaid="S" value="0" constant="true"/>
      <parameter id="k1" metaid="k1" value="0.1" constant="true"/>
      <parameter id="k2" metaid="k2" value="0.1" constant="true"/>
      <parameter id="k3" metaid="k3" value="0.1" constant="true"/>
      <parameter id="k4" metaid="k4" value="0.1" constant="true"/>
    </listOfParameters>
    <listOfReactions>
      <reaction id="R1" metaid="R1" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="A" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="B" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction1">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k1 </ci>
              <ci> A </ci>
              <ci> S </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R2" metaid="R2" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="B" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="A" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction2">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k2 </ci>
              <ci> B </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R3" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="C" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="D" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <listOfModifiers>
          <modifierSpeciesReference species="B"/>
        </listOfModifiers>
        <kineticLaw>
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k3 </ci>
              <ci> C </ci>
              <ci> B </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R4" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="D" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="C" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw>
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k4 </ci>
              <ci> D </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
    </listOfReactions>
  </model>
</sbml>
"""

# create an empty RDF object
rdf = RDF()

editor = rdf.to_editor(sbml, generate_new_metaids=True, sbml_semantic_extraction=False)

sbml_with_metaids = editor.get_xml()
print(sbml_with_metaids)

with editor.new_physical_entity() as physical_entity:
    physical_entity \
        .about("SpeciesA", eUriType.MODEL_URI) \
        .identity("uniprot/PD12345") \
        .is_part_of("fma:1234") \
        .has_property( "EntityProperty", eUriType.LOCAL_URI, "opb:OPB_12345")


print(rdf)


Output
[16:33:33 +01:00]  warning : XML parser warning: Unsupported version '1.1' 
<?xml version="1.1" encoding="UTF-8"?>
<!-- Created by libAntimony version v2.12.0 with libSBML version 5.18.1. -->
<sbml xmlns="http://www.sbml.org/sbml/level3/version1/core" level="3" version="1">
  <model metaid="SimpleRegulation" id="SimpleRegulation">
    <listOfCompartments>
      <compartment id="cell" spatialDimensions="3" size="1" constant="true" metaid="compartment0000"/>
    </listOfCompartments>
    <listOfSpecies>
      <species id="A" metaid="SpeciesA" compartment="cell" initialConcentration="10" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="B" metaid="SpeciesB" compartment="cell" initialConcentration="0" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="C" metaid="SpeciesC" compartment="cell" initialConcentration="10" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="D" metaid="SpeciesD" compartment="cell" initialConcentration="0" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
    </listOfSpecies>
    <listOfParameters>
      <parameter id="S" metaid="S" value="0" constant="true"/>
      <parameter id="k1" metaid="k1" value="0.1" constant="true"/>
      <parameter id="k2" metaid="k2" value="0.1" constant="true"/>
      <parameter id="k3" metaid="k3" value="0.1" constant="true"/>
      <parameter id="k4" metaid="k4" value="0.1" constant="true"/>
    </listOfParameters>
    <listOfReactions>
      <reaction id="R1" metaid="R1" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="A" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="B" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction1">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k1 </ci>
              <ci> A </ci>
              <ci> S </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R2" metaid="R2" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="B" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="A" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction2">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k2 </ci>
              <ci> B </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R3" reversible="false" fast="false" metaid="reaction0000">
        <listOfReactants>
          <speciesReference species="C" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="D" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <listOfModifiers>
          <modifierSpeciesReference species="B"/>
        </listOfModifiers>
        <kineticLaw metaid="kineticLaw0000">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k3 </ci>
              <ci> C </ci>
              <ci> B </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R4" reversible="false" fast="false" metaid="reaction0001">
        <listOfReactants>
          <speciesReference species="D" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="C" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="kineticLaw0001">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k4 </ci>
              <ci> D </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
    </listOfReactions>
  </model>
</sbml>

@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqbiol: <http://biomodels.net/biology-qualifiers/> .
@prefix OMEXlib: <http://omex-library.org/> .
@prefix local: <http://omex-library.org/NewOmex.omex/NewModel.rdf#> .

local:EntityProperty
    bqbiol:isPropertyOf <http://omex-library.org/NewOmex.omex/NewModel.xml#SpeciesA> ;
    bqbiol:isVersionOf <https://identifiers.org/opb:OPB_12345> .

<http://omex-library.org/NewOmex.omex/NewModel.xml#SpeciesA>
    bqbiol:is <https://identifiers.org/uniprot/PD12345> ;
    bqbiol:isPartOf <https://identifiers.org/fma:1234> .


Annotations on a Physical Process

Demonstrate the creation of a physical process composite annotation in Python
from pyomexmeta import RDF, eUriType

sbml = """<?xml version="1.1" encoding="UTF-8"?>
<sbml xmlns="http://www.sbml.org/sbml/level3/version1/core" level="3" version="1">
  <model metaid="SimpleRegulation" id="SimpleRegulation">
    <listOfCompartments>
      <compartment id="cell" spatialDimensions="3" size="1" constant="true"/>
    </listOfCompartments>
    <listOfSpecies>
      <species id="A" metaid="A" compartment="cell" initialConcentration="10" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="B" metaid="B" compartment="cell" initialConcentration="0" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="C" metaid="C" compartment="cell" initialConcentration="10" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
      <species id="D" metaid="D" compartment="cell" initialConcentration="0" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
    </listOfSpecies>
    <listOfParameters>
      <parameter id="S"  metaid="S"  value="0" constant="true"/>
      <parameter id="k1" metaid="k1" value="0.1" constant="true"/>
      <parameter id="k2" metaid="k2" value="0.1" constant="true"/>
      <parameter id="k3" metaid="k3" value="0.1" constant="true"/>
      <parameter id="k4" metaid="k4" value="0.1" constant="true"/>
    </listOfParameters>
    <listOfReactions>
      <reaction id="R1" metaid="R1" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="A" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="B" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction1">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k1 </ci>
              <ci> A </ci>
              <ci> S </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R2" metaid="R2" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="B" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="A" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction2">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k2 </ci>
              <ci> B </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R3" metaid="R3" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="C" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="D" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <listOfModifiers>
          <modifierSpeciesReference species="B"/>
        </listOfModifiers>
        <kineticLaw metaid="KineticLawForReaction3">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k3 </ci>
              <ci> C </ci>
              <ci> B </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id="R4" metaid="R4" reversible="false" fast="false">
        <listOfReactants>
          <speciesReference species="D" stoichiometry="1" constant="true"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species="C" stoichiometry="1" constant="true"/>
        </listOfProducts>
        <kineticLaw metaid="KineticLawForReaction4">
          <math xmlns="http://www.w3.org/1998/Math/MathML">
            <apply>
              <times/>
              <ci> k4 </ci>
              <ci> D </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
    </listOfReactions>
  </model>
</sbml>
"""

# create an empty RDF object
rdf = RDF()

editor = rdf.to_editor(sbml, generate_new_metaids=False, sbml_semantic_extraction=False)

# For when generate_new_metaids=True
# sbml_with_metaids = editor.get_xml()
# print(sbml_with_metaids)

# annotate reaction R1as the chemical concentration
# flow rate OPB:OPB_00593 of Smad3 to phosphorylated smad3
with editor.new_physical_process() as physical_process:
    physical_process \
        .about("R1", eUriType.MODEL_URI) \
        .is_version_of("GO/GO12345") \
        .add_source(physical_entity_reference="SpeciesA", uri_type=eUriType.MODEL_URI, multiplier=1) \
        .add_sink(physical_entity_reference="SpeciesB", uri_type=eUriType.MODEL_URI, multiplier=1) \
        .add_mediator(physical_entity_reference="SpeciesC", uri_type=eUriType.MODEL_URI) \
        .has_property(is_version_of="OPB:OPB_00593")


print(rdf)

Output
[16:33:33 +01:00]  warning : XML parser warning: Unsupported version '1.1' 
@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqbiol: <http://biomodels.net/biology-qualifiers/> .
@prefix semsim: <http://bime.uw.edu/semsim/> .
@prefix OMEXlib: <http://omex-library.org/> .
@prefix local: <http://omex-library.org/NewOmex.omex/NewModel.rdf#> .

local:MediatorParticipant0000
    semsim:hasPhysicalEntityReference <http://omex-library.org/NewOmex.omex/NewModel.xml#SpeciesC> .

local:ProcessProperty0000
    bqbiol:isPropertyOf <http://omex-library.org/NewOmex.omex/NewModel.xml#R1> ;
    bqbiol:isVersionOf <https://identifiers.org/OPB:OPB_00593> .

local:SinkParticipant0000
    semsim:hasMultiplier "1"^^rdf:double ;
    semsim:hasPhysicalEntityReference <http://omex-library.org/NewOmex.omex/NewModel.xml#SpeciesB> .

local:SourceParticipant0000
    semsim:hasMultiplier "1"^^rdf:double ;
    semsim:hasPhysicalEntityReference <http://omex-library.org/NewOmex.omex/NewModel.xml#SpeciesA> .

<http://omex-library.org/NewOmex.omex/NewModel.xml#R1>
    semsim:hasMediatorParticipant local:MediatorParticipant0000 ;
    semsim:hasSinkParticipant local:SinkParticipant0000 ;
    semsim:hasSourceParticipant local:SourceParticipant0000 ;
    bqbiol:isVersionOf <https://identifiers.org/GO/GO12345> .


Annotations on an energy differential

Demonstrate the creation of a energy differential composite annotation in Python
from pyomexmeta import RDF, eUriType

sbml = """<?xml version=\"1.0\" encoding=\"UTF-8\"?>
<sbml xmlns=\"http://www.sbml.org/sbml/level3/version1/core\" level=\"3\" version=\"1\">
  <model metaid=\"EnergyDiff\" id=\"EnergyDiff\">
    <listOfCompartments>
      <compartment id=\"cytoplasm\" metaid=\"cytoplasm\" spatialDimensions=\"3\" size=\"1\" constant=\"true\"/>
      <compartment id=\"extracellular\" metaid=\"extracellular\" spatialDimensions=\"3\" size=\"1\" constant=\"true\"/>
    </listOfCompartments>
    <listOfSpecies>
      <species id=\"Ca_ex\" metaid=\"Ca_ex\" compartment=\"extracellular\" initialConcentration=\"1\" hasOnlySubstanceUnits=\"false\" boundaryCondition=\"false\" constant=\"false\"/>
      <species id=\"Ca_cyt\" metaid=\"Ca_cyt\" compartment=\"cytoplasm\" initialConcentration=\"0.3\" hasOnlySubstanceUnits=\"false\" boundaryCondition=\"false\" constant=\"false\"/>
    </listOfSpecies>
    <listOfParameters>
      <parameter id=\"k1\" metaid=\"k1\" value=\"0.1\" constant=\"true\"/>
      <parameter id=\"k2\" metaid=\"k2\" value=\"0.1\" constant=\"true\"/>
      <parameter id=\"N_Ca\" metaid=\"NernstPotential\"/>
    </listOfParameters>
    <listOfReactions>
      <reaction id=\"NernstReversalPotential_in\" metaid=\"NernstReversalPotential_in\" reversible=\"false\" fast=\"false\">
        <listOfReactants>
          <speciesReference species=\"Ca_ex\" stoichiometry=\"1\" constant=\"true\"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species=\"Ca_cyt\" stoichiometry=\"1\" constant=\"true\"/>
        </listOfProducts>
        <kineticLaw>
          <math xmlns=\"http://www.w3.org/1998/Math/MathML\">
            <apply>
              <times/>
              <ci> cytoplasm </ci>
              <ci> k1 </ci>
              <ci> Ca_ex </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
      <reaction id=\"NernstReversalPotential_out\" metaid=\"NernstReversalPotential_out\" reversible=\"false\" fast=\"false\">
        <listOfReactants>
          <speciesReference species=\"Ca_cyt\" stoichiometry=\"1\" constant=\"true\"/>
        </listOfReactants>
        <listOfProducts>
          <speciesReference species=\"Ca_ex\" stoichiometry=\"1\" constant=\"true\"/>
        </listOfProducts>
        <kineticLaw>
          <math xmlns=\"http://www.w3.org/1998/Math/MathML\">
            <apply>
              <times/>
              <ci> extracellular </ci>
              <ci> k2 </ci>
              <ci> Ca_cyt </ci>
            </apply>
          </math>
        </kineticLaw>
      </reaction>
    </listOfReactions>
  </model>
</sbml>
"""

# create an empty RDF object
rdf = RDF()

editor = rdf.to_editor(sbml, generate_new_metaids=True, sbml_semantic_extraction=False)

# Property: Osmotic Pressure
with editor.new_energy_diff() as energy_diff:
    # about(): this could be named anything provided it is unique to the rdf graph because it is a LOCAL_URI
    energy_diff.about(about="MembranePotential", uri_type=eUriType.LOCAL_URI)
    energy_diff.add_source(physical_entity_reference="Ca_ex", uri_type=eUriType.MODEL_URI)
    energy_diff.add_sink(physical_entity_reference="Ca_cyt",uri_type= eUriType.MODEL_URI)
    # for energy differentials we use the full signature of hasProperty.
    # The arguments are as follows:
    #  - NernstPotential: the "about" value that is used for the property associated with the energy differential
    #  - MODEL_URI: the uri type for the previous argument. In SBML this is a <parameter> and therefore a MODEL_URI
    #      while in CellML this is a construct local to the rdf document (so LOCAL_URI)
    #  - OPB:OPB_01581: the OPB term associated with the property (isVersionOf)
    energy_diff.has_property(property_about="NernstPotential", about_uri_type=eUriType.MODEL_URI, is_version_of="OPB:OPB_01581")

print(rdf)

Output
@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqbiol: <http://biomodels.net/biology-qualifiers/> .
@prefix semsim: <http://bime.uw.edu/semsim/> .
@prefix OMEXlib: <http://omex-library.org/> .
@prefix local: <http://omex-library.org/NewOmex.omex/NewModel.rdf#> .

local:MembranePotential
    semsim:hasSinkParticipant local:SinkParticipant0000 ;
    semsim:hasSourceParticipant local:SourceParticipant0000 .

local:SinkParticipant0000
    semsim:hasPhysicalEntityReference <http://omex-library.org/NewOmex.omex/NewModel.xml#Ca_cyt> .

local:SourceParticipant0000
    semsim:hasPhysicalEntityReference <http://omex-library.org/NewOmex.omex/NewModel.xml#Ca_ex> .

<http://omex-library.org/NewOmex.omex/NewModel.xml#NernstPotential>
    bqbiol:isPropertyOf local:MembranePotential ;
    bqbiol:isVersionOf <https://identifiers.org/OPB:OPB_01581> .


Composite Annotations: CellML Examples

Annotations on a Physical Entity

Demonstrate the creation of a physical entity composite annotation in Python
from pyomexmeta import RDF, eUriType

cellml = """<?xml version=\"1.0\" encoding=\"UTF-8\"?>
<model xmlns=\"http://www.cellml.org/cellml/1.1#\" xmlns:cmeta=\"http://www.cellml.org/metadata/1.0#\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" xmlns:rdf=\"http://www.w3.org/1999/02/22-rdf-syntax-ns#\" xmlns:bqs=\"http://www.cellml.org/bqs/1.0#\" xmlns:semsim=\"http://bime.uw.edu/semsim/#\" xmlns:dc=\"http://purl.org/dc/terms/\" xmlns:vCard=\"http://www.w3.org/2001/vcard-rdf/3.0#\" name=\"annotation_examples\" cmeta:id=\"annExamples\">
  <component name=\"main\">
    <variable cmeta:id=\"main.Volume\" initial_value=\"100\" name=\"Volume\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.MembraneVoltage\" initial_value=\"-80\" name=\"MembraneVoltage\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.ReactionRate\" initial_value=\"1\" name=\"ReactionRate\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity1\" initial_value=\"1\" name=\"entity1\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity2\" initial_value=\"1\" name=\"entity2\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity3\" initial_value=\"1\" name=\"entity3\" units=\"dimensionless\" />
  </component>
</model>
"""

# create an empty RDF object
rdf = RDF()


editor = rdf.to_editor(cellml, generate_new_metaids=True, sbml_semantic_extraction=False)

sbml_with_metaids = editor.get_xml()
print(sbml_with_metaids)

with editor.new_physical_entity() as physical_entity:
    physical_entity \
        .about("entity0", eUriType.LOCAL_URI) \
        .identity("fma:9570") \
        .is_part_of("fma:18228") \
        .has_property("main.Volume", eUriType.MODEL_URI, "opb:OPB_00154")

print(rdf)


Output
<?xml version="1.0" encoding="UTF-8"?>
<model xmlns="http://www.cellml.org/cellml/1.1#" xmlns:cmeta="http://www.cellml.org/metadata/1.0#" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#" xmlns:bqs="http://www.cellml.org/bqs/1.0#" xmlns:semsim="http://bime.uw.edu/semsim/#" xmlns:dc="http://purl.org/dc/terms/" xmlns:vCard="http://www.w3.org/2001/vcard-rdf/3.0#" name="annotation_examples" cmeta:id="annExamples">
  <component name="main" cmeta:id="component0000">
    <variable cmeta:id="main.Volume" initial_value="100" name="Volume" units="dimensionless"/>
    <variable cmeta:id="main.MembraneVoltage" initial_value="-80" name="MembraneVoltage" units="dimensionless"/>
    <variable cmeta:id="main.ReactionRate" initial_value="1" name="ReactionRate" units="dimensionless"/>
    <variable cmeta:id="main.entity1" initial_value="1" name="entity1" units="dimensionless"/>
    <variable cmeta:id="main.entity2" initial_value="1" name="entity2" units="dimensionless"/>
    <variable cmeta:id="main.entity3" initial_value="1" name="entity3" units="dimensionless"/>
  </component>
</model>

@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqbiol: <http://biomodels.net/biology-qualifiers/> .
@prefix OMEXlib: <http://omex-library.org/> .
@prefix local: <http://omex-library.org/NewOmex.omex/NewModel.rdf#> .

local:entity0
    bqbiol:is <https://identifiers.org/fma:9570> ;
    bqbiol:isPartOf <https://identifiers.org/fma:18228> .

<http://omex-library.org/NewOmex.omex/NewModel.xml#main.Volume>
    bqbiol:isPropertyOf local:entity0 ;
    bqbiol:isVersionOf <https://identifiers.org/opb:OPB_00154> .


Annotations on a Physical Process

Demonstrate the creation of a physical process composite annotation in Python
from pyomexmeta import RDF, eUriType

cellml = """<?xml version=\"1.0\" encoding=\"UTF-8\"?>
<model xmlns=\"http://www.cellml.org/cellml/1.1#\" xmlns:cmeta=\"http://www.cellml.org/metadata/1.0#\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" xmlns:rdf=\"http://www.w3.org/1999/02/22-rdf-syntax-ns#\" xmlns:bqs=\"http://www.cellml.org/bqs/1.0#\" xmlns:semsim=\"http://bime.uw.edu/semsim/#\" xmlns:dc=\"http://purl.org/dc/terms/\" xmlns:vCard=\"http://www.w3.org/2001/vcard-rdf/3.0#\" name=\"annotation_examples\" cmeta:id=\"annExamples\">
  <component name=\"main\">
    <variable cmeta:id=\"main.Volume\" initial_value=\"100\" name=\"Volume\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.MembraneVoltage\" initial_value=\"-80\" name=\"MembraneVoltage\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.ReactionRate\" initial_value=\"1\" name=\"ReactionRate\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity1\" initial_value=\"1\" name=\"entity1\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity2\" initial_value=\"1\" name=\"entity2\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity3\" initial_value=\"1\" name=\"entity3\" units=\"dimensionless\" />
  </component>
</model>
"""
# create an empty RDF object
rdf = RDF()

editor = rdf.to_editor(cellml, generate_new_metaids=True, sbml_semantic_extraction=False)
sbml_with_metaids = editor.get_xml()


with editor.new_physical_process() as physical_process:
    physical_process.about("Process", eUriType.LOCAL_URI) \
        .add_source("entity1", eUriType.LOCAL_URI, 1) \
        .add_sink("entity2", eUriType.LOCAL_URI, 1) \
        .add_mediator("entity3", eUriType.LOCAL_URI) \
        .has_property("main.ReactionRate", eUriType.MODEL_URI, "opb:OPB_00592")

print(rdf)

Output
@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqbiol: <http://biomodels.net/biology-qualifiers/> .
@prefix semsim: <http://bime.uw.edu/semsim/> .
@prefix OMEXlib: <http://omex-library.org/> .
@prefix local: <http://omex-library.org/NewOmex.omex/NewModel.rdf#> .

local:MediatorParticipant0000
    semsim:hasPhysicalEntityReference local:entity3 .

local:Process
    semsim:hasMediatorParticipant local:MediatorParticipant0000 ;
    semsim:hasSinkParticipant local:SinkParticipant0000 ;
    semsim:hasSourceParticipant local:SourceParticipant0000 .

local:SinkParticipant0000
    semsim:hasMultiplier "1"^^rdf:double ;
    semsim:hasPhysicalEntityReference local:entity2 .

local:SourceParticipant0000
    semsim:hasMultiplier "1"^^rdf:double ;
    semsim:hasPhysicalEntityReference local:entity1 .

<http://omex-library.org/NewOmex.omex/NewModel.xml#main.ReactionRate>
    bqbiol:isPropertyOf local:Process ;
    bqbiol:isVersionOf <https://identifiers.org/opb:OPB_00592> .


Annotations on an energy differential

Demonstrate the creation of a energy differential composite annotation in Python
from pyomexmeta import RDF, eUriType

cellml = """<?xml version=\"1.0\" encoding=\"UTF-8\"?>
<model xmlns=\"http://www.cellml.org/cellml/1.1#\" xmlns:cmeta=\"http://www.cellml.org/metadata/1.0#\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" xmlns:rdf=\"http://www.w3.org/1999/02/22-rdf-syntax-ns#\" xmlns:bqs=\"http://www.cellml.org/bqs/1.0#\" xmlns:semsim=\"http://bime.uw.edu/semsim/#\" xmlns:dc=\"http://purl.org/dc/terms/\" xmlns:vCard=\"http://www.w3.org/2001/vcard-rdf/3.0#\" name=\"annotation_examples\" cmeta:id=\"annExamples\">
  <component name=\"main\">
    <variable cmeta:id=\"main.Volume\" initial_value=\"100\" name=\"Volume\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.MembraneVoltage\" initial_value=\"-80\" name=\"MembraneVoltage\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.ReactionRate\" initial_value=\"1\" name=\"ReactionRate\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity1\" initial_value=\"1\" name=\"entity1\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity2\" initial_value=\"1\" name=\"entity2\" units=\"dimensionless\" />
    <variable cmeta:id=\"main.entity3\" initial_value=\"1\" name=\"entity3\" units=\"dimensionless\" />
  </component>
</model>
"""

# create an empty RDF object
rdf = RDF()

editor = rdf.to_editor(cellml, generate_new_metaids=True, sbml_semantic_extraction=False)

with editor.new_physical_entity() as physical_entity:
    physical_entity.about("entity0", eUriType.LOCAL_URI) \
        .identity("fma:9570") \
        .is_part_of("fma:18228") \
        .has_property("main.Volume", eUriType.MODEL_URI, "opb:OPB_00154")


print(rdf)

Output
@prefix rdf: <http://www.w3.org/1999/02/22-rdf-syntax-ns#> .
@prefix bqbiol: <http://biomodels.net/biology-qualifiers/> .
@prefix OMEXlib: <http://omex-library.org/> .
@prefix local: <http://omex-library.org/NewOmex.omex/NewModel.rdf#> .

local:entity0
    bqbiol:is <https://identifiers.org/fma:9570> ;
    bqbiol:isPartOf <https://identifiers.org/fma:18228> .

<http://omex-library.org/NewOmex.omex/NewModel.xml#main.Volume>
    bqbiol:isPropertyOf local:entity0 ;
    bqbiol:isVersionOf <https://identifiers.org/opb:OPB_00154> .