<jupyter_start><jupyter_text>Details of defining a model$$\def\n{\mathbf{n}}\def\x{\mathbf{x}}\def\N{\mathbb{\mathbb{N}}}\def\X{\mathbb{X}}\def\NX{\mathbb{\N_0^\X}}\def\C{\mathcal{C}}\def\Jc{\mathcal{J}_c}\def\DM{\Delta M_{c,j}}\newcommand\diff{\mathop{}\!\mathrm{d}}\def\Xc{\mathbf{X}_c}\def\Yc{\mathbf{Y}_c}\newcommand{\muset}[1]{\dot{\{}1\dot{\}}}$$This notebook provides a more detailed view of how `TransitionClass`es of a population model are defined.It describes the objects hidden behind syntactic sugar when writing something like`TransitionClass([x] + [y] -to> [x+y], 'k_F')`This will hopefully help deciphering error messages when making mistakes in the "short-hand" notation.<jupyter_code># initialize sympy printing (for latex output)
from sympy import init_printing
init_printing()

# import functions and classes for compartment models
from compartor import *<jupyter_output><empty_output><jupyter_text>BasicsThe left-hand-side and right-hand-site of transitions are expressions involving zero or more `Compartments`.To define a compartment $[\x]$, first create a content variable $\x \in \N_0^D$, then define a compartment $[\x]$ with content $\x$.<jupyter_code>x = Content('x')
Compartment(x)<jupyter_output><empty_output><jupyter_text>Content variables are $D$-dimensional. In equations, `x[d]` can be used to denote the copy number of chemical species `d`.Next, we can build a stoichiometric equation using the `Transition` class, which takes two parameter: the reactants expression and the products expression. For empty reactants or products, use `EmptySet`.The following defines the "Exit" transition from the nested birth-death process example.<jupyter_code>transition_E = Transition(Compartment(x), EmptySet, name='E')
display(transition_E)<jupyter_output><empty_output><jupyter_text>`name="E"` is an optional parameter, defining the subscript used for the transition class. (It is only used for displaying equations, e.g., the $h_E$ over the arrow in the above example.)To define a propensity for the transition class, we have to define a content-independent rate constant $k_c$ and a reactant tuning function $g_c$.<jupyter_code>k_E = Constant('k_E')
g_E = 1<jupyter_output><empty_output><jupyter_text>To specify the transition class with its propensity, use the constructor `TransitionClass(transition_c, k_c, g_c)`.<jupyter_code>Exit = TransitionClass(transition_E, k_E, g_E)
display(Exit)<jupyter_output><empty_output><jupyter_text>The optional 4th argument of `TransitionClass` specifies the outcome distribution $\pi_c$.This must be defined, if content variables occur in products that do not occur in reactants.The type of $\pi_c$ is `OutcomeDistribution`, which is a class comprising* an expression or symbol to use for displaying $\pi_c$ in compound expressions* a function `conditional_expectation` that takes an expression over reactants $\Xc$, and returns it's expectation over all product compartment variables.There are generators for several predefined outcome distributions. If nothing is specified, as in the above "Exit" transition example, `OutcomeDistribution.Identity()` is used. This has symbol $1$ and `conditional_expectation` returns the expression over reactants $\Xc$ as is. This is appropriate when $\Yc = \emptyset$ or, more precisely, if all content variables occurring in products already occur in reactants.Other generators include* `Poisson()` which specifies a poisson distribution over a scalar, and* `Uniform()` which specifies a uniform distribution over a scalar.For example, `Poisson()` can be used to define an "Intake" transition class as follows:<jupyter_code>from sympy import Symbol
y = Content('y')
transition_I = Transition(EmptySet, Compartment(y), name='I')
k_I = Constant('k_I')
g_I = 1
pi_I = OutcomeDistribution.Poisson(
    Symbol("\pi_{Poiss}(y; \lambda)"),
    y[0],
    Symbol("\lambda", positive=True))
