Skip to content

Commit

Permalink
Merge pull request #229 from SciML/no-parameters-no-worry
Browse files Browse the repository at this point in the history
Fix for ODEs with no parameters
  • Loading branch information
sumiya11 authored Nov 6, 2023
2 parents 790a655 + a26a525 commit d1c2260
Show file tree
Hide file tree
Showing 7 changed files with 73 additions and 7 deletions.
8 changes: 4 additions & 4 deletions Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -29,23 +29,23 @@ AbstractAlgebra = "0.13, 0.18, 0.19, 0.20, 0.21, 0.22, 0.23, 0.24, 0.25, 0.26, 0
BenchmarkTools = "1"
Combinatorics = "1"
DataStructures = "0.18"
Dates = "1.6, 1.7"
Groebner = "0.4, 0.5"
IterTools = "1"
LinearAlgebra = "1.6, 1.7"
Logging = "1.6, 1.7"
MacroTools = "0.5"
ModelingToolkit = "7, 8"
Nemo = "0.24, 0.25, 0.26, 0.27, 0.28, 0.29, 0.30, 0.31, 0.32, 0.33, 0.34, 0.35, 0.36, 0.37"
ParamPunPam = "0.2"
PrecompileTools = "1"
Primes = "0.5"
Random = "1.6, 1.7"
SpecialFunctions = "1, 2"
SymbolicUtils = "1"
Symbolics = "5"
TestSetExtensions = "2"
julia = "1.6, 1.7"
Dates = "1.6, 1.7"
Random = "1.6, 1.7"
Logging = "1.6, 1.7"
LinearAlgebra = "1.6, 1.7"

[extras]
CPUSummary = "2a0fbf3d-bb9c-48f3-b0a9-814d99fd7ab9"
Expand Down
2 changes: 1 addition & 1 deletion docs/src/tutorials/identifiable_functions.md
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ find_identifiable_functions(LLW1987, with_states = true)
```

By default, `find_identifiable_functions` tries to simplify the output functions as much as possible, and it has `simplify` keyword responsible for
the degree of simplification. The default value is `:standard` but one could use `:string` to try to simplify further
the degree of simplification. The default value is `:standard` but one could use `:strong` to try to simplify further
(at the expense of heavier computation) or use `:weak` to simplify less (but compute faster).

[^1]: > Y. Lecourtier, F. Lamnabhi-Lagarrigue, and E. Walter, [*A method to prove that nonlinear models can be unidentifiable*](https://doi.org/10.1109/CDC.1987.272467), Proceedings of the 26th Conference on Decision and Control, December 1987, 2144-2145;
2 changes: 1 addition & 1 deletion src/StructuralIdentifiability.jl
Original file line number Diff line number Diff line change
Expand Up @@ -97,7 +97,7 @@ The function returns a dictionary from the functions to check to their identifia
"""
function assess_identifiability(
ode::ODE{P};
funcs_to_check::Array{<:RingElem, 1} = Array{RingElem, 1}(),
funcs_to_check = Vector(),
p::Float64 = 0.99,
) where {P <: MPolyElem{fmpq}}
p_glob = 1 - (1 - p) * 0.9
Expand Down
34 changes: 33 additions & 1 deletion src/global_identifiability.jl
Original file line number Diff line number Diff line change
@@ -1,3 +1,17 @@
"""
extract_identifiable_functions_raw(io_equations, ode, known, with_states)
Takes as input input-output equations, the corresponding ode, a list of functions assumed to be known
and a flag `with_states`.
Extracts generators of the field of identifiable functions (with or without states) withous
any simplifications.
Returns a tuple consisting of
- a dictionary with at most two keys: `:no_states` and `:with_states`. The respective values
are identifiable functions containing or not the state variables
- a polynomial ring containing all returned identifiable functions
(parameters or parameters + states)
"""
function extract_identifiable_functions_raw(
io_equations::Dict{P, P},
ode::ODE{P},
Expand Down Expand Up @@ -79,6 +93,10 @@ This function takes the following optional arguments:
- `with_states`: Also report the identifiabile functions in states. Default is
`false`. If this is `true`, the identifiable functions involving parameters only
will be simplified
The function returns a tuple containing the following:
- a list of identifiable functions (as pairs [num, denum])
- the ring containing all these functuons (either parameters only of with states)
"""
function initial_identifiable_functions(
ode::ODE{T};
Expand Down Expand Up @@ -113,7 +131,7 @@ function initial_identifiable_functions(
_runtime_logger[:rank_time] = rank_times

if any([dim != rk + 1 for (dim, rk) in zip(dims, wranks)])
@warn "One of the Wronskians has corank greater than one, so the results of the algorithm will be valid only for multiexperiment identifiability. If you still would like to assess single-experiment identifiability, we recommend using SIAN (https://github.com/alexeyovchinnikov/SIAN-Julia)"
@warn "One of the Wronskians has corank greater than one, so the results of the algorithm will be valid only for multiexperiment identifiability. If you still would like to assess single-experiment identifiability, we recommend using SIAN (https://github.com/alexeyovchinnikov/SIAN-Julia)"
end
end

Expand Down Expand Up @@ -159,7 +177,18 @@ function initial_identifiable_functions(
end

# ------------------------------------------------------------------------------
"""
check_identifiability(ode, funcs_to_check; known, p, var_change_policy)
Input:
- `ode` - the ODE model
- `funcs_to_check` - the functions to check identifiability for
- `known` - a list of functions in states which are assumed to be known and generic
- `p` - probability of correctness
- `var_change` - a policy for variable change (`:default`, `:yes`, `:no`), affects only the runtime
Output: a list L of booleans with L[i] being the identifiability status of the i-th function to check
"""
function check_identifiability(
ode::ODE{P},
funcs_to_check::Array{<:Any, 1};
Expand All @@ -176,6 +205,9 @@ function check_identifiability(
break
end
end
if !states_needed && isempty(ode.parameters)
return [true for _ in funcs_to_check]
end

half_p = 0.5 + p / 2
identifiable_functions_raw, bring = initial_identifiable_functions(
Expand Down
10 changes: 10 additions & 0 deletions src/identifiable_functions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -56,6 +56,16 @@ function find_identifiable_functions(
_runtime_logger[:id_npoints_interpolation] = 0
_runtime_logger[:id_beautifulization] = 0.0
runtime_start = time_ns()
if isempty(ode.parameters) && !with_states
@warn """
There are no parameters in the given ODE, thus no identifiabile
functions.
Use `find_identifiable_functions` with keyword `with_states=true` to
compute the functions with the ODE states included."""
bring = parent(ode)
id_funcs = [one(bring)]
return id_funcs
end
half_p = 0.5 + p / 2
id_funcs, bring = initial_identifiable_functions(
ode,
Expand Down
20 changes: 20 additions & 0 deletions test/identifiability.jl
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,26 @@
),
)

#--------------------------------------------------------------------------
# No parameters no worry

ode = @ODEmodel(x1'(t) = x1, x2'(t) = x2, y(t) = x1 + x2(t))
funcs_to_test = [x1, x2, x1 + x2]
correct = [:nonidentifiable, :nonidentifiable, :globally]
push!(
test_cases,
Dict(
:ode => ode,
:funcs => funcs_to_test,
:correct => Dict(funcs_to_test .=> correct),
),
)

# Also test when `funcs_to_test` is empty!
funcs_to_test = Vector{typeof(x1)}()
correct = Dict(x1 => :nonidentifiable, x2 => :nonidentifiable)
push!(test_cases, Dict(:ode => ode, :funcs => funcs_to_test, :correct => correct))

#--------------------------------------------------------------------------

ode = @ODEmodel(
Expand Down
4 changes: 4 additions & 0 deletions test/identifiable_functions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -935,6 +935,10 @@ ode = StructuralIdentifiability.@ODEmodel(
ident_funcs = [k3, k1 // k7, k5 // k2, k6 // k2, k7 * EpoR_A]
push!(test_cases, (ode = ode, ident_funcs = ident_funcs))

ode = @ODEmodel(x1'(t) = x1, x2'(t) = x2, y(t) = x1 + x2(t))
ident_funcs = [x1 + x2]
push!(test_cases, (ode = ode, ident_funcs = ident_funcs, with_states = true))

# TODO: verify that Maple returns the same
@testset "Identifiable functions of parameters" begin
p = 0.99
Expand Down

0 comments on commit d1c2260

Please sign in to comment.