Skip to content

Commit

Permalink
update more docs with internal references
Browse files Browse the repository at this point in the history
  • Loading branch information
TorkelE committed May 28, 2024
1 parent aa9735c commit 0d13a40
Show file tree
Hide file tree
Showing 8 changed files with 46 additions and 42 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@ Catalyst models are created through the `@reaction_network` *macro*. For more in

The `@reaction_network` command is followed by the `begin` keyword, which is followed by one line for each *reaction* of the model. Each reaction consists of a *reaction rate*, followed by the reaction itself. The reaction contains a set of *substrates* and a set of *products* (what is consumed and produced by the reaction, respectively). These are separated by a `-->` arrow. Finally, the model ends with the `end` keyword.

Here, we create a simple *birth-death* model, where a single species ($X$) is created at rate $b$, and degraded at rate $d$. The model is stored in the variable `rn`.
Here, we create a simple [*birth-death* model](@ref basic_CRN_library_bd), where a single species ($X$) is created at rate $b$, and degraded at rate $d$. The model is stored in the variable `rn`.
```@example ex2
rn = @reaction_network begin
b, 0 --> X
Expand Down Expand Up @@ -136,7 +136,7 @@ Pkg.add("JumpProcesses")
using JumpProcesses
```

This time, we will declare a so-called [SIR model for an infectious disease](https://en.wikipedia.org/wiki/Compartmental_models_in_epidemiology#The_SIR_model). Note that even if this model does not describe a set of chemical reactions, it can be modelled using the same framework. The model consists of 3 species:
This time, we will declare a so-called [SIR model for an infectious disease](@ref basic_CRN_library_sir). Note that even if this model does not describe a set of chemical reactions, it can be modelled using the same framework. The model consists of 3 species:
* $S$, the amount of *susceptible* individuals.
* $I$, the amount of *infected* individuals.
* $R$, the amount of *recovered* (or *removed*) individuals.
Expand Down
17 changes: 9 additions & 8 deletions docs/src/model_creation/chemistry_related_functionality.md
Original file line number Diff line number Diff line change
@@ -1,9 +1,10 @@
# [Chemistry-related functionality](@id chemistry_functionality)

While Catalyst has primarily been designed around the modelling of biological systems, reaction network models are also common across chemistry. This section describes two types of functionality, that while of general interest, should be especially useful in the modelling of chemical systems.
While Catalyst has primarily been designed around the modelling of biological systems, reaction network models are also common in chemistry. This section describes two types of functionality, that while of general interest, should be especially useful in the modelling of chemical systems.
- The `@compound` option, which enables the user to designate that a specific species is composed of certain subspecies.
- The `balance_reaction` function, which enables the user to balance a reaction so the same number of components occur on both sides.


## Modelling with compound species

### Creating compound species programmatically
Expand All @@ -17,7 +18,7 @@ Next, we create the `CO2` compound species:
```@example chem1
@compound CO2 ~ C + 2O
```
Here, the compound is the first argument to the macro, followed by its component (with the left-hand and right-hand sides separated by a `~` sign). While non-compound species (such as `C` and `O`) have their independent variable (in this case `t`) designated, independent variables are generally not designated for compounds (these are instead directly inferred from their components). Components with non-unit stoichiometries have this value written before the component (generally, the rules for designating the components of a compound are identical to those of designating the substrates or products of a reaction). The created compound, `CO2`, is also a species, and can be used wherever e.g. `C` can be used:
Here, the compound is the first argument to the macro, followed by its component (with the left-hand and right-hand sides separated by a `~` sign). While non-compound species (such as `C` and `O`) have their independent variable (in this case `t`) designated, independent variables are generally not designated for compounds (these are instead directly inferred from their components). Components with non-unitary stoichiometries have this value written before the component (generally, the rules for designating the components of a compound are identical to those of designating the substrates or products of a reaction). The created compound, `CO2`, is also a species, and can be used wherever e.g. `C` can be used:
```@example chem1
isspecies(CO2)
```
Expand All @@ -32,7 +33,7 @@ Alternatively, we can retrieve the components and their stoichiometric coefficie
```@example chem1
component_coefficients(CO2)
```
Finally, it is possible to check whether a species is a compound or not using the `iscompound` function:
Finally, it is possible to check whether a species is a compound using the `iscompound` function:
```@example chem1
iscompound(CO2)
```
Expand Down Expand Up @@ -68,7 +69,7 @@ end
```
When creating compound species using the DSL, it is important to note that *every component must be known to the system as a species, either by being declared using the `@species` or `@compound` options, or by appearing in a reaction*. E.g. the following is not valid
```julia
rn = @reaction_network begin``
rn = @reaction_network begin
@compounds begin
C2O ~ C + 2O
H2O ~ 2H + O
Expand All @@ -77,7 +78,7 @@ rn = @reaction_network begin``
(k1,k2), H2O+ CO2 <--> H2CO3
end
```
as the components `C`, `H`, and `O` are not declared as a species anywhere. Please also note that only `@compounds` can be used as an option in the DSL, not `@compound`.
as the components `C`, `H`, and `O` are not declared as species anywhere. Please also note that only `@compounds` can be used as an option in the DSL, not `@compound`.

### Designating metadata and default values for compounds
Just like for normal species, it is possible to designate metadata and default values for compounds. Metadata is provided after the compound name, but separated from it by a `,`:
Expand All @@ -92,10 +93,10 @@ If both default values and meta data are provided, the metadata is provided afte
```@example chem1
@compound (CO2 = 2.0, [unit="mol"]) ~ C + 2O
```
In all of these cases, the side to the left of the `~` must be enclosed within `()`.
In all of these cases, the left-hand side must be enclosed within `()`.

### Compounds with multiple independent variables
While we generally do not need to specify independent variables for compound, if the components (together) have more than one independent variable, this have to be done:
While we generally do not need to specify independent variables for compound, if the components (together) have more than one independent variable, this *must be done*:
```@example chem1
t = default_t()
@variables s
Expand All @@ -113,7 +114,7 @@ Here, the reaction rate (`k`) is not involved in the reaction balancing. We use
```@example chem1
balance_reaction(rx)
```
which correctly finds the (rather trivial) solution `C + 2O --> CO2`. Here we note that `balance_reaction` actually returns a vector. The reason is that the reaction balancing problem may have several solutions. Typically, there is only a single solution (in which case this is the vector's only element). No, or an infinite number of, solutions is also possible depending on the given reaction.
which correctly finds the (rather trivial) solution `C + 2O --> CO2`. Here we note that `balance_reaction` actually returns a vector. The reason is that, in some cases, the reaction balancing problem does not have a single obvious solution. Typically, a single solution is the obvious candidate (in which case this is the vector's only element). However, when this is not the case, the vector instead contain several reactions (from which a balanced reaction cab be generated).

Let us consider a more elaborate example, the reaction between ammonia (NH₃) and oxygen (O₂) to form nitrogen monoxide (NO) and water (H₂O). Let us first create the components and the unbalanced reaction:
```@example chem2
Expand Down
32 changes: 16 additions & 16 deletions docs/src/model_creation/dsl_advanced.md
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ using Catalyst
```

## [Explicit specification of network species and parameters](@id dsl_advanced_options_declaring_species_and_parameters)
[Previously](@ref ref), we mentioned that the DSL automatically determines which symbols correspond to species and which to parameters. This is done by designating everything that appears as either a substrate or a product as a species, and all remaining quantities as parameters (i.e. those only appearing within rates or [stoichiometric constants](@ref ref)). Sometimes, one might want to manually override this default behaviour for a given symbol. I.e. consider the following model, where the conversion of a protein `P` from its inactive form (`Pᵢ`) to its active form (`Pₐ`) is catalysed by an enzyme (`E`). Using the most natural description:
Previously, we mentioned that the DSL automatically determines which symbols correspond to species and which to parameters. This is done by designating everything that appears as either a substrate or a product as a species, and all remaining quantities as parameters (i.e. those only appearing within rates or [stoichiometric constants](@ref dsl_description_stoichiometries_parameters)). Sometimes, one might want to manually override this default behaviour for a given symbol. I.e. consider the following model, where the conversion of a protein `P` from its inactive form (`Pᵢ`) to its active form (`Pₐ`) is catalysed by an enzyme (`E`). Using the most natural description:
```@example dsl_advanced_explicit_definitions
catalysis_sys = @reaction_network begin
k*E, Pᵢ --> Pₐ
Expand Down Expand Up @@ -74,7 +74,7 @@ end
parameters(dimerisation)
```
!!! danger
Generally, Catalyst and the SciML ecosystem *do not* guarantee that parameter and species order are preserved throughout various operations on a model. Writing programs that depend on these orders is *strongly discouraged*. There are, however, some legacy packages which still depend on order (one example is provided [here](@ref ref)). In these situations, this might be useful. However, in these cases, it is recommended that the user is extra wary, and also checks the order manually.
Generally, Catalyst and the SciML ecosystem *do not* guarantee that parameter and species order are preserved throughout various operations on a model. Writing programs that depend on these orders is *strongly discouraged*. There are, however, some legacy packages which still depend on order (one example can be found [here](@ref optimization_parameter_fitting_basics)). In these situations, this might be useful. However, in these cases, it is recommended that the user is extra wary, and also checks the order manually.

!!! note
The syntax of the `@species` and `@parameters` options is identical to that used by the `@species` and `@parameters` macros [used in programmatic modelling in Catalyst](@ref programmatic_CRN_construction) (for e.g. designating metadata or initial conditions). Hence, if one has learnt how to specify species/parameters using either approach, that knowledge can be transferred to the other one.
Expand Down Expand Up @@ -125,7 +125,7 @@ rn = @reaction_network begin
end
tspan = (0.0, 10.0)
p = [:p => 1.0, :D => 0.2]
p = [:p => 1.0, :d => 0.2]
oprob = ODEProblem(rn, u0, tspan, p)
sol = solve(oprob)
plot(sol)
Expand Down Expand Up @@ -171,18 +171,18 @@ two_state_system = @reaction_network begin
(ka,kD), Xi <--> Xa
end
```
A metadata can be given to only a subset of a system's species/parameters, and a quantity can be given several metadata entries. To give several metadata, separate each by a `,`. Here we only provide a description for `kA`, for which we also provide a [bounds metadata](@ref https://docs.sciml.ai/ModelingToolkit/dev/basics/Variable_metadata/#Bounds),
A metadata can be given to only a subset of a system's species/parameters, and a quantity can be given several metadata entries. To give several metadata, separate each by a `,`. Here we only provide a description for `kA`, for which we also provide a [`bounds` metadata](@ref https://docs.sciml.ai/ModelingToolkit/dev/basics/Variable_metadata/#Bounds),
```@example dsl_advanced_metadata
two_state_system = @reaction_network begin
@parameters kA [description="X's activation rate", bound=(0.01,10.0)]
@parameters kA [description="X's activation rate", bounds=(0.01,10.0)]
(ka,kD), Xi <--> Xa
end
```

It is possible to add both default values and metadata to a parameter/species. In this case, first provide the default value, and next the metadata. I.e. to in the above example set $kA$'s default value to $1.0$ we use
```@example dsl_advanced_metadata
two_state_system = @reaction_network begin
@parameters kA=1.0 [description="X's activation rate", bound=(0.01,10.0)]
@parameters kA=1.0 [description="X's activation rate", bounds=(0.01,10.0)]
(ka,kD), Xi <--> Xa
end
```
Expand All @@ -191,10 +191,10 @@ When designating metadata for species/parameters in `begin ... end` blocks the s
```@example dsl_advanced_metadata
two_state_system = @reaction_network begin
@parameters begin
kA, [description="X's activation rate", bound=(0.01,10.0)]
kD = 1.0, [description="X's deactivation rate"]
kA, [description="X's activation rate", bounds=(0.01,10.0)]
kD = 1.0, [description="X's deactivation rate"]
end
(ka,kD), Xi <--> Xa
(kA,kD), Xi <--> Xa
end
```

Expand Down Expand Up @@ -229,7 +229,7 @@ end
A common use-case for constant species is when modelling systems where some species are present in such surplus that their amounts the reactions' effect on it is negligible. A system which is commonly modelled this way is the [Brusselator](https://en.wikipedia.org/wiki/Brusselator).

### [Designating parameter types](@id dsl_advanced_options_parameter_types)
Sometimes it is desired to designate that a parameter should have a specific [type](@ref ref). When supplying this parameter's value to e.g. an `ODEProblem`, that parameter will then be restricted to that specific type. Designating a type is done by appending the parameter with `::` followed by its type. E.g. in the following example we specify that the parameter `n` (the number of `X` molecules in the `Xn` polymer) must be an integer (`Int64`)
Sometimes it is desired to designate that a parameter should have a specific [type](https://docs.julialang.org/en/v1/manual/types/). When supplying this parameter's value to e.g. an `ODEProblem`, that parameter will then be restricted to that specific type. Designating a type is done by appending the parameter with `::` followed by its type. E.g. in the following example we specify that the parameter `n` (the number of `X` molecules in the `Xn` polymer) must be an integer (`Int64`)
```@example dsl_advanced_parameter_types
using Catalyst # hide
polymerisation_network = @reaction_network begin
Expand All @@ -238,7 +238,7 @@ polymerisation_network = @reaction_network begin
end
nothing # hide
```
Generally, when simulating models with mixed parameter types, it is recommended to [declare parameter values as tuples, rather than vectors](@ref ref), e.g.:
Generally, when simulating models with mixed parameter types, it is recommended to [declare parameter values as tuples, rather than vectors](@ref simulation_intro_ODEs_input_forms), e.g.:
```@example dsl_advanced_parameter_types
ps = (:kB => 0.2, :kD => 1.0, :n => 2)
nothing # hide
Expand All @@ -254,7 +254,7 @@ nothing # hide
```

### [Vector-valued species or parameters](@id dsl_advanced_options_vector_variables)
Sometimes, one wishes to declare a large number of similar parameters or species. This can be done by *creating them as vectors*. E.g. below we create a [two-state system](@ref ref). However, instead of declaring `X1` and `X2` (and `k1` and `k2`) as separate entities, we declare them as vectors:
Sometimes, one wishes to declare a large number of similar parameters or species. This can be done by *creating them as vectors*. E.g. below we create a [two-state system](@ref basic_CRN_library_two_states). However, instead of declaring `X1` and `X2` (and `k1` and `k2`) as separate entities, we declare them as vectors:
```@example dsl_advanced_vector_variables
using Catalyst # hide
two_state_model = @reaction_network begin
Expand Down Expand Up @@ -316,7 +316,7 @@ end
rn1 == rn2
```

Setting model names is primarily useful for [hierarchical modelling](@ref ref), where network names are appended to the display names of subnetworks' species and parameters.
Setting model names is primarily useful for [hierarchical modelling](@ref compositional_modeling), where network names are appended to the display names of subnetworks' species and parameters.

## [Creating observables](@id dsl_advanced_options_observables)
Sometimes one might want to use observable variables. These are variables with values that can be computed directly from a system's state (rather than having their values implicitly given by reactions or equations). Observables can be designated using the `@observables` option. Here, the `@observables` option is followed by a `begin ... end` block with one line for each observable. Each line first gives the observable, followed by a `~` (*not* a `=`!), followed by an expression describing how to compute it.
Expand Down Expand Up @@ -350,11 +350,11 @@ sol[:Xtot]
to get a vector with `Xtot`'s value throughout the simulation. We can also use
```@example dsl_advanced_observables
using Plots
plot(sol; idxs = [:Xtot, :Ytot])
plot(sol; idxs = :Xtot)
```
to plot the observables (rather than the species).

Observables can be defined using complicated expressions containing species, parameters, and [variables](@ref ref) (but not other observables). In the following example (which uses a [parametric stoichiometry](@ref ref)) `X` polymerises to form a complex `Xn` containing `n` copies of `X`. Here, we create an observable describing the total number of `X` molecules in the system:
Observables can be defined using complicated expressions containing species, parameters, and [variables](@ref ref) (but not other observables). In the following example (which uses a [parametric stoichiometry](@ref dsl_description_stoichiometries_parameters)) `X` polymerises to form a complex `Xn` containing `n` copies of `X`. Here, we create an observable describing the total number of `X` molecules in the system:
```@example dsl_advanced_observables
rn = @reaction_network begin
@observables Xtot ~ X + n*Xn
Expand All @@ -378,7 +378,7 @@ Observables are by default considered [variables](@ref ref) (not species). To de
```@example dsl_advanced_observables
rn = @reaction_network begin
@species Xtot(t)
@observables Xtot ~ X + n*XnXY
@observables Xtot ~ X + n*Xn
(kB,kD), n*X <--> Xn
end
nothing # hide
Expand Down
Loading

0 comments on commit 0d13a40

Please sign in to comment.