Intake = TransitionClass(transition_I, k_I, g_I, pi_I)
display(Intake)<jupyter_output><empty_output><jupyter_text>When displaying a `TransitionClass`, the assembled propensity function is shown, where $w(\n;\Xc)$ is derived from the reactants in the stoichiometric equation.The function `display_propensity_details()` displays the constituent expressions $k_c$, $g_c$, and $\pi_c$.<jupyter_code>display_propensity_details(Exit)<jupyter_output><empty_output><jupyter_text>To specify a reaction network, collect the transition classes in a list.The function `display_transition_classes()` displays the equations and propensities of the network.<jupyter_code>transitions = [Exit, Intake]
display_transition_classes(transitions)<jupyter_output><empty_output><jupyter_text>Chemical EventsTo specify chemical reactions, we need the `ContentChange()` class. It is constructed with $D$ (number of species) integer arguments, specifying the change (to a compartment variable) in each species.For example, the "birth" transition from the nested birth-death process example can be defined as<jupyter_code>x = Content('x')
Transition(Compartment(x), Compartment(x + ContentChange(1)), name='b')<jupyter_output><empty_output><jupyter_text>Note, that `Compartment(...)` in the product encloses `x + ContentChange(1)`, so the `+` acts on the compartment content, rather than specifying a set of compartments.Besides addition, also multiplication can be used inside compartment content equations, e.g.<jupyter_code>Transition(Compartment(x), Compartment(2*x + ContentChange(1)))<jupyter_output><empty_output><jupyter_text>Compartment EventsTo specify compartment events, simply construct (multi-)sets of compartments using `+`.For example, a "Coagulation" transition can be specified as<jupyter_code>x = Content('x')
y = Content('y')
Transition(Compartment(x) + Compartment(y), Compartment(x + y), name='C')<jupyter_output><empty_output><jupyter_text>In principle, also multiplication can be used, e.g., `2*Compartment(x)` means "two compartments of content 'x'".But the moment derivation mechanics is not ready for that yet...  Example: nested birth-death processFinally, as a complete example, here is the definition of nested birth-death process.<jupyter_code>y = Content('y')
x = Content('x')

# Intake
transition_I = Transition(EmptySet, Compartment(y), name='I')
k_I = Constant('k_I')
g_I = 1
pi_I = OutcomeDistribution.Poisson(Symbol("\pi_{Poiss}(y; \lambda)"), y[0], Symbol("\lambda"))
Intake = TransitionClass(transition_I, k_I, g_I, pi_I)

# Exit
transition_E = Transition(Compartment(x), EmptySet, name='E')
k_E = Constant('k_E')
g_E = 1
Exit = TransitionClass(transition_E, k_E, g_E)

# birth
transition_b = Transition(Compartment(x), Compartment(x + ContentChange(1)), name='b')
k_b = Constant('k_b')
g_b = 1
Birth = TransitionClass(transition_b, k_b, g_b)

# death
transition_d = Transition(Compartment(x), Compartment(x + ContentChange(-1)), name='d')
k_d = Constant('k_d')
g_d = x[0] # TODO x should be enough here, in case D=1.
Death = TransitionClass(transition_d, k_d, g_d)

transitions = [Intake, Exit, Birth, Death]
display_transition_classes(transitions)<jupyter_output><empty_output><jupyter_text>The short-hand definition of the same model is given below. The two are 100% equivalent, they generate the same objects in the `transitions` list.<jupyter_code>x = Content('x')
y = Content('y')

pi_I = OutcomeDistribution.Poisson(Symbol("\pi_{Poiss}(y; \lambda)"), y[0], Symbol("\lambda"))

Intake = TransitionClass( {}  -to> [y],   'k_I', pi=pi_I, name='I')
Exit   = TransitionClass( [x] -to> {},    'k_E',          name='E')
Birth  = TransitionClass( [x] -to> [x+1], 'k_b',          name='b')
Death  = TransitionClass( [x] -to> [x-1], 'k_d', x[0],    name='d')

transitions = [Intake, Exit, Birth, Death]
display_transition_classes(transitions)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12042761463163318, 'fraction_numerical': 0.0037351880473982485, 'mean_word_length': 6.62770137524558, 'pattern_counts': {'":': 0, '<': 53, '<?xml version=': 0, '>': 58, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '650602', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 2438, 'n_tokens_neox': 2327, 'n_words': 765, 'path': '(EXTRA) Model definition details.ipynb', 'repo_name': 'zechnerlab/Compartor', 'size': '7764'}