diff --git a/Project.toml b/Project.toml index c3d1e59..2cd9c41 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "DiffEqParamEstim" uuid = "1130ab10-4a5a-5621-a13d-e4788d82bd4c" -version = "2.4.0" +version = "2.5.0" authors = ["Chris Rackauckas "] [deps] @@ -21,13 +21,13 @@ StatsAPI = "82ae8749-77ed-4fe6-ae5f-f523153014b0" Calculus = "0.5" CommonSolve = "0.2.6" Dierckx = "0.4, 0.5" -DiffEqBase = "6" +DiffEqBase = "6, 7" Distributions = "0.25" ForwardDiff = "0.10" PenaltyFunctions = "0.1, 0.2, 0.3" PreallocationTools = "0.2, 0.3, 0.4, 1.0" RecursiveArrayTools = "1.0, 2.0, 3, 4" -SciMLBase = "1.69, 2, 3.1" +SciMLBase = "1.69, 2, 3" SciMLSensitivity = "7" Statistics = "1.10" StatsAPI = "1.8.0" diff --git a/src/cost_functions.jl b/src/cost_functions.jl index 8edfd88..e329893 100644 --- a/src/cost_functions.jl +++ b/src/cost_functions.jl @@ -40,7 +40,7 @@ function (f::L2Loss)(sol::DiffEqBase.AbstractNoTimeSolution) dudt = f.dudt if sol isa DiffEqBase.AbstractEnsembleSolution - failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol) + failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol.u) else failure = !SciMLBase.successful_retcode(sol.retcode) end @@ -48,16 +48,17 @@ function (f::L2Loss)(sol::DiffEqBase.AbstractNoTimeSolution) sumsq = 0.0 + solu = sol.u if weight === nothing - @inbounds for i in 1:length(sol) - sumsq += (data[i] - sol[i])^2 + @inbounds for i in 1:length(solu) + sumsq += (data[i] - solu[i])^2 end else - @inbounds for i in 1:length(sol) + @inbounds for i in 1:length(solu) if weight isa Real - sumsq = sumsq + ((data[i] - sol[i])^2) * weight + sumsq = sumsq + ((data[i] - solu[i])^2) * weight else - sumsq = sumsq + ((data[i] - sol[i])^2) * weight[i] + sumsq = sumsq + ((data[i] - solu[i])^2) * weight[i] end end end @@ -72,7 +73,7 @@ function (f::L2Loss)(sol::SciMLBase.AbstractSciMLSolution) dudt = f.dudt if sol isa DiffEqBase.AbstractEnsembleSolution - failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol) + failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol.u) else failure = !SciMLBase.successful_retcode(sol.retcode) end @@ -80,8 +81,9 @@ function (f::L2Loss)(sol::SciMLBase.AbstractSciMLSolution) sumsq = 0.0 + nsteps = length(sol.u) if weight === nothing - @inbounds for i in 1:length(sol) + @inbounds for i in 1:nsteps for j in 1:length(sol.u[i]) sumsq += (data[j, i] - sol[j, i])^2 end @@ -108,7 +110,7 @@ function (f::L2Loss)(sol::SciMLBase.AbstractSciMLSolution) end end else - @inbounds for i in 1:length(sol) + @inbounds for i in 1:nsteps if weight isa Real for j in 1:length(sol.u[i]) sumsq = sumsq + ((data[j, i] - sol[j, i])^2) * weight @@ -195,7 +197,7 @@ end function (f::LogLikeLoss)(sol::SciMLBase.AbstractSciMLSolution) distributions = f.data_distributions if sol isa DiffEqBase.AbstractEnsembleSolution - failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol) + failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol.u) else failure = !SciMLBase.successful_retcode(sol.retcode) end @@ -244,7 +246,7 @@ end function (f::LogLikeLoss)(sol::DiffEqBase.AbstractEnsembleSolution) distributions = f.data_distributions - failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol) + failure = any(!SciMLBase.successful_retcode(s.retcode) for s in sol.u) failure && return Inf ll = 0.0 if eltype(distributions) <: UnivariateDistribution @@ -252,7 +254,7 @@ function (f::LogLikeLoss)(sol::DiffEqBase.AbstractEnsembleSolution) # i is the number of time points # j is the size of the system # corresponds to distributions[i,j] - vals = [s[i, j] for s in sol] + vals = [s[i, j] for s in sol.u] ll -= loglikelihood(distributions[i, j], vals) end else @@ -260,7 +262,7 @@ function (f::LogLikeLoss)(sol::DiffEqBase.AbstractEnsembleSolution) # i is the number of time points # j is the size of the system # corresponds to distributions[i,j] - vals = [s[i, j] for i in 1:length(sol.u[1].u[1]), s in sol] + vals = [s[i, j] for i in 1:length(sol.u[1].u[1]), s in sol.u] ll -= loglikelihood(distributions[j], vals) end end @@ -271,12 +273,12 @@ function (f::LogLikeLoss)(sol::DiffEqBase.AbstractEnsembleSolution) if eltype(distributions) <: UnivariateDistribution for j in 2:length(f.t), i in 1:length(sol.u[1].u[1]) - vals = [s[i, j] - s[i, j - 1] for s in sol] + vals = [s[i, j] - s[i, j - 1] for s in sol.u] fdll -= logpdf(distributions[j - 1, i], vals)[1] end else for j in 2:length(f.t) - vals = [s[i, j] - s[i, j - 1] for i in 1:length(sol.u[1].u[1]), s in sol] + vals = [s[i, j] - s[i, j - 1] for i in 1:length(sol.u[1].u[1]), s in sol.u] fdll -= logpdf(distributions[j - 1], vals)[1] end end diff --git a/test/tests_on_odes/regularization_test.jl b/test/tests_on_odes/regularization_test.jl index 54f2e4b..039672f 100644 --- a/test/tests_on_odes/regularization_test.jl +++ b/test/tests_on_odes/regularization_test.jl @@ -36,8 +36,8 @@ result = solve(optprob, Optim.BFGS()) optprob = Optimization.OptimizationProblem(cost_function_2, [1.2, 2.7]) result = solve(optprob, Optim.BFGS()) -@test result.minimizer ≈ [1.5; 3.0] atol = 3.0e-1 +@test result.u ≈ [1.5; 3.0] atol = 3.0e-1 optprob = Optimization.OptimizationProblem(cost_function_3, [1.3, 0.8, 2.8, 1.2]) result = solve(optprob, Optim.BFGS()) -@test result.minimizer ≈ [1.5; 1.0; 3.0; 1.0] atol = 5.0e-1 +@test result.u ≈ [1.5; 1.0; 3.0; 1.0] atol = 5.0e-1