Definitions and Evaluations of Rules in SimBiology Models
Rules are mathematical expressions that allow you to define or modify model quantities, namely compartment capacity, species amount, or parameter value.
Rules can take the form of initial assignments, assignments during the course of a simulation (repeated assignments), algebraic relationships, or differential equations (rate rules). Details of each type of rule are described next.
An initial assignment rule lets you specify the initial value of a model quantity as a numeric value or as a function of other model quantities. It is evaluated once at the beginning of a simulation.
An initial assignment rule is expressed as
Expression, and the rule is specified as the
Expression. For example, you could write an initial
assignment rule to set the initial amount of
species2 to be
species2 = k *
k is a constant parameter.
A repeated assignment rule lets you specify the value of a quantity as a numeric value or as a function of other quantities repeatedly during the simulation. It is evaluated at every time step, which is determined by the solver during the simulation process.
A repeated assignment rule is expressed as
Expression, and the rule is specified as the
Expression. For example, to repeatedly evaluate the total
species amount by summing up the species in different compartments, you could enter:
xTotal = c1.X + c2.X, where
xTotal is a
c2 are the name
of compartments where species x resides.
An algebraic rule lets you specify mathematical constraints on one or more model quantities that must hold during a simulation. It is evaluated continuously during a simulation.
An algebraic rule takes the form
0 = Expression, and the rule
is specified as the
Expression. For example, if you have a mass
conservation equation such as
species_total = species1 +
species2, write the corresponding algebraic rule as
species2 - species_total.
However, repeated assignment rules are mathematically equivalent to algebraic rules, but result in exact solutions instead of approximated solutions. Therefore, it is recommended that you use repeated assignment rules instead of algebraic rules whenever possible. Use algebraic rules only when:
You cannot analytically solve the equations to get a closed-form solution. For example, there is no closed-form solution for
x^4 + ax^3 + bx^2 + cx + k = 0whereas the closed-form solution for
kx – c = 0is
x = c/k.
You have multiple equations with multiple unknowns, and they could be inconvenient to solve.
If you use an algebraic rule, rate rule, or repeated assignment to vary the value
of a parameter or compartment during the simulation, make sure the
ConstantValue property of the
ConstantCapacity of the compartment
is set to
Repeated Assignment vs. Algebraic Rules
Repeated assignment rules are mathematically equivalent to algebraic rules, but result in exact solutions. However, algebraic rules are solved numerically, and the accuracy depends on the error tolerances specified in the simulation settings. In addition, there are several advantages to repeated assignment rules such as better computational performance, more accurate results since no rules have to be solved numerically (hence no approximations), and sensitivity analysis support.
If you can analytically solve for a variable, use a repeated assignment rule instead of an algebraic rule.
In repeated assignment rules, the constrained variable is explicitly defined as the left-hand side, whereas in algebraic rules it is inferred from the degrees of freedom in the system of equations. See also Considerations When Imposing Constraints.
A rate rule represents a differential equation and lets you specify the time derivative of a model quantity. It is evaluated continuously during a simulation.
A rate rule is represented as , which is expressed in SimBiology as
Expression. For example, if you have a differential equation for
species x, , write the rate rule as:
x = k * (y +
For more examples, see Rate Rule Examples.
Evaluation Order of Rules
At the start of the simulation (that is, at simulation time = 0), SimBiology® evaluates the initial assignment and repeated assignment rules as a set of simultaneous constraints. SimBiology treats the rules as a unified system of constraints and automatically reorders and evaluates them. The order in which the rules appear in the model has no effect on the simulation results.
If a quantity is being modified by an assignment rule, the rule replaces initial
value properties, such as
Value. Similarly, a
variant altering such quantities has no effect because the value is superseded by
the assignment rules.
SimBiology throws an error if the model has circular dependencies in the initial assignment and repeated assignment rules. In other words, initial assignments and repeated assignments cannot have a variable that is explicitly or implicitly referenced on both the left- and right-hand sides of the equation.
For instance, you cannot create circular sets of assignments such as
b + 1 and
b = a + 1, where a and
b are explicitly referenced on both sides of the equation. An
example of an implicit reference is when an assignment rule references a species in
concentration. In this case, the compartment that contains the species is implicitly
You might observe different simulation results with respect to initial
assignments for previous releases of SimBiology (R2017a or earlier). To recover
the same simulation results at time = 0, as in R2017a or earlier, use the
updateInitialAssignments function in the command line. If you
are using the SimBiology app, right-click the model from the
Project Workspace and select Remove Order
Conservation of Amounts During Simulation
During a simulation (that is, at simulation time > 0), SimBiology conserves species amounts rather than concentrations if there are any changes to the volume of a compartment where the species reside. In other words, if you have a repeated assignment rule or an event that changes the volume, then you see the effect of conservation of species amounts at time > 0.
However, at the beginning of a simulation (that is, at simulation time = 0), the concept of amount conservation does not apply because there are no changes before time = 0. Only one set of initial conditions exists and SimBiology uses the conditions at the start of the simulation. Specifically, at time = 0, SimBiology:
Initializes variables for species, compartments, and parameters using the corresponding
Updates the values by replacing them with the corresponding alternate values from variants, if any.
Updates the values by evaluating initial assignment and repeated assignment rules as a set of simultaneous constraints. Therefore, the assignment rules replace initial values if model quantities are being modified by such rules or variants.
In previous releases (R2017a or earlier), if a repeated assignment changed a
compartment volume, SimBiology used the compartment capacity to determine the
initial amount and conserved it when the compartment volume changed at time = 0.
In R2017b or later, SimBiology uses the
property of the species as the initial condition at time = 0. Consider the
m = sbiomodel('m1') v = addcompartment(m,'v',10,'ConstantCapacity',0,'CapacityUnit','liter') p = addparameter(m,'p','ValueUnit','liter') r = addrule(m,'v = 100 * p','repeatedAssignment') s = addspecies(v,'s',50,'InitialAmountUnit','milligram/liter')
50 milligram/liter * 10 liter = 500 milligram, and then applied the repeated assignment rule
v = 100 liter. So, the concentration of s was then calculated and reported as
500 milligram/100 liter = 5 milligram/literat time = 0.
In R2017b or later, SimBiology uses the
property of species s, and reports the initial amount of
Writing Rule Expressions
Use MATLAB® syntax to write a mathematical expression for a rule. Note that no semicolon or comma is needed at the end of a rule expression. If your algebraic, repeated assignment, or rate rule expression is not continuous or differentiable, see Using Events to Address Discontinuities in Rule and Reaction Rate Expressions before simulating your model.
Considerations When Imposing Constraints
Suppose that you have a species
y whose amount is determined by
y = m * x - c. In SimBiology, the algebraic rule to
describe this constraint is written as
m * x - c - y. If you want
to use this rule to determine the value of
c must be
variables or constants whose values are known or determined by other equations.
Therefore, you must ensure that the system of equations is not overconstrained or
underconstrained. For instance, if you have more equations than unknowns, then the
system is overconstrained. Conversely, if you have more unknowns than the equations,
then the system is underconstrained.
Rate Rule Examples
The following examples show how to create rate rules for different applications.
Create a Rate Rule for a Constant Rate of Change
This example shows how to increase the amount or concentration of a species by a constant value using the zero-order rate rule. For example, suppose species
x increases by a constant rate
k. The rate of change is:
Set the initial amount of species
x to 2, and the value of parameter
k to 1. Use the following commands to set up a SimBiology model accordingly and simulate it.
m = sbiomodel('m'); c = addcompartment(m,'comp'); s = addspecies(m,'x','InitialAmount',2); p = addparameter(m,'k','Value',1); r = addrule(m,'x = k','RuleType','rate'); [t,sd,species] = sbiosimulate(m); plot(t,sd); legend(species) xlabel('Time'); ylabel('Species Amount');
Alternatively, you could model a constant increase in a species using the Mass Action reaction
x with the forward rate constant
clear m = sbiomodel('m'); c = addcompartment(m,'comp'); s = addspecies(m,'x','InitialAmount',2); r = addreaction(m,'null -> x'); kl = addkineticlaw(r,'MassAction'); p = addparameter(kl,'k','Value',1); kl.ParameterVariableNames = 'k'; [t,sd,species] = sbiosimulate(m); plot(t,sd); legend(species) xlabel('Time'); ylabel('Species Amount');
Create a Rate Rule for an Exponential Rate of Change
This example shows how to change the amount of a species similar to a first-order reaction using the first-order rate rule. For example, suppose the species
x decays exponentially. The rate of change of species
The analytical solution is:
where is the amount of species at time t, and is the initial amount. Use the following commands to set up a SimBiology model accordingly and simulate it.
m = sbiomodel('m'); c = addcompartment(m,'comp'); s = addspecies(m,'x','InitialAmount',2); p = addparameter(m,'k','Value',1); r = addrule(m,'x = -k * x','RuleType','rate'); [t,sd,species] = sbiosimulate(m); plot(t,sd); legend(species); xlabel('Time'); ylabel('Species Amount');
If the amount of a species
x is determined by a rate rule and
x is also in a reaction,
x must have its
BoundaryCondition property set to
true. For example, with a reaction
a -> x and a rate rule , set the
BoundaryCondition property of species
true so that a differential rate term is not created from the reaction. The amount of
x is determined solely by a differential rate term from the rate rule. If the
BoundaryCondition property is set to
false, you will get the following error message such as
Invalid rule variable 'x' in rate rule or reaction.
Create a Rate Rule to Define a Differential Rate Equation
Many mathematical models in the literature are described with differential
rate equations for the species. You could manually convert the equations to
reactions, or you could enter the equations as rate rules. For example, you
could enter the following differential rate equation for a species
as a rate rule in SimBiology:
C = vi - (vd*X*C)/(Kc + C) - kd*C
Create a Rate Rule for the Rate of Change That Is Determined by Another Species
This example shows how to create a rate rule where a species from one reaction can determine the rate of another reaction if it is in the second reaction rate equation. Similarly, a species from a reaction can determine the rate of another species if it is in the rate rule that defines that other species. Suppose you have a SimBiology model with three species (
c), one reaction (
a -> b), and two parameters (
k2). The rate equation is defined as , and rate rule is . The solution for the species in the reaction are:
Since the rate rule is dependent on the reaction, . The solution is:
Enter the following commands to set up a SimBiology model accordingly and simulate it.
m = sbiomodel('m'); c = addcompartment(m,'comp'); s1 = addspecies(m,'a','InitialAmount',10,'InitialAmountUnits','mole'); s2 = addspecies(m,'b','InitialAmount',0,'InitialAmountUnits','mole'); s3 = addspecies(m,'c','InitialAmount',5,'InitialAmountUnits','mole'); rxn = addreaction(m,'a -> b'); kl = addkineticlaw(rxn,'MassAction'); p1 = addparameter(kl,'k1','Value',1,'ValueUnits','1/second'); rule = addrule(m,'c = k2 * a','RuleType','rate'); kl.ParameterVariableNames = 'k1'; p2 = addparameter(m,'k2','Value',1,'ValueUnits','1/second'); [t,sd,species] = sbiosimulate(m); plot(t,sd); legend(species); xlabel('Time'); ylabel('Species Amount');