diff --git a/.gitignore b/.gitignore index ca42cbdc0..4bd358137 100644 --- a/.gitignore +++ b/.gitignore @@ -5,5 +5,6 @@ Manifest.toml .DS_Store /docs/build/* /docs/src/libs/*/examples +/docs/src/assets/Project.toml hall_of_fame.* -/tmp/* \ No newline at end of file +/tmp/* diff --git a/Project.toml b/Project.toml index 8b02e9c4a..af1f38232 100644 --- a/Project.toml +++ b/Project.toml @@ -17,6 +17,7 @@ QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" Reexport = "189a3867-3050-52da-a836-e630ba90ab69" +SciMLStructures = "53ae85a6-f571-4167-b2af-e1d143709226" Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" @@ -35,6 +36,8 @@ ProgressMeter = "1.6" QuadGK = "2.4" RecipesBase = "1" Reexport = "1.0" +OrdinaryDiffEqTsit5 = "1" +SciMLStructures = "1" Setfield = "1" Statistics = "1" StatsBase = "0.32.0, 0.33, 0.34" @@ -43,10 +46,11 @@ Symbolics = "5.30.1, 6" julia = "1.10" [extras] +OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a" Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Test", "Random", "SafeTestsets", "Pkg"] +test = ["Test", "Random", "SafeTestsets", "Pkg", "OrdinaryDiffEqTsit5"] diff --git a/docs/make.jl b/docs/make.jl index 383b10d87..58b068531 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -12,6 +12,7 @@ dev_subpkg("DataDrivenDMD") dev_subpkg("DataDrivenSparse") dev_subpkg("DataDrivenSR") + using Documenter using DataDrivenDiffEq using DataDrivenDMD @@ -21,8 +22,8 @@ using DataDrivenSR using StatsBase using Literate -cp("./docs/Manifest.toml", "./docs/src/assets/Manifest.toml", force = true) -cp("./docs/Project.toml", "./docs/src/assets/Project.toml", force = true) +cp(joinpath(@__DIR__, "Manifest.toml"), joinpath(@__DIR__, "src", "assets", "Manifest.toml"), force = true) +cp(joinpath(@__DIR__, "Project.toml"), joinpath(@__DIR__, "src", "assets", "Project.toml"), force = true) ENV["GKSwstype"] = "100" @@ -73,9 +74,12 @@ makedocs(sitename = "DataDrivenDiffEq.jl", modules = [DataDrivenDiffEq, DataDrivenDMD, DataDrivenSparse, DataDrivenSR], clean = true, doctest = false, linkcheck = true, warnonly = [:missing_docs, :cross_references], - linkcheck_ignore = ["http://cwrowley.princeton.edu/papers/Hemati-2017a.pdf", + linkcheck_ignore = [ + "http://cwrowley.princeton.edu/papers/Hemati-2017a.pdf", "https://royalsocietypublishing.org/doi/10.1098/rspa.2020.0279", - "https://www.pnas.org/doi/10.1073/pnas.1517384113"], + "https://www.pnas.org/doi/10.1073/pnas.1517384113", + "https://link.springer.com/article/10.1007/s00332-015-9258-5", + ], format = Documenter.HTML(assets = ["assets/favicon.ico"], canonical = "https://docs.sciml.ai/DataDrivenDiffEq/stable/"), pages = pages) diff --git a/docs/src/assets/Project.toml b/docs/src/assets/Project.toml deleted file mode 100644 index c7a92b3bd..000000000 --- a/docs/src/assets/Project.toml +++ /dev/null @@ -1,16 +0,0 @@ -[deps] -DataDrivenDMD = "3c9adf31-5280-42ff-b439-b71cc6b07807" -DataDrivenDiffEq = "2445eb08-9709-466a-b3fc-47e12bd697a2" -DataDrivenSR = "7fed8a53-d475-4873-af3a-ba53cceea094" -DataDrivenSparse = "5b588203-7d8b-4fab-a537-c31a7f73f46b" -Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" -LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" -Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306" -ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" -OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" -Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" -StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3" -StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" - -[compat] -Documenter = "0.27" diff --git a/docs/src/citations.md b/docs/src/citations.md index cb1d187b7..5574347a3 100644 --- a/docs/src/citations.md +++ b/docs/src/citations.md @@ -15,16 +15,19 @@ If you are using `DataDrivenDiffEq.jl` for research, please cite } ``` -If you are using the [SymbolicRegression.jl](https://docs.sciml.ai/SymbolicRegression/stable/) API, please cite +If you are using the [SymbolicRegression.jl](https://ai.damtp.cam.ac.uk/symbolicregression/dev/) API, please cite ```bibtex -@software{pysr, - author = {Miles Cranmer}, - title = {PySR: Fast \& Parallelized Symbolic Regression in Python/Julia}, - month = sep, - year = 2020, - publisher = {Zenodo}, - doi = {10.5281/zenodo.4041459}, - url = {http://doi.org/10.5281/zenodo.4041459} +@misc{cranmerInterpretableMachineLearning2023, + title = {Interpretable {Machine} {Learning} for {Science} with {PySR} and {SymbolicRegression}.jl}, + url = {http://arxiv.org/abs/2305.01582}, + doi = {10.48550/arXiv.2305.01582}, + urldate = {2023-07-17}, + publisher = {arXiv}, + author = {Cranmer, Miles}, + month = may, + year = {2023}, + note = {arXiv:2305.01582 [astro-ph, physics:physics]}, + keywords = {Astrophysics - Instrumentation and Methods for Astrophysics, Computer Science - Machine Learning, Computer Science - Neural and Evolutionary Computing, Computer Science - Symbolic Computation, Physics - Data Analysis, Statistics and Probability}, } ``` diff --git a/docs/src/index.md b/docs/src/index.md index 73e752247..639558981 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -94,7 +94,7 @@ Pkg.add("DataDrivenSR") ## Contributing - Please refer to the - [SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages](https://github.com/SciML/ColPrac/blob/master/README.md) + [SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages](https://github.com/SciML/ColPrac) for guidance on PRs, issues, and other matters relating to contributing to SciML. - See the [SciML Style Guide](https://github.com/SciML/SciMLStyle) for common coding practices and other style decisions. diff --git a/docs/src/libs/datadrivendmd/example_04.jl b/docs/src/libs/datadrivendmd/example_04.jl index 8add48ee2..a06eac202 100644 --- a/docs/src/libs/datadrivendmd/example_04.jl +++ b/docs/src/libs/datadrivendmd/example_04.jl @@ -1,6 +1,6 @@ # # [Nonlinear Time Continuous System](@id nonlinear_continuos) # -# Similarly, we can use the [Extended Dynamic Mode Decomposition](https://link.springer.com/article/10.1007/s00332-015-9258-5) via a nonlinear [`Basis`](@ref) of observables. Here, we will look at a rather [famous example](https://arxiv.org/pdf/1510.03007.pdf) with a finite dimensional solution. +# Similarly, we can use the [Extended Dynamic Mode Decomposition](https://link.springer.com/article/10.1007/s00332-015-9258-5) via a nonlinear [`Basis`](@ref) of observables. Here, we will look at a rather [famous example](https://arxiv.org/pdf/1510.03007) with a finite dimensional solution. using DataDrivenDiffEq using LinearAlgebra diff --git a/docs/src/libs/datadrivensparse/example_03.jl b/docs/src/libs/datadrivensparse/example_03.jl index 75b67b0a5..027bf010a 100644 --- a/docs/src/libs/datadrivensparse/example_03.jl +++ b/docs/src/libs/datadrivensparse/example_03.jl @@ -1,7 +1,7 @@ # # [Implicit Nonlinear Dynamics : Michaelis-Menten](@id michaelis_menten) # # What if you want to estimate an implicitly defined system of the form ``f(u_t, u, p, t) = 0``? -# The solution : Implicit Sparse Identification. This method was originally described in [this paper](http://ieeexplore.ieee.org/document/7809160/), and currently there exist [robust algorithms](https://royalsocietypublishing.org/doi/10.1098/rspa.2020.0279) to identify these systems. +# The solution : Implicit Sparse Identification. This method was originally described in [this paper](https://ieeexplore.ieee.org/document/7809160/), and currently there exist [robust algorithms](https://royalsocietypublishing.org/doi/10.1098/rspa.2020.0279) to identify these systems. # We will focus on [Michaelis-Menten Kinetics](https://en.wikipedia.org/wiki/Michaelis%E2%80%93Menten_kinetics). As before, we will define the [`DataDrivenProblem`](@ref) and the [`Basis`](@ref) containing possible candidate functions for our [sparse regression](@ref sparse_algorithms). # Let's generate some data! We will use two experiments starting from different initial conditions. diff --git a/docs/src/libs/datadrivensparse/example_04.jl b/docs/src/libs/datadrivensparse/example_04.jl index 031b8903b..48dcbe887 100644 --- a/docs/src/libs/datadrivensparse/example_04.jl +++ b/docs/src/libs/datadrivensparse/example_04.jl @@ -5,6 +5,7 @@ using DataDrivenDiffEq using ModelingToolkit +using ModelingToolkit: t_nounits as t, D_nounits as D using OrdinaryDiffEq using LinearAlgebra using DataDrivenSparse @@ -12,31 +13,25 @@ using DataDrivenSparse #md gr() using Test #src -@parameters begin - t - α = 1.0 - β = 1.3 - γ = 2.0 - δ = 0.5 +@mtkmodel Autoregulation begin + @parameters begin + α = 1.0 + β = 1.3 + γ = 2.0 + δ = 0.5 + end + @variables begin + (x(t))[1:2] = [20.0; 12.0] + end + @equations begin + D(x[1]) ~ α / (1 + x[2]) - β * x[1] + D(x[2]) ~ γ / (1 + x[1]) - δ * x[2] + end end -@variables begin - x[1:2](t) = [20.0; 12.0] -end - -x = collect(x) -D = Differential(t) - -eqs = [D(x[1]) ~ α / (1 + x[2]) - β * x[1]; - D(x[2]) ~ γ / (1 + x[1]) - δ * x[2]] - -sys = ODESystem(eqs, t, x, [α, β, γ, δ], name = :Autoregulation) - -x0 = [x[1] => 20.0; x[2] => 12.0] - +@mtkbuild sys = Autoregulation() tspan = (0.0, 5.0) - -de_problem = ODEProblem{true, SciMLBase.NoSpecialize}(sys, x0, tspan) +de_problem = ODEProblem{true, SciMLBase.NoSpecialize}(sys, [], tspan, []) de_solution = solve(de_problem, Tsit5(), saveat = 0.005); #md plot(de_solution) @@ -44,6 +39,7 @@ de_solution = solve(de_problem, Tsit5(), saveat = 0.005); dd_prob = DataDrivenProblem(de_solution) +x = sys.x eqs = [polynomial_basis(x, 4); D.(x); x .* D(x[1]); x .* D(x[2])] basis = Basis(eqs, x, independent_variable = t, implicits = D.(x)) diff --git a/docs/src/libs/datadrivensr/example_01.jl b/docs/src/libs/datadrivensr/example_01.jl index 7073fc654..35ea046c6 100644 --- a/docs/src/libs/datadrivensr/example_01.jl +++ b/docs/src/libs/datadrivensr/example_01.jl @@ -5,7 +5,7 @@ # Hence, the performance might differ and depends strongly on the hyperparameters of the optimization. # This example might not recover the groundtruth, but is showing off the use within `DataDrivenDiffEq.jl`. # -# DataDrivenDiffEq offers an interface to [`SymbolicRegression.jl`](https://docs.sciml.ai/SymbolicRegression/stable/) to infer more complex functions. +# DataDrivenDiffEq offers an interface to [`SymbolicRegression.jl`](https://ai.damtp.cam.ac.uk/symbolicregression/dev/) to infer more complex functions. using DataDrivenDiffEq using LinearAlgebra diff --git a/lib/DataDrivenDMD/src/algorithms.jl b/lib/DataDrivenDMD/src/algorithms.jl index 25592ad82..9548abe8b 100644 --- a/lib/DataDrivenDMD/src/algorithms.jl +++ b/lib/DataDrivenDMD/src/algorithms.jl @@ -189,7 +189,7 @@ $(TYPEDEF) Approximates the Koopman operator `K` via the forward-backward DMD. -It is assumed that `K = sqrt(K₁*inv(K₂))`, where `K₁` is the approximation via forward and `K₂` via [DMDSVD](@ref). Based on [this paper](https://arxiv.org/pdf/1507.02264.pdf). +It is assumed that `K = sqrt(K₁*inv(K₂))`, where `K₁` is the approximation via forward and `K₂` via [DMDSVD](@ref). Based on [this paper](https://arxiv.org/pdf/1507.02264). If `truncation` ∈ (0, 1) is given, the singular value decomposition is reduced to include only entries bigger than `truncation*maximum(Σ)`. If `truncation` is an integer, the reduced SVD up to `truncation` is used for computation. diff --git a/lib/DataDrivenSparse/src/algorithms/STLSQ.jl b/lib/DataDrivenSparse/src/algorithms/STLSQ.jl index 262447cd5..6c0f80d61 100644 --- a/lib/DataDrivenSparse/src/algorithms/STLSQ.jl +++ b/lib/DataDrivenSparse/src/algorithms/STLSQ.jl @@ -2,7 +2,7 @@ $(TYPEDEF) `STLSQ` is taken from the [original paper on SINDY](https://www.pnas.org/doi/10.1073/pnas.1517384113) and implements a sequentially thresholded least squares iteration. `λ` is the threshold of the iteration. -It is based upon [this Matlab implementation](https://github.com/eurika-kaiser/SINDY-MPC/blob/master/utils/sparsifyDynamics.m). +It is based upon [this Matlab implementation](https://github.com/eurika-kaiser/SINDY-MPC/blob/e1dfd9908b2b56af303ee9fb30a133aced4fd757/utils/sparsifyDynamics.m). It solves the following problem ```math \\argmin_{x} \\frac{1}{2} \\| Ax-b\\|_2 + \\rho \\|x\\|_2 diff --git a/src/DataDrivenDiffEq.jl b/src/DataDrivenDiffEq.jl index 7378a10c0..8e353f629 100644 --- a/src/DataDrivenDiffEq.jl +++ b/src/DataDrivenDiffEq.jl @@ -14,6 +14,7 @@ using Setfield @reexport using ModelingToolkit using ModelingToolkit: AbstractSystem, AbstractTimeDependentSystem +using SciMLStructures: SciMLStructures as SS using SymbolicUtils: operation, arguments, iscall, issym using Symbolics using Symbolics: scalarize, variable, value diff --git a/src/problem/type.jl b/src/problem/type.jl index 9dc6cbd28..39244e474 100644 --- a/src/problem/type.jl +++ b/src/problem/type.jl @@ -152,10 +152,15 @@ function DataDrivenProblem(X::AbstractMatrix; DX::AbstractMatrix = Array{eltype(X)}(undef, 0, 0), Y::AbstractMatrix = Array{eltype(X)}(undef, 0, 0), U::F = Array{eltype(X)}(undef, 0, 0), - p::AbstractVector = Array{eltype(X)}(undef, 0), + p::Union{AbstractVector, MTKParameters} = Array{eltype(X)}(undef, 0), probtype = nothing, kwargs...) where {F <: Union{AbstractMatrix, Function}} - return DataDrivenProblem(probtype, X, t, DX, Y, U, p; kwargs...) + if SS.isscimlstructure(p) + _p, _, _ = SS.canonicalize(SS.Tunable(), p) + else + _p = p + end + return DataDrivenProblem(probtype, X, t, DX, Y, U, _p; kwargs...) end function Base.summary(io::IO, x::DataDrivenProblem{N, C, P}) where {N, C, P} @@ -313,7 +318,7 @@ function ModelingToolkit.get_dvs(p::ABSTRACT_DIRECT_PROB, i = :, j = :) end function ModelingToolkit.get_dvs(p::ABSTRACT_DISCRETE_PROB, i = :, j = :) - ModelingToolkit.states(p, i, j) + states(p, i, j) end function ModelingToolkit.observed(p::AbstractDataDrivenProblem, i = :, j = :) @@ -327,7 +332,7 @@ function ModelingToolkit.controls(p::AbstractDataDrivenProblem, i = :, j = :) end function Base.getindex(p::AbstractDataDrivenProblem, i = :, j = :) - return (ModelingToolkit.states(p, i, j), + return (states(p, i, j), ModelingToolkit.parameters(p), ModelingToolkit.independent_variable(p, j), ModelingToolkit.controls(p, i, j)) diff --git a/src/utils/collocation.jl b/src/utils/collocation.jl index 70daa5fe8..9448deacc 100644 --- a/src/utils/collocation.jl +++ b/src/utils/collocation.jl @@ -142,7 +142,7 @@ $(SIGNATURES) Unified interface for collocation techniques. The input can either be a `CollocationKernel` (see list below) or a wrapped `InterpolationMethod` from -[DataInterpolations.jl](https://github.com/PumasAI/DataInterpolations.jl). +[DataInterpolations.jl](https://github.com/SciML/DataInterpolations.jl). Computes a non-parametrically smoothed estimate of `u'` and `u` given the `data`, where each column is a snapshot of the timeseries at @@ -156,7 +156,7 @@ u′,u,t = collocate_data(data,tpoints,interp) ``` # Collocation Kernels -See [this paper](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2631937/) for more information. +See [this paper](https://pmc.ncbi.nlm.nih.gov/articles/PMC2631937/) for more information. + EpanechnikovKernel + UniformKernel + TriangularKernel @@ -170,7 +170,7 @@ See [this paper](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2631937/) for more + SilvermanKernel # Interpolation Methods -See [DataInterpolations.jl](https://github.com/PumasAI/DataInterpolations.jl) for more information. +See [DataInterpolations.jl](https://github.com/SciML/DataInterpolations.jl) for more information. + ConstantInterpolation + LinearInterpolation + QuadraticInterpolation diff --git a/test/problem/problem.jl b/test/problem/problem.jl index bbe48caec..2611ce796 100644 --- a/test/problem/problem.jl +++ b/test/problem/problem.jl @@ -162,3 +162,31 @@ end prob2 = (X = X, t = t, Y = Y)) @test_throws ArgumentError ContinuousDataset(wrong_data) end + +@testset "DataDrivenProblem from ODEProblem solution" begin + using OrdinaryDiffEqTsit5 + using ModelingToolkit: t_nounits as time, D_nounits as D + + @mtkmodel Autoregulation begin + @parameters begin + α = 1.0 + β = 1.3 + γ = 2.0 + δ = 0.5 + end + @variables begin + (x(time))[1:2] = [20.0; 12.0] + end + @equations begin + D(x[1]) ~ α / (1 + x[2]) - β * x[1] + D(x[2]) ~ γ / (1 + x[1]) - δ * x[2] + end + end + + @mtkbuild sys = Autoregulation() + tspan = (0.0, 5.0) + de_problem = ODEProblem{true, SciMLBase.NoSpecialize}(sys, [], tspan, []) + de_solution = solve(de_problem, Tsit5(), saveat = 0.005) + prob = DataDrivenProblem(de_solution) + @test is_valid(prob) +end