diff --git a/mira/dkg/model.py b/mira/dkg/model.py index f0f989ee8..eb885ae62 100644 --- a/mira/dkg/model.py +++ b/mira/dkg/model.py @@ -21,11 +21,12 @@ from fastapi.responses import FileResponse from pydantic import BaseModel, Field -from mira.examples.sir import sir_bilayer, sir +from mira.examples.sir import sir_bilayer, sir, sir_parameterized_init from mira.metamodel import ( NaturalConversion, Template, ControlledConversion, stratify, Concept, ModelComparisonGraphdata, TemplateModelDelta, TemplateModel, Parameter, simplify_rate_laws, aggregate_parameters, + counts_to_dimensionless ) from mira.modeling import Model from mira.modeling.askenet.petrinet import AskeNetPetriNetModel, ModelSpecification @@ -87,6 +88,9 @@ #: PetriNetModel json example petrinet_json = PetriNetModel(Model(sir)).to_pydantic() askenet_petrinet_json = AskeNetPetriNetModel(Model(sir)).to_pydantic() +askenet_petrinet_json_units_values = AskeNetPetriNetModel( + Model(sir_parameterized_init) +).to_pydantic() @model_blueprint.post( @@ -221,6 +225,58 @@ def model_stratification( return template_model +@model_blueprint.post( + "/counts_to_dimensionless_mira", + response_model=TemplateModel, + tags=["modeling"] +) +def dimension_transform( + query: Dict[str, Any] = Body( + ..., + example={ + "model": sir_parameterized_init, + "counts_unit": "person", + "norm_factor": 1e5, + }, + ) +): + """Convert all entity concentrations to dimensionless units""" + # convert to template model + tm_json = query.pop("model") + tm = TemplateModel.from_json(tm_json) + # The concepts should have their units' expressions as sympy.Expr, + # currently they are strings + tm_dimless = counts_to_dimensionless(tm=tm, **query) + return tm_dimless + + +@model_blueprint.post( + "/counts_to_dimensionless_amr", + response_model=ModelSpecification, + tags=["modeling"] +) +def dimension_transform( + query: Dict[str, Any] = Body( + ..., + example={ + "model": askenet_petrinet_json_units_values, + "counts_units": "persons", + "norm_factor": 1e5, + }, + ) +): + """Convert all entity concentrations to dimensionless units""" + # convert to template model + amr_json = query.pop("model") + tm = template_model_from_askenet_json(amr_json) + + # Create a dimensionless model + dimless_model = counts_to_dimensionless(tm=tm, **query) + + # Transform back to askenet model + return AskeNetPetriNetModel(Model(dimless_model)).to_pydantic() + + @model_blueprint.get( "/biomodels/{model_id}", response_model=TemplateModel, diff --git a/mira/examples/sir.py b/mira/examples/sir.py index 377d17651..0fff0cad0 100644 --- a/mira/examples/sir.py +++ b/mira/examples/sir.py @@ -6,7 +6,7 @@ from mira.metamodel import ControlledConversion, NaturalConversion, \ GroupedControlledConversion, TemplateModel, Initial, Parameter, \ - safe_parse_expr + safe_parse_expr, Unit from .concepts import susceptible, infected, recovered, infected_symptomatic, \ infected_asymptomatic @@ -16,6 +16,7 @@ "sir_2_city", "sir_bilayer", "sir_parameterized", + "sir_parameterized_init", "svir", ] @@ -136,3 +137,22 @@ infection_asymptomatic, ], ) + +# SIR Parameterized Model with initial values and units, used as example in +# docs and tests +sir_parameterized_init = _d(sir_parameterized) +sir_init_val_norm = 1e5 +for template in sir_parameterized_init.templates: + for concept in template.get_concepts(): + concept.units = Unit(expression=sympy.Symbol('person')) +sir_parameterized_init.initials['susceptible_population'].value = \ + sir_init_val_norm - 1 +sir_parameterized_init.initials['infected_population'].value = 1 +sir_parameterized_init.initials['immune_population'].value = 0 + +sir_parameterized_init.parameters['beta'].units = \ + Unit(expression=1 / (sympy.Symbol('person') * sympy.Symbol('day'))) +old_beta = sir_parameterized_init.parameters['beta'].value + +for initial in sir_parameterized_init.initials.values(): + initial.concept.units = Unit(expression=sympy.Symbol('person')) diff --git a/mira/metamodel/__init__.py b/mira/metamodel/__init__.py index ecdfe1b52..ef67b5aae 100644 --- a/mira/metamodel/__init__.py +++ b/mira/metamodel/__init__.py @@ -6,4 +6,5 @@ from .schema import * from .search import * from .ops import * +from .units import * from .utils import * diff --git a/mira/metamodel/ops.py b/mira/metamodel/ops.py index 9ad683a5d..36fddd71a 100644 --- a/mira/metamodel/ops.py +++ b/mira/metamodel/ops.py @@ -9,11 +9,16 @@ from .template_model import TemplateModel, Initial, Parameter from .templates import * +from .units import Unit, dimensionless_units +from .utils import SympyExprStr + __all__ = [ "stratify", "simplify_rate_laws", - "aggregate_parameters" + "aggregate_parameters", + "get_term_roles", + "counts_to_dimensionless" ] @@ -393,3 +398,69 @@ def get_term_roles(term, template, parameters): else: term_roles['other'].append(symbol.name) return dict(term_roles) + + +def counts_to_dimensionless(tm: TemplateModel, + counts_unit: str, + norm_factor: float): + """Convert all entity concentrations to dimensionless units. + + Parameters + ---------- + tm : + A template model. + counts_unit : + The unit of the counts. + norm_factor : + The normalization factor to convert counts to concentration. + + Returns + ------- + : + A template model with all entity concentrations converted to + dimensionless units. + """ + # Make a deepcopy up front so we don't change the original template model + tm = deepcopy(tm) + # Make a symbol of the counts unit for calculations + counts_unit_symbol = sympy.Symbol(counts_unit) + + initials_normalized = set() + # First we normalize concepts and their initials + for template in tm.templates: + # Since concepts can be distributed across templates, we have to go + # template by template + for concept in template.get_concepts(): + if concept.units: + # We figure out what the exponent of the counts unit is + # if it appears in the units of the concept + (coeff, exponent) = \ + concept.units.expression.args[0].as_coeff_exponent(counts_unit_symbol) + # If the exponent is other than zero then normalization is needed + if exponent: + concept.units.expression = \ + SympyExprStr(concept.units.expression.args[0] / + (counts_unit_symbol ** exponent)) + # We not try to see if there is a corresponding initial condition + # for the concept and if so, we normalize it as well + if concept.name in tm.initials and concept.name not in initials_normalized: + init = tm.initials[concept.name] + if init.value is not None: + init.value /= (norm_factor ** exponent) + if init.concept.units: + init.concept.units.expression = \ + SympyExprStr(init.concept.units.expression.args[0] / + (counts_unit_symbol ** exponent)) + initials_normalized.add(concept.name) + # Now we do the same for parameters + for p_name, p in tm.parameters.items(): + if p.units: + (coeff, exponent) = \ + p.units.expression.args[0].as_coeff_exponent(counts_unit_symbol) + if exponent: + p.units.expression = \ + SympyExprStr(p.units.expression.args[0] / + (counts_unit_symbol ** exponent)) + p.value /= (norm_factor ** exponent) + + return tm diff --git a/mira/metamodel/template_model.py b/mira/metamodel/template_model.py index a7c91023f..950247a79 100644 --- a/mira/metamodel/template_model.py +++ b/mira/metamodel/template_model.py @@ -3,14 +3,15 @@ import datetime import sys -from typing import List, Dict, Set, Optional, Mapping, Tuple +from typing import List, Dict, Set, Optional, Mapping, Tuple, Any import networkx as nx import sympy from pydantic import BaseModel, Field from .templates import * -from .utils import safe_parse_expr +from .units import Unit +from .utils import safe_parse_expr, SympyExprStr class Initial(BaseModel): @@ -19,6 +20,23 @@ class Initial(BaseModel): concept: Concept value: float + class Config: + arbitrary_types_allowed = True + json_encoders = { + SympyExprStr: lambda e: str(e), + } + json_decoders = { + SympyExprStr: lambda e: sympy.parse_expr(e) + } + + @classmethod + def from_json(cls, data: Dict[str, Any]) -> "Initial": + value = data.pop('value') + concept_json = data.pop('concept') + # Get Concept + concept = Concept.from_json(concept_json) + return cls(concept=concept, value=value) + class Distribution(BaseModel): """A distribution of values for a parameter.""" @@ -363,24 +381,28 @@ def from_json(cls, data) -> "TemplateModel": for concept in template.get_concepts() } - initials = { - name: ( - Initial( + initials = {} + for name, value in data.get('initials', {}).items(): + if isinstance(value, float): + # If the data is just a float, upgrade it to + # a :class:`Initial` instance + initials[name] = Initial( concept=concepts[name], value=value, ) - # If the data is just a float, upgrade it to - # a :class:`Initial` instance - if isinstance(value, float) + else: # If the data is not a float, assume it's JSON - # for a :class:`Initial` instance - else value - ) - for name, value in data.get('initials', {}).items() + # for a :class:`Initial` instance and parse it to Initial + initials[name] = Initial.from_json(value) + + # Handle parameters + parameters = { + par_key: Parameter.from_json(par_dict) + for par_key, par_dict in data.get('parameters', {}).items() } return cls(templates=templates, - parameters=data.get('parameters', {}), + parameters=parameters, initials=initials, annotations=data.get('annotations')) diff --git a/mira/metamodel/templates.py b/mira/metamodel/templates.py index 963393f2e..cb389046b 100644 --- a/mira/metamodel/templates.py +++ b/mira/metamodel/templates.py @@ -7,7 +7,6 @@ "Concept", "Template", "Provenance", - "Unit", "ControlledConversion", "ControlledProduction", "ControlledDegradation", @@ -18,14 +17,11 @@ "GroupedControlledProduction", "GroupedControlledDegradation", "SpecifiedTemplate", - "SympyExprStr", "templates_equal", "context_refinement", - "UNIT_SYMBOLS" ] import logging -import os import sys from collections import ChainMap from itertools import product @@ -52,7 +48,9 @@ except ImportError: from typing_extensions import Annotated -from .utils import safe_parse_expr +from .units import Unit, UNIT_SYMBOLS +from .utils import safe_parse_expr, SympyExprStr + IS_EQUAL = "is_equal" REFINEMENT_OF = "refinement_of" @@ -81,44 +79,6 @@ class Config(BaseModel): ) -class SympyExprStr(sympy.Expr): - @classmethod - def __get_validators__(cls): - yield cls.validate - - @classmethod - def validate(cls, v): - if isinstance(v, cls): - return v - return cls(v) - - @classmethod - def __modify_schema__(cls, field_schema): - field_schema.update(type="string", example="2*x") - - def __str__(self): - return super().__str__()[len(self.__class__.__name__)+1:-1] - - def __repr__(self): - return str(self) - - -class Unit(BaseModel): - """A unit of measurement.""" - class Config: - arbitrary_types_allowed = True - json_encoders = { - SympyExprStr: lambda e: str(e), - } - json_decoders = { - SympyExprStr: lambda e: safe_parse_expr(e) - } - - expression: SympyExprStr = Field( - description="The expression for the unit." - ) - - class Concept(BaseModel): """A concept is specified by its identifier(s), name, and - optionally - its context. @@ -142,6 +102,15 @@ class Concept(BaseModel): ) _base_name: str = pydantic.PrivateAttr(None) + class Config: + arbitrary_types_allowed = True + json_encoders = { + SympyExprStr: lambda e: str(e), + } + json_decoders = { + SympyExprStr: lambda e: sympy.parse_expr(e) + } + def with_context(self, do_rename=False, **context) -> "Concept": """Return this concept with extra context. @@ -297,6 +266,16 @@ def refinement_of( return ontological_refinement and contextual_refinement + @classmethod + def from_json(cls, data) -> "Concept": + # Handle Units + if isinstance(data, Concept): + return data + elif data.get('units'): + data['units'] = Unit.from_json(data['units']) + + return cls(**data) + class Template(BaseModel): """The Template is a parent class for model processes""" @@ -330,6 +309,12 @@ def from_json(cls, data, rate_symbols=None) -> "Template": rate = safe_parse_expr(rate_str, local_dict=rate_symbols) else: rate = None + + # Handle concepts + for concept_key in stmt_cls.concept_keys: + if concept_key in data: + data[concept_key] = Concept.from_json(data[concept_key]) + return stmt_cls(**{k: v for k, v in data.items() if k not in {'rate_law', 'type'}}, rate_law=rate) @@ -726,7 +711,6 @@ def with_context(self, do_rename=False, **context) -> "GroupedControlledProducti ) - class ControlledProduction(Template): """Specifies a process of production controlled by one controller""" @@ -1082,17 +1066,3 @@ def has_controller(template: Template, controller: Concept) -> bool: return template.controller == controller else: raise NotImplementedError - - -def load_units(): - path = os.path.join(os.path.dirname(os.path.abspath(__file__)), - os.pardir, 'dkg', 'resources', 'unit_names.tsv') - with open(path, 'r') as fh: - units = {} - for line in fh.readlines(): - symbol = line.strip() - units[symbol] = sympy.Symbol(symbol) - return units - - -UNIT_SYMBOLS = load_units() diff --git a/mira/metamodel/units.py b/mira/metamodel/units.py new file mode 100644 index 000000000..ffa30a051 --- /dev/null +++ b/mira/metamodel/units.py @@ -0,0 +1,66 @@ +__all__ = [ + 'Unit', + 'person_units', + 'day_units', + 'per_day_units', + 'dimensionless_units', + 'per_day_per_person_units', + 'UNIT_SYMBOLS' +] + +import os +from typing import Dict, Any + +import sympy +from pydantic import BaseModel, Field +from .utils import SympyExprStr + + +def load_units(): + path = os.path.join(os.path.dirname(os.path.abspath(__file__)), + os.pardir, 'dkg', 'resources', 'unit_names.tsv') + with open(path, 'r') as fh: + units = {} + for line in fh.readlines(): + symbol = line.strip() + units[symbol] = sympy.Symbol(symbol) + return units + + +UNIT_SYMBOLS = load_units() + + +class Unit(BaseModel): + """A unit of measurement.""" + class Config: + arbitrary_types_allowed = True + json_encoders = { + SympyExprStr: lambda e: str(e), + } + json_decoders = { + SympyExprStr: lambda e: sympy.parse_expr(e) + } + + expression: SympyExprStr = Field( + description="The expression for the unit." + ) + + @classmethod + def from_json(cls, data: Dict[str, Any]) -> "Unit": + expr_str = data.get('expression') + if expr_str: + data['expression'] = sympy.parse_expr( + expr_str, local_dict=UNIT_SYMBOLS + ) + + assert data.get('expression') is None or not isinstance( + data['expression'], str + ) + return cls(**data) + + +person_units = Unit(expression=sympy.Symbol('person')) +day_units = Unit(expression=sympy.Symbol('day')) +per_day_units = Unit(expression=1/sympy.Symbol('day')) +dimensionless_units = Unit(expression=sympy.Integer('1')) +per_day_per_person_units = Unit(expression=1/(sympy.Symbol('day')*sympy.Symbol('person'))) diff --git a/mira/metamodel/utils.py b/mira/metamodel/utils.py index 038fddc71..f8ff6d440 100644 --- a/mira/metamodel/utils.py +++ b/mira/metamodel/utils.py @@ -1,5 +1,5 @@ __all__ = ['get_parseable_expression', 'revert_parseable_expression', - 'safe_parse_expr'] + 'safe_parse_expr', 'SympyExprStr'] import sympy @@ -20,3 +20,25 @@ def safe_parse_expr(s: str, local_dict=None) -> sympy.Expr: local_dict={get_parseable_expression(k): v for k, v in local_dict.items()} if local_dict else None) + + +class SympyExprStr(sympy.Expr): + @classmethod + def __get_validators__(cls): + yield cls.validate + + @classmethod + def validate(cls, v): + if isinstance(v, cls): + return v + return cls(v) + + @classmethod + def __modify_schema__(cls, field_schema): + field_schema.update(type="string", example="2*x") + + def __str__(self): + return super().__str__()[len(self.__class__.__name__)+1:-1] + + def __repr__(self): + return str(self) diff --git a/mira/modeling/askenet/petrinet.py b/mira/modeling/askenet/petrinet.py index eccf91e3b..22978beaf 100644 --- a/mira/modeling/askenet/petrinet.py +++ b/mira/modeling/askenet/petrinet.py @@ -88,7 +88,7 @@ def __init__(self, model: Model): # } initial = var.data.get('initial_value') if initial is not None: - if isinstance(initial, float): + if isinstance(initial, (float, int)): initial = safe_parse_expr(str(initial)) initial_data = { 'target': name, diff --git a/mira/sources/askenet/petrinet.py b/mira/sources/askenet/petrinet.py index 75e3e0a06..5f25feab2 100644 --- a/mira/sources/askenet/petrinet.py +++ b/mira/sources/askenet/petrinet.py @@ -121,8 +121,10 @@ def template_model_from_askenet_json(model_json) -> TemplateModel: except TypeError: continue - initial = Initial(concept=concepts[initial_state['target']], - value=initial_val) + initial = Initial( + concept=concepts[initial_state['target']].copy(deep=True), + value=initial_val + ) initials[initial.concept.name] = initial # We get observables from the semantics @@ -188,9 +190,9 @@ def template_model_from_askenet_json(model_json) -> TemplateModel: both = set(inputs) & set(outputs) # We can now get the appropriate concepts for each group - input_concepts = [concepts[i] for i in inputs] - output_concepts = [concepts[i] for i in outputs] - controller_concepts = [concepts[i] for i in controllers] + input_concepts = [concepts[i].copy(deep=True) for i in inputs] + output_concepts = [concepts[i].copy(deep=True) for i in outputs] + controller_concepts = [concepts[i].copy(deep=True) for i in controllers] templates.extend(transition_to_templates(rate_obj, input_concepts, @@ -253,9 +255,13 @@ def parameter_to_mira(parameter): """Return a MIRA parameter from a parameter""" distr = Distribution(**parameter['distribution']) \ if parameter.get('distribution') else None - return Parameter(name=parameter['id'], - value=parameter.get('value'), - distribution=distr) + data = { + "name": parameter['id'], + "value": parameter.get('value'), + "distribution": distr, + "units": parameter.get('units') + } + return Parameter.from_json(data) def transition_to_templates(transition_rate, input_concepts, output_concepts, diff --git a/tests/test_model_api.py b/tests/test_model_api.py index 67314b6c2..ecf44f686 100644 --- a/tests/test_model_api.py +++ b/tests/test_model_api.py @@ -6,15 +6,17 @@ from pathlib import Path from typing import List, Union +import sympy from fastapi import FastAPI from fastapi.testclient import TestClient -from mira.examples.sir import sir_parameterized, sir +from mira.examples.sir import sir_parameterized, sir, \ + sir_parameterized_init, sir_init_val_norm from mira.dkg.model import model_blueprint, ModelComparisonResponse from mira.dkg.api import RelationQuery from mira.dkg.web_client import is_ontological_child_web, get_relations_web from mira.metamodel import Concept, ControlledConversion, NaturalConversion, \ - TemplateModel, Distribution + TemplateModel, Distribution, Unit from mira.metamodel.ops import stratify from mira.metamodel.templates import SympyExprStr from mira.metamodel.comparison import TemplateModelComparison, \ @@ -451,3 +453,80 @@ def test_n_way_comparison(self): sorted_json_str(local_response.dict()), sorted_json_str(resp_model.dict()), ) + + def test_counts_to_dimensionless_mira(self): + # Test counts_to_dimensionless + old_beta = sir_parameterized_init.parameters['beta'].value + + response = self.client.post( + "/api/counts_to_dimensionless_mira", + json={ + "model": json.loads(sir_parameterized_init.json()), + "counts_unit": "person", + "norm_factor": sir_init_val_norm, + }, + ) + self.assertEqual(200, response.status_code) + + tm_dimless_json = response.json() + tm_dimless = TemplateModel.from_json(tm_dimless_json) + + for template in tm_dimless.templates: + for concept in template.get_concepts(): + assert concept.units.expression.args[0].equals(1), \ + concept.units + + assert tm_dimless.parameters['beta'].units.expression.args[0].equals( + 1 / sympy.Symbol('day')) + assert tm_dimless.parameters['beta'].value == \ + old_beta * sir_init_val_norm + + assert tm_dimless.initials['susceptible_population'].value == \ + (sir_init_val_norm - 1) / sir_init_val_norm + assert tm_dimless.initials['infected_population'].value == \ + 1 / sir_init_val_norm + assert tm_dimless.initials['immune_population'].value == 0 + + for initial in tm_dimless.initials.values(): + assert initial.concept.units.expression.args[0].equals(1) + + def test_counts_to_dimensionless_amr(self): + # Same test as test_counts_to_dimensionless_mira but with sending + # the model as an askenetpetrinet instead of a mira model + norm = sir_init_val_norm + old_beta = sir_parameterized_init.parameters['beta'].value + + # transform to askenetpetrinet + amr_json = AskeNetPetriNetModel(Model(sir_parameterized_init)).to_json() + + # Post to /api/counts_to_dimensionless_amr + response = self.client.post( + "/api/counts_to_dimensionless_amr", + json={ + "model": amr_json, + "counts_unit": "person", + "norm_factor": norm, + }, + ) + self.assertEqual(200, response.status_code) + + # transform json > amr > + amr_dimless_json = response.json() + tm_dimless = template_model_from_askenet_json(amr_dimless_json) + + for template in tm_dimless.templates: + for concept in template.get_concepts(): + assert concept.units.expression.args[0].equals(1), \ + concept.units + + assert tm_dimless.parameters['beta'].units.expression.args[0].equals( + 1 / sympy.Symbol('day')) + assert tm_dimless.parameters['beta'].value == old_beta * norm + + assert tm_dimless.initials['susceptible_population'].value \ + == (norm - 1) / norm + assert tm_dimless.initials['infected_population'].value == 1 / norm + assert tm_dimless.initials['immune_population'].value == 0 + + for initial in tm_dimless.initials.values(): + assert initial.concept.units.expression.args[0].equals(1) diff --git a/tests/test_ops.py b/tests/test_ops.py index 9eb1bfa9e..f968077e3 100644 --- a/tests/test_ops.py +++ b/tests/test_ops.py @@ -16,7 +16,7 @@ TemplateModel, safe_parse_expr ) -from mira.metamodel.ops import stratify, simplify_rate_law +from mira.metamodel.ops import stratify, simplify_rate_law, counts_to_dimensionless from mira.examples.sir import cities, sir, sir_2_city, sir_parameterized from mira.examples.concepts import infected, susceptible from mira.examples.chime import sviivr @@ -320,3 +320,38 @@ def _make_template(rate_law): (1 - _s('alpha')) * _s('S') * _s('A')) assert templates[1].rate_law.args[0].equals( (1 - _s('alpha')) * _s('beta') * _s('S') * _s('B')) + + +def test_counts_to_dimensionless(): + """Test that counts are converted to dimensionless.""" + from mira.metamodel import Unit + tm = _d(sir_parameterized) + + for template in tm.templates: + for concept in template.get_concepts(): + concept.units = Unit(expression=sympy.Symbol('person')) + tm.initials['susceptible_population'].value = 1e5-1 + tm.initials['infected_population'].value = 1 + tm.initials['immune_population'].value = 0 + + tm.parameters['beta'].units = \ + Unit(expression=1/(sympy.Symbol('person')*sympy.Symbol('day'))) + old_beta = tm.parameters['beta'].value + + for initial in tm.initials.values(): + initial.concept.units = Unit(expression=sympy.Symbol('person')) + + tm = counts_to_dimensionless(tm, 'person', 1e5) + for template in tm.templates: + for concept in template.get_concepts(): + assert concept.units.expression.args[0].equals(1), concept.units + + assert tm.parameters['beta'].units.expression.args[0].equals(1/sympy.Symbol('day')) + assert tm.parameters['beta'].value == old_beta*1e5 + + assert tm.initials['susceptible_population'].value == (1e5-1)/1e5 + assert tm.initials['infected_population'].value == 1/1e5 + assert tm.initials['immune_population'].value == 0 + + for initial in tm.initials.values(): + assert initial.concept.units.expression.args[0].equals(1)