|
2 | 2 |
|
3 | 3 | # For a reaction system, list of measured quantities and known parameters, generate a StructuralIdentifiability compatible ODE.
|
4 | 4 | """
|
5 |
| - si_ode(rs::ReactionSystem; measured_quantities=observed(rs), known_p = Num[], ignore_no_measured_warn=false) |
| 5 | +make_si_ode(rs::ReactionSystem; measured_quantities=observed(rs), known_p = Num[], ignore_no_measured_warn=false) |
6 | 6 |
|
7 | 7 | Creates a ODE system of the form used within the StructuralIdentifiability.jl package. The output system is compatible with all StructuralIdentifiability functions.
|
8 | 8 |
|
9 | 9 | Arguments:
|
10 | 10 | - `rs::ReactionSystem`; The reaction system we wish to convert to an ODE.
|
11 |
| -- `measured_quantities=observed(rs)`: The quantities of the system we can measure. May either be equations (e.g. `x1 + x2`), or single species (e.g. the symbolic `x`, `rs.s`, or the symbol `:x`). Defaults to the system's observables. |
| 11 | +- `measured_quantities=[]`: The quantities of the system we can measure. May either be equations (e.g. `x1 + x2`), or single species (e.g. the symbolic `x`, `rs.s`, or the symbol `:x`). Defaults to the system's observables. |
12 | 12 | - `known_p = Num[]`: List of parameters which values are known.
|
13 |
| -- `ignore_no_measured_warn=false`: If set to `true`, no warning is provided when the `measured_quantities` vector is empty. |
| 13 | +- `ignore_no_measured_warn = false`: If set to `true`, no warning is provided when the `measured_quantities` vector is empty. |
| 14 | +- `remove_conserved = true`: Whether to eliminate conservation laws when computing the ode (this can reduce runtime of identifiability analysis significantly). |
14 | 15 |
|
15 | 16 | Example:
|
16 | 17 | ```julia
|
| 18 | +using Catalyst, StructuralIdentifiability |
17 | 19 | rs = @reaction_network begin
|
18 | 20 | (p,d), 0 <--> X
|
19 | 21 | end
|
20 |
| -si_ode(rs; measured_quantities = [:X], known_p = [:p]) |
| 22 | +make_si_ode(rs; measured_quantities = [:X], known_p = [:p]) |
| 23 | +``` |
21 | 24 |
|
22 | 25 | Notes:
|
23 |
| -This function is part of the StructuralIdentifiability.jl extension. StructuralIdentifiability.jl must be imported to access it. |
| 26 | +- This function is part of the StructuralIdentifiability.jl extension. StructuralIdentifiability.jl must be imported to access it. |
| 27 | +- `measured_quantities` and `known_p` input may also be symbolic (e.g. measured_quantities = [rs.X]) |
24 | 28 | ```
|
25 | 29 | """
|
26 | 30 | function Catalyst.make_si_ode(rs::ReactionSystem; measured_quantities = [], known_p = [],
|
27 |
| - ignore_no_measured_warn=false, remove_conserved = true) |
| 31 | + ignore_no_measured_warn = false, remove_conserved = true) |
28 | 32 | # Creates a MTK ODESystem, and a list of measured quantities (there are equations).
|
29 | 33 | # Gives these to SI to create an SI ode model of its preferred form.
|
30 | 34 | osys, conseqs, _ = make_osys(rs; remove_conserved)
|
@@ -111,12 +115,12 @@ function make_measured_quantities(rs::ReactionSystem, measured_quantities::Vecto
|
111 | 115 | # Appends the known parameters to the measured_quantities vector. Converts any Symbols to symbolics.
|
112 | 116 | mqiterator = Iterators.flatten((measured_quantities, known_p))
|
113 | 117 | mqs = [(q isa Symbol) ? Catalyst._symbol_to_var(rs, q) : q for q in mqiterator]
|
114 |
| - mqs = vector_subs(measured_quantities, conseqs) |
| 118 | + mqs = vector_subs(mqs, conseqs) |
115 | 119 |
|
116 | 120 | # Creates one internal observation variable for each measured quantity (`___internal_observables`).
|
117 | 121 | # Creates a vector of equations, setting each measured quantity equal to one observation variable.
|
118 |
| - @variables t (___internal_observables(t))[1:length(measured_quantities)] |
119 |
| - return Equation[(q isa Equation) ? q : (___internal_observables[i] ~ q) for (i,q) in enumerate(measured_quantities)] |
| 122 | + @variables t (___internal_observables(Catalyst.get_iv(rs)))[1:length(mqs)] |
| 123 | + return Equation[(q isa Equation) ? q : (___internal_observables[i] ~ q) for (i,q) in enumerate(mqs)] |
120 | 124 | end
|
121 | 125 |
|
122 | 126 | # Creates the functions that we wish to check for identifiability.
|
|
135 | 139 | function make_output(out, funcs_to_check, conseqs)
|
136 | 140 | funcs_to_check = vector_subs(funcs_to_check, conseqs)
|
137 | 141 | out = Dict(zip(vector_subs(keys(out), conseqs), values(out)))
|
138 |
| - out = sort(out; by = x -> findfirst(isequal(x, ftc) for ftc in funcs_to_check)) |
139 |
| - return out |
| 142 | + sortdict = Dict(ftc => i for (i,ftc) in enumerate(funcs_to_check)) |
| 143 | + return sort(out; by = x -> sortdict[x]) |
140 | 144 | end
|
141 | 145 |
|
142 | 146 | # For a vector of expressions and a conservation law, substitutes the law into every equation.
|
|
0 commit comments