From ab437dde78583f9d26e84340b605b56352528924 Mon Sep 17 00:00:00 2001 From: "Documenter.jl" Date: Thu, 25 Jan 2024 14:33:27 +0000 Subject: [PATCH] build based on 190a585 --- previews/PR266/.documenter-siteinfo.json | 2 +- previews/PR266/assets/Manifest.toml | 30 ++++++------ previews/PR266/assets/Project.toml | 2 +- previews/PR266/export/export/index.html | 2 +- .../identifiability/index.html | 4 +- previews/PR266/index.html | 20 ++++---- previews/PR266/input/input/index.html | 12 +++-- .../PR266/ioequations/ioequations/index.html | 2 +- previews/PR266/search_index.js | 2 +- .../PR266/tutorials/creating_ode/index.html | 6 +-- .../PR266/tutorials/discrete_time/index.html | 49 +++++++++++++++---- .../tutorials/identifiability/index.html | 4 +- .../identifiable_functions/index.html | 4 +- .../tutorials/reparametrization/index.html | 4 +- previews/PR266/utils/elimination/index.html | 2 +- .../utils/global_identifiability/index.html | 4 +- .../utils/local_identifiability/index.html | 2 +- previews/PR266/utils/ode/index.html | 4 +- .../PR266/utils/power_series_utils/index.html | 2 +- previews/PR266/utils/primality/index.html | 2 +- .../PR266/utils/reparametrization/index.html | 4 +- previews/PR266/utils/util/index.html | 4 +- previews/PR266/utils/wronskian/index.html | 2 +- 23 files changed, 103 insertions(+), 66 deletions(-) diff --git a/previews/PR266/.documenter-siteinfo.json b/previews/PR266/.documenter-siteinfo.json index 6ded9cad..155f2c4b 100644 --- a/previews/PR266/.documenter-siteinfo.json +++ b/previews/PR266/.documenter-siteinfo.json @@ -1 +1 @@ -{"documenter":{"julia_version":"1.10.0","generation_timestamp":"2024-01-21T19:37:22","documenter_version":"1.2.1"}} \ No newline at end of file +{"documenter":{"julia_version":"1.10.0","generation_timestamp":"2024-01-25T14:33:19","documenter_version":"1.2.1"}} \ No newline at end of file diff --git a/previews/PR266/assets/Manifest.toml b/previews/PR266/assets/Manifest.toml index f5f68c45..13c3cc63 100644 --- a/previews/PR266/assets/Manifest.toml +++ b/previews/PR266/assets/Manifest.toml @@ -2,7 +2,7 @@ julia_version = "1.10.0" manifest_format = "2.0" -project_hash = "88a670dcfb901d5c3347fe8f16ba784750a1af5b" +project_hash = "79a87c8b6351eb6b4461d22a8b16148b365aecfc" [[deps.ADTypes]] git-tree-sha1 = "41c37aa88889c171f1300ceac1313c06e891d245" @@ -144,9 +144,9 @@ version = "0.5.1" [[deps.ChainRulesCore]] deps = ["Compat", "LinearAlgebra"] -git-tree-sha1 = "c1deebd76f7a443d527fc0430d5758b8b2112ed8" +git-tree-sha1 = "0d12ee16b3f62e4e33c3277773730a5b21a74152" uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" -version = "1.19.1" +version = "1.20.0" weakdeps = ["SparseArrays"] [deps.ChainRulesCore.extensions] @@ -228,9 +228,9 @@ uuid = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f" version = "4.1.1" [[deps.DataAPI]] -git-tree-sha1 = "8da84edb865b0b5b0100c0666a9bc9a0b71c553c" +git-tree-sha1 = "abe83f3a2f1b857aac70ef8b269080af17764bbe" uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a" -version = "1.15.0" +version = "1.16.0" [[deps.DataStructures]] deps = ["Compat", "InteractiveUtils", "OrderedCollections"] @@ -943,10 +943,10 @@ uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" version = "1.2.0" [[deps.NonlinearSolve]] -deps = ["ADTypes", "ArrayInterface", "ConcreteStructs", "DiffEqBase", "EnumX", "FastBroadcast", "FastClosures", "FiniteDiff", "ForwardDiff", "LazyArrays", "LineSearches", "LinearAlgebra", "LinearSolve", "MaybeInplace", "PrecompileTools", "Printf", "RecursiveArrayTools", "Reexport", "SciMLBase", "SciMLOperators", "SimpleNonlinearSolve", "SparseArrays", "SparseDiffTools", "StaticArrays", "UnPack"] -git-tree-sha1 = "72b036b728461272ae1b1c3f7096cb4c319d8793" +deps = ["ADTypes", "ArrayInterface", "ConcreteStructs", "DiffEqBase", "FastBroadcast", "FastClosures", "FiniteDiff", "ForwardDiff", "LazyArrays", "LineSearches", "LinearAlgebra", "LinearSolve", "MaybeInplace", "PrecompileTools", "Preferences", "Printf", "RecursiveArrayTools", "Reexport", "SciMLBase", "SimpleNonlinearSolve", "SparseArrays", "SparseDiffTools", "StaticArraysCore", "TimerOutputs"] +git-tree-sha1 = "78bdd3a4a62865cf43c53d63783b0cbfddcdbbe6" uuid = "8913a72c-1f9b-4ce2-8d82-65094dcecaec" -version = "3.4.0" +version = "3.5.0" [deps.NonlinearSolve.extensions] NonlinearSolveBandedMatricesExt = "BandedMatrices" @@ -1162,14 +1162,15 @@ version = "1.3.4" [[deps.RecursiveArrayTools]] deps = ["Adapt", "ArrayInterface", "DocStringExtensions", "GPUArraysCore", "IteratorInterfaceExtensions", "LinearAlgebra", "RecipesBase", "SparseArrays", "StaticArraysCore", "Statistics", "SymbolicIndexingInterface", "Tables"] -git-tree-sha1 = "dd7fc1923fde0cc6cdff451352d17924b0704ca1" +git-tree-sha1 = "16f1bb9de02b8bce31a7b2495345532901214cae" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" -version = "3.5.4" +version = "3.6.2" [deps.RecursiveArrayTools.extensions] RecursiveArrayToolsFastBroadcastExt = "FastBroadcast" RecursiveArrayToolsMeasurementsExt = "Measurements" RecursiveArrayToolsMonteCarloMeasurementsExt = "MonteCarloMeasurements" + RecursiveArrayToolsReverseDiffExt = ["ReverseDiff", "Zygote"] RecursiveArrayToolsTrackerExt = "Tracker" RecursiveArrayToolsZygoteExt = "Zygote" @@ -1177,6 +1178,7 @@ version = "3.5.4" FastBroadcast = "7034ab61-46d4-4ed7-9d0f-46aef9175898" Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" @@ -1238,9 +1240,9 @@ version = "0.6.42" [[deps.SciMLBase]] deps = ["ADTypes", "ArrayInterface", "CommonSolve", "ConstructionBase", "Distributed", "DocStringExtensions", "EnumX", "FillArrays", "FunctionWrappersWrappers", "IteratorInterfaceExtensions", "LinearAlgebra", "Logging", "Markdown", "PrecompileTools", "Preferences", "Printf", "RecipesBase", "RecursiveArrayTools", "Reexport", "RuntimeGeneratedFunctions", "SciMLOperators", "StaticArraysCore", "Statistics", "SymbolicIndexingInterface", "Tables", "TruncatedStacktraces"] -git-tree-sha1 = "de41474ac529bf81598e064587421cc5ebc28fa0" +git-tree-sha1 = "d91985cfda7d730a885d7dbc89889e8184b72802" uuid = "0bca4576-84f4-4d90-8ffe-ffa030f20462" -version = "2.20.0" +version = "2.21.0" [deps.SciMLBase.extensions] SciMLBaseChainRulesCoreExt = "ChainRulesCore" @@ -1447,9 +1449,9 @@ uuid = "bea87d4a-7f5b-5778-9afe-8cc45184846c" version = "7.2.1+1" [[deps.SymbolicIndexingInterface]] -git-tree-sha1 = "74502f408d99fc217a9d7cd901d9ffe45af892b1" +git-tree-sha1 = "34573fc920adfd457c5be704098d0168e4f20e54" uuid = "2efcf032-c050-4f8e-a9bb-153293bab1f5" -version = "0.3.3" +version = "0.3.4" [[deps.SymbolicUtils]] deps = ["AbstractTrees", "Bijections", "ChainRulesCore", "Combinatorics", "ConstructionBase", "DataStructures", "DocStringExtensions", "DynamicPolynomials", "IfElse", "LabelledArrays", "LinearAlgebra", "MultivariatePolynomials", "NaNMath", "Setfield", "SparseArrays", "SpecialFunctions", "StaticArrays", "SymbolicIndexingInterface", "TimerOutputs", "Unityper"] diff --git a/previews/PR266/assets/Project.toml b/previews/PR266/assets/Project.toml index 21e7b181..e449e2a5 100644 --- a/previews/PR266/assets/Project.toml +++ b/previews/PR266/assets/Project.toml @@ -7,5 +7,5 @@ StructuralIdentifiability = "220ca800-aa68-49bb-acd8-6037fa93a544" [compat] BenchmarkTools = "1.3" Documenter = "0.27, 1" -ModelingToolkit = "8.34" +ModelingToolkit = "8.74" StructuralIdentifiability = "0.5" diff --git a/previews/PR266/export/export/index.html b/previews/PR266/export/export/index.html index 9b79df17..1d8b727f 100644 --- a/previews/PR266/export/export/index.html +++ b/previews/PR266/export/export/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Exporting to Other Systems

Here we put some helpful utilities to export your code to other identifiability software.

StructuralIdentifiability.print_for_mapleFunction
print_for_maple(ode, package)

Prints the ODE in the format accepted by maple packages

  • SIAN (https://github.com/pogudingleb/SIAN) if package=:SIAN
  • DifferentialAlgebra if package=:DifferentialAlgebra
  • DifferentialThomas if package=:DifferentialThomas
source
+

Exporting to Other Systems

Here we put some helpful utilities to export your code to other identifiability software.

StructuralIdentifiability.print_for_mapleFunction
print_for_maple(ode, package)

Prints the ODE in the format accepted by maple packages

  • SIAN (https://github.com/pogudingleb/SIAN) if package=:SIAN
  • DifferentialAlgebra if package=:DifferentialAlgebra
  • DifferentialThomas if package=:DifferentialThomas
source
diff --git a/previews/PR266/identifiability/identifiability/index.html b/previews/PR266/identifiability/identifiability/index.html index 485d732f..a5a1d6b1 100644 --- a/previews/PR266/identifiability/identifiability/index.html +++ b/previews/PR266/identifiability/identifiability/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Functions to Assess Identifiability

Assessing All Types of Identifiability

StructuralIdentifiability.assess_identifiabilityFunction
assess_identifiability(ode; funcs_to_check = [], prob_threshold=0.99, loglevel=Logging.Info)

Input:

  • ode - the ODE model
  • funcs_to_check - list of functions to check identifiability for; if empty, all parameters and states are taken
  • prob_threshold - probability of correctness.
  • loglevel - the minimal level of log messages to display (Logging.Info by default)

Assesses identifiability of a given ODE model. The result is guaranteed to be correct with the probability at least prob_threshold. The function returns an (ordered) dictionary from the functions to check to their identifiability properties (one of :nonidentifiable, :locally, :globally).

source

Assessing Local Identifiability

StructuralIdentifiability.assess_local_identifiabilityFunction
assess_local_identifiability(ode::ODE{P}; funcs_to_check::Array{<: Any, 1}, prob_threshold::Float64=0.99, type=:SE, loglevel=Logging.Info) where P <: MPolyRingElem{Nemo.QQFieldElem}

Checks the local identifiability/observability of the functions in funcs_to_check. The result is correct with probability at least prob_threshold.

Call this function if you have a specific collection of parameters of which you would like to check local identifiability.

type can be either :SE (single-experiment identifiability) or :ME (multi-experiment identifiability). If the type is :ME, states are not allowed to appear in the funcs_to_check.

source

Finding Identifiable Functions

StructuralIdentifiability.find_identifiable_functionsFunction
find_identifiable_functions(ode::ODE; options...)

Finds all functions of parameters/states that are identifiable in the given ODE system.

Options

This functions takes the following optional arguments:

  • with_states: When true, also reports the identifiabile functions in the ODE states. Default is false.
  • simplify: The extent to which the output functions are simplified. Stronger simplification may require more time. Possible options are:
    • :standard: Default simplification.
    • :weak: Weak simplification. This option is the fastest, but the output functions can be quite complex.
    • :strong: Strong simplification. This option is the slowest, but the output
    functions are nice and simple.
    • :absent: No simplification.
  • prob_threshold: A float in the range from 0 to 1, the probability of correctness. Default is 0.99.
  • seed: The rng seed. Default value is 42.
  • loglevel - the minimal level of log messages to display (Logging.Info by default)

Example

using StructuralIdentifiability
+

Functions to Assess Identifiability

Assessing All Types of Identifiability

StructuralIdentifiability.assess_identifiabilityFunction
assess_identifiability(ode; funcs_to_check = [], prob_threshold=0.99, loglevel=Logging.Info)

Input:

  • ode - the ODE model
  • funcs_to_check - list of functions to check identifiability for; if empty, all parameters and states are taken
  • prob_threshold - probability of correctness.
  • loglevel - the minimal level of log messages to display (Logging.Info by default)

Assesses identifiability of a given ODE model. The result is guaranteed to be correct with the probability at least prob_threshold. The function returns an (ordered) dictionary from the functions to check to their identifiability properties (one of :nonidentifiable, :locally, :globally).

source

Assessing Local Identifiability

StructuralIdentifiability.assess_local_identifiabilityFunction
assess_local_identifiability(ode::ODE{P}; funcs_to_check::Array{<: Any, 1}, prob_threshold::Float64=0.99, type=:SE, loglevel=Logging.Info) where P <: MPolyRingElem{Nemo.QQFieldElem}

Checks the local identifiability/observability of the functions in funcs_to_check. The result is correct with probability at least prob_threshold.

Call this function if you have a specific collection of parameters of which you would like to check local identifiability.

type can be either :SE (single-experiment identifiability) or :ME (multi-experiment identifiability). If the type is :ME, states are not allowed to appear in the funcs_to_check.

source
assess_local_identifiability(dds::DDS{P}; funcs_to_check::Array{<: Any, 1}, known_ic, prob_threshold::Float64=0.99, loglevel=Logging.Info) where P <: MPolyRingElem{Nemo.QQFieldElem}

Checks the local identifiability/observability of the functions in funcs_to_check. The result is correct with probability at least prob_threshold. A list of quantities can be provided as known_ic for which the initial conditions can be assumed to be known and generic.

source

Finding Identifiable Functions

StructuralIdentifiability.find_identifiable_functionsFunction
find_identifiable_functions(ode::ODE; options...)

Finds all functions of parameters/states that are identifiable in the given ODE system.

Options

This functions takes the following optional arguments:

  • with_states: When true, also reports the identifiabile functions in the ODE states. Default is false.
  • simplify: The extent to which the output functions are simplified. Stronger simplification may require more time. Possible options are:
    • :standard: Default simplification.
    • :weak: Weak simplification. This option is the fastest, but the output functions can be quite complex.
    • :strong: Strong simplification. This option is the slowest, but the output
    functions are nice and simple.
    • :absent: No simplification.
  • prob_threshold: A float in the range from 0 to 1, the probability of correctness. Default is 0.99.
  • seed: The rng seed. Default value is 42.
  • loglevel - the minimal level of log messages to display (Logging.Info by default)

Example

using StructuralIdentifiability
 
 ode = @ODEmodel(
     x0'(t) = -(a01 + a21) * x0(t) + a12 * x1(t),
@@ -16,4 +16,4 @@
 # prints
 3-element Vector{AbstractAlgebra.Generic.Frac{Nemo.QQMPolyRingElem}}:
  a12 + a01 + a21
- a12*a01
source
+ a12*a01
source
diff --git a/previews/PR266/index.html b/previews/PR266/index.html index a4ee19ec..ac99aa7d 100644 --- a/previews/PR266/index.html +++ b/previews/PR266/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

StructuralIdentifiability.jl

StructuralIdentifiability.jl is a comprehensive toolbox for assessing identifiability parameters.

Installation

To install StructuralIdentifiability.jl, use the Julia package manager:

using Pkg
+

StructuralIdentifiability.jl

StructuralIdentifiability.jl is a comprehensive toolbox for assessing identifiability parameters.

Installation

To install StructuralIdentifiability.jl, use the Julia package manager:

using Pkg
 Pkg.add("StructuralIdentifiability")

Citation

@article{structidjl,
   author  = {Dong, R. and Goodbrake, C. and Harrington, H. and Pogudin G.},
   title   = {Differential Elimination for Dynamical Models via Projections with Applications to Structural Identifiability},
@@ -30,7 +30,7 @@
   Threads: 1 on 4 virtual cores
A more complete overview of all dependencies and their versions is also provided.
Status `~/work/StructuralIdentifiability.jl/StructuralIdentifiability.jl/docs/Manifest.toml`
   [47edcb42] ADTypes v0.2.6
   [a4c015fc] ANSIColoredPrinters v0.0.1
-  [c3fe647b] AbstractAlgebra v0.35.3
+⌅ [c3fe647b] AbstractAlgebra v0.35.3
   [1520ce14] AbstractTrees v0.4.4
 ⌅ [79e6a3ab] Adapt v3.7.2
   [ec485272] ArnoldiMethod v0.2.0
@@ -43,7 +43,7 @@
   [2a0fbf3d] CPUSummary v0.2.4
   [00ebfdb7] CSTParser v3.4.0
   [49dc2e85] Calculus v0.5.1
-  [d360d2e6] ChainRulesCore v1.19.1
+  [d360d2e6] ChainRulesCore v1.20.0
   [fb6a15b2] CloseOpenIntervals v0.1.12
   [861a8166] Combinatorics v1.0.2
   [a80b9123] CommonMark v0.8.12
@@ -55,7 +55,7 @@
   [187b0558] ConstructionBase v1.5.4
   [adafc99b] CpuId v0.3.1
   [a8cc5b0e] Crayons v4.1.1
-  [9a962f9c] DataAPI v1.15.0
+  [9a962f9c] DataAPI v1.16.0
   [864edb3b] DataStructures v0.18.16
   [e2d170a0] DataValueInterfaces v1.0.0
   [2b5f629d] DiffEqBase v6.146.0
@@ -132,8 +132,8 @@
   [d41bc354] NLSolversBase v7.8.3
   [2774e3e8] NLsolve v4.5.1
   [77ba4419] NaNMath v1.0.2
-  [2edaba10] Nemo v0.39.1
-  [8913a72c] NonlinearSolve v3.4.0
+⌅ [2edaba10] Nemo v0.39.1
+  [8913a72c] NonlinearSolve v3.5.0
   [6fe1bfb0] OffsetArrays v1.13.0
   [bac558e1] OrderedCollections v1.6.3
   [1dea7af3] OrdinaryDiffEq v6.69.0
@@ -155,7 +155,7 @@
   [fb686558] RandomExtensions v0.4.4
   [e6cf234a] RandomNumbers v1.5.3
   [3cdcf5f2] RecipesBase v1.3.4
-  [731186ca] RecursiveArrayTools v3.5.4
+  [731186ca] RecursiveArrayTools v3.6.2
   [f2c3362d] RecursiveFactorization v0.2.21
   [189a3867] Reexport v1.2.2
   [2792f1a3] RegistryInstances v0.1.0
@@ -164,7 +164,7 @@
   [7e49a35a] RuntimeGeneratedFunctions v0.5.12
   [94e857df] SIMDTypes v0.1.0
   [476501e8] SLEEFPirates v0.6.42
-  [0bca4576] SciMLBase v2.20.0
+  [0bca4576] SciMLBase v2.21.0
   [c0aeaf25] SciMLOperators v0.3.7
   [efcf1570] Setfield v1.1.1
   [727e6d20] SimpleNonlinearSolve v1.3.1
@@ -184,7 +184,7 @@
   [7792a7ef] StrideArraysCore v0.5.2
   [892a3eda] StringManipulation v0.3.4
   [220ca800] StructuralIdentifiability v0.5.1 `~/work/StructuralIdentifiability.jl/StructuralIdentifiability.jl`
-  [2efcf032] SymbolicIndexingInterface v0.3.3
+  [2efcf032] SymbolicIndexingInterface v0.3.4
   [d1185830] SymbolicUtils v1.5.0
   [0c5d862f] Symbolics v5.16.1
   [3783bdb8] TableTraits v1.0.1
@@ -268,4 +268,4 @@
   [3f19e933] p7zip_jll v17.4.0+2
 Info Packages marked with ⌃ and ⌅ have new versions available. Those with ⌃ may be upgradable, but those with ⌅ are restricted by compatibility constraints from upgrading. To see why use `status --outdated -m`
You can also download the manifest file and the -project file.
+project file.
diff --git a/previews/PR266/input/input/index.html b/previews/PR266/input/input/index.html index f41986e7..0ad9d235 100644 --- a/previews/PR266/input/input/index.html +++ b/previews/PR266/input/input/index.html @@ -1,12 +1,18 @@ -Parsing input ODE system · StructuralIdentifiability.jl

Parsing input ODE system

StructuralIdentifiability.@ODEmodelMacro
macro ODEmodel

Macro for creating an ODE from a list of equations. It also injects all variables into the global scope.

Example

Creating a simple ODE:

using StructuralIdentifiability
+

Parsing input system

StructuralIdentifiability.@ODEmodelMacro
macro ODEmodel

Macro for creating an ODE from a list of equations. It also injects all variables into the global scope.

Example

Creating a simple ODE:

using StructuralIdentifiability
 
 ode = @ODEmodel(
     x1'(t) = a * x1(t) + u(t),
     x2'(t) = b * x2(t) + c*x1(t)*x2(t),
     y(t) = x1(t)
-)

Here,

  • x1, x2 are state variables
  • y is an output variable
  • u is an input variable
  • a, b, c are time-indepdendent parameters
source
StructuralIdentifiability.ODEType

The main structure that represents input ODE system.

Stores information about states (x_vars), outputs (y_vars), inputs (u_vars), parameters (parameters) and the equations.

This structure is constructed via @ODEmodel macro.

source
StructuralIdentifiability.set_parameter_valuesFunction
set_parameter_values(ode, param_values)

Input:

  • ode - an ODE as above
  • param_values - values for (possibly, some of) the parameters as dictionary parameter => value

Output:

  • new ode with the parameters in param_values plugged with the given numbers
source

Create Compartmental Model

StructuralIdentifiability.linear_compartment_modelFunction
linear_compartment_model(graph, inputs, outputs, leaks)

Input: defines a linear compartment model with nodes numbered from 1 to n by

  • graph - and array of integer arrays representing the adjacency lists of the graph
  • inputs - array of input nodes
  • outputs - array of output nodes
  • leaks - array of sink nodes

Output:

  • the corresponding ODE system in the notation of https://doi.org/10.1007/s11538-015-0098-0
source
+)

Here,

  • x1, x2 are state variables
  • y is an output variable
  • u is an input variable
  • a, b, c are time-indepdendent parameters
source
StructuralIdentifiability.ODEType

The main structure that represents input ODE system.

Stores information about states (x_vars), outputs (y_vars), inputs (u_vars), parameters (parameters) and the equations.

This structure is constructed via @ODEmodel macro.

source
StructuralIdentifiability.set_parameter_valuesFunction
set_parameter_values(ode, param_values)

Input:

  • ode - an ODE as above
  • param_values - values for (possibly, some of) the parameters as dictionary parameter => value

Output:

  • new ode with the parameters in param_values plugged with the given numbers
source

Create Compartmental Model

StructuralIdentifiability.linear_compartment_modelFunction
linear_compartment_model(graph, inputs, outputs, leaks)

Input: defines a linear compartment model with nodes numbered from 1 to n by

  • graph - and array of integer arrays representing the adjacency lists of the graph
  • inputs - array of input nodes
  • outputs - array of output nodes
  • leaks - array of sink nodes

Output:

  • the corresponding ODE system in the notation of https://doi.org/10.1007/s11538-015-0098-0
source

Discrete-time systems

StructuralIdentifiability.@DDSmodelMacro
macro DDSmodel

Macro for creating a DDS (discrete dynamical system) from a list of equations. It also injects all variables into the global scope.

Example

Creating a simple DDS:

using StructuralIdentifiability
+
+ode = @ODEmodel(
+    x1(t + 1) = a * x1(t) + u(t),
+    x2(t + 1) = b * x2(t) + c*x1(t)*x2(t),
+    y(t) = x1(t)
+)

Here,

  • x1, x2 are state variables
  • y is an output variable
  • u is an input variable
  • a, b, c are time-indepdendent parameters
source
diff --git a/previews/PR266/ioequations/ioequations/index.html b/previews/PR266/ioequations/ioequations/index.html index 277d6edc..a432895c 100644 --- a/previews/PR266/ioequations/ioequations/index.html +++ b/previews/PR266/ioequations/ioequations/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Finding Input-Output Equations

StructuralIdentifiability.find_ioequationsFunction
find_ioequations(ode, [var_change_policy=:default])

Finds the input-output equations of an ODE system Input:

  • ode - the ODE system
  • var_change_policy - whether to perform automatic variable change, can be one of :default, :yes, :no
  • loglevel - logging level (default: Logging.Info)

Output:

  • a dictionary from “leaders” to the corresponding input-output equations; if an extra projection is needed, it will be the value corresponding to rand_proj_var
source

Reducing with respect to Input-Output Equations

StructuralIdentifiability.PBRepresentationType

The structure for storing a projection-based representation of differential ideal (see Section 2.3 https://arxiv.org/abs/2111.00991). Contains the following fields:

  • y_names - the names of the variables with finite order in the profile (typically, outputs)
  • u_names - the names of the variables with infinite order in the profile (typically, inputs)
  • param_names - the names of the parameters
  • profile - the profile of the PB-representation (see Definition 2.13) as a dict from y_names with finite orders to the orders
  • projections - the corresponding projections (see Definition 2.15) as a dict from y_names to the projections
source
StructuralIdentifiability.pseudodivisionFunction
pseudodivision(f, g, x)

Computes the result of pseudodivision of f by g as univariate polynomials in x Input:

  • f - the polynomial to be divided
  • g - the polynomial to divide by
  • x - the variable for the division

Output: the pseudoremainder of f divided by g w.r.t. x

source
StructuralIdentifiability.diffreduceFunction
diffreduce(diffpoly, pbr)

Computes the result of differential reduction of a differential polynomial diffpoly with respect to the charset defined by a PB-representation pbr Input:

  • diffpoly - a polynomial representing a differential polynomial to be reduced
  • pbr - a projection-based representation

Output: the result of differential reduction of diffpoly by pbr considered as a characteristic set (see Remark 2.20 in the paper)

source
+

Finding Input-Output Equations

StructuralIdentifiability.find_ioequationsFunction
find_ioequations(ode, [var_change_policy=:default])

Finds the input-output equations of an ODE system Input:

  • ode - the ODE system
  • var_change_policy - whether to perform automatic variable change, can be one of :default, :yes, :no
  • loglevel - logging level (default: Logging.Info)

Output:

  • a dictionary from “leaders” to the corresponding input-output equations; if an extra projection is needed, it will be the value corresponding to rand_proj_var
source

Reducing with respect to Input-Output Equations

StructuralIdentifiability.PBRepresentationType

The structure for storing a projection-based representation of differential ideal (see Section 2.3 https://arxiv.org/abs/2111.00991). Contains the following fields:

  • y_names - the names of the variables with finite order in the profile (typically, outputs)
  • u_names - the names of the variables with infinite order in the profile (typically, inputs)
  • param_names - the names of the parameters
  • profile - the profile of the PB-representation (see Definition 2.13) as a dict from y_names with finite orders to the orders
  • projections - the corresponding projections (see Definition 2.15) as a dict from y_names to the projections
source
StructuralIdentifiability.pseudodivisionFunction
pseudodivision(f, g, x)

Computes the result of pseudodivision of f by g as univariate polynomials in x Input:

  • f - the polynomial to be divided
  • g - the polynomial to divide by
  • x - the variable for the division

Output: the pseudoremainder of f divided by g w.r.t. x

source
StructuralIdentifiability.diffreduceFunction
diffreduce(diffpoly, pbr)

Computes the result of differential reduction of a differential polynomial diffpoly with respect to the charset defined by a PB-representation pbr Input:

  • diffpoly - a polynomial representing a differential polynomial to be reduced
  • pbr - a projection-based representation

Output: the result of differential reduction of diffpoly by pbr considered as a characteristic set (see Remark 2.20 in the paper)

source
diff --git a/previews/PR266/search_index.js b/previews/PR266/search_index.js index c8f4bacb..a5fcddc6 100644 --- a/previews/PR266/search_index.js +++ b/previews/PR266/search_index.js @@ -1,3 +1,3 @@ var documenterSearchIndex = {"docs": -[{"location":"tutorials/reparametrization/#Reparametrizations","page":"Reparametrizations","title":"Reparametrizations","text":"","category":"section"},{"location":"tutorials/reparametrization/#Overview","page":"Reparametrizations","title":"Overview","text":"","category":"section"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Once one has found that not all parameters and/or states of the model at hand are identifiable, one natural desire is to reparametrize the model into a one with better identifiability properties. StructuralIdentifiability offers such a functionality via the function reparametrize_global. It takes as input an ODE model and produces its transformation into another model with the same input-output behaviour but with the states and parameters being globally identifiable. Note that, in general, such a transformation may not exist in the class of rational models, so sometimes the function returns an ODE not on the whole affine space but on a manifold.","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"More precisely, the function returns a dictionary with three keys:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":":new_vars is a dictionary which maps the new parameters and new states into the formulas expressing them in terms of the original parameters and states;\n:new_ode is the ODE satisfied by these new states (and the expression of the output in terms of the new states);\n:implicit_relations is a list of algebraic relations between the new states and parameters. Being nonempty, this is exactly the list of equations defining the manifold, on which the new ODE model is defined. In many interesting cases, however, this list is empty meaning that the new ODE is a standard rational ODE model.","category":"page"},{"location":"tutorials/reparametrization/#Example:-SEUIR-model","page":"Reparametrizations","title":"Example: SEUIR model","text":"","category":"section"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Consider a SEUIR epidemiological model from[1]:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"begincases\nS(t) = -beta frac(U(t) + I(t))S(t)N\nE(t) = beta frac(U(t) + I(t))S(t)N - gamma E(t)\nU(t) = (1 - alpha) gamma E(t) - delta U(t)\nI(t) = alpha gamma E(t) - delta I(t)\nR(t) = delta(U(t) + I(t))\ny(t) = I(t)\nendcases","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"In this model S is, as usually, the number of susceptible people, E is the number of people exposed to virus but not yet infected (as in a simple SEIR model[1]), and I and U correspond to number of infected people who report the infection and who do not, respectively. We define the model but omit R compartment since it does not affect the output dynamics:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n S'(t) = -b * (U(t) + I(t)) * S(t) / N,\n E'(t) = b * (U(t) + I(t)) * S(t) / N - g * E(t),\n U'(t) = (1 - a) * g * E(t) - d * U(t),\n I'(t) = a * g * E(t) - d * I(t),\n y(t) = I(t)\n)","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Majority of the states and parameters are not identifiable in this case:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"assess_identifiability(ode)","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Let us attempt to reparametrize the model, and print new variables:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"reparam = reparametrize_global(ode)\n@assert isempty(reparam[:implicit_relations]) # checking that the result is an ODE on the whole space, not on a manifold\nreparam[:new_vars]","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"In these new variables and parameters, the original ODE can be rewritten as follows:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"reparam[:new_ode]","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"In order to analyze this result, let us give more interpretable names to the new variables and parameters:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"I = I widetildeE = alpha E widetildeS = alpha widetildeI = alpha (I + U) gamma = gammadelta = deltawidetildebeta = fracbetaalpha N","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Then the reparametrize system becomes","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"begincases\nwidetildeS(t) = -widetildebeta widetildeS(t) widetildeI(t)\nwidetildeE(t) = widetildebeta widetildeS(t) widetildeI(t) - gamma widetildeE(t)\nwidetildeI(t) = -delta widetildeI(t) + gammawidetildeE(t)\nI(t) = gammawidetildeE(t) - delta I(t)\ny(t) = I(t)\nendcases","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"This reparametrization not only reduces the dimension of the parameter space from 5 to 3 but reveals interesting structural properties of the model:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"The first three equations form a self-contained model which is equivalent to a simple SEIR model, so the model gets \"decomposed\";\nNew variables widetildeS, widetildeE, widetildeI are obtained from S, E, and I by scaling by alpha which is the ratio of people who report being infected. One can interpret this as there is a part of population who would report infection and the other part who would not. Ultimately, we can model only the ones who would as this is mainly they who contribute to the output.","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Finally, we can check that the new model is indeed globally identifiable:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"assess_identifiability(reparam[:new_ode])","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"[1]: T. Sauer, T. Berry, D. Ebeigbe, M. Norton, A. Whalen, S. Schiff, Identifiability of infection model parameters early in an epidemic, SIAM Journal on Control and Optimization, 2022;","category":"page"},{"location":"utils/ode/#Functions-to-work-with-the-ODE-structure","page":"ODE Tools","title":"Functions to work with the ODE structure","text":"","category":"section"},{"location":"utils/ode/","page":"ODE Tools","title":"ODE Tools","text":"Modules = [StructuralIdentifiability]\nPages = [\"ODE.jl\", \"submodels.jl\"]\nOrder = [:function]","category":"page"},{"location":"utils/ode/#StructuralIdentifiability._get_var-Tuple{Any}","page":"ODE Tools","title":"StructuralIdentifiability._get_var","text":"For an expression of the form f'(t) or f(t) returns (f, true) and (f, false), resp\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability._reduce_mod_p-Tuple{Nemo.QQMPolyRingElem, Int64}","page":"ODE Tools","title":"StructuralIdentifiability._reduce_mod_p","text":"_reduce_mod_p(f, p)\n\nReduces a polynomial/rational function over Q modulo p\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.power_series_solution-Union{Tuple{P}, Tuple{T}, Tuple{ODE{P}, Dict{P, T}, Dict{P, T}, Dict{P, Vector{T}}, Int64}} where {T<:AbstractAlgebra.FieldElem, P<:AbstractAlgebra.MPolyRingElem{T}}","page":"ODE Tools","title":"StructuralIdentifiability.power_series_solution","text":"power_series_solution(ode, param_values, initial_conditions, input_values, prec)\n\nInput:\n\node - an ode to solve\nparam_values - parameter values, must be a dictionary mapping parameter to a value\ninitial_conditions - initial conditions of ode, must be a dictionary mapping state variable to a value\ninput_values - power series for the inputs presented as a dictionary variable => list of coefficients\nprec - the precision of solutions\n\nOutput:\n\ncomputes a power series solution with precision prec presented as a dictionary variable => corresponding coordinate of the solution\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.reduce_ode_mod_p-Tuple{ODE{<:AbstractAlgebra.MPolyRingElem{Nemo.QQFieldElem}}, Int64}","page":"ODE Tools","title":"StructuralIdentifiability.reduce_ode_mod_p","text":"reduce_ode_mod_p(ode, p)\n\nInput: ode is an ODE over QQ, p is a prime number Output: the reduction mod p, throws an exception if p divides one of the denominators\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.set_parameter_values-Union{Tuple{P}, Tuple{T}, Tuple{ODE{P}, Dict{P, T}}} where {T<:AbstractAlgebra.FieldElem, P<:AbstractAlgebra.MPolyRingElem{T}}","page":"ODE Tools","title":"StructuralIdentifiability.set_parameter_values","text":"set_parameter_values(ode, param_values)\n\nInput:\n\node - an ODE as above\nparam_values - values for (possibly, some of) the parameters as dictionary parameter => value\n\nOutput:\n\nnew ode with the parameters in param_values plugged with the given numbers\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.find_submodels-Union{Tuple{ODE{P}}, Tuple{P}} where P<:AbstractAlgebra.MPolyRingElem","page":"ODE Tools","title":"StructuralIdentifiability.find_submodels","text":"find_submodels(ode)\n\nThe function calculates and returns all valid submodels given a system of ODEs.\n\nInput:\n\node - an ODEs system to be studied\n\nOutput: \n\nA list of submodels represented as ode objects\n\nExample:\n\n>ode = @ODEmodel(x1'(t) = x1(t)^2, \n x2'(t) = x1(t) * x2(t), \n y1(t) = x1(t), \n y2(t) = x2(t))\n>find_submodels(ode)\n ODE{QQMPolyRingElem}[\n \n x1'(t) = a(t)*x2(t)^2 + x1(t)\n y1(t) = x1(t)\n ]\n\n\n\n\n\n","category":"method"},{"location":"utils/primality/#Primality-Checks","page":"Primality Checks","title":"Primality Checks","text":"","category":"section"},{"location":"utils/primality/","page":"Primality Checks","title":"Primality Checks","text":"Pages=[\"primality.md\"]","category":"page"},{"location":"utils/primality/","page":"Primality Checks","title":"Primality Checks","text":"StructuralIdentifiability.check_primality","category":"page"},{"location":"utils/primality/#StructuralIdentifiability.check_primality","page":"Primality Checks","title":"StructuralIdentifiability.check_primality","text":"check_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem}, extra_relations::Array{QQMPolyRingElem, 1})\n\nThe function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.\n\nThe extra_relations allows adding more polynomials to the generators (not affecting the saturation).\n\n\n\n\n\ncheck_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem})\n\nThe function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.\n\n\n\n\n\n","category":"function"},{"location":"tutorials/identifiable_functions/#Globally-Identifiable-Functions","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"","category":"section"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"In addition to assessing identifiabuility of given functions of parameters or states, StructuralIdentifiability.jl provides the function find_identifiable_functions which finds a set of identifiable functions such that any other identifiable function can be expressed using them. This allows to find out what actually is identifiable and what does non-identifiability in the model at hand looks like.","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"For example, consider the following model[1].","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"using StructuralIdentifiability # hide\nLLW1987 = @ODEmodel(\n x1'(t) = -p1 * x1(t) + p2 * u(t),\n x2'(t) = -p3 * x2(t) + p4 * u(t),\n x3'(t) = -(p1 + p3) * x3(t) + (p4 * x1(t) + p2 * x2(t)) * u(t),\n y1(t) = x3(t)\n)","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"Several decades ago, this model was introduced to demonstrate the presence of nontrivial unidentifiability in nonlinear systems of ODEs. Nowadays, we can automatically find the identifiable combinations of parameters:","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"using StructuralIdentifiability # hide\nfind_identifiable_functions(LLW1987)","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"From these expressions, we see that the values of p1 and p3 are not identifiable but an unordered pair of numbers {p1, p3} is uniquely determined since p1 + p3 and p1 * p3 are known. Furthermore, we see that, for fixed input and output, p2 and p4 can take infinitely many values but knowing one of them, we would also be able to determine the other.","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"Moreover, we can find generators of all identifiable functions in parameters and states:","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"find_identifiable_functions(LLW1987, with_states = true)","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"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 :strong to try to simplify further (at the expense of heavier computation) or use :weak to simplify less (but compute faster).","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"As assess_identifiability and assess_local_identifiability, find_identifiable_functions accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"[1]: Y. Lecourtier, F. Lamnabhi-Lagarrigue, and E. Walter, A method to prove that nonlinear models can be unidentifiable, Proceedings of the 26th Conference on Decision and Control, December 1987, 2144-2145;","category":"page"},{"location":"tutorials/identifiability/#Identifiability-of-Differential-Models-(Local-and-Global)","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"","category":"section"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Recall that we consider ODE models in the state-space form","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"begincases\nmathbfx(t) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"where mathbfx(t) mathbfy(t), and mathbfu(t) are time-dependent states, outputs, and inputs, respectively, and mathbfp are scalar parameters. We will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs. Typically, two types of identifiability are distinguished","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"local identifiability: the value can be recovered up to finitely many options;\nglobal identifiability: the value can be recovered uniquely.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Note that in the case of states, term observability it typically used. We will use identifiability for both states and parameters for brevity and uniformity. While the notion of global identifiability is much more precise, assessing local identifiability is typically much faster, and can be performed for the models whose global identifiability analysis is out of reach.","category":"page"},{"location":"tutorials/identifiability/#Local-identifiability","page":"Identifiability of Differential Models (Local and Global)","title":"Local identifiability","text":"","category":"section"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"We consider the Lotka-Volterra model:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"begincases\nx_1(t) = a x_1(t) - b x_1(t) x_2(t) + u(t)\nx_2(t) = -c x_2(t) + d x_1(t) x_2(t)\ny(t) = x_1(t)\nendcases","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"The local identifiability of all parameters and states in this model can be assessed as follows","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = a * x1(t) - b * x1(t) * x2(t) + u(t),\n x2'(t) = -c * x2(t) + d * x1(t) * x2(t),\n y(t) = x1(t)\n)\n\nassess_local_identifiability(ode)","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"We see that x_1(t) is locally identifiable (no surprises, this is an output), a c and d are identifiable as well. On the other hand, x_2(t) and b are nonidentifiable. This can be explained by the following scaling symmetry","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"x_2(t) to lambda x_2(t) quad b to fracblambda","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"which preserves input and output for every nonzero lambda. The algorithm behind this call is the one due to Sedoglavic[1].","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Function assess_local_identifiability has several optional parameters","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"funcs_to_check a list of specific functions of parameters and states to check identifiability for (see an example below). If not provided, the identifiability is assessed for all parameters and states.\nprob_threshold (default 099, i.e. 99%) is the probability of correctness. The algorithm can, in theory, produce wrong result, but the probability that it is correct is guaranteed to be at least prob_threshold. However, the probability bounds we use are quite conservative, so the actual probability of correctness is likely to be much higher.\ntype (default :SE). By default, the algorithm checks the standard single-experiment identifiability. If one sets type = :ME, then the algorithm checks multi-experiment identifiability, that is, identifiability from several experiments with independent initial conditions (the algorithm from [2] is used).\nloglevel (default Logging.Info). The minimal level of logging messages to be displayed. Available options: Logging.Debug, Logging.Info, Logging.Warn, and Logging.Error.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Note that the scaling symmetry given above suggests that b x_2(t) may in fact be identifiable. This can be checked using funcs_to_check parameter:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"assess_local_identifiability(ode, funcs_to_check = [b * x2])","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Indeed!","category":"page"},{"location":"tutorials/identifiability/#Global-identifiability","page":"Identifiability of Differential Models (Local and Global)","title":"Global identifiability","text":"","category":"section"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"One can obtain more refined information about a model using assess_identifiability function. We will showcase it using the Goodwin oscillator model [3].","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = -b * x1(t) + 1 / (c + x4(t)),\n x2'(t) = alpha * x1(t) - beta * x2(t),\n x3'(t) = gama * x2(t) - delta * x3(t),\n x4'(t) = sigma * x4(t) * (gama * x2(t) - delta * x3(t)) / x3(t),\n y(t) = x1(t)\n)\n\nassess_identifiability(ode)","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"As a result, each parameter/state is assigned one of the labels :globally (globally identifiable), :locally (locally but not globally identifiable), or :nonidentifiable (not identifiable, even locally). The algorithm behind this computation follows [4].","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Similarly to assess_local_identifiability, this function has optional parameters:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"funcs_to_check a list of specific functions of parameters and states to check identifiability for (see an example below). If not provided, the identifiability is assessed for all parameters and states. Note that the computations for states may be more involved than for the parameters, so one may want to call the function with funcs_to_check = ode.parameters if the call assess_identifiability(ode) takes too long.\nprob_threshold (default 099, i.e. 99%) is the probability of correctness. Same story as above: the probability estimates are very conservative, so the actual error probability is much lower than 1%. Also, currently, the probability of correctness does not include the probability of correctness of the modular reconstruction for Groebner bases. This probability is ensured by an additional check modulo a large prime, and can be neglected for practical purposes.\nloglevel (default Logging.Info). The minimal level of logging messages to be displayed. Available options: Logging.Debug, Logging.Info, Logging.Warn, and Logging.Error.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Using funcs_to_check parameter, one can further inverstigate the nature of the lack of identifiability in the model at hand. For example, for the Goodwin oscillator, we can check if beta + delta and beta * delta are identifiable:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"assess_identifiability(ode, funcs_to_check = [beta + delta, beta * delta])","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"And we see that they indeed are. This means, in particular, that the reason why beta and delta are not identifiable is because their values can be exchanged. One may wonder how could we guess these functions beta + delta, beta * delta. In fact, they can be just computed using find_identifiable_functions function as we will explain in the next tutorial. Stay tuned!","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[1]: A. Sedoglavic, A probabilistic algorithm to test local algebraic observability in polynomial time, Journal of Symbolic Computation, 2002.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[2]: A. Ovchinnikov, A. Pillay, G. Pogudin, T. Scanlon, Multi-experiment Parameter Identifiability of ODEs and Model Theory, SIAM Journal on Applied Algebra and Geometry, 2022.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[3]: D. Gonze, P. Ruoff, The Goodwin Oscillator and its Legacy, Acta Biotheoretica, 2020.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[4]: R. Dong, C. Goodbrake, H. Harrington, G. Pogudin, Differential elimination for dynamical models via projections with applications to structural identifiability, SIAM Journal on Applied Algebra and Geometry, 2023.","category":"page"},{"location":"export/export/#Exporting-to-Other-Systems","page":"Exporting to Other Systems","title":"Exporting to Other Systems","text":"","category":"section"},{"location":"export/export/","page":"Exporting to Other Systems","title":"Exporting to Other Systems","text":"Here we put some helpful utilities to export your code to other identifiability software.","category":"page"},{"location":"export/export/","page":"Exporting to Other Systems","title":"Exporting to Other Systems","text":"print_for_maple\nprint_for_DAISY\nprint_for_GenSSI\nprint_for_COMBOS","category":"page"},{"location":"export/export/#StructuralIdentifiability.print_for_maple","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_maple","text":"print_for_maple(ode, package)\n\nPrints the ODE in the format accepted by maple packages\n\nSIAN (https://github.com/pogudingleb/SIAN) if package=:SIAN\nDifferentialAlgebra if package=:DifferentialAlgebra\nDifferentialThomas if package=:DifferentialThomas\n\n\n\n\n\n","category":"function"},{"location":"export/export/#StructuralIdentifiability.print_for_DAISY","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_DAISY","text":"print_for_DAISY(ode)\n\nPrints the ODE in the format accepted by DAISY (https://daisy.dei.unipd.it/)\n\n\n\n\n\n","category":"function"},{"location":"export/export/#StructuralIdentifiability.print_for_GenSSI","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_GenSSI","text":"print_for_GenSSI(ode)\n\nPrints the ODE in the format accepted by GenSSI 2.0 (https://github.com/genssi-developer/GenSSI)\n\n\n\n\n\n","category":"function"},{"location":"export/export/#StructuralIdentifiability.print_for_COMBOS","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_COMBOS","text":"print_for_COMBOS(ode)\n\nPrints the ODE in the format accepted by COMBOS (http://biocyb1.cs.ucla.edu/combos/)\n\n\n\n\n\n","category":"function"},{"location":"utils/util/#Other-Helpful-Functions","page":"Other Utilities","title":"Other Helpful Functions","text":"","category":"section"},{"location":"utils/util/","page":"Other Utilities","title":"Other Utilities","text":"Pages=[\"util.md\"]","category":"page"},{"location":"utils/util/","page":"Other Utilities","title":"Other Utilities","text":"Modules = [StructuralIdentifiability]\nPages = [\"util.jl\"]","category":"page"},{"location":"utils/util/#StructuralIdentifiability.compare_rational_func_by-Union{Tuple{By}, Tuple{Any, Any, By}, Tuple{Any, Any, By, Any}} where By","page":"Other Utilities","title":"StructuralIdentifiability.compare_rational_func_by","text":"compare_rational_func_by(f, g, by)\n\nReturns \n\n-1 if f < g,\n0 if f = g, and \n1 if f > g.\n\nFunctions' numerators and denominators are compared using by.\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.decompose_derivative-Tuple{String, Array{String}}","page":"Other Utilities","title":"StructuralIdentifiability.decompose_derivative","text":"decompose_derivative(varname, prefixes)\n\nDetermines if it is possible to represent the varname as a_number where a is an element of prefixes If yes, returns a pair (a, number), otherwise nothing\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.dennums_to_fractions-Union{Tuple{Array{Vector{T}, 1}}, Tuple{T}} where T","page":"Other Utilities","title":"StructuralIdentifiability.dennums_to_fractions","text":"dennums_to_fractions(dennums)\n\nReturns the field generators represented by fractions.\n\nInput: an array of arrays of polynomials, as in [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]\n\nOutput: an array of fractions [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.eval_at_dict-Union{Tuple{P}, Tuple{P, Dict{P, <:AbstractAlgebra.RingElem}}} where P<:AbstractAlgebra.MPolyRingElem","page":"Other Utilities","title":"StructuralIdentifiability.eval_at_dict","text":"eval_at_dict(f, d)\n\nEvaluates a polynomial/rational function on a dictionary of type var => val and missing values are replaced with zeroes\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.extract_coefficients-Union{Tuple{P}, Tuple{P, Vector{P}}} where P<:AbstractAlgebra.MPolyRingElem","page":"Other Utilities","title":"StructuralIdentifiability.extract_coefficients","text":"extract_coefficients(poly, variables)\n\nInput:\n\npoly - multivariate polynomial\nvariables - a list of variables from the generators of the ring of p\n\nOutput:\n\ndictionary with keys being tuples of length lenght(variables) and values being polynomials in the variables other than those which are the coefficients at the corresponding monomials (in a smaller polynomial ring)\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.fractions_to_dennums-Tuple{Any}","page":"Other Utilities","title":"StructuralIdentifiability.fractions_to_dennums","text":"fractions_to_dennums(fractions)\n\nReturns the field generators represented by lists of denominators and numerators.\n\nInput: an array of fractions, as in [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]\n\nOutput: an array of arrays of polynomials, [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.gen_tag_name","page":"Other Utilities","title":"StructuralIdentifiability.gen_tag_name","text":"gen_tag_name(base; stop_words)\ngen_tag_names(n, base; stop_words)\n\nGenerates a string which will not collide with the words in stop_words.\n\nArguments\n\nn: Generates a sequence of unique strings of length n\nbase: A string or a vector of strings, the base for the generated sequence\nstop_words: A vector of strings, stop words\n\n\n\n\n\n","category":"function"},{"location":"utils/util/#StructuralIdentifiability.make_substitution-Union{Tuple{P}, NTuple{4, P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Other Utilities","title":"StructuralIdentifiability.make_substitution","text":"make_substitution(f, var_sub, val_numer, val_denom)\n\nSubstitute a variable in a polynomial with an expression\n\nInput:\n\nf - the polynomial\nvar_sub - the variable to be substituted\nvar_numer - numerator of the substitution expression\nvar_denom - denominator of the substitution expression\n\nOutput:\n\npolynomial - result of substitution\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.parent_ring_change-Tuple{AbstractAlgebra.MPolyRingElem, AbstractAlgebra.MPolyRing}","page":"Other Utilities","title":"StructuralIdentifiability.parent_ring_change","text":"parent_ring_change(poly, new_ring)\n\nConverts a polynomial to a different polynomial ring Input\n\npoly - a polynomial to be converted\nnew_ring - a polynomial ring such that every variable name appearing in poly appears among the generators\n\nOutput:\n\na polynomial in new_ring “equal” to poly\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.switch_ring-Tuple{AbstractAlgebra.MPolyRingElem, AbstractAlgebra.MPolyRing}","page":"Other Utilities","title":"StructuralIdentifiability.switch_ring","text":"switch_ring(v, ring)\n\nFor a variable v, returns a variable in ring with the same name\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.uncertain_factorization-Tuple{AbstractAlgebra.MPolyRingElem{Nemo.QQFieldElem}}","page":"Other Utilities","title":"StructuralIdentifiability.uncertain_factorization","text":"uncertain_factorization(f)\n\nInput:\n\nf - polynomial with rational coefficients\n\nOutput:\n\nlist of pairs (div, certainty) where\ndiv's are divisors of f such that f is their product with certain powers\nif certainty is true, div is Q-irreducible\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#Wronskian-Tools","page":"Wronskian Tools","title":"Wronskian Tools","text":"","category":"section"},{"location":"utils/wronskian/","page":"Wronskian Tools","title":"Wronskian Tools","text":"Modules = [StructuralIdentifiability]\nPages = [\"wronskian.jl\"]","category":"page"},{"location":"utils/wronskian/#StructuralIdentifiability.get_max_below-Tuple{StructuralIdentifiability.ExpVectTrie, Vector{Int64}}","page":"Wronskian Tools","title":"StructuralIdentifiability.get_max_below","text":"get_max_below(t, vect)\n\nInput:\n\nt - a trie with exponent vectors\nvect - yet another exponent vector\n\nOutput:\n\na pair (d, v) where v is a vector in the trie which is componentwise ≤ vect and the difference d is as small as possible\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#StructuralIdentifiability.massive_eval-Tuple{Any, Any}","page":"Wronskian Tools","title":"StructuralIdentifiability.massive_eval","text":"massive_eval(polys, eval_dict)\n\nInput:\n\npolys - a list of polynomials\neval_dict - dictionary from variables to the values. Missing values are treated as zeroes\n\nOutput:\n\na list of values of the polynomials\n\nEvaluates a list of polynomials at a point. Assumes that multiplications are relatively expensive (like in truncated power series) so all the monomials are precomputed first and the values of monomials of lower degree are cached and used to compute the values of the monomials of higher degree\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#StructuralIdentifiability.monomial_compress-Tuple{Any, ODE}","page":"Wronskian Tools","title":"StructuralIdentifiability.monomial_compress","text":"monomial_compress(io_equation, ode)\n\nCompresses an input-output equation for the rank computation Input:\n\nio_equation - input-output equation\node - the corresponding ODE model\n\nOutput:\n\npair (coeffs, terms) such that:\nsum of coeffs[i] * terms[i] = io_equation\ncoeffs involve only parameters, terms involve only inputs and outputs\nlength of the representation is the smallest possible\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#StructuralIdentifiability.wronskian-Union{Tuple{P}, Tuple{Dict{P, P}, ODE{P}}} where P<:AbstractAlgebra.MPolyRingElem","page":"Wronskian Tools","title":"StructuralIdentifiability.wronskian","text":"wronskian(io_equations, ode)\n\nInput:\n\nio_equations - a set of io-equations in the form of the Dict as returned by find_ioequations\node - the ODE object\n\nOutput:\n\na list of Wronskians evaluated at a point modulo prime\n\nComputes the Wronskians of io_equations\n\n\n\n\n\n","category":"method"},{"location":"utils/global_identifiability/#Global-Identifiability-Tools","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"","category":"section"},{"location":"utils/global_identifiability/","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"Pages=[\"global_identifiability.md\"]","category":"page"},{"location":"utils/global_identifiability/","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"CurrentModule=StructuralIdentifiability","category":"page"},{"location":"utils/global_identifiability/","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"StructuralIdentifiability.RationalFunctionField\nStructuralIdentifiability.field_contains\nStructuralIdentifiability.get_degree_and_coeffsize","category":"page"},{"location":"utils/global_identifiability/#StructuralIdentifiability.RationalFunctionField","page":"Global Identifiability Tools","title":"StructuralIdentifiability.RationalFunctionField","text":"RationalFunctionField\n\nA subfield of the field of rational functions over the rationals.\n\nExample\n\nusing Nemo\nusing StructuralIdentifiability: RationalFunctionField\n\nR, (x, y, z) = QQ[\"x\", \"y\", \"z\"]\n\n# Constructs a subfield generated by x / y, y / z\nrff = RationalFunctionField([x // y, y // z])\n\n# Constructs a subfield generated by y / x, 1 / x, z / y\nrff = RationalFunctionField([[x, y, R(1)], [y, z]])\n\n\n\n\n\n","category":"type"},{"location":"utils/global_identifiability/#StructuralIdentifiability.field_contains","page":"Global Identifiability Tools","title":"StructuralIdentifiability.field_contains","text":"field_contains(field, ratfuncs, prob_threshold)\n\nChecks whether given rational function field field contains given rational functions ratfuncs (represented as a list of lists). The result is correct with probability at least prob_threshold\n\nInputs:\n\nfield - a rational function field\nratfuncs - a list of lists of polynomials. Each of the lists, say, [f1, ..., fn], defines generators f2/f1, ..., fn/f1.\nprob_threshold real number from (0, 1)\n\nOutput:\n\na list L[i] of bools of length length(rat_funcs) such that L[i] is true iff the i-th function belongs to field\n\n\n\n\n\n","category":"function"},{"location":"utils/global_identifiability/#StructuralIdentifiability.get_degree_and_coeffsize","page":"Global Identifiability Tools","title":"StructuralIdentifiability.get_degree_and_coeffsize","text":"get_degree_and_coeffsize(f)\n\nfor f being a polynomial/rational function over rationals (QQ) returns a tuple (degree, max_coef_size)\n\n\n\n\n\n","category":"function"},{"location":"utils/reparametrization/#Model-reparametrization","page":"Model reparametrization","title":"Model reparametrization","text":"","category":"section"},{"location":"utils/reparametrization/","page":"Model reparametrization","title":"Model reparametrization","text":"StructuralIdentifiability.reparametrize_global","category":"page"},{"location":"utils/reparametrization/#StructuralIdentifiability.reparametrize_global","page":"Model reparametrization","title":"StructuralIdentifiability.reparametrize_global","text":"reparametrize_global(ode, options...)\n\nFinds a reparametrization of ode in terms of globally identifiabile functions.\n\nReturns a tuple (new_ode, new_vars, implicit_relations), such that:\n\nnew_ode is the reparametrized ODE system\nnew_vars is a mapping from the new variables to the original ones\nrelations is the array of implicit algebraic relations among new_vars. Usually, relations is empty.\n\nOptions\n\nThe function accepts the following optional arguments.\n\nseed: A float in the range from 0 to 1, random seed (default is seed = 42). \nprob_threshold: The probability of correctness (default is prob_threshold = 0.99).\n\nExample\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = a * x1(t) - b*x1(t)*x2(t),\n x2'(t) = -c * x2(t) + d*x1(t)*x2(t),\n y(t) = x1(t)\n)\n\nnew_ode, new_vars, relations = reparametrize_global(ode)\n\nThen, we have the following:\n\n# new_ode\nX2'(t) = X1(t)*X2(t)*a2 - X2(t)*a1\nX1'(t) = -X1(t)*X2(t) + X1(t)*a3\ny1(t) = X1(t)\n\n# new_vars\nDict{Nemo.QQMPolyRingElem, AbstractAlgebra.Generic.Frac{Nemo.QQMPolyRingElem}} with 6 entries:\n X2 => b*x2\n y1 => y\n X1 => x1\n a2 => d\n a3 => a\n a1 => c\n\nNotice that the new_ode is fully identifiabile, and has 1 less parameter compared to the original one.\n\n\n\n\n\n","category":"function"},{"location":"utils/power_series_utils/#Power-Series-Utilities","page":"Power Series Tools","title":"Power Series Utilities","text":"","category":"section"},{"location":"utils/power_series_utils/","page":"Power Series Tools","title":"Power Series Tools","text":"Pages =[\"power_series_utils.md\"]","category":"page"},{"location":"utils/power_series_utils/","page":"Power Series Tools","title":"Power Series Tools","text":"Modules = [StructuralIdentifiability]\nPages = [\"power_series_utils.jl\"]","category":"page"},{"location":"utils/power_series_utils/#StructuralIdentifiability._matrix_inv_newton_iteration-Union{Tuple{T}, Tuple{AbstractAlgebra.MatElem{T}, AbstractAlgebra.MatElem{T}}} where T<:(AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.FieldElem})","page":"Power Series Tools","title":"StructuralIdentifiability._matrix_inv_newton_iteration","text":"_matrix_inv_newton_iteration(M, Minv)\n\nPerforms a single step of Newton iteration for inverting M with Minv being a partial result\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_diff-Tuple{AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.RingElem}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_diff","text":"ps_diff(ps)\n\nInput:\n\nps - (absolute capped) univariate power series\n\nOutput:\n\nthe derivative of ps\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_integrate-Tuple{AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.FieldElem}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_integrate","text":"ps_integrate(ps)\n\nInput:\n\nps - (absolute capped) univariate power series\n\nOutput:\n\nthe integral of ps without constant term\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_homlinear_de-Union{Tuple{T}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}, Int64}} where T<:AbstractAlgebra.FieldElem","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_homlinear_de","text":"ps_matrix_homlinear_de(A, Y0, prec)\n\nInput:\n\nA - a square matrix with entries in a univariate power series ring\nY0 - a square invertible matrix over the base field\n\nOutput:\n\nmatrix Y such that Y' = AY up to precision of A - 1 and Y(0) = Y0\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_inv","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_inv","text":"ps_matrix_inv(M, prec)\n\nInput:\n\nM - a square matrix with entries in a univariate power series ring it is assumed that M(0) is invertible and all entries having the same precision\nprec - an integer, precision, if -1 then defaults to precision of M\n\nOutput:\n\nthe inverse of M computed up to prec\n\n\n\n\n\n","category":"function"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_linear_de-Union{Tuple{T}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}, Int64}} where T<:AbstractAlgebra.FieldElem","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_linear_de","text":"ps_matrix_linear_de(A, B, Y0, prec)\n\nInput:\n\nA, B - square matrices with entries in a univariate power series ring\nY0 - a matrix over the base field with the rows number the same as A\n\nOutput:\n\nmatrix Y such that Y' = AY + B up to precision of A - 1 and Y(0) = Y0\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_log-Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.FieldElem}}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_log","text":"ps_matrix_log(M)\n\nInput:\n\nM - a square matrix with entries in a univariate power series ring it is assumed that M(0) is the identity\n\nOutput:\n\nthe natural log of M\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_ode_solution-Union{Tuple{P}, Tuple{T}, Tuple{Vector{P}, Dict{P, T}, Dict{P, Vector{T}}, Int64}} where {T<:AbstractAlgebra.FieldElem, P<:AbstractAlgebra.MPolyRingElem{T}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_ode_solution","text":"ps_ode_solution(equations, ic, inputs, prec)\n\nInput:\n\nequations - a system of the form A(x u mu)x - B(x u mu) = 0, where A is a generically nonsingular square matrix. Assumption: A is nonzero at zero\nic - initial conditions for x's (dictionary)\ninputs - power series for inputs represented as arrays (dictionary)\nprec - precision of the solution\n\nOutput:\n\npower series solution of the system\n\n\n\n\n\n","category":"method"},{"location":"identifiability/identifiability/#Functions-to-Assess-Identifiability","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"","category":"section"},{"location":"identifiability/identifiability/#Assessing-All-Types-of-Identifiability","page":"Functions to Assess Identifiability","title":"Assessing All Types of Identifiability","text":"","category":"section"},{"location":"identifiability/identifiability/","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"assess_identifiability","category":"page"},{"location":"identifiability/identifiability/#StructuralIdentifiability.assess_identifiability","page":"Functions to Assess Identifiability","title":"StructuralIdentifiability.assess_identifiability","text":"assess_identifiability(ode; funcs_to_check = [], prob_threshold=0.99, loglevel=Logging.Info)\n\nInput:\n\node - the ODE model\nfuncs_to_check - list of functions to check identifiability for; if empty, all parameters and states are taken\nprob_threshold - probability of correctness.\nloglevel - the minimal level of log messages to display (Logging.Info by default)\n\nAssesses identifiability of a given ODE model. The result is guaranteed to be correct with the probability at least prob_threshold. The function returns an (ordered) dictionary from the functions to check to their identifiability properties (one of :nonidentifiable, :locally, :globally).\n\n\n\n\n\n","category":"function"},{"location":"identifiability/identifiability/#Assessing-Local-Identifiability","page":"Functions to Assess Identifiability","title":"Assessing Local Identifiability","text":"","category":"section"},{"location":"identifiability/identifiability/","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"assess_local_identifiability","category":"page"},{"location":"identifiability/identifiability/#StructuralIdentifiability.assess_local_identifiability","page":"Functions to Assess Identifiability","title":"StructuralIdentifiability.assess_local_identifiability","text":"assess_local_identifiability(ode::ODE{P}; funcs_to_check::Array{<: Any, 1}, prob_threshold::Float64=0.99, type=:SE, loglevel=Logging.Info) where P <: MPolyRingElem{Nemo.QQFieldElem}\n\nChecks the local identifiability/observability of the functions in funcs_to_check. The result is correct with probability at least prob_threshold.\n\nCall this function if you have a specific collection of parameters of which you would like to check local identifiability.\n\ntype can be either :SE (single-experiment identifiability) or :ME (multi-experiment identifiability). If the type is :ME, states are not allowed to appear in the funcs_to_check.\n\n\n\n\n\n","category":"function"},{"location":"identifiability/identifiability/#Finding-Identifiable-Functions","page":"Functions to Assess Identifiability","title":"Finding Identifiable Functions","text":"","category":"section"},{"location":"identifiability/identifiability/","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"find_identifiable_functions","category":"page"},{"location":"identifiability/identifiability/#StructuralIdentifiability.find_identifiable_functions","page":"Functions to Assess Identifiability","title":"StructuralIdentifiability.find_identifiable_functions","text":"find_identifiable_functions(ode::ODE; options...)\n\nFinds all functions of parameters/states that are identifiable in the given ODE system.\n\nOptions\n\nThis functions takes the following optional arguments:\n\nwith_states: When true, also reports the identifiabile functions in the ODE states. Default is false.\nsimplify: The extent to which the output functions are simplified. Stronger simplification may require more time. Possible options are:\n:standard: Default simplification.\n:weak: Weak simplification. This option is the fastest, but the output functions can be quite complex.\n:strong: Strong simplification. This option is the slowest, but the output\nfunctions are nice and simple.\n:absent: No simplification.\nprob_threshold: A float in the range from 0 to 1, the probability of correctness. Default is 0.99.\nseed: The rng seed. Default value is 42.\nloglevel - the minimal level of log messages to display (Logging.Info by default)\n\nExample\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x0'(t) = -(a01 + a21) * x0(t) + a12 * x1(t),\n x1'(t) = a21 * x0(t) - a12 * x1(t),\n y(t) = x0(t)\n)\n\nfind_identifiable_functions(ode)\n\n# prints\n3-element Vector{AbstractAlgebra.Generic.Frac{Nemo.QQMPolyRingElem}}:\n a12 + a01 + a21\n a12*a01\n\n\n\n\n\n","category":"function"},{"location":"utils/elimination/#Elimination","page":"Elimination","title":"Elimination","text":"","category":"section"},{"location":"utils/elimination/","page":"Elimination","title":"Elimination","text":"Pages=[\"elimination.md\"]","category":"page"},{"location":"utils/elimination/","page":"Elimination","title":"Elimination","text":"Modules = [StructuralIdentifiability]\nPages = [\"elimination.jl\"]","category":"page"},{"location":"utils/elimination/#StructuralIdentifiability.Bezout_matrix-Union{Tuple{P}, Tuple{P, P, P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Elimination","title":"StructuralIdentifiability.Bezout_matrix","text":"Bezout_matrix(f, g, var_elim)\n\nCompute the Bezout matrix of two polynomials f, g with respect to var_elim\n\nInputs:\n\nf - first polynomial\ng - second polynomial\nvar_elim - variable, of which f and g are considered as polynomials\n\nOutput:\n\nM::MatrixElem - The Bezout matrix\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.Sylvester_matrix-Union{Tuple{P}, Tuple{P, P, P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Elimination","title":"StructuralIdentifiability.Sylvester_matrix","text":"Sylvester_matrix(f, g, var_elim)\n\nCompute the Sylvester matrix of two polynomials f, g with respect to var_elim Inputs:\n\nf - first polynomial\ng - second polynomial\nvar_elim - variable, of which f and g are considered as polynomials\n\nOutput:\n\nM::MatrixElem - The Sylvester matrix\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.choose-Union{Tuple{P}, Tuple{Vector{P}, Any}} where P<:(AbstractAlgebra.MPolyRingElem{<:AbstractAlgebra.FieldElem})","page":"Elimination","title":"StructuralIdentifiability.choose","text":"choose(polys, generic_point_generator)\n\nInput:\n\npolys - an array of distinct irreducible polynomials in the same ring\ngeneric_point_generator - a generic point generator as described above for one of polys\n\nOutput:\n\nthe polynomial that vanishes at the generic_point_generator\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.eliminate_var-Union{Tuple{P}, Tuple{P, P, P, Any}} where P<:(AbstractAlgebra.MPolyRingElem{<:AbstractAlgebra.FieldElem})","page":"Elimination","title":"StructuralIdentifiability.eliminate_var","text":"eliminate_var(f, g, var_elim, generic_point_generator)\n\nEliminate a variable from a pair of polynomials\n\nInput:\n\nf and g - polynomials\nvar_elim - variable to be eliminated\ngeneric_point_generator - a generic point generator object for the factor of the resultant of f and g of interest\n\nOutput:\n\npolynomial - the desired factor of the resultant of f and g\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.simplify_matrix-Union{Tuple{AbstractAlgebra.MatElem{P}}, Tuple{P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Elimination","title":"StructuralIdentifiability.simplify_matrix","text":"simplify_matrix(M)\n\nEliminate GCD of entries of every row and column\n\nInput:\n\nM::MatrixElem - matrix to be simplified\n\nOutput:\n\nM::MatrixElem - Simplified matrix\nextra_factors::Vector{AbstractAlgebra.MPolyRingElem} - array of GCDs eliminated from M.\n\n\n\n\n\n","category":"method"},{"location":"tutorials/discrete_time/#Identifiability-of-Discrete-Time-Models-(Local)","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"","category":"section"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Now we consider a discrete-time model in the state-space form","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"begincases\nDeltamathbfx(t) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases quad textwhere quad Delta mathbfx(t) = mathbfx(t + 1) - mathbfx(t)","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"and mathbfx(t) mathbfy(t), and mathbfu(t) are time-dependent states, outputs, and inputs, respectively, and mathbfp are scalar parameters. As in the ODE case, we will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs (in the generic case, see Definition 3 in [1] for details). Again, we will distinguish two types of identifiability","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"local identifiability: the value can be recovered up to finitely many options;\nglobal identifiability: the value can be recovered uniquely.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Currently, StructuralIdentifiability.jl allows to assess only local identifiability for discrete-time models, and below we will describe how this can be done. As a running example, we will use the following discrete version of the SIR model:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"begincases\nDelta S(t) = S(t) - beta S(t) I(t)\nDelta I(t) = I(t) + beta S(t) I(t) - alpha I(t)\nDelta R(t) = R(t) + alpha I(t)\ny(t) = I(t)\nendcases","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"where the observable is I, the number of infected people. We start with creating a system as a DiscreteSystem from ModelingToolkit:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"using ModelingToolkit\nusing StructuralIdentifiability\n\n@parameters α β\n@variables t S(t) I(t) R(t) y(t)\nD = Difference(t; dt = 1.0)\n\neqs = [D(S) ~ S - β * S * I, D(I) ~ I + β * S * I - α * I, D(R) ~ R + α * I]\n@named sir = DiscreteSystem(eqs)","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Once the model is defined, we can assess identifiability by providing the formula for the observable:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(sir; measured_quantities = [y ~ I])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"For each parameter or state, the value in the returned dictionary is true (1) if the parameter is locally identifiable and false (0) otherwise. We see that R(t) is not identifiable, which makes sense: it does not affect the dynamics of the observable in any way.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"In principle, it is not required to give a name to the observable, so one can write this shorter","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(sir; measured_quantities = [I])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"The assess_local_identifiability function has three important keyword arguments:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"funcs_to_check is a list of functions for which one want to assess identifiability, for example, the following code will check if β * S is locally identifiable.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(sir; measured_quantities = [I], funcs_to_check = [β * S])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"prob_threshold is the probability of correctness (default value 0.99, i.e., 99%). The underlying algorithm is a Monte-Carlo algorithm, so in principle it may produce incorrect result but the probability of correctness of the returned result is guaranteed to be at least prob_threshold (in fact, the employed bounds are quite conservative, so in practice incorrect result is almost never produced).\nknown_ic is a list of the states for which initial conditions are known. In this case, the identifiability results will be valid not at any time point t but only at t = 0.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"As other main functions in the package, assess_local_identifiability accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"The implementation is based on a version of the observability rank criterion and will be described in a forthcoming paper.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"[1]: S. Nõmm, C. Moog, Identifiability of discrete-time nonlinear systems, IFAC Proceedings Volumes, 2004.","category":"page"},{"location":"tutorials/creating_ode/#Creating-ODE-System","page":"Creating ODE System","title":"Creating ODE System","text":"","category":"section"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Most of the algorithms in StructuralIdentifiability.jl take as input system of ordinary differential equations (ODEs) in the state space form, that is:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"begincases\nmathbfx(t) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"which involves","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"a vector mathbfx(t) of the state variables of the system,\na vector mathbfu(t) of extermal inputs,\na vector mathbfp of scalar parameters,\na vector mathbfy(t) of outputs (i.e., observations),\nand vectors of rational functions mathbff and mathbfg (for discussion of the non-rational case, see this issue).","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"In the standard setup, inputs and outputs are assumed to be known, and the goal is to assess identifiability of parameters and/or states from the input-output data. In the case of states, this property is also called observability.","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"There are two ways to define such a system to be processed using StructuralIdentifiability.jl. We will demonstrate them using the following example system (Wright's population model of two mutualist species with control[1]):","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"begincases\nx_1(t) = r_1 x_1(t)(1 - c_1 x_1(t)) + fracbeta_1 x_1(t)x_2(t)chi_1 + x_2(t) + u(t)\nx_2(t) = r_2 x_2(t)(1 - c_2 x_2(t)) + fracbeta_2 x_1(t)x_2(t)chi_2 + x_1(t)\ny(t) = x_1(t)\nendcases","category":"page"},{"location":"tutorials/creating_ode/#Defining-the-model-using-@ODEmodel-macro","page":"Creating ODE System","title":"Defining the model using @ODEmodel macro","text":"","category":"section"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"One way to define the model is to use the @ODEmodel macro provided by the StructuralIdentifiability.jl package.","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) =\n r1 * x1(t) * (1 - c1 * x1(t)) + beta1 * x1(t) * x2(t) / (chi1 + x2(t)) + u(t),\n x2'(t) = r2 * x2(t) * (1 - c2 * x2(t)) + beta2 * x1(t) * x2(t) / (chi2 + x1(t)),\n y(t) = x1(t)\n)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Then one can, for example, assess identifiability of the parameters and states by","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"assess_identifiability(ode)","category":"page"},{"location":"tutorials/creating_ode/#Defining-using-ModelingToolkit","page":"Creating ODE System","title":"Defining using ModelingToolkit","text":"","category":"section"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Alternatively, one can use ModelingToolkit: encode the equations for the states as ODESystem and specify the outputs separately. In order to do this, we first introduce all functions and scalars:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"using StructuralIdentifiability, ModelingToolkit\n\n@parameters r1, r2, c1, c2, beta1, beta2, chi1, chi2\n@variables t, x1(t), x2(t), y(t), u(t)\n\nD = Differential(t)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"And then defined the system and separately the outputs:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"eqs = [\n D(x1) ~ r1 * x1 * (1 - c1 * x1) + beta1 * x1 * x2 / (chi1 + x2) + u,\n D(x2) ~ r2 * x2 * (1 - c2 * x2) + beta2 * x1 * x2 / (chi2 + x1),\n]\n\nmeasured_quantities = [y ~ x1]\n\node_mtk = ODESystem(eqs, t, name = :mutualist)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Then, for example, the identifiability of parameters and states can be assessed as follows:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"assess_identifiability(ode_mtk, measured_quantities = measured_quantities)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"[1]: D. H. Wright, A Simple, Stable Model of Mutualism Incorporating Handling Time, The American Naturalist, 1989, 134(4).","category":"page"},{"location":"ioequations/ioequations/#Finding-Input-Output-Equations","page":"Input-Output Equation tools","title":"Finding Input-Output Equations","text":"","category":"section"},{"location":"ioequations/ioequations/","page":"Input-Output Equation tools","title":"Input-Output Equation tools","text":"find_ioequations","category":"page"},{"location":"ioequations/ioequations/#StructuralIdentifiability.find_ioequations","page":"Input-Output Equation tools","title":"StructuralIdentifiability.find_ioequations","text":"find_ioequations(ode, [var_change_policy=:default])\n\nFinds the input-output equations of an ODE system Input:\n\node - the ODE system\nvar_change_policy - whether to perform automatic variable change, can be one of :default, :yes, :no\nloglevel - logging level (default: Logging.Info)\n\nOutput:\n\na dictionary from “leaders” to the corresponding input-output equations; if an extra projection is needed, it will be the value corresponding to rand_proj_var\n\n\n\n\n\n","category":"function"},{"location":"ioequations/ioequations/#Reducing-with-respect-to-Input-Output-Equations","page":"Input-Output Equation tools","title":"Reducing with respect to Input-Output Equations","text":"","category":"section"},{"location":"ioequations/ioequations/","page":"Input-Output Equation tools","title":"Input-Output Equation tools","text":"PBRepresentation\npseudodivision\ndiffreduce\nio_switch!","category":"page"},{"location":"ioequations/ioequations/#StructuralIdentifiability.PBRepresentation","page":"Input-Output Equation tools","title":"StructuralIdentifiability.PBRepresentation","text":"The structure for storing a projection-based representation of differential ideal (see Section 2.3 https://arxiv.org/abs/2111.00991). Contains the following fields:\n\ny_names - the names of the variables with finite order in the profile (typically, outputs)\nu_names - the names of the variables with infinite order in the profile (typically, inputs)\nparam_names - the names of the parameters\nprofile - the profile of the PB-representation (see Definition 2.13) as a dict from y_names with finite orders to the orders\nprojections - the corresponding projections (see Definition 2.15) as a dict from y_names to the projections\n\n\n\n\n\n","category":"type"},{"location":"ioequations/ioequations/#StructuralIdentifiability.pseudodivision","page":"Input-Output Equation tools","title":"StructuralIdentifiability.pseudodivision","text":"pseudodivision(f, g, x)\n\nComputes the result of pseudodivision of f by g as univariate polynomials in x Input:\n\nf - the polynomial to be divided\ng - the polynomial to divide by\nx - the variable for the division\n\nOutput: the pseudoremainder of f divided by g w.r.t. x\n\n\n\n\n\n","category":"function"},{"location":"ioequations/ioequations/#StructuralIdentifiability.diffreduce","page":"Input-Output Equation tools","title":"StructuralIdentifiability.diffreduce","text":"diffreduce(diffpoly, pbr)\n\nComputes the result of differential reduction of a differential polynomial diffpoly with respect to the charset defined by a PB-representation pbr Input:\n\ndiffpoly - a polynomial representing a differential polynomial to be reduced\npbr - a projection-based representation\n\nOutput: the result of differential reduction of diffpoly by pbr considered as a characteristic set (see Remark 2.20 in the paper)\n\n\n\n\n\n","category":"function"},{"location":"ioequations/ioequations/#StructuralIdentifiability.io_switch!","page":"Input-Output Equation tools","title":"StructuralIdentifiability.io_switch!","text":"io_switch(pbr)\n\nIn a single-output pb-representation pbr makes the leading variable to be the first of the inputs\n\n\n\n\n\n","category":"function"},{"location":"utils/local_identifiability/#Local-Identifiability-Tools","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"","category":"section"},{"location":"utils/local_identifiability/","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"Pages=[\"local_identifiability.md\"]","category":"page"},{"location":"utils/local_identifiability/","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"CurrentModule=StructuralIdentifiability","category":"page"},{"location":"utils/local_identifiability/","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"StructuralIdentifiability.differentiate_solution\nStructuralIdentifiability.differentiate_output","category":"page"},{"location":"utils/local_identifiability/#StructuralIdentifiability.differentiate_solution","page":"Local Identifiability Tools","title":"StructuralIdentifiability.differentiate_solution","text":"differentiate_solution(ode, params, ic, inputs, prec)\n\nInput:\n\nthe same as for power_series_solutions\n\nOutput:\n\na tuple consisting of the power series solution and a dictionary of the form (u, v) => power series, where u is a state variable v is a state or parameter, and the power series is the partial derivative of the function u w.r.t. v evaluated at the solution\n\n\n\n\n\n","category":"function"},{"location":"utils/local_identifiability/#StructuralIdentifiability.differentiate_output","page":"Local Identifiability Tools","title":"StructuralIdentifiability.differentiate_output","text":"differentiate_output(ode, params, ic, inputs, prec)\n\nSimilar to differentiate_solution but computes partial derivatives of prescribed outputs returns a dictionary of the form y_function => Dict(var => dy/dvar) where dy/dvar is the derivative of y_function with respect to var.\n\n\n\n\n\n","category":"function"},{"location":"#StructuralIdentifiability.jl","page":"Home","title":"StructuralIdentifiability.jl","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"StructuralIdentifiability.jl is a comprehensive toolbox for assessing identifiability parameters.","category":"page"},{"location":"#Installation","page":"Home","title":"Installation","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"To install StructuralIdentifiability.jl, use the Julia package manager:","category":"page"},{"location":"","page":"Home","title":"Home","text":"using Pkg\nPkg.add(\"StructuralIdentifiability\")","category":"page"},{"location":"#Citation","page":"Home","title":"Citation","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"@article{structidjl,\n author = {Dong, R. and Goodbrake, C. and Harrington, H. and Pogudin G.},\n title = {Differential Elimination for Dynamical Models via Projections with Applications to Structural Identifiability},\n journal = {SIAM Journal on Applied Algebra and Geometry},\n url = {https://doi.org/10.1137/22M1469067},\n year = {2023}\n volume = {7},\n number = {1},\n pages = {194-235}\n}","category":"page"},{"location":"#Feature-Summary","page":"Home","title":"Feature Summary","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"StructuralIdentifiability.jl can assess local and global identifiability of ODE models. In addition to these straightforward identifiability queries on individual parameters, the package can distinguish between single- and multi-experiment identifiability.","category":"page"},{"location":"#Feature-List","page":"Home","title":"Feature List","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"Local identifiability checks\nGlobal identifiability checks\nAssessment of identifiable functions of parameters and states\nModel reparametrization (experimental)","category":"page"},{"location":"#Contributing","page":"Home","title":"Contributing","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"Please refer to the SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages for guidance on PRs, issues, and other matters relating to contributing to StructuralIdentifiability.\nThere are a few community forums:\nThe #diffeq-bridged channel in the Julia Slack\nJuliaDiffEq on Gitter\nOn the Julia Discourse forums\nSee also SciML Community page","category":"page"},{"location":"#Reproducibility","page":"Home","title":"Reproducibility","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"
The documentation of this SciML package was built using these direct dependencies,","category":"page"},{"location":"","page":"Home","title":"Home","text":"using Pkg # hide\nPkg.status() # hide","category":"page"},{"location":"","page":"Home","title":"Home","text":"
","category":"page"},{"location":"","page":"Home","title":"Home","text":"
and using this machine and Julia version.","category":"page"},{"location":"","page":"Home","title":"Home","text":"using InteractiveUtils # hide\nversioninfo() # hide","category":"page"},{"location":"","page":"Home","title":"Home","text":"
","category":"page"},{"location":"","page":"Home","title":"Home","text":"
A more complete overview of all dependencies and their versions is also provided.","category":"page"},{"location":"","page":"Home","title":"Home","text":"using Pkg # hide\nPkg.status(; mode = PKGMODE_MANIFEST) # hide","category":"page"},{"location":"","page":"Home","title":"Home","text":"
","category":"page"},{"location":"","page":"Home","title":"Home","text":"You can also download the \nmanifest file and the\nproject file.","category":"page"},{"location":"input/input/#Parsing-input-ODE-system","page":"Parsing input ODE system","title":"Parsing input ODE system","text":"","category":"section"},{"location":"input/input/","page":"Parsing input ODE system","title":"Parsing input ODE system","text":"@ODEmodel(ex::Expr...)\nODE\nset_parameter_values","category":"page"},{"location":"input/input/#StructuralIdentifiability.@ODEmodel-Tuple{Vararg{Expr}}","page":"Parsing input ODE system","title":"StructuralIdentifiability.@ODEmodel","text":"macro ODEmodel\n\nMacro for creating an ODE from a list of equations. It also injects all variables into the global scope.\n\nExample\n\nCreating a simple ODE:\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = a * x1(t) + u(t),\n x2'(t) = b * x2(t) + c*x1(t)*x2(t),\n y(t) = x1(t)\n)\n\nHere,\n\nx1, x2 are state variables\ny is an output variable\nu is an input variable\na, b, c are time-indepdendent parameters\n\n\n\n\n\n","category":"macro"},{"location":"input/input/#StructuralIdentifiability.ODE","page":"Parsing input ODE system","title":"StructuralIdentifiability.ODE","text":"The main structure that represents input ODE system.\n\nStores information about states (x_vars), outputs (y_vars), inputs (u_vars), parameters (parameters) and the equations.\n\nThis structure is constructed via @ODEmodel macro.\n\n\n\n\n\n","category":"type"},{"location":"input/input/#StructuralIdentifiability.set_parameter_values","page":"Parsing input ODE system","title":"StructuralIdentifiability.set_parameter_values","text":"set_parameter_values(ode, param_values)\n\nInput:\n\node - an ODE as above\nparam_values - values for (possibly, some of) the parameters as dictionary parameter => value\n\nOutput:\n\nnew ode with the parameters in param_values plugged with the given numbers\n\n\n\n\n\n","category":"function"},{"location":"input/input/#Create-Compartmental-Model","page":"Parsing input ODE system","title":"Create Compartmental Model","text":"","category":"section"},{"location":"input/input/","page":"Parsing input ODE system","title":"Parsing input ODE system","text":"linear_compartment_model","category":"page"},{"location":"input/input/#StructuralIdentifiability.linear_compartment_model","page":"Parsing input ODE system","title":"StructuralIdentifiability.linear_compartment_model","text":"linear_compartment_model(graph, inputs, outputs, leaks)\n\nInput: defines a linear compartment model with nodes numbered from 1 to n by\n\ngraph - and array of integer arrays representing the adjacency lists of the graph\ninputs - array of input nodes\noutputs - array of output nodes\nleaks - array of sink nodes\n\nOutput:\n\nthe corresponding ODE system in the notation of https://doi.org/10.1007/s11538-015-0098-0\n\n\n\n\n\n","category":"function"}] +[{"location":"tutorials/reparametrization/#Reparametrizations","page":"Reparametrizations","title":"Reparametrizations","text":"","category":"section"},{"location":"tutorials/reparametrization/#Overview","page":"Reparametrizations","title":"Overview","text":"","category":"section"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Once one has found that not all parameters and/or states of the model at hand are identifiable, one natural desire is to reparametrize the model into a one with better identifiability properties. StructuralIdentifiability offers such a functionality via the function reparametrize_global. It takes as input an ODE model and produces its transformation into another model with the same input-output behaviour but with the states and parameters being globally identifiable. Note that, in general, such a transformation may not exist in the class of rational models, so sometimes the function returns an ODE not on the whole affine space but on a manifold.","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"More precisely, the function returns a dictionary with three keys:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":":new_vars is a dictionary which maps the new parameters and new states into the formulas expressing them in terms of the original parameters and states;\n:new_ode is the ODE satisfied by these new states (and the expression of the output in terms of the new states);\n:implicit_relations is a list of algebraic relations between the new states and parameters. Being nonempty, this is exactly the list of equations defining the manifold, on which the new ODE model is defined. In many interesting cases, however, this list is empty meaning that the new ODE is a standard rational ODE model.","category":"page"},{"location":"tutorials/reparametrization/#Example:-SEUIR-model","page":"Reparametrizations","title":"Example: SEUIR model","text":"","category":"section"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Consider a SEUIR epidemiological model from[1]:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"begincases\nS(t) = -beta frac(U(t) + I(t))S(t)N\nE(t) = beta frac(U(t) + I(t))S(t)N - gamma E(t)\nU(t) = (1 - alpha) gamma E(t) - delta U(t)\nI(t) = alpha gamma E(t) - delta I(t)\nR(t) = delta(U(t) + I(t))\ny(t) = I(t)\nendcases","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"In this model S is, as usually, the number of susceptible people, E is the number of people exposed to virus but not yet infected (as in a simple SEIR model[1]), and I and U correspond to number of infected people who report the infection and who do not, respectively. We define the model but omit R compartment since it does not affect the output dynamics:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n S'(t) = -b * (U(t) + I(t)) * S(t) / N,\n E'(t) = b * (U(t) + I(t)) * S(t) / N - g * E(t),\n U'(t) = (1 - a) * g * E(t) - d * U(t),\n I'(t) = a * g * E(t) - d * I(t),\n y(t) = I(t)\n)","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Majority of the states and parameters are not identifiable in this case:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"assess_identifiability(ode)","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Let us attempt to reparametrize the model, and print new variables:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"reparam = reparametrize_global(ode)\n@assert isempty(reparam[:implicit_relations]) # checking that the result is an ODE on the whole space, not on a manifold\nreparam[:new_vars]","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"In these new variables and parameters, the original ODE can be rewritten as follows:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"reparam[:new_ode]","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"In order to analyze this result, let us give more interpretable names to the new variables and parameters:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"I = I widetildeE = alpha E widetildeS = alpha widetildeI = alpha (I + U) gamma = gammadelta = deltawidetildebeta = fracbetaalpha N","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Then the reparametrize system becomes","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"begincases\nwidetildeS(t) = -widetildebeta widetildeS(t) widetildeI(t)\nwidetildeE(t) = widetildebeta widetildeS(t) widetildeI(t) - gamma widetildeE(t)\nwidetildeI(t) = -delta widetildeI(t) + gammawidetildeE(t)\nI(t) = gammawidetildeE(t) - delta I(t)\ny(t) = I(t)\nendcases","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"This reparametrization not only reduces the dimension of the parameter space from 5 to 3 but reveals interesting structural properties of the model:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"The first three equations form a self-contained model which is equivalent to a simple SEIR model, so the model gets \"decomposed\";\nNew variables widetildeS, widetildeE, widetildeI are obtained from S, E, and I by scaling by alpha which is the ratio of people who report being infected. One can interpret this as there is a part of population who would report infection and the other part who would not. Ultimately, we can model only the ones who would as this is mainly they who contribute to the output.","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"Finally, we can check that the new model is indeed globally identifiable:","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"assess_identifiability(reparam[:new_ode])","category":"page"},{"location":"tutorials/reparametrization/","page":"Reparametrizations","title":"Reparametrizations","text":"[1]: T. Sauer, T. Berry, D. Ebeigbe, M. Norton, A. Whalen, S. Schiff, Identifiability of infection model parameters early in an epidemic, SIAM Journal on Control and Optimization, 2022;","category":"page"},{"location":"utils/ode/#Functions-to-work-with-the-ODE-structure","page":"ODE Tools","title":"Functions to work with the ODE structure","text":"","category":"section"},{"location":"utils/ode/","page":"ODE Tools","title":"ODE Tools","text":"Modules = [StructuralIdentifiability]\nPages = [\"ODE.jl\", \"submodels.jl\"]\nOrder = [:function]","category":"page"},{"location":"utils/ode/#StructuralIdentifiability._reduce_mod_p-Tuple{Nemo.QQMPolyRingElem, Int64}","page":"ODE Tools","title":"StructuralIdentifiability._reduce_mod_p","text":"_reduce_mod_p(f, p)\n\nReduces a polynomial/rational function over Q modulo p\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.power_series_solution-Union{Tuple{P}, Tuple{T}, Tuple{ODE{P}, Dict{P, T}, Dict{P, T}, Dict{P, Vector{T}}, Int64}} where {T<:AbstractAlgebra.FieldElem, P<:AbstractAlgebra.MPolyRingElem{T}}","page":"ODE Tools","title":"StructuralIdentifiability.power_series_solution","text":"power_series_solution(ode, param_values, initial_conditions, input_values, prec)\n\nInput:\n\node - an ode to solve\nparam_values - parameter values, must be a dictionary mapping parameter to a value\ninitial_conditions - initial conditions of ode, must be a dictionary mapping state variable to a value\ninput_values - power series for the inputs presented as a dictionary variable => list of coefficients\nprec - the precision of solutions\n\nOutput:\n\ncomputes a power series solution with precision prec presented as a dictionary variable => corresponding coordinate of the solution\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.reduce_ode_mod_p-Tuple{ODE{<:AbstractAlgebra.MPolyRingElem{Nemo.QQFieldElem}}, Int64}","page":"ODE Tools","title":"StructuralIdentifiability.reduce_ode_mod_p","text":"reduce_ode_mod_p(ode, p)\n\nInput: ode is an ODE over QQ, p is a prime number Output: the reduction mod p, throws an exception if p divides one of the denominators\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.set_parameter_values-Union{Tuple{P}, Tuple{T}, Tuple{ODE{P}, Dict{P, T}}} where {T<:AbstractAlgebra.FieldElem, P<:AbstractAlgebra.MPolyRingElem{T}}","page":"ODE Tools","title":"StructuralIdentifiability.set_parameter_values","text":"set_parameter_values(ode, param_values)\n\nInput:\n\node - an ODE as above\nparam_values - values for (possibly, some of) the parameters as dictionary parameter => value\n\nOutput:\n\nnew ode with the parameters in param_values plugged with the given numbers\n\n\n\n\n\n","category":"method"},{"location":"utils/ode/#StructuralIdentifiability.find_submodels-Union{Tuple{ODE{P}}, Tuple{P}} where P<:AbstractAlgebra.MPolyRingElem","page":"ODE Tools","title":"StructuralIdentifiability.find_submodels","text":"find_submodels(ode)\n\nThe function calculates and returns all valid submodels given a system of ODEs.\n\nInput:\n\node - an ODEs system to be studied\n\nOutput: \n\nA list of submodels represented as ode objects\n\nExample:\n\n>ode = @ODEmodel(x1'(t) = x1(t)^2, \n x2'(t) = x1(t) * x2(t), \n y1(t) = x1(t), \n y2(t) = x2(t))\n>find_submodels(ode)\n ODE{QQMPolyRingElem}[\n \n x1'(t) = a(t)*x2(t)^2 + x1(t)\n y1(t) = x1(t)\n ]\n\n\n\n\n\n","category":"method"},{"location":"utils/primality/#Primality-Checks","page":"Primality Checks","title":"Primality Checks","text":"","category":"section"},{"location":"utils/primality/","page":"Primality Checks","title":"Primality Checks","text":"Pages=[\"primality.md\"]","category":"page"},{"location":"utils/primality/","page":"Primality Checks","title":"Primality Checks","text":"StructuralIdentifiability.check_primality","category":"page"},{"location":"utils/primality/#StructuralIdentifiability.check_primality","page":"Primality Checks","title":"StructuralIdentifiability.check_primality","text":"check_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem}, extra_relations::Array{QQMPolyRingElem, 1})\n\nThe function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.\n\nThe extra_relations allows adding more polynomials to the generators (not affecting the saturation).\n\n\n\n\n\ncheck_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem})\n\nThe function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.\n\n\n\n\n\n","category":"function"},{"location":"tutorials/identifiable_functions/#Globally-Identifiable-Functions","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"","category":"section"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"In addition to assessing identifiabuility of given functions of parameters or states, StructuralIdentifiability.jl provides the function find_identifiable_functions which finds a set of identifiable functions such that any other identifiable function can be expressed using them. This allows to find out what actually is identifiable and what does non-identifiability in the model at hand looks like.","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"For example, consider the following model[1].","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"using StructuralIdentifiability # hide\nLLW1987 = @ODEmodel(\n x1'(t) = -p1 * x1(t) + p2 * u(t),\n x2'(t) = -p3 * x2(t) + p4 * u(t),\n x3'(t) = -(p1 + p3) * x3(t) + (p4 * x1(t) + p2 * x2(t)) * u(t),\n y1(t) = x3(t)\n)","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"Several decades ago, this model was introduced to demonstrate the presence of nontrivial unidentifiability in nonlinear systems of ODEs. Nowadays, we can automatically find the identifiable combinations of parameters:","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"using StructuralIdentifiability # hide\nfind_identifiable_functions(LLW1987)","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"From these expressions, we see that the values of p1 and p3 are not identifiable but an unordered pair of numbers {p1, p3} is uniquely determined since p1 + p3 and p1 * p3 are known. Furthermore, we see that, for fixed input and output, p2 and p4 can take infinitely many values but knowing one of them, we would also be able to determine the other.","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"Moreover, we can find generators of all identifiable functions in parameters and states:","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"find_identifiable_functions(LLW1987, with_states = true)","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"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 :strong to try to simplify further (at the expense of heavier computation) or use :weak to simplify less (but compute faster).","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"As assess_identifiability and assess_local_identifiability, find_identifiable_functions accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.","category":"page"},{"location":"tutorials/identifiable_functions/","page":"Globally Identifiable Functions","title":"Globally Identifiable Functions","text":"[1]: Y. Lecourtier, F. Lamnabhi-Lagarrigue, and E. Walter, A method to prove that nonlinear models can be unidentifiable, Proceedings of the 26th Conference on Decision and Control, December 1987, 2144-2145;","category":"page"},{"location":"tutorials/identifiability/#Identifiability-of-Differential-Models-(Local-and-Global)","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"","category":"section"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Recall that we consider ODE models in the state-space form","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"begincases\nmathbfx(t) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"where mathbfx(t) mathbfy(t), and mathbfu(t) are time-dependent states, outputs, and inputs, respectively, and mathbfp are scalar parameters. We will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs. Typically, two types of identifiability are distinguished","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"local identifiability: the value can be recovered up to finitely many options;\nglobal identifiability: the value can be recovered uniquely.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Note that in the case of states, term observability it typically used. We will use identifiability for both states and parameters for brevity and uniformity. While the notion of global identifiability is much more precise, assessing local identifiability is typically much faster, and can be performed for the models whose global identifiability analysis is out of reach.","category":"page"},{"location":"tutorials/identifiability/#Local-identifiability","page":"Identifiability of Differential Models (Local and Global)","title":"Local identifiability","text":"","category":"section"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"We consider the Lotka-Volterra model:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"begincases\nx_1(t) = a x_1(t) - b x_1(t) x_2(t) + u(t)\nx_2(t) = -c x_2(t) + d x_1(t) x_2(t)\ny(t) = x_1(t)\nendcases","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"The local identifiability of all parameters and states in this model can be assessed as follows","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = a * x1(t) - b * x1(t) * x2(t) + u(t),\n x2'(t) = -c * x2(t) + d * x1(t) * x2(t),\n y(t) = x1(t)\n)\n\nassess_local_identifiability(ode)","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"We see that x_1(t) is locally identifiable (no surprises, this is an output), a c and d are identifiable as well. On the other hand, x_2(t) and b are nonidentifiable. This can be explained by the following scaling symmetry","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"x_2(t) to lambda x_2(t) quad b to fracblambda","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"which preserves input and output for every nonzero lambda. The algorithm behind this call is the one due to Sedoglavic[1].","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Function assess_local_identifiability has several optional parameters","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"funcs_to_check a list of specific functions of parameters and states to check identifiability for (see an example below). If not provided, the identifiability is assessed for all parameters and states.\nprob_threshold (default 099, i.e. 99%) is the probability of correctness. The algorithm can, in theory, produce wrong result, but the probability that it is correct is guaranteed to be at least prob_threshold. However, the probability bounds we use are quite conservative, so the actual probability of correctness is likely to be much higher.\ntype (default :SE). By default, the algorithm checks the standard single-experiment identifiability. If one sets type = :ME, then the algorithm checks multi-experiment identifiability, that is, identifiability from several experiments with independent initial conditions (the algorithm from [2] is used).\nloglevel (default Logging.Info). The minimal level of logging messages to be displayed. Available options: Logging.Debug, Logging.Info, Logging.Warn, and Logging.Error.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Note that the scaling symmetry given above suggests that b x_2(t) may in fact be identifiable. This can be checked using funcs_to_check parameter:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"assess_local_identifiability(ode, funcs_to_check = [b * x2])","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Indeed!","category":"page"},{"location":"tutorials/identifiability/#Global-identifiability","page":"Identifiability of Differential Models (Local and Global)","title":"Global identifiability","text":"","category":"section"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"One can obtain more refined information about a model using assess_identifiability function. We will showcase it using the Goodwin oscillator model [3].","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = -b * x1(t) + 1 / (c + x4(t)),\n x2'(t) = alpha * x1(t) - beta * x2(t),\n x3'(t) = gama * x2(t) - delta * x3(t),\n x4'(t) = sigma * x4(t) * (gama * x2(t) - delta * x3(t)) / x3(t),\n y(t) = x1(t)\n)\n\nassess_identifiability(ode)","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"As a result, each parameter/state is assigned one of the labels :globally (globally identifiable), :locally (locally but not globally identifiable), or :nonidentifiable (not identifiable, even locally). The algorithm behind this computation follows [4].","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Similarly to assess_local_identifiability, this function has optional parameters:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"funcs_to_check a list of specific functions of parameters and states to check identifiability for (see an example below). If not provided, the identifiability is assessed for all parameters and states. Note that the computations for states may be more involved than for the parameters, so one may want to call the function with funcs_to_check = ode.parameters if the call assess_identifiability(ode) takes too long.\nprob_threshold (default 099, i.e. 99%) is the probability of correctness. Same story as above: the probability estimates are very conservative, so the actual error probability is much lower than 1%. Also, currently, the probability of correctness does not include the probability of correctness of the modular reconstruction for Groebner bases. This probability is ensured by an additional check modulo a large prime, and can be neglected for practical purposes.\nloglevel (default Logging.Info). The minimal level of logging messages to be displayed. Available options: Logging.Debug, Logging.Info, Logging.Warn, and Logging.Error.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"Using funcs_to_check parameter, one can further inverstigate the nature of the lack of identifiability in the model at hand. For example, for the Goodwin oscillator, we can check if beta + delta and beta * delta are identifiable:","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"assess_identifiability(ode, funcs_to_check = [beta + delta, beta * delta])","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"And we see that they indeed are. This means, in particular, that the reason why beta and delta are not identifiable is because their values can be exchanged. One may wonder how could we guess these functions beta + delta, beta * delta. In fact, they can be just computed using find_identifiable_functions function as we will explain in the next tutorial. Stay tuned!","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[1]: A. Sedoglavic, A probabilistic algorithm to test local algebraic observability in polynomial time, Journal of Symbolic Computation, 2002.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[2]: A. Ovchinnikov, A. Pillay, G. Pogudin, T. Scanlon, Multi-experiment Parameter Identifiability of ODEs and Model Theory, SIAM Journal on Applied Algebra and Geometry, 2022.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[3]: D. Gonze, P. Ruoff, The Goodwin Oscillator and its Legacy, Acta Biotheoretica, 2020.","category":"page"},{"location":"tutorials/identifiability/","page":"Identifiability of Differential Models (Local and Global)","title":"Identifiability of Differential Models (Local and Global)","text":"[4]: R. Dong, C. Goodbrake, H. Harrington, G. Pogudin, Differential elimination for dynamical models via projections with applications to structural identifiability, SIAM Journal on Applied Algebra and Geometry, 2023.","category":"page"},{"location":"export/export/#Exporting-to-Other-Systems","page":"Exporting to Other Systems","title":"Exporting to Other Systems","text":"","category":"section"},{"location":"export/export/","page":"Exporting to Other Systems","title":"Exporting to Other Systems","text":"Here we put some helpful utilities to export your code to other identifiability software.","category":"page"},{"location":"export/export/","page":"Exporting to Other Systems","title":"Exporting to Other Systems","text":"print_for_maple\nprint_for_DAISY\nprint_for_GenSSI\nprint_for_COMBOS","category":"page"},{"location":"export/export/#StructuralIdentifiability.print_for_maple","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_maple","text":"print_for_maple(ode, package)\n\nPrints the ODE in the format accepted by maple packages\n\nSIAN (https://github.com/pogudingleb/SIAN) if package=:SIAN\nDifferentialAlgebra if package=:DifferentialAlgebra\nDifferentialThomas if package=:DifferentialThomas\n\n\n\n\n\n","category":"function"},{"location":"export/export/#StructuralIdentifiability.print_for_DAISY","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_DAISY","text":"print_for_DAISY(ode)\n\nPrints the ODE in the format accepted by DAISY (https://daisy.dei.unipd.it/)\n\n\n\n\n\n","category":"function"},{"location":"export/export/#StructuralIdentifiability.print_for_GenSSI","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_GenSSI","text":"print_for_GenSSI(ode)\n\nPrints the ODE in the format accepted by GenSSI 2.0 (https://github.com/genssi-developer/GenSSI)\n\n\n\n\n\n","category":"function"},{"location":"export/export/#StructuralIdentifiability.print_for_COMBOS","page":"Exporting to Other Systems","title":"StructuralIdentifiability.print_for_COMBOS","text":"print_for_COMBOS(ode)\n\nPrints the ODE in the format accepted by COMBOS (http://biocyb1.cs.ucla.edu/combos/)\n\n\n\n\n\n","category":"function"},{"location":"utils/util/#Other-Helpful-Functions","page":"Other Utilities","title":"Other Helpful Functions","text":"","category":"section"},{"location":"utils/util/","page":"Other Utilities","title":"Other Utilities","text":"Pages=[\"util.md\"]","category":"page"},{"location":"utils/util/","page":"Other Utilities","title":"Other Utilities","text":"Modules = [StructuralIdentifiability]\nPages = [\"util.jl\"]","category":"page"},{"location":"utils/util/#StructuralIdentifiability.compare_rational_func_by-Union{Tuple{By}, Tuple{Any, Any, By}, Tuple{Any, Any, By, Any}} where By","page":"Other Utilities","title":"StructuralIdentifiability.compare_rational_func_by","text":"compare_rational_func_by(f, g, by)\n\nReturns \n\n-1 if f < g,\n0 if f = g, and \n1 if f > g.\n\nFunctions' numerators and denominators are compared using by.\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.decompose_derivative-Tuple{String, Array{String}}","page":"Other Utilities","title":"StructuralIdentifiability.decompose_derivative","text":"decompose_derivative(varname, prefixes)\n\nDetermines if it is possible to represent the varname as a_number where a is an element of prefixes If yes, returns a pair (a, number), otherwise nothing\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.dennums_to_fractions-Union{Tuple{Array{Vector{T}, 1}}, Tuple{T}} where T","page":"Other Utilities","title":"StructuralIdentifiability.dennums_to_fractions","text":"dennums_to_fractions(dennums)\n\nReturns the field generators represented by fractions.\n\nInput: an array of arrays of polynomials, as in [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]\n\nOutput: an array of fractions [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.eval_at_dict-Union{Tuple{P}, Tuple{P, Dict{P, <:AbstractAlgebra.RingElem}}} where P<:AbstractAlgebra.MPolyRingElem","page":"Other Utilities","title":"StructuralIdentifiability.eval_at_dict","text":"eval_at_dict(f, d)\n\nEvaluates a polynomial/rational function on a dictionary of type var => val and missing values are replaced with zeroes\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.extract_coefficients-Union{Tuple{P}, Tuple{P, Vector{P}}} where P<:AbstractAlgebra.MPolyRingElem","page":"Other Utilities","title":"StructuralIdentifiability.extract_coefficients","text":"extract_coefficients(poly, variables)\n\nInput:\n\npoly - multivariate polynomial\nvariables - a list of variables from the generators of the ring of p\n\nOutput:\n\ndictionary with keys being tuples of length lenght(variables) and values being polynomials in the variables other than those which are the coefficients at the corresponding monomials (in a smaller polynomial ring)\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.fractions_to_dennums-Tuple{Any}","page":"Other Utilities","title":"StructuralIdentifiability.fractions_to_dennums","text":"fractions_to_dennums(fractions)\n\nReturns the field generators represented by lists of denominators and numerators.\n\nInput: an array of fractions, as in [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]\n\nOutput: an array of arrays of polynomials, [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.gen_tag_name","page":"Other Utilities","title":"StructuralIdentifiability.gen_tag_name","text":"gen_tag_name(base; stop_words)\ngen_tag_names(n, base; stop_words)\n\nGenerates a string which will not collide with the words in stop_words.\n\nArguments\n\nn: Generates a sequence of unique strings of length n\nbase: A string or a vector of strings, the base for the generated sequence\nstop_words: A vector of strings, stop words\n\n\n\n\n\n","category":"function"},{"location":"utils/util/#StructuralIdentifiability.make_substitution-Union{Tuple{P}, NTuple{4, P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Other Utilities","title":"StructuralIdentifiability.make_substitution","text":"make_substitution(f, var_sub, val_numer, val_denom)\n\nSubstitute a variable in a polynomial with an expression\n\nInput:\n\nf - the polynomial\nvar_sub - the variable to be substituted\nvar_numer - numerator of the substitution expression\nvar_denom - denominator of the substitution expression\n\nOutput:\n\npolynomial - result of substitution\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.parent_ring_change-Tuple{AbstractAlgebra.MPolyRingElem, AbstractAlgebra.MPolyRing}","page":"Other Utilities","title":"StructuralIdentifiability.parent_ring_change","text":"parent_ring_change(poly, new_ring)\n\nConverts a polynomial to a different polynomial ring Input\n\npoly - a polynomial to be converted\nnew_ring - a polynomial ring such that every variable name appearing in poly appears among the generators\n\nOutput:\n\na polynomial in new_ring “equal” to poly\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.switch_ring-Tuple{AbstractAlgebra.MPolyRingElem, AbstractAlgebra.MPolyRing}","page":"Other Utilities","title":"StructuralIdentifiability.switch_ring","text":"switch_ring(v, ring)\n\nFor a variable v, returns a variable in ring with the same name\n\n\n\n\n\n","category":"method"},{"location":"utils/util/#StructuralIdentifiability.uncertain_factorization-Tuple{AbstractAlgebra.MPolyRingElem{Nemo.QQFieldElem}}","page":"Other Utilities","title":"StructuralIdentifiability.uncertain_factorization","text":"uncertain_factorization(f)\n\nInput:\n\nf - polynomial with rational coefficients\n\nOutput:\n\nlist of pairs (div, certainty) where\ndiv's are divisors of f such that f is their product with certain powers\nif certainty is true, div is Q-irreducible\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#Wronskian-Tools","page":"Wronskian Tools","title":"Wronskian Tools","text":"","category":"section"},{"location":"utils/wronskian/","page":"Wronskian Tools","title":"Wronskian Tools","text":"Modules = [StructuralIdentifiability]\nPages = [\"wronskian.jl\"]","category":"page"},{"location":"utils/wronskian/#StructuralIdentifiability.get_max_below-Tuple{StructuralIdentifiability.ExpVectTrie, Vector{Int64}}","page":"Wronskian Tools","title":"StructuralIdentifiability.get_max_below","text":"get_max_below(t, vect)\n\nInput:\n\nt - a trie with exponent vectors\nvect - yet another exponent vector\n\nOutput:\n\na pair (d, v) where v is a vector in the trie which is componentwise ≤ vect and the difference d is as small as possible\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#StructuralIdentifiability.massive_eval-Tuple{Any, Any}","page":"Wronskian Tools","title":"StructuralIdentifiability.massive_eval","text":"massive_eval(polys, eval_dict)\n\nInput:\n\npolys - a list of polynomials\neval_dict - dictionary from variables to the values. Missing values are treated as zeroes\n\nOutput:\n\na list of values of the polynomials\n\nEvaluates a list of polynomials at a point. Assumes that multiplications are relatively expensive (like in truncated power series) so all the monomials are precomputed first and the values of monomials of lower degree are cached and used to compute the values of the monomials of higher degree\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#StructuralIdentifiability.monomial_compress-Tuple{Any, ODE}","page":"Wronskian Tools","title":"StructuralIdentifiability.monomial_compress","text":"monomial_compress(io_equation, ode)\n\nCompresses an input-output equation for the rank computation Input:\n\nio_equation - input-output equation\node - the corresponding ODE model\n\nOutput:\n\npair (coeffs, terms) such that:\nsum of coeffs[i] * terms[i] = io_equation\ncoeffs involve only parameters, terms involve only inputs and outputs\nlength of the representation is the smallest possible\n\n\n\n\n\n","category":"method"},{"location":"utils/wronskian/#StructuralIdentifiability.wronskian-Union{Tuple{P}, Tuple{Dict{P, P}, ODE{P}}} where P<:AbstractAlgebra.MPolyRingElem","page":"Wronskian Tools","title":"StructuralIdentifiability.wronskian","text":"wronskian(io_equations, ode)\n\nInput:\n\nio_equations - a set of io-equations in the form of the Dict as returned by find_ioequations\node - the ODE object\n\nOutput:\n\na list of Wronskians evaluated at a point modulo prime\n\nComputes the Wronskians of io_equations\n\n\n\n\n\n","category":"method"},{"location":"utils/global_identifiability/#Global-Identifiability-Tools","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"","category":"section"},{"location":"utils/global_identifiability/","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"Pages=[\"global_identifiability.md\"]","category":"page"},{"location":"utils/global_identifiability/","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"CurrentModule=StructuralIdentifiability","category":"page"},{"location":"utils/global_identifiability/","page":"Global Identifiability Tools","title":"Global Identifiability Tools","text":"StructuralIdentifiability.RationalFunctionField\nStructuralIdentifiability.field_contains\nStructuralIdentifiability.get_degree_and_coeffsize","category":"page"},{"location":"utils/global_identifiability/#StructuralIdentifiability.RationalFunctionField","page":"Global Identifiability Tools","title":"StructuralIdentifiability.RationalFunctionField","text":"RationalFunctionField\n\nA subfield of the field of rational functions over the rationals.\n\nExample\n\nusing Nemo\nusing StructuralIdentifiability: RationalFunctionField\n\nR, (x, y, z) = QQ[\"x\", \"y\", \"z\"]\n\n# Constructs a subfield generated by x / y, y / z\nrff = RationalFunctionField([x // y, y // z])\n\n# Constructs a subfield generated by y / x, 1 / x, z / y\nrff = RationalFunctionField([[x, y, R(1)], [y, z]])\n\n\n\n\n\n","category":"type"},{"location":"utils/global_identifiability/#StructuralIdentifiability.field_contains","page":"Global Identifiability Tools","title":"StructuralIdentifiability.field_contains","text":"field_contains(field, ratfuncs, prob_threshold)\n\nChecks whether given rational function field field contains given rational functions ratfuncs (represented as a list of lists). The result is correct with probability at least prob_threshold\n\nInputs:\n\nfield - a rational function field\nratfuncs - a list of lists of polynomials. Each of the lists, say, [f1, ..., fn], defines generators f2/f1, ..., fn/f1.\nprob_threshold real number from (0, 1)\n\nOutput:\n\na list L[i] of bools of length length(rat_funcs) such that L[i] is true iff the i-th function belongs to field\n\n\n\n\n\n","category":"function"},{"location":"utils/global_identifiability/#StructuralIdentifiability.get_degree_and_coeffsize","page":"Global Identifiability Tools","title":"StructuralIdentifiability.get_degree_and_coeffsize","text":"get_degree_and_coeffsize(f)\n\nfor f being a polynomial/rational function over rationals (QQ) returns a tuple (degree, max_coef_size)\n\n\n\n\n\n","category":"function"},{"location":"utils/reparametrization/#Model-reparametrization","page":"Model reparametrization","title":"Model reparametrization","text":"","category":"section"},{"location":"utils/reparametrization/","page":"Model reparametrization","title":"Model reparametrization","text":"StructuralIdentifiability.reparametrize_global","category":"page"},{"location":"utils/reparametrization/#StructuralIdentifiability.reparametrize_global","page":"Model reparametrization","title":"StructuralIdentifiability.reparametrize_global","text":"reparametrize_global(ode, options...)\n\nFinds a reparametrization of ode in terms of globally identifiabile functions.\n\nReturns a tuple (new_ode, new_vars, implicit_relations), such that:\n\nnew_ode is the reparametrized ODE system\nnew_vars is a mapping from the new variables to the original ones\nrelations is the array of implicit algebraic relations among new_vars. Usually, relations is empty.\n\nOptions\n\nThe function accepts the following optional arguments.\n\nseed: A float in the range from 0 to 1, random seed (default is seed = 42). \nprob_threshold: The probability of correctness (default is prob_threshold = 0.99).\n\nExample\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = a * x1(t) - b*x1(t)*x2(t),\n x2'(t) = -c * x2(t) + d*x1(t)*x2(t),\n y(t) = x1(t)\n)\n\nnew_ode, new_vars, relations = reparametrize_global(ode)\n\nThen, we have the following:\n\n# new_ode\nX2'(t) = X1(t)*X2(t)*a2 - X2(t)*a1\nX1'(t) = -X1(t)*X2(t) + X1(t)*a3\ny1(t) = X1(t)\n\n# new_vars\nDict{Nemo.QQMPolyRingElem, AbstractAlgebra.Generic.Frac{Nemo.QQMPolyRingElem}} with 6 entries:\n X2 => b*x2\n y1 => y\n X1 => x1\n a2 => d\n a3 => a\n a1 => c\n\nNotice that the new_ode is fully identifiabile, and has 1 less parameter compared to the original one.\n\n\n\n\n\n","category":"function"},{"location":"utils/power_series_utils/#Power-Series-Utilities","page":"Power Series Tools","title":"Power Series Utilities","text":"","category":"section"},{"location":"utils/power_series_utils/","page":"Power Series Tools","title":"Power Series Tools","text":"Pages =[\"power_series_utils.md\"]","category":"page"},{"location":"utils/power_series_utils/","page":"Power Series Tools","title":"Power Series Tools","text":"Modules = [StructuralIdentifiability]\nPages = [\"power_series_utils.jl\"]","category":"page"},{"location":"utils/power_series_utils/#StructuralIdentifiability._matrix_inv_newton_iteration-Union{Tuple{T}, Tuple{AbstractAlgebra.MatElem{T}, AbstractAlgebra.MatElem{T}}} where T<:(AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.FieldElem})","page":"Power Series Tools","title":"StructuralIdentifiability._matrix_inv_newton_iteration","text":"_matrix_inv_newton_iteration(M, Minv)\n\nPerforms a single step of Newton iteration for inverting M with Minv being a partial result\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_diff-Tuple{AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.RingElem}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_diff","text":"ps_diff(ps)\n\nInput:\n\nps - (absolute capped) univariate power series\n\nOutput:\n\nthe derivative of ps\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_integrate-Tuple{AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.FieldElem}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_integrate","text":"ps_integrate(ps)\n\nInput:\n\nps - (absolute capped) univariate power series\n\nOutput:\n\nthe integral of ps without constant term\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_homlinear_de-Union{Tuple{T}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}, Int64}} where T<:AbstractAlgebra.FieldElem","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_homlinear_de","text":"ps_matrix_homlinear_de(A, Y0, prec)\n\nInput:\n\nA - a square matrix with entries in a univariate power series ring\nY0 - a square invertible matrix over the base field\n\nOutput:\n\nmatrix Y such that Y' = AY up to precision of A - 1 and Y(0) = Y0\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_inv","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_inv","text":"ps_matrix_inv(M, prec)\n\nInput:\n\nM - a square matrix with entries in a univariate power series ring it is assumed that M(0) is invertible and all entries having the same precision\nprec - an integer, precision, if -1 then defaults to precision of M\n\nOutput:\n\nthe inverse of M computed up to prec\n\n\n\n\n\n","category":"function"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_linear_de-Union{Tuple{T}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}}, Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{T}}, AbstractAlgebra.MatElem{<:T}, Int64}} where T<:AbstractAlgebra.FieldElem","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_linear_de","text":"ps_matrix_linear_de(A, B, Y0, prec)\n\nInput:\n\nA, B - square matrices with entries in a univariate power series ring\nY0 - a matrix over the base field with the rows number the same as A\n\nOutput:\n\nmatrix Y such that Y' = AY + B up to precision of A - 1 and Y(0) = Y0\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_matrix_log-Tuple{AbstractAlgebra.MatElem{<:AbstractAlgebra.AbsPowerSeriesRingElem{<:AbstractAlgebra.FieldElem}}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_matrix_log","text":"ps_matrix_log(M)\n\nInput:\n\nM - a square matrix with entries in a univariate power series ring it is assumed that M(0) is the identity\n\nOutput:\n\nthe natural log of M\n\n\n\n\n\n","category":"method"},{"location":"utils/power_series_utils/#StructuralIdentifiability.ps_ode_solution-Union{Tuple{P}, Tuple{T}, Tuple{Vector{P}, Dict{P, T}, Dict{P, Vector{T}}, Int64}} where {T<:AbstractAlgebra.FieldElem, P<:AbstractAlgebra.MPolyRingElem{T}}","page":"Power Series Tools","title":"StructuralIdentifiability.ps_ode_solution","text":"ps_ode_solution(equations, ic, inputs, prec)\n\nInput:\n\nequations - a system of the form A(x u mu)x - B(x u mu) = 0, where A is a generically nonsingular square matrix. Assumption: A is nonzero at zero\nic - initial conditions for x's (dictionary)\ninputs - power series for inputs represented as arrays (dictionary)\nprec - precision of the solution\n\nOutput:\n\npower series solution of the system\n\n\n\n\n\n","category":"method"},{"location":"identifiability/identifiability/#Functions-to-Assess-Identifiability","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"","category":"section"},{"location":"identifiability/identifiability/#Assessing-All-Types-of-Identifiability","page":"Functions to Assess Identifiability","title":"Assessing All Types of Identifiability","text":"","category":"section"},{"location":"identifiability/identifiability/","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"assess_identifiability","category":"page"},{"location":"identifiability/identifiability/#StructuralIdentifiability.assess_identifiability","page":"Functions to Assess Identifiability","title":"StructuralIdentifiability.assess_identifiability","text":"assess_identifiability(ode; funcs_to_check = [], prob_threshold=0.99, loglevel=Logging.Info)\n\nInput:\n\node - the ODE model\nfuncs_to_check - list of functions to check identifiability for; if empty, all parameters and states are taken\nprob_threshold - probability of correctness.\nloglevel - the minimal level of log messages to display (Logging.Info by default)\n\nAssesses identifiability of a given ODE model. The result is guaranteed to be correct with the probability at least prob_threshold. The function returns an (ordered) dictionary from the functions to check to their identifiability properties (one of :nonidentifiable, :locally, :globally).\n\n\n\n\n\n","category":"function"},{"location":"identifiability/identifiability/#Assessing-Local-Identifiability","page":"Functions to Assess Identifiability","title":"Assessing Local Identifiability","text":"","category":"section"},{"location":"identifiability/identifiability/","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"assess_local_identifiability","category":"page"},{"location":"identifiability/identifiability/#StructuralIdentifiability.assess_local_identifiability","page":"Functions to Assess Identifiability","title":"StructuralIdentifiability.assess_local_identifiability","text":"assess_local_identifiability(ode::ODE{P}; funcs_to_check::Array{<: Any, 1}, prob_threshold::Float64=0.99, type=:SE, loglevel=Logging.Info) where P <: MPolyRingElem{Nemo.QQFieldElem}\n\nChecks the local identifiability/observability of the functions in funcs_to_check. The result is correct with probability at least prob_threshold.\n\nCall this function if you have a specific collection of parameters of which you would like to check local identifiability.\n\ntype can be either :SE (single-experiment identifiability) or :ME (multi-experiment identifiability). If the type is :ME, states are not allowed to appear in the funcs_to_check.\n\n\n\n\n\nassess_local_identifiability(dds::DDS{P}; funcs_to_check::Array{<: Any, 1}, known_ic, prob_threshold::Float64=0.99, loglevel=Logging.Info) where P <: MPolyRingElem{Nemo.QQFieldElem}\n\nChecks the local identifiability/observability of the functions in funcs_to_check. The result is correct with probability at least prob_threshold. A list of quantities can be provided as known_ic for which the initial conditions can be assumed to be known and generic.\n\n\n\n\n\n","category":"function"},{"location":"identifiability/identifiability/#Finding-Identifiable-Functions","page":"Functions to Assess Identifiability","title":"Finding Identifiable Functions","text":"","category":"section"},{"location":"identifiability/identifiability/","page":"Functions to Assess Identifiability","title":"Functions to Assess Identifiability","text":"find_identifiable_functions","category":"page"},{"location":"identifiability/identifiability/#StructuralIdentifiability.find_identifiable_functions","page":"Functions to Assess Identifiability","title":"StructuralIdentifiability.find_identifiable_functions","text":"find_identifiable_functions(ode::ODE; options...)\n\nFinds all functions of parameters/states that are identifiable in the given ODE system.\n\nOptions\n\nThis functions takes the following optional arguments:\n\nwith_states: When true, also reports the identifiabile functions in the ODE states. Default is false.\nsimplify: The extent to which the output functions are simplified. Stronger simplification may require more time. Possible options are:\n:standard: Default simplification.\n:weak: Weak simplification. This option is the fastest, but the output functions can be quite complex.\n:strong: Strong simplification. This option is the slowest, but the output\nfunctions are nice and simple.\n:absent: No simplification.\nprob_threshold: A float in the range from 0 to 1, the probability of correctness. Default is 0.99.\nseed: The rng seed. Default value is 42.\nloglevel - the minimal level of log messages to display (Logging.Info by default)\n\nExample\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x0'(t) = -(a01 + a21) * x0(t) + a12 * x1(t),\n x1'(t) = a21 * x0(t) - a12 * x1(t),\n y(t) = x0(t)\n)\n\nfind_identifiable_functions(ode)\n\n# prints\n3-element Vector{AbstractAlgebra.Generic.Frac{Nemo.QQMPolyRingElem}}:\n a12 + a01 + a21\n a12*a01\n\n\n\n\n\n","category":"function"},{"location":"utils/elimination/#Elimination","page":"Elimination","title":"Elimination","text":"","category":"section"},{"location":"utils/elimination/","page":"Elimination","title":"Elimination","text":"Pages=[\"elimination.md\"]","category":"page"},{"location":"utils/elimination/","page":"Elimination","title":"Elimination","text":"Modules = [StructuralIdentifiability]\nPages = [\"elimination.jl\"]","category":"page"},{"location":"utils/elimination/#StructuralIdentifiability.Bezout_matrix-Union{Tuple{P}, Tuple{P, P, P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Elimination","title":"StructuralIdentifiability.Bezout_matrix","text":"Bezout_matrix(f, g, var_elim)\n\nCompute the Bezout matrix of two polynomials f, g with respect to var_elim\n\nInputs:\n\nf - first polynomial\ng - second polynomial\nvar_elim - variable, of which f and g are considered as polynomials\n\nOutput:\n\nM::MatrixElem - The Bezout matrix\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.Sylvester_matrix-Union{Tuple{P}, Tuple{P, P, P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Elimination","title":"StructuralIdentifiability.Sylvester_matrix","text":"Sylvester_matrix(f, g, var_elim)\n\nCompute the Sylvester matrix of two polynomials f, g with respect to var_elim Inputs:\n\nf - first polynomial\ng - second polynomial\nvar_elim - variable, of which f and g are considered as polynomials\n\nOutput:\n\nM::MatrixElem - The Sylvester matrix\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.choose-Union{Tuple{P}, Tuple{Vector{P}, Any}} where P<:(AbstractAlgebra.MPolyRingElem{<:AbstractAlgebra.FieldElem})","page":"Elimination","title":"StructuralIdentifiability.choose","text":"choose(polys, generic_point_generator)\n\nInput:\n\npolys - an array of distinct irreducible polynomials in the same ring\ngeneric_point_generator - a generic point generator as described above for one of polys\n\nOutput:\n\nthe polynomial that vanishes at the generic_point_generator\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.eliminate_var-Union{Tuple{P}, Tuple{P, P, P, Any}} where P<:(AbstractAlgebra.MPolyRingElem{<:AbstractAlgebra.FieldElem})","page":"Elimination","title":"StructuralIdentifiability.eliminate_var","text":"eliminate_var(f, g, var_elim, generic_point_generator)\n\nEliminate a variable from a pair of polynomials\n\nInput:\n\nf and g - polynomials\nvar_elim - variable to be eliminated\ngeneric_point_generator - a generic point generator object for the factor of the resultant of f and g of interest\n\nOutput:\n\npolynomial - the desired factor of the resultant of f and g\n\n\n\n\n\n","category":"method"},{"location":"utils/elimination/#StructuralIdentifiability.simplify_matrix-Union{Tuple{AbstractAlgebra.MatElem{P}}, Tuple{P}} where P<:AbstractAlgebra.MPolyRingElem","page":"Elimination","title":"StructuralIdentifiability.simplify_matrix","text":"simplify_matrix(M)\n\nEliminate GCD of entries of every row and column\n\nInput:\n\nM::MatrixElem - matrix to be simplified\n\nOutput:\n\nM::MatrixElem - Simplified matrix\nextra_factors::Vector{AbstractAlgebra.MPolyRingElem} - array of GCDs eliminated from M.\n\n\n\n\n\n","category":"method"},{"location":"tutorials/discrete_time/#Identifiability-of-Discrete-Time-Models-(Local)","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"","category":"section"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Now we consider a discrete-time model in the state-space form. Such a model is typically written either in terms of shift:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"begincases\nmathbfx(t + 1) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"or in terms of difference","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"begincases\nDeltamathbfx(t) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases quad textwhere quad Delta mathbfx(t) = mathbfx(t + 1) - mathbfx(t)","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"In both cases,mathbfx(t) mathbfy(t), and mathbfu(t) are time-dependent states, outputs, and inputs, respectively, and mathbfp are scalar parameters. As in the ODE case, we will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs (in the generic case, see Definition 3 in [1] for details). Again, we will distinguish two types of identifiability","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"local identifiability: the value can be recovered up to finitely many options;\nglobal identifiability: the value can be recovered uniquely.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Currently, StructuralIdentifiability.jl allows to assess only local identifiability for discrete-time models, and below we will describe how this can be done. As a running example, we will use the following discrete version of the SIR model:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"\nbegincases\nS(t + 1) = S(t) - beta S(t) I(t)\nI(t + 1) = I(t) + beta S(t) I(t) - alpha I(t)\nR(t + 1) = R(t) + alpha I(t)\ny(t) = I(t)\nendcases\nquad textor\nquad\nbegincases\nDelta S(t) = -beta S(t) I(t)\nDelta I(t) = beta S(t) I(t) - alpha I(t)\nDelta R(t) = alpha I(t)\ny(t) = I(t)\nendcases","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"where the observable is I, the number of infected people. The native way to define such a model in StructuralIdentifiability is to use @DDSmodel macro which uses the shift notation:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"using StructuralIdentifiability\n\ndds = @DDSmodel(\n S(t + 1) = S(t) - β * S(t) * I(t),\n I(t + 1) = I(t) + β * S(t) * I(t) - α * I(t),\n R(t + 1) = R(t) + α * I(t),\n y(t) = I(t)\n)","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Then local identifiability can be assessed using assess_local_identifiability function:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(dds)","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"For each parameter or state, the value in the returned dictionary is true (1) if the parameter is locally identifiable and false (0) otherwise. We see that R(t) is not identifiable, which makes sense: it does not affect the dynamics of the observable in any way.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"The assess_local_identifiability function has three important keyword arguments:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"funcs_to_check is a list of functions for which one want to assess identifiability, for example, the following code will check if β * S is locally identifiable.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(dds; funcs_to_check = [β * S])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"prob_threshold is the probability of correctness (default value 0.99, i.e., 99%). The underlying algorithm is a Monte-Carlo algorithm, so in principle it may produce incorrect result but the probability of correctness of the returned result is guaranteed to be at least prob_threshold (in fact, the employed bounds are quite conservative, so in practice incorrect result is almost never produced).\nknown_ic is a list of the states for which initial conditions are known. In this case, the identifiability results will be valid not at any time point t but only at t = 0.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"As other main functions in the package, assess_local_identifiability accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"If one loads ModelingToolkit (and thus the ModelingToolkitExt extension), one can use DiscreteSystem from ModelingToolkit to describe the input model (now in terms of difference!):","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"using ModelingToolkit\nusing StructuralIdentifiability\n\n@parameters α β\n@variables t S(t) I(t) R(t) y(t)\nD = Difference(t; dt = 1.0)\n\neqs = [D(S) ~ S - β * S * I, D(I) ~ I + β * S * I - α * I, D(R) ~ R + α * I]\n@named sir = DiscreteSystem(eqs)","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"Then the same computation can be carried out with the models defined this way:","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(sir; measured_quantities = [y ~ I])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"In principle, it is not required to give a name to the observable, so one can write this shorter","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(sir; measured_quantities = [I])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"The same example but with specified functions to check","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"assess_local_identifiability(sir; measured_quantities = [I], funcs_to_check = [β * S])","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"The implementation is based on a version of the observability rank criterion and will be described in a forthcoming paper.","category":"page"},{"location":"tutorials/discrete_time/","page":"Identifiability of Discrete-Time Models (Local)","title":"Identifiability of Discrete-Time Models (Local)","text":"[1]: S. Nõmm, C. Moog, Identifiability of discrete-time nonlinear systems, IFAC Proceedings Volumes, 2004.","category":"page"},{"location":"tutorials/creating_ode/#Creating-ODE-System","page":"Creating ODE System","title":"Creating ODE System","text":"","category":"section"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Most of the algorithms in StructuralIdentifiability.jl take as input system of ordinary differential equations (ODEs) in the state space form, that is:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"begincases\nmathbfx(t) = mathbff(mathbfx(t) mathbfp mathbfu(t))\nmathbfy(t) = mathbfg(mathbfx(t) mathbfp mathbfu(t))\nendcases","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"which involves","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"a vector mathbfx(t) of the state variables of the system,\na vector mathbfu(t) of extermal inputs,\na vector mathbfp of scalar parameters,\na vector mathbfy(t) of outputs (i.e., observations),\nand vectors of rational functions mathbff and mathbfg (for discussion of the non-rational case, see this issue).","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"In the standard setup, inputs and outputs are assumed to be known, and the goal is to assess identifiability of parameters and/or states from the input-output data. In the case of states, this property is also called observability.","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"There are two ways to define such a system to be processed using StructuralIdentifiability.jl. We will demonstrate them using the following example system (Wright's population model of two mutualist species with control[1]):","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"begincases\nx_1(t) = r_1 x_1(t)(1 - c_1 x_1(t)) + fracbeta_1 x_1(t)x_2(t)chi_1 + x_2(t) + u(t)\nx_2(t) = r_2 x_2(t)(1 - c_2 x_2(t)) + fracbeta_2 x_1(t)x_2(t)chi_2 + x_1(t)\ny(t) = x_1(t)\nendcases","category":"page"},{"location":"tutorials/creating_ode/#Defining-the-model-using-@ODEmodel-macro","page":"Creating ODE System","title":"Defining the model using @ODEmodel macro","text":"","category":"section"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"One way to define the model is to use the @ODEmodel macro provided by the StructuralIdentifiability.jl package.","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"using StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) =\n r1 * x1(t) * (1 - c1 * x1(t)) + beta1 * x1(t) * x2(t) / (chi1 + x2(t)) + u(t),\n x2'(t) = r2 * x2(t) * (1 - c2 * x2(t)) + beta2 * x1(t) * x2(t) / (chi2 + x1(t)),\n y(t) = x1(t)\n)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Then one can, for example, assess identifiability of the parameters and states by","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"assess_identifiability(ode)","category":"page"},{"location":"tutorials/creating_ode/#Defining-using-ModelingToolkit","page":"Creating ODE System","title":"Defining using ModelingToolkit","text":"","category":"section"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"StructuralIdentifiability has an extension ModelingToolkitExt which allows to use ODESystem from ModelingToolkit to descibe a model. The extension is loaded automatically once ModelingToolkit is loaded via using ModelingToolkit. In this case, one should encode the equations for the states as ODESystem and specify the outputs separately. In order to do this, we first introduce all functions and scalars:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"using StructuralIdentifiability, ModelingToolkit\n\n@parameters r1, r2, c1, c2, beta1, beta2, chi1, chi2\n@variables t, x1(t), x2(t), y(t), u(t)\n\nD = Differential(t)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"And then defined the system and separately the outputs:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"eqs = [\n D(x1) ~ r1 * x1 * (1 - c1 * x1) + beta1 * x1 * x2 / (chi1 + x2) + u,\n D(x2) ~ r2 * x2 * (1 - c2 * x2) + beta2 * x1 * x2 / (chi2 + x1),\n]\n\nmeasured_quantities = [y ~ x1]\n\node_mtk = ODESystem(eqs, t, name = :mutualist)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"Then, for example, the identifiability of parameters and states can be assessed as follows:","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"assess_identifiability(ode_mtk, measured_quantities = measured_quantities)","category":"page"},{"location":"tutorials/creating_ode/","page":"Creating ODE System","title":"Creating ODE System","text":"[1]: D. H. Wright, A Simple, Stable Model of Mutualism Incorporating Handling Time, The American Naturalist, 1989, 134(4).","category":"page"},{"location":"ioequations/ioequations/#Finding-Input-Output-Equations","page":"Input-Output Equation tools","title":"Finding Input-Output Equations","text":"","category":"section"},{"location":"ioequations/ioequations/","page":"Input-Output Equation tools","title":"Input-Output Equation tools","text":"find_ioequations","category":"page"},{"location":"ioequations/ioequations/#StructuralIdentifiability.find_ioequations","page":"Input-Output Equation tools","title":"StructuralIdentifiability.find_ioequations","text":"find_ioequations(ode, [var_change_policy=:default])\n\nFinds the input-output equations of an ODE system Input:\n\node - the ODE system\nvar_change_policy - whether to perform automatic variable change, can be one of :default, :yes, :no\nloglevel - logging level (default: Logging.Info)\n\nOutput:\n\na dictionary from “leaders” to the corresponding input-output equations; if an extra projection is needed, it will be the value corresponding to rand_proj_var\n\n\n\n\n\n","category":"function"},{"location":"ioequations/ioequations/#Reducing-with-respect-to-Input-Output-Equations","page":"Input-Output Equation tools","title":"Reducing with respect to Input-Output Equations","text":"","category":"section"},{"location":"ioequations/ioequations/","page":"Input-Output Equation tools","title":"Input-Output Equation tools","text":"PBRepresentation\npseudodivision\ndiffreduce\nio_switch!","category":"page"},{"location":"ioequations/ioequations/#StructuralIdentifiability.PBRepresentation","page":"Input-Output Equation tools","title":"StructuralIdentifiability.PBRepresentation","text":"The structure for storing a projection-based representation of differential ideal (see Section 2.3 https://arxiv.org/abs/2111.00991). Contains the following fields:\n\ny_names - the names of the variables with finite order in the profile (typically, outputs)\nu_names - the names of the variables with infinite order in the profile (typically, inputs)\nparam_names - the names of the parameters\nprofile - the profile of the PB-representation (see Definition 2.13) as a dict from y_names with finite orders to the orders\nprojections - the corresponding projections (see Definition 2.15) as a dict from y_names to the projections\n\n\n\n\n\n","category":"type"},{"location":"ioequations/ioequations/#StructuralIdentifiability.pseudodivision","page":"Input-Output Equation tools","title":"StructuralIdentifiability.pseudodivision","text":"pseudodivision(f, g, x)\n\nComputes the result of pseudodivision of f by g as univariate polynomials in x Input:\n\nf - the polynomial to be divided\ng - the polynomial to divide by\nx - the variable for the division\n\nOutput: the pseudoremainder of f divided by g w.r.t. x\n\n\n\n\n\n","category":"function"},{"location":"ioequations/ioequations/#StructuralIdentifiability.diffreduce","page":"Input-Output Equation tools","title":"StructuralIdentifiability.diffreduce","text":"diffreduce(diffpoly, pbr)\n\nComputes the result of differential reduction of a differential polynomial diffpoly with respect to the charset defined by a PB-representation pbr Input:\n\ndiffpoly - a polynomial representing a differential polynomial to be reduced\npbr - a projection-based representation\n\nOutput: the result of differential reduction of diffpoly by pbr considered as a characteristic set (see Remark 2.20 in the paper)\n\n\n\n\n\n","category":"function"},{"location":"ioequations/ioequations/#StructuralIdentifiability.io_switch!","page":"Input-Output Equation tools","title":"StructuralIdentifiability.io_switch!","text":"io_switch(pbr)\n\nIn a single-output pb-representation pbr makes the leading variable to be the first of the inputs\n\n\n\n\n\n","category":"function"},{"location":"utils/local_identifiability/#Local-Identifiability-Tools","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"","category":"section"},{"location":"utils/local_identifiability/","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"Pages=[\"local_identifiability.md\"]","category":"page"},{"location":"utils/local_identifiability/","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"CurrentModule=StructuralIdentifiability","category":"page"},{"location":"utils/local_identifiability/","page":"Local Identifiability Tools","title":"Local Identifiability Tools","text":"StructuralIdentifiability.differentiate_solution\nStructuralIdentifiability.differentiate_output","category":"page"},{"location":"utils/local_identifiability/#StructuralIdentifiability.differentiate_solution","page":"Local Identifiability Tools","title":"StructuralIdentifiability.differentiate_solution","text":"differentiate_solution(ode, params, ic, inputs, prec)\n\nInput:\n\nthe same as for power_series_solutions\n\nOutput:\n\na tuple consisting of the power series solution and a dictionary of the form (u, v) => power series, where u is a state variable v is a state or parameter, and the power series is the partial derivative of the function u w.r.t. v evaluated at the solution\n\n\n\n\n\n","category":"function"},{"location":"utils/local_identifiability/#StructuralIdentifiability.differentiate_output","page":"Local Identifiability Tools","title":"StructuralIdentifiability.differentiate_output","text":"differentiate_output(ode, params, ic, inputs, prec)\n\nSimilar to differentiate_solution but computes partial derivatives of prescribed outputs returns a dictionary of the form y_function => Dict(var => dy/dvar) where dy/dvar is the derivative of y_function with respect to var.\n\n\n\n\n\n","category":"function"},{"location":"#StructuralIdentifiability.jl","page":"Home","title":"StructuralIdentifiability.jl","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"StructuralIdentifiability.jl is a comprehensive toolbox for assessing identifiability parameters.","category":"page"},{"location":"#Installation","page":"Home","title":"Installation","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"To install StructuralIdentifiability.jl, use the Julia package manager:","category":"page"},{"location":"","page":"Home","title":"Home","text":"using Pkg\nPkg.add(\"StructuralIdentifiability\")","category":"page"},{"location":"#Citation","page":"Home","title":"Citation","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"@article{structidjl,\n author = {Dong, R. and Goodbrake, C. and Harrington, H. and Pogudin G.},\n title = {Differential Elimination for Dynamical Models via Projections with Applications to Structural Identifiability},\n journal = {SIAM Journal on Applied Algebra and Geometry},\n url = {https://doi.org/10.1137/22M1469067},\n year = {2023}\n volume = {7},\n number = {1},\n pages = {194-235}\n}","category":"page"},{"location":"#Feature-Summary","page":"Home","title":"Feature Summary","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"StructuralIdentifiability.jl can assess local and global identifiability of ODE models. In addition to these straightforward identifiability queries on individual parameters, the package can distinguish between single- and multi-experiment identifiability.","category":"page"},{"location":"#Feature-List","page":"Home","title":"Feature List","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"Local identifiability checks\nGlobal identifiability checks\nAssessment of identifiable functions of parameters and states\nModel reparametrization (experimental)","category":"page"},{"location":"#Contributing","page":"Home","title":"Contributing","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"Please refer to the SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages for guidance on PRs, issues, and other matters relating to contributing to StructuralIdentifiability.\nThere are a few community forums:\nThe #diffeq-bridged channel in the Julia Slack\nJuliaDiffEq on Gitter\nOn the Julia Discourse forums\nSee also SciML Community page","category":"page"},{"location":"#Reproducibility","page":"Home","title":"Reproducibility","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"
The documentation of this SciML package was built using these direct dependencies,","category":"page"},{"location":"","page":"Home","title":"Home","text":"using Pkg # hide\nPkg.status() # hide","category":"page"},{"location":"","page":"Home","title":"Home","text":"
","category":"page"},{"location":"","page":"Home","title":"Home","text":"
and using this machine and Julia version.","category":"page"},{"location":"","page":"Home","title":"Home","text":"using InteractiveUtils # hide\nversioninfo() # hide","category":"page"},{"location":"","page":"Home","title":"Home","text":"
","category":"page"},{"location":"","page":"Home","title":"Home","text":"
A more complete overview of all dependencies and their versions is also provided.","category":"page"},{"location":"","page":"Home","title":"Home","text":"using Pkg # hide\nPkg.status(; mode = PKGMODE_MANIFEST) # hide","category":"page"},{"location":"","page":"Home","title":"Home","text":"
","category":"page"},{"location":"","page":"Home","title":"Home","text":"You can also download the \nmanifest file and the\nproject file.","category":"page"},{"location":"input/input/#Parsing-input-system","page":"Parsing input system","title":"Parsing input system","text":"","category":"section"},{"location":"input/input/","page":"Parsing input system","title":"Parsing input system","text":"@ODEmodel(ex::Expr...)\nODE\nset_parameter_values","category":"page"},{"location":"input/input/#StructuralIdentifiability.@ODEmodel-Tuple{Vararg{Expr}}","page":"Parsing input system","title":"StructuralIdentifiability.@ODEmodel","text":"macro ODEmodel\n\nMacro for creating an ODE from a list of equations. It also injects all variables into the global scope.\n\nExample\n\nCreating a simple ODE:\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x1'(t) = a * x1(t) + u(t),\n x2'(t) = b * x2(t) + c*x1(t)*x2(t),\n y(t) = x1(t)\n)\n\nHere,\n\nx1, x2 are state variables\ny is an output variable\nu is an input variable\na, b, c are time-indepdendent parameters\n\n\n\n\n\n","category":"macro"},{"location":"input/input/#StructuralIdentifiability.ODE","page":"Parsing input system","title":"StructuralIdentifiability.ODE","text":"The main structure that represents input ODE system.\n\nStores information about states (x_vars), outputs (y_vars), inputs (u_vars), parameters (parameters) and the equations.\n\nThis structure is constructed via @ODEmodel macro.\n\n\n\n\n\n","category":"type"},{"location":"input/input/#StructuralIdentifiability.set_parameter_values","page":"Parsing input system","title":"StructuralIdentifiability.set_parameter_values","text":"set_parameter_values(ode, param_values)\n\nInput:\n\node - an ODE as above\nparam_values - values for (possibly, some of) the parameters as dictionary parameter => value\n\nOutput:\n\nnew ode with the parameters in param_values plugged with the given numbers\n\n\n\n\n\n","category":"function"},{"location":"input/input/#Create-Compartmental-Model","page":"Parsing input system","title":"Create Compartmental Model","text":"","category":"section"},{"location":"input/input/","page":"Parsing input system","title":"Parsing input system","text":"linear_compartment_model","category":"page"},{"location":"input/input/#StructuralIdentifiability.linear_compartment_model","page":"Parsing input system","title":"StructuralIdentifiability.linear_compartment_model","text":"linear_compartment_model(graph, inputs, outputs, leaks)\n\nInput: defines a linear compartment model with nodes numbered from 1 to n by\n\ngraph - and array of integer arrays representing the adjacency lists of the graph\ninputs - array of input nodes\noutputs - array of output nodes\nleaks - array of sink nodes\n\nOutput:\n\nthe corresponding ODE system in the notation of https://doi.org/10.1007/s11538-015-0098-0\n\n\n\n\n\n","category":"function"},{"location":"input/input/#Discrete-time-systems","page":"Parsing input system","title":"Discrete-time systems","text":"","category":"section"},{"location":"input/input/","page":"Parsing input system","title":"Parsing input system","text":"@DDSmodel","category":"page"},{"location":"input/input/#StructuralIdentifiability.@DDSmodel","page":"Parsing input system","title":"StructuralIdentifiability.@DDSmodel","text":"macro DDSmodel\n\nMacro for creating a DDS (discrete dynamical system) from a list of equations. It also injects all variables into the global scope.\n\nExample\n\nCreating a simple DDS:\n\nusing StructuralIdentifiability\n\node = @ODEmodel(\n x1(t + 1) = a * x1(t) + u(t),\n x2(t + 1) = b * x2(t) + c*x1(t)*x2(t),\n y(t) = x1(t)\n)\n\nHere,\n\nx1, x2 are state variables\ny is an output variable\nu is an input variable\na, b, c are time-indepdendent parameters\n\n\n\n\n\n","category":"macro"}] } diff --git a/previews/PR266/tutorials/creating_ode/index.html b/previews/PR266/tutorials/creating_ode/index.html index 9b146c31..a392b39a 100644 --- a/previews/PR266/tutorials/creating_ode/index.html +++ b/previews/PR266/tutorials/creating_ode/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Creating ODE System

Most of the algorithms in StructuralIdentifiability.jl take as input system of ordinary differential equations (ODEs) in the state space form, that is:

\[\begin{cases} +

Creating ODE System

Most of the algorithms in StructuralIdentifiability.jl take as input system of ordinary differential equations (ODEs) in the state space form, that is:

\[\begin{cases} \mathbf{x}'(t) = \mathbf{f}(\mathbf{x}(t), \mathbf{p}, \mathbf{u}(t)),\\ \mathbf{y}(t) = \mathbf{g}(\mathbf{x}(t), \mathbf{p}, \mathbf{u(t)}), \end{cases}\]

which involves

  • a vector $\mathbf{x}(t)$ of the state variables of the system,

  • a vector $\mathbf{u}(t)$ of extermal inputs,

  • a vector $\mathbf{p}$ of scalar parameters,

  • a vector $\mathbf{y}(t)$ of outputs (i.e., observations),

  • and vectors of rational functions $\mathbf{f}$ and $\mathbf{g}$ (for discussion of the non-rational case, see this issue).

In the standard setup, inputs and outputs are assumed to be known, and the goal is to assess identifiability of parameters and/or states from the input-output data. In the case of states, this property is also called observability.

There are two ways to define such a system to be processed using StructuralIdentifiability.jl. We will demonstrate them using the following example system (Wright's population model of two mutualist species with control[1]):

\[\begin{cases} @@ -30,7 +30,7 @@ chi1 => :nonidentifiable chi2 => :globally r1 => :globally - r2 => :globally

Defining using ModelingToolkit

Alternatively, one can use ModelingToolkit: encode the equations for the states as ODESystem and specify the outputs separately. In order to do this, we first introduce all functions and scalars:

using StructuralIdentifiability, ModelingToolkit
+  r2    => :globally

Defining using ModelingToolkit

StructuralIdentifiability has an extension ModelingToolkitExt which allows to use ODESystem from ModelingToolkit to descibe a model. The extension is loaded automatically once ModelingToolkit is loaded via using ModelingToolkit. In this case, one should encode the equations for the states as ODESystem and specify the outputs separately. In order to do this, we first introduce all functions and scalars:

using StructuralIdentifiability, ModelingToolkit
 
 @parameters r1, r2, c1, c2, beta1, beta2, chi1, chi2
 @variables t, x1(t), x2(t), y(t), u(t)
@@ -57,4 +57,4 @@
   c2    => :nonidentifiable
   chi2  => :globally
   r2    => :globally
-  beta2 => :globally
+ beta2 => :globally
diff --git a/previews/PR266/tutorials/discrete_time/index.html b/previews/PR266/tutorials/discrete_time/index.html index 06b18586..0fa55d89 100644 --- a/previews/PR266/tutorials/discrete_time/index.html +++ b/previews/PR266/tutorials/discrete_time/index.html @@ -3,15 +3,44 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Identifiability of Discrete-Time Models (Local)

Now we consider a discrete-time model in the state-space form

\[\begin{cases} +

Identifiability of Discrete-Time Models (Local)

Now we consider a discrete-time model in the state-space form. Such a model is typically written either in terms of shift:

\[\begin{cases} +\mathbf{x}(t + 1) = \mathbf{f}(\mathbf{x}(t), \mathbf{p}, \mathbf{u}(t)),\\ +\mathbf{y}(t) = \mathbf{g}(\mathbf{x}(t), \mathbf{p}, \mathbf{u(t)}), +\end{cases}\]

or in terms of difference

\[\begin{cases} \Delta\mathbf{x}(t) = \mathbf{f}(\mathbf{x}(t), \mathbf{p}, \mathbf{u}(t)),\\ \mathbf{y}(t) = \mathbf{g}(\mathbf{x}(t), \mathbf{p}, \mathbf{u(t)}), -\end{cases} \quad \text{where} \quad \Delta \mathbf{x}(t) := \mathbf{x}(t + 1) - \mathbf{x}(t)\]

and $\mathbf{x}(t), \mathbf{y}(t)$, and $\mathbf{u}(t)$ are time-dependent states, outputs, and inputs, respectively, and $\mathbf{p}$ are scalar parameters. As in the ODE case, we will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs (in the generic case, see Definition 3 in [1] for details). Again, we will distinguish two types of identifiability

  • local identifiability: the value can be recovered up to finitely many options;

  • global identifiability: the value can be recovered uniquely.

Currently, StructuralIdentifiability.jl allows to assess only local identifiability for discrete-time models, and below we will describe how this can be done. As a running example, we will use the following discrete version of the SIR model:

\[\begin{cases} -\Delta S(t) = S(t) - \beta S(t) I(t),\\ -\Delta I(t) = I(t) + \beta S(t) I(t) - \alpha I(t),\\ -\Delta R(t) = R(t) + \alpha I(t),\\ +\end{cases} \quad \text{where} \quad \Delta \mathbf{x}(t) := \mathbf{x}(t + 1) - \mathbf{x}(t).\]

In both cases,$\mathbf{x}(t), \mathbf{y}(t)$, and $\mathbf{u}(t)$ are time-dependent states, outputs, and inputs, respectively, and $\mathbf{p}$ are scalar parameters. As in the ODE case, we will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs (in the generic case, see Definition 3 in [1] for details). Again, we will distinguish two types of identifiability

  • local identifiability: the value can be recovered up to finitely many options;

  • global identifiability: the value can be recovered uniquely.

Currently, StructuralIdentifiability.jl allows to assess only local identifiability for discrete-time models, and below we will describe how this can be done. As a running example, we will use the following discrete version of the SIR model:

\[ +\begin{cases} +S(t + 1) = S(t) - \beta S(t) I(t),\\ +I(t + 1) = I(t) + \beta S(t) I(t) - \alpha I(t),\\ +R(t + 1) = R(t) + \alpha I(t),\\ +y(t) = I(t), +\end{cases} +\quad \text{or} +\quad +\begin{cases} +\Delta S(t) = -\beta S(t) I(t),\\ +\Delta I(t) = \beta S(t) I(t) - \alpha I(t),\\ +\Delta R(t) = \alpha I(t),\\ y(t) = I(t), -\end{cases}\]

where the observable is I, the number of infected people. We start with creating a system as a DiscreteSystem from ModelingToolkit:

using ModelingToolkit
+\end{cases}\]

where the observable is I, the number of infected people. The native way to define such a model in StructuralIdentifiability is to use @DDSmodel macro which uses the shift notation:

using StructuralIdentifiability
+
+dds = @DDSmodel(
+    S(t + 1) = S(t) - β * S(t) * I(t),
+    I(t + 1) = I(t) + β * S(t) * I(t) - α * I(t),
+    R(t + 1) = R(t) + α * I(t),
+    y(t) = I(t)
+)
S(t + 1) = -S(t)*I(t)*β + S(t)
+I(t + 1) = S(t)*I(t)*β - I(t)*α + I(t)
+R(t + 1) = I(t)*α + R(t)
+y(t) = I(t)
+

Then local identifiability can be assessed using assess_local_identifiability function:

assess_local_identifiability(dds)
OrderedCollections.OrderedDict{Any, Bool} with 5 entries:
+  α    => 1
+  β    => 1
+  S(t) => 1
+  I(t) => 1
+  R(t) => 0

For each parameter or state, the value in the returned dictionary is true (1) if the parameter is locally identifiable and false (0) otherwise. We see that R(t) is not identifiable, which makes sense: it does not affect the dynamics of the observable in any way.

The assess_local_identifiability function has three important keyword arguments:

  • funcs_to_check is a list of functions for which one want to assess identifiability, for example, the following code will check if β * S is locally identifiable.
assess_local_identifiability(dds; funcs_to_check = [β * S])
OrderedCollections.OrderedDict{Any, Bool} with 1 entry:
+  S(t)*β => 1
  • prob_threshold is the probability of correctness (default value 0.99, i.e., 99%). The underlying algorithm is a Monte-Carlo algorithm, so in principle it may produce incorrect result but the probability of correctness of the returned result is guaranteed to be at least prob_threshold (in fact, the employed bounds are quite conservative, so in practice incorrect result is almost never produced).

  • known_ic is a list of the states for which initial conditions are known. In this case, the identifiability results will be valid not at any time point t but only at t = 0.

As other main functions in the package, assess_local_identifiability accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.

If one loads ModelingToolkit (and thus the ModelingToolkitExt extension), one can use DiscreteSystem from ModelingToolkit to describe the input model (now in terms of difference!):

using ModelingToolkit
 using StructuralIdentifiability
 
 @parameters α β
@@ -24,15 +53,15 @@
 \mathrm{Difference(t; dt=1.0, update=false)}\left( I\left( t \right) \right) =& I\left( t \right) - I\left( t \right) \alpha + S\left( t \right) I\left( t \right) \beta \\
 \mathrm{Difference(t; dt=1.0, update=false)}\left( R\left( t \right) \right) =& R\left( t \right) + I\left( t \right) \alpha
 \end{align}
- \]

Once the model is defined, we can assess identifiability by providing the formula for the observable:

assess_local_identifiability(sir; measured_quantities = [y ~ I])
OrderedCollections.OrderedDict{SymbolicUtils.BasicSymbolic{Real}, Bool} with 5 entries:
+ \]

Then the same computation can be carried out with the models defined this way:

assess_local_identifiability(sir; measured_quantities = [y ~ I])
OrderedCollections.OrderedDict{SymbolicUtils.BasicSymbolic{Real}, Bool} with 5 entries:
   S(t) => 1
   I(t) => 1
   R(t) => 0
   β    => 1
-  α    => 1

For each parameter or state, the value in the returned dictionary is true (1) if the parameter is locally identifiable and false (0) otherwise. We see that R(t) is not identifiable, which makes sense: it does not affect the dynamics of the observable in any way.

In principle, it is not required to give a name to the observable, so one can write this shorter

assess_local_identifiability(sir; measured_quantities = [I])
OrderedCollections.OrderedDict{SymbolicUtils.BasicSymbolic{Real}, Bool} with 5 entries:
+  α    => 1

In principle, it is not required to give a name to the observable, so one can write this shorter

assess_local_identifiability(sir; measured_quantities = [I])
OrderedCollections.OrderedDict{SymbolicUtils.BasicSymbolic{Real}, Bool} with 5 entries:
   S(t) => 1
   I(t) => 1
   R(t) => 0
   β    => 1
-  α    => 1

The assess_local_identifiability function has three important keyword arguments:

  • funcs_to_check is a list of functions for which one want to assess identifiability, for example, the following code will check if β * S is locally identifiable.
assess_local_identifiability(sir; measured_quantities = [I], funcs_to_check = [β * S])
OrderedCollections.OrderedDict{Symbolics.Num, Bool} with 1 entry:
-  S(t)*β => 1
  • prob_threshold is the probability of correctness (default value 0.99, i.e., 99%). The underlying algorithm is a Monte-Carlo algorithm, so in principle it may produce incorrect result but the probability of correctness of the returned result is guaranteed to be at least prob_threshold (in fact, the employed bounds are quite conservative, so in practice incorrect result is almost never produced).

  • known_ic is a list of the states for which initial conditions are known. In this case, the identifiability results will be valid not at any time point t but only at t = 0.

As other main functions in the package, assess_local_identifiability accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.

The implementation is based on a version of the observability rank criterion and will be described in a forthcoming paper.

+ α => 1

The same example but with specified functions to check

assess_local_identifiability(sir; measured_quantities = [I], funcs_to_check = [β * S])
OrderedCollections.OrderedDict{Symbolics.Num, Bool} with 1 entry:
+  S(t)*β => 1

The implementation is based on a version of the observability rank criterion and will be described in a forthcoming paper.

diff --git a/previews/PR266/tutorials/identifiability/index.html b/previews/PR266/tutorials/identifiability/index.html index 6feb8322..dafc2887 100644 --- a/previews/PR266/tutorials/identifiability/index.html +++ b/previews/PR266/tutorials/identifiability/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Identifiability of Differential Models (Local and Global)

Recall that we consider ODE models in the state-space form

\[\begin{cases} +

Identifiability of Differential Models (Local and Global)

Recall that we consider ODE models in the state-space form

\[\begin{cases} \mathbf{x}'(t) = \mathbf{f}(\mathbf{x}(t), \mathbf{p}, \mathbf{u}(t)),\\ \mathbf{y}(t) = \mathbf{g}(\mathbf{x}(t), \mathbf{p}, \mathbf{u(t)}), \end{cases}\]

where $\mathbf{x}(t), \mathbf{y}(t)$, and $\mathbf{u}(t)$ are time-dependent states, outputs, and inputs, respectively, and $\mathbf{p}$ are scalar parameters. We will call that a parameter or a states (or a function of them) is identifiable if its value can be recovered from time series for inputs and outputs. Typically, two types of identifiability are distinguished

  • local identifiability: the value can be recovered up to finitely many options;

  • global identifiability: the value can be recovered uniquely.

Note that in the case of states, term observability it typically used. We will use identifiability for both states and parameters for brevity and uniformity. While the notion of global identifiability is much more precise, assessing local identifiability is typically much faster, and can be performed for the models whose global identifiability analysis is out of reach.

Local identifiability

We consider the Lotka-Volterra model:

\[\begin{cases} @@ -48,4 +48,4 @@ gama => :nonidentifiable sigma => :globally

As a result, each parameter/state is assigned one of the labels :globally (globally identifiable), :locally (locally but not globally identifiable), or :nonidentifiable (not identifiable, even locally). The algorithm behind this computation follows [4].

Similarly to assess_local_identifiability, this function has optional parameters:

  • funcs_to_check a list of specific functions of parameters and states to check identifiability for (see an example below). If not provided, the identifiability is assessed for all parameters and states. Note that the computations for states may be more involved than for the parameters, so one may want to call the function with funcs_to_check = ode.parameters if the call assess_identifiability(ode) takes too long.

  • prob_threshold (default $0.99$, i.e. 99%) is the probability of correctness. Same story as above: the probability estimates are very conservative, so the actual error probability is much lower than 1%. Also, currently, the probability of correctness does not include the probability of correctness of the modular reconstruction for Groebner bases. This probability is ensured by an additional check modulo a large prime, and can be neglected for practical purposes.

  • loglevel (default Logging.Info). The minimal level of logging messages to be displayed. Available options: Logging.Debug, Logging.Info, Logging.Warn, and Logging.Error.

Using funcs_to_check parameter, one can further inverstigate the nature of the lack of identifiability in the model at hand. For example, for the Goodwin oscillator, we can check if beta + delta and beta * delta are identifiable:

assess_identifiability(ode, funcs_to_check = [beta + delta, beta * delta])
OrderedCollections.OrderedDict{Any, Symbol} with 2 entries:
   beta + delta => :globally
-  beta*delta   => :globally

And we see that they indeed are. This means, in particular, that the reason why beta and delta are not identifiable is because their values can be exchanged. One may wonder how could we guess these functions beta + delta, beta * delta. In fact, they can be just computed using find_identifiable_functions function as we will explain in the next tutorial. Stay tuned!

+ beta*delta => :globally

And we see that they indeed are. This means, in particular, that the reason why beta and delta are not identifiable is because their values can be exchanged. One may wonder how could we guess these functions beta + delta, beta * delta. In fact, they can be just computed using find_identifiable_functions function as we will explain in the next tutorial. Stay tuned!

diff --git a/previews/PR266/tutorials/identifiable_functions/index.html b/previews/PR266/tutorials/identifiable_functions/index.html index 55aa312b..f559b1df 100644 --- a/previews/PR266/tutorials/identifiable_functions/index.html +++ b/previews/PR266/tutorials/identifiable_functions/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Globally Identifiable Functions

In addition to assessing identifiabuility of given functions of parameters or states, StructuralIdentifiability.jl provides the function find_identifiable_functions which finds a set of identifiable functions such that any other identifiable function can be expressed using them. This allows to find out what actually is identifiable and what does non-identifiability in the model at hand looks like.

For example, consider the following model[1].

LLW1987 = @ODEmodel(
+

Globally Identifiable Functions

In addition to assessing identifiabuility of given functions of parameters or states, StructuralIdentifiability.jl provides the function find_identifiable_functions which finds a set of identifiable functions such that any other identifiable function can be expressed using them. This allows to find out what actually is identifiable and what does non-identifiability in the model at hand looks like.

For example, consider the following model[1].

LLW1987 = @ODEmodel(
     x1'(t) = -p1 * x1(t) + p2 * u(t),
     x2'(t) = -p3 * x2(t) + p4 * u(t),
     x3'(t) = -(p1 + p3) * x3(t) + (p4 * x1(t) + p2 * x2(t)) * u(t),
@@ -22,4 +22,4 @@
  p1*p3
  p1 + p3
  x1(t)*p4 + x2(t)*p2
- (-p1 + p3)//(x1(t)*p4 - x2(t)*p2)

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 :strong to try to simplify further (at the expense of heavier computation) or use :weak to simplify less (but compute faster).

As assess_identifiability and assess_local_identifiability, find_identifiable_functions accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.

+ (-p1 + p3)//(x1(t)*p4 - x2(t)*p2)

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 :strong to try to simplify further (at the expense of heavier computation) or use :weak to simplify less (but compute faster).

As assess_identifiability and assess_local_identifiability, find_identifiable_functions accepts an optional parameter loglevel (default: Logging.Info) to adjust the verbosity of logging.

diff --git a/previews/PR266/tutorials/reparametrization/index.html b/previews/PR266/tutorials/reparametrization/index.html index 3f7e5883..62d8702e 100644 --- a/previews/PR266/tutorials/reparametrization/index.html +++ b/previews/PR266/tutorials/reparametrization/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Reparametrizations

Overview

Once one has found that not all parameters and/or states of the model at hand are identifiable, one natural desire is to reparametrize the model into a one with better identifiability properties. StructuralIdentifiability offers such a functionality via the function reparametrize_global. It takes as input an ODE model and produces its transformation into another model with the same input-output behaviour but with the states and parameters being globally identifiable. Note that, in general, such a transformation may not exist in the class of rational models, so sometimes the function returns an ODE not on the whole affine space but on a manifold.

More precisely, the function returns a dictionary with three keys:

  • :new_vars is a dictionary which maps the new parameters and new states into the formulas expressing them in terms of the original parameters and states;

  • :new_ode is the ODE satisfied by these new states (and the expression of the output in terms of the new states);

  • :implicit_relations is a list of algebraic relations between the new states and parameters. Being nonempty, this is exactly the list of equations defining the manifold, on which the new ODE model is defined. In many interesting cases, however, this list is empty meaning that the new ODE is a standard rational ODE model.

Example: SEUIR model

Consider a SEUIR epidemiological model from[1]:

\[\begin{cases} +

Reparametrizations

Overview

Once one has found that not all parameters and/or states of the model at hand are identifiable, one natural desire is to reparametrize the model into a one with better identifiability properties. StructuralIdentifiability offers such a functionality via the function reparametrize_global. It takes as input an ODE model and produces its transformation into another model with the same input-output behaviour but with the states and parameters being globally identifiable. Note that, in general, such a transformation may not exist in the class of rational models, so sometimes the function returns an ODE not on the whole affine space but on a manifold.

More precisely, the function returns a dictionary with three keys:

  • :new_vars is a dictionary which maps the new parameters and new states into the formulas expressing them in terms of the original parameters and states;

  • :new_ode is the ODE satisfied by these new states (and the expression of the output in terms of the new states);

  • :implicit_relations is a list of algebraic relations between the new states and parameters. Being nonempty, this is exactly the list of equations defining the manifold, on which the new ODE model is defined. In many interesting cases, however, this list is empty meaning that the new ODE is a standard rational ODE model.

Example: SEUIR model

Consider a SEUIR epidemiological model from[1]:

\[\begin{cases} S(t)' = -\beta \frac{(U(t) + I(t))S(t)}{N},\\ E(t)' = \beta \frac{(U(t) + I(t))S(t)}{N} - \gamma E(t),\\ U(t)' = (1 - \alpha) \gamma E(t) - \delta U(t),\\ @@ -60,4 +60,4 @@ X4(t) => :globally a1 => :globally a2 => :globally - a3 => :globally

+ a3 => :globally
diff --git a/previews/PR266/utils/elimination/index.html b/previews/PR266/utils/elimination/index.html index c1c82898..437561f1 100644 --- a/previews/PR266/utils/elimination/index.html +++ b/previews/PR266/utils/elimination/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Elimination

StructuralIdentifiability.Bezout_matrixMethod
Bezout_matrix(f, g, var_elim)

Compute the Bezout matrix of two polynomials f, g with respect to var_elim

Inputs:

  • f - first polynomial
  • g - second polynomial
  • var_elim - variable, of which f and g are considered as polynomials

Output:

  • M::MatrixElem - The Bezout matrix
source
StructuralIdentifiability.Sylvester_matrixMethod
Sylvester_matrix(f, g, var_elim)

Compute the Sylvester matrix of two polynomials f, g with respect to var_elim Inputs:

  • f - first polynomial
  • g - second polynomial
  • var_elim - variable, of which f and g are considered as polynomials

Output:

  • M::MatrixElem - The Sylvester matrix
source
StructuralIdentifiability.chooseMethod
choose(polys, generic_point_generator)

Input:

  • polys - an array of distinct irreducible polynomials in the same ring
  • generic_point_generator - a generic point generator as described above for one of polys

Output:

  • the polynomial that vanishes at the generic_point_generator
source
StructuralIdentifiability.eliminate_varMethod
eliminate_var(f, g, var_elim, generic_point_generator)

Eliminate a variable from a pair of polynomials

Input:

  • f and g - polynomials
  • var_elim - variable to be eliminated
  • generic_point_generator - a generic point generator object for the factor of the resultant of f and g of interest

Output:

  • polynomial - the desired factor of the resultant of f and g
source
StructuralIdentifiability.simplify_matrixMethod
simplify_matrix(M)

Eliminate GCD of entries of every row and column

Input:

  • M::MatrixElem - matrix to be simplified

Output:

  • M::MatrixElem - Simplified matrix
  • extra_factors::Vector{AbstractAlgebra.MPolyRingElem} - array of GCDs eliminated from M.
source
+

Elimination

StructuralIdentifiability.Bezout_matrixMethod
Bezout_matrix(f, g, var_elim)

Compute the Bezout matrix of two polynomials f, g with respect to var_elim

Inputs:

  • f - first polynomial
  • g - second polynomial
  • var_elim - variable, of which f and g are considered as polynomials

Output:

  • M::MatrixElem - The Bezout matrix
source
StructuralIdentifiability.Sylvester_matrixMethod
Sylvester_matrix(f, g, var_elim)

Compute the Sylvester matrix of two polynomials f, g with respect to var_elim Inputs:

  • f - first polynomial
  • g - second polynomial
  • var_elim - variable, of which f and g are considered as polynomials

Output:

  • M::MatrixElem - The Sylvester matrix
source
StructuralIdentifiability.chooseMethod
choose(polys, generic_point_generator)

Input:

  • polys - an array of distinct irreducible polynomials in the same ring
  • generic_point_generator - a generic point generator as described above for one of polys

Output:

  • the polynomial that vanishes at the generic_point_generator
source
StructuralIdentifiability.eliminate_varMethod
eliminate_var(f, g, var_elim, generic_point_generator)

Eliminate a variable from a pair of polynomials

Input:

  • f and g - polynomials
  • var_elim - variable to be eliminated
  • generic_point_generator - a generic point generator object for the factor of the resultant of f and g of interest

Output:

  • polynomial - the desired factor of the resultant of f and g
source
StructuralIdentifiability.simplify_matrixMethod
simplify_matrix(M)

Eliminate GCD of entries of every row and column

Input:

  • M::MatrixElem - matrix to be simplified

Output:

  • M::MatrixElem - Simplified matrix
  • extra_factors::Vector{AbstractAlgebra.MPolyRingElem} - array of GCDs eliminated from M.
source
diff --git a/previews/PR266/utils/global_identifiability/index.html b/previews/PR266/utils/global_identifiability/index.html index 7ede406f..6d22255c 100644 --- a/previews/PR266/utils/global_identifiability/index.html +++ b/previews/PR266/utils/global_identifiability/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Global Identifiability Tools

StructuralIdentifiability.RationalFunctionFieldType
RationalFunctionField

A subfield of the field of rational functions over the rationals.

Example

using Nemo
+

Global Identifiability Tools

StructuralIdentifiability.RationalFunctionFieldType
RationalFunctionField

A subfield of the field of rational functions over the rationals.

Example

using Nemo
 using StructuralIdentifiability: RationalFunctionField
 
 R, (x, y, z) = QQ["x", "y", "z"]
@@ -12,4 +12,4 @@
 rff = RationalFunctionField([x // y, y // z])
 
 # Constructs a subfield generated by y / x, 1 / x, z / y
-rff = RationalFunctionField([[x, y, R(1)], [y, z]])
source
StructuralIdentifiability.field_containsFunction
field_contains(field, ratfuncs, prob_threshold)

Checks whether given rational function field field contains given rational functions ratfuncs (represented as a list of lists). The result is correct with probability at least prob_threshold

Inputs:

  • field - a rational function field
  • ratfuncs - a list of lists of polynomials. Each of the lists, say, [f1, ..., fn], defines generators f2/f1, ..., fn/f1.
  • prob_threshold real number from (0, 1)

Output:

  • a list L[i] of bools of length length(rat_funcs) such that L[i] is true iff the i-th function belongs to field
source
+rff = RationalFunctionField([[x, y, R(1)], [y, z]])
source
StructuralIdentifiability.field_containsFunction
field_contains(field, ratfuncs, prob_threshold)

Checks whether given rational function field field contains given rational functions ratfuncs (represented as a list of lists). The result is correct with probability at least prob_threshold

Inputs:

  • field - a rational function field
  • ratfuncs - a list of lists of polynomials. Each of the lists, say, [f1, ..., fn], defines generators f2/f1, ..., fn/f1.
  • prob_threshold real number from (0, 1)

Output:

  • a list L[i] of bools of length length(rat_funcs) such that L[i] is true iff the i-th function belongs to field
source
diff --git a/previews/PR266/utils/local_identifiability/index.html b/previews/PR266/utils/local_identifiability/index.html index f9d8a6d7..619e0743 100644 --- a/previews/PR266/utils/local_identifiability/index.html +++ b/previews/PR266/utils/local_identifiability/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Local Identifiability Tools

StructuralIdentifiability.differentiate_solutionFunction
differentiate_solution(ode, params, ic, inputs, prec)

Input:

  • the same as for power_series_solutions

Output:

  • a tuple consisting of the power series solution and a dictionary of the form (u, v) => power series, where u is a state variable v is a state or parameter, and the power series is the partial derivative of the function u w.r.t. v evaluated at the solution
source
StructuralIdentifiability.differentiate_outputFunction
differentiate_output(ode, params, ic, inputs, prec)

Similar to differentiate_solution but computes partial derivatives of prescribed outputs returns a dictionary of the form y_function => Dict(var => dy/dvar) where dy/dvar is the derivative of y_function with respect to var.

source
+

Local Identifiability Tools

StructuralIdentifiability.differentiate_solutionFunction
differentiate_solution(ode, params, ic, inputs, prec)

Input:

  • the same as for power_series_solutions

Output:

  • a tuple consisting of the power series solution and a dictionary of the form (u, v) => power series, where u is a state variable v is a state or parameter, and the power series is the partial derivative of the function u w.r.t. v evaluated at the solution
source
StructuralIdentifiability.differentiate_outputFunction
differentiate_output(ode, params, ic, inputs, prec)

Similar to differentiate_solution but computes partial derivatives of prescribed outputs returns a dictionary of the form y_function => Dict(var => dy/dvar) where dy/dvar is the derivative of y_function with respect to var.

source
diff --git a/previews/PR266/utils/ode/index.html b/previews/PR266/utils/ode/index.html index b814ead6..20ee6821 100644 --- a/previews/PR266/utils/ode/index.html +++ b/previews/PR266/utils/ode/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Functions to work with the ODE structure

StructuralIdentifiability.power_series_solutionMethod
power_series_solution(ode, param_values, initial_conditions, input_values, prec)

Input:

  • ode - an ode to solve
  • param_values - parameter values, must be a dictionary mapping parameter to a value
  • initial_conditions - initial conditions of ode, must be a dictionary mapping state variable to a value
  • input_values - power series for the inputs presented as a dictionary variable => list of coefficients
  • prec - the precision of solutions

Output:

  • computes a power series solution with precision prec presented as a dictionary variable => corresponding coordinate of the solution
source
StructuralIdentifiability.set_parameter_valuesMethod
set_parameter_values(ode, param_values)

Input:

  • ode - an ODE as above
  • param_values - values for (possibly, some of) the parameters as dictionary parameter => value

Output:

  • new ode with the parameters in param_values plugged with the given numbers
source
StructuralIdentifiability.find_submodelsMethod
find_submodels(ode)

The function calculates and returns all valid submodels given a system of ODEs.

Input:

  • ode - an ODEs system to be studied

Output:

  • A list of submodels represented as ode objects

Example:

>ode = @ODEmodel(x1'(t) = x1(t)^2, 
+

Functions to work with the ODE structure

StructuralIdentifiability.power_series_solutionMethod
power_series_solution(ode, param_values, initial_conditions, input_values, prec)

Input:

  • ode - an ode to solve
  • param_values - parameter values, must be a dictionary mapping parameter to a value
  • initial_conditions - initial conditions of ode, must be a dictionary mapping state variable to a value
  • input_values - power series for the inputs presented as a dictionary variable => list of coefficients
  • prec - the precision of solutions

Output:

  • computes a power series solution with precision prec presented as a dictionary variable => corresponding coordinate of the solution
source
StructuralIdentifiability.set_parameter_valuesMethod
set_parameter_values(ode, param_values)

Input:

  • ode - an ODE as above
  • param_values - values for (possibly, some of) the parameters as dictionary parameter => value

Output:

  • new ode with the parameters in param_values plugged with the given numbers
source
StructuralIdentifiability.find_submodelsMethod
find_submodels(ode)

The function calculates and returns all valid submodels given a system of ODEs.

Input:

  • ode - an ODEs system to be studied

Output:

  • A list of submodels represented as ode objects

Example:

>ode = @ODEmodel(x1'(t) = x1(t)^2, 
                  x2'(t) = x1(t) * x2(t), 
                  y1(t) = x1(t), 
                  y2(t) = x2(t))
@@ -12,4 +12,4 @@
         
         x1'(t) = a(t)*x2(t)^2 + x1(t)
         y1(t) = x1(t)
-    ]
source
+ ]
source
diff --git a/previews/PR266/utils/power_series_utils/index.html b/previews/PR266/utils/power_series_utils/index.html index 6958ecae..d500b978 100644 --- a/previews/PR266/utils/power_series_utils/index.html +++ b/previews/PR266/utils/power_series_utils/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Power Series Utilities

StructuralIdentifiability.ps_matrix_homlinear_deMethod
ps_matrix_homlinear_de(A, Y0, prec)

Input:

  • A - a square matrix with entries in a univariate power series ring
  • Y0 - a square invertible matrix over the base field

Output:

  • matrix Y such that Y' = AY up to precision of A - 1 and Y(0) = Y0
source
StructuralIdentifiability.ps_matrix_invFunction
ps_matrix_inv(M, prec)

Input:

  • M - a square matrix with entries in a univariate power series ring it is assumed that M(0) is invertible and all entries having the same precision
  • prec - an integer, precision, if -1 then defaults to precision of M

Output:

  • the inverse of M computed up to prec
source
StructuralIdentifiability.ps_matrix_linear_deMethod
ps_matrix_linear_de(A, B, Y0, prec)

Input:

  • A, B - square matrices with entries in a univariate power series ring
  • Y0 - a matrix over the base field with the rows number the same as A

Output:

  • matrix Y such that Y' = AY + B up to precision of A - 1 and Y(0) = Y0
source
StructuralIdentifiability.ps_ode_solutionMethod
ps_ode_solution(equations, ic, inputs, prec)

Input:

  • equations - a system of the form $A(x, u, mu)x' - B(x, u, mu) = 0$, where A is a generically nonsingular square matrix. Assumption: A is nonzero at zero
  • ic - initial conditions for x's (dictionary)
  • inputs - power series for inputs represented as arrays (dictionary)
  • prec - precision of the solution

Output:

  • power series solution of the system
source
+

Power Series Utilities

StructuralIdentifiability.ps_matrix_homlinear_deMethod
ps_matrix_homlinear_de(A, Y0, prec)

Input:

  • A - a square matrix with entries in a univariate power series ring
  • Y0 - a square invertible matrix over the base field

Output:

  • matrix Y such that Y' = AY up to precision of A - 1 and Y(0) = Y0
source
StructuralIdentifiability.ps_matrix_invFunction
ps_matrix_inv(M, prec)

Input:

  • M - a square matrix with entries in a univariate power series ring it is assumed that M(0) is invertible and all entries having the same precision
  • prec - an integer, precision, if -1 then defaults to precision of M

Output:

  • the inverse of M computed up to prec
source
StructuralIdentifiability.ps_matrix_linear_deMethod
ps_matrix_linear_de(A, B, Y0, prec)

Input:

  • A, B - square matrices with entries in a univariate power series ring
  • Y0 - a matrix over the base field with the rows number the same as A

Output:

  • matrix Y such that Y' = AY + B up to precision of A - 1 and Y(0) = Y0
source
StructuralIdentifiability.ps_ode_solutionMethod
ps_ode_solution(equations, ic, inputs, prec)

Input:

  • equations - a system of the form $A(x, u, mu)x' - B(x, u, mu) = 0$, where A is a generically nonsingular square matrix. Assumption: A is nonzero at zero
  • ic - initial conditions for x's (dictionary)
  • inputs - power series for inputs represented as arrays (dictionary)
  • prec - precision of the solution

Output:

  • power series solution of the system
source
diff --git a/previews/PR266/utils/primality/index.html b/previews/PR266/utils/primality/index.html index 776d87eb..671f2306 100644 --- a/previews/PR266/utils/primality/index.html +++ b/previews/PR266/utils/primality/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Primality Checks

StructuralIdentifiability.check_primalityFunction
check_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem}, extra_relations::Array{QQMPolyRingElem, 1})

The function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.

The extra_relations allows adding more polynomials to the generators (not affecting the saturation).

source
check_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem})

The function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.

source
+

Primality Checks

StructuralIdentifiability.check_primalityFunction
check_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem}, extra_relations::Array{QQMPolyRingElem, 1})

The function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.

The extra_relations allows adding more polynomials to the generators (not affecting the saturation).

source
check_primality(polys::Dict{QQMPolyRingElem, QQMPolyRingElem})

The function checks if the ideal generated by the polynomials and saturated at the leading coefficient with respect to the corresponding variables is prime over rationals.

source
diff --git a/previews/PR266/utils/reparametrization/index.html b/previews/PR266/utils/reparametrization/index.html index da1f88f8..dcbd226e 100644 --- a/previews/PR266/utils/reparametrization/index.html +++ b/previews/PR266/utils/reparametrization/index.html @@ -3,7 +3,7 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Model reparametrization

StructuralIdentifiability.reparametrize_globalFunction
reparametrize_global(ode, options...)

Finds a reparametrization of ode in terms of globally identifiabile functions.

Returns a tuple (new_ode, new_vars, implicit_relations), such that:

  • new_ode is the reparametrized ODE system
  • new_vars is a mapping from the new variables to the original ones
  • relations is the array of implicit algebraic relations among new_vars. Usually, relations is empty.

Options

The function accepts the following optional arguments.

  • seed: A float in the range from 0 to 1, random seed (default is seed = 42).
  • prob_threshold: The probability of correctness (default is prob_threshold = 0.99).

Example

using StructuralIdentifiability
+

Model reparametrization

StructuralIdentifiability.reparametrize_globalFunction
reparametrize_global(ode, options...)

Finds a reparametrization of ode in terms of globally identifiabile functions.

Returns a tuple (new_ode, new_vars, implicit_relations), such that:

  • new_ode is the reparametrized ODE system
  • new_vars is a mapping from the new variables to the original ones
  • relations is the array of implicit algebraic relations among new_vars. Usually, relations is empty.

Options

The function accepts the following optional arguments.

  • seed: A float in the range from 0 to 1, random seed (default is seed = 42).
  • prob_threshold: The probability of correctness (default is prob_threshold = 0.99).

Example

using StructuralIdentifiability
 
 ode = @ODEmodel(
     x1'(t) = a * x1(t) - b*x1(t)*x2(t),
@@ -23,4 +23,4 @@
   X1 => x1
   a2 => d
   a3 => a
-  a1 => c

Notice that the new_ode is fully identifiabile, and has 1 less parameter compared to the original one.

source
+ a1 => c

Notice that the new_ode is fully identifiabile, and has 1 less parameter compared to the original one.

source
diff --git a/previews/PR266/utils/util/index.html b/previews/PR266/utils/util/index.html index ea2d831e..d5b82641 100644 --- a/previews/PR266/utils/util/index.html +++ b/previews/PR266/utils/util/index.html @@ -3,5 +3,5 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Other Helpful Functions

StructuralIdentifiability.dennums_to_fractionsMethod
dennums_to_fractions(dennums)

Returns the field generators represented by fractions.

Input: an array of arrays of polynomials, as in [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]

Output: an array of fractions [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]

source
StructuralIdentifiability.extract_coefficientsMethod
extract_coefficients(poly, variables)

Input:

  • poly - multivariate polynomial
  • variables - a list of variables from the generators of the ring of p

Output:

  • dictionary with keys being tuples of length lenght(variables) and values being polynomials in the variables other than those which are the coefficients at the corresponding monomials (in a smaller polynomial ring)
source
StructuralIdentifiability.fractions_to_dennumsMethod
fractions_to_dennums(fractions)

Returns the field generators represented by lists of denominators and numerators.

Input: an array of fractions, as in [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]

Output: an array of arrays of polynomials, [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]

source
StructuralIdentifiability.gen_tag_nameFunction
gen_tag_name(base; stop_words)
-gen_tag_names(n, base; stop_words)

Generates a string which will not collide with the words in stop_words.

Arguments

  • n: Generates a sequence of unique strings of length n
  • base: A string or a vector of strings, the base for the generated sequence
  • stop_words: A vector of strings, stop words
source
StructuralIdentifiability.make_substitutionMethod
make_substitution(f, var_sub, val_numer, val_denom)

Substitute a variable in a polynomial with an expression

Input:

  • f - the polynomial
  • var_sub - the variable to be substituted
  • var_numer - numerator of the substitution expression
  • var_denom - denominator of the substitution expression

Output:

  • polynomial - result of substitution
source
StructuralIdentifiability.parent_ring_changeMethod
parent_ring_change(poly, new_ring)

Converts a polynomial to a different polynomial ring Input

  • poly - a polynomial to be converted
  • new_ring - a polynomial ring such that every variable name appearing in poly appears among the generators

Output:

  • a polynomial in new_ring “equal” to poly
source
StructuralIdentifiability.uncertain_factorizationMethod
uncertain_factorization(f)

Input:

  • f - polynomial with rational coefficients

Output:

  • list of pairs (div, certainty) where
    • div's are divisors of f such that f is their product with certain powers
    • if certainty is true, div is $Q$-irreducible
source
+

Other Helpful Functions

StructuralIdentifiability.dennums_to_fractionsMethod
dennums_to_fractions(dennums)

Returns the field generators represented by fractions.

Input: an array of arrays of polynomials, as in [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]

Output: an array of fractions [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]

source
StructuralIdentifiability.extract_coefficientsMethod
extract_coefficients(poly, variables)

Input:

  • poly - multivariate polynomial
  • variables - a list of variables from the generators of the ring of p

Output:

  • dictionary with keys being tuples of length lenght(variables) and values being polynomials in the variables other than those which are the coefficients at the corresponding monomials (in a smaller polynomial ring)
source
StructuralIdentifiability.fractions_to_dennumsMethod
fractions_to_dennums(fractions)

Returns the field generators represented by lists of denominators and numerators.

Input: an array of fractions, as in [f2/f1, f3/f1, ..., g2/g1, g3/g1, ...]

Output: an array of arrays of polynomials, [[f1, f2, f3, ...], [g1, g2, g3, ...], ...]

source
StructuralIdentifiability.gen_tag_nameFunction
gen_tag_name(base; stop_words)
+gen_tag_names(n, base; stop_words)

Generates a string which will not collide with the words in stop_words.

Arguments

  • n: Generates a sequence of unique strings of length n
  • base: A string or a vector of strings, the base for the generated sequence
  • stop_words: A vector of strings, stop words
source
StructuralIdentifiability.make_substitutionMethod
make_substitution(f, var_sub, val_numer, val_denom)

Substitute a variable in a polynomial with an expression

Input:

  • f - the polynomial
  • var_sub - the variable to be substituted
  • var_numer - numerator of the substitution expression
  • var_denom - denominator of the substitution expression

Output:

  • polynomial - result of substitution
source
StructuralIdentifiability.parent_ring_changeMethod
parent_ring_change(poly, new_ring)

Converts a polynomial to a different polynomial ring Input

  • poly - a polynomial to be converted
  • new_ring - a polynomial ring such that every variable name appearing in poly appears among the generators

Output:

  • a polynomial in new_ring “equal” to poly
source
StructuralIdentifiability.uncertain_factorizationMethod
uncertain_factorization(f)

Input:

  • f - polynomial with rational coefficients

Output:

  • list of pairs (div, certainty) where
    • div's are divisors of f such that f is their product with certain powers
    • if certainty is true, div is $Q$-irreducible
source
diff --git a/previews/PR266/utils/wronskian/index.html b/previews/PR266/utils/wronskian/index.html index 416bd463..82d1560a 100644 --- a/previews/PR266/utils/wronskian/index.html +++ b/previews/PR266/utils/wronskian/index.html @@ -3,4 +3,4 @@ function gtag(){dataLayer.push(arguments);} gtag('js', new Date()); gtag('config', 'UA-90474609-3', {'page_path': location.pathname + location.search + location.hash}); -

Wronskian Tools

StructuralIdentifiability.get_max_belowMethod
get_max_below(t, vect)

Input:

  • t - a trie with exponent vectors
  • vect - yet another exponent vector

Output:

  • a pair (d, v) where v is a vector in the trie which is componentwise ≤ vect and the difference d is as small as possible
source
StructuralIdentifiability.massive_evalMethod
massive_eval(polys, eval_dict)

Input:

  • polys - a list of polynomials
  • eval_dict - dictionary from variables to the values. Missing values are treated as zeroes

Output:

  • a list of values of the polynomials

Evaluates a list of polynomials at a point. Assumes that multiplications are relatively expensive (like in truncated power series) so all the monomials are precomputed first and the values of monomials of lower degree are cached and used to compute the values of the monomials of higher degree

source
StructuralIdentifiability.monomial_compressMethod
monomial_compress(io_equation, ode)

Compresses an input-output equation for the rank computation Input:

  • io_equation - input-output equation
  • ode - the corresponding ODE model

Output:

  • pair (coeffs, terms) such that:
    • sum of coeffs[i] * terms[i] = io_equation
    • coeffs involve only parameters, terms involve only inputs and outputs
    • length of the representation is the smallest possible
source
StructuralIdentifiability.wronskianMethod
wronskian(io_equations, ode)

Input:

  • io_equations - a set of io-equations in the form of the Dict as returned by find_ioequations
  • ode - the ODE object

Output:

  • a list of Wronskians evaluated at a point modulo prime

Computes the Wronskians of io_equations

source
+

Wronskian Tools

StructuralIdentifiability.get_max_belowMethod
get_max_below(t, vect)

Input:

  • t - a trie with exponent vectors
  • vect - yet another exponent vector

Output:

  • a pair (d, v) where v is a vector in the trie which is componentwise ≤ vect and the difference d is as small as possible
source
StructuralIdentifiability.massive_evalMethod
massive_eval(polys, eval_dict)

Input:

  • polys - a list of polynomials
  • eval_dict - dictionary from variables to the values. Missing values are treated as zeroes

Output:

  • a list of values of the polynomials

Evaluates a list of polynomials at a point. Assumes that multiplications are relatively expensive (like in truncated power series) so all the monomials are precomputed first and the values of monomials of lower degree are cached and used to compute the values of the monomials of higher degree

source
StructuralIdentifiability.monomial_compressMethod
monomial_compress(io_equation, ode)

Compresses an input-output equation for the rank computation Input:

  • io_equation - input-output equation
  • ode - the corresponding ODE model

Output:

  • pair (coeffs, terms) such that:
    • sum of coeffs[i] * terms[i] = io_equation
    • coeffs involve only parameters, terms involve only inputs and outputs
    • length of the representation is the smallest possible
source
StructuralIdentifiability.wronskianMethod
wronskian(io_equations, ode)

Input:

  • io_equations - a set of io-equations in the form of the Dict as returned by find_ioequations
  • ode - the ODE object

Output:

  • a list of Wronskians evaluated at a point modulo prime

Computes the Wronskians of io_equations

source