Commit f7c87a98 authored by Javier's avatar Javier

LsqFit -> Optim

fit_routine uses Optim
LsqFit deps removed.
parent 1e233c9a
...@@ -153,12 +153,6 @@ version = "1.0.2" ...@@ -153,12 +153,6 @@ version = "1.0.2"
deps = ["Random", "Serialization", "Sockets"] deps = ["Random", "Serialization", "Sockets"]
uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"
[[Distributions]]
deps = ["FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "StaticArrays", "Statistics", "StatsBase", "StatsFuns"]
git-tree-sha1 = "501c11d708917ca09ce357bed163dbaf0f30229f"
uuid = "31c24e10-a181-5473-b8eb-7969acd0382f"
version = "0.23.12"
[[DocStringExtensions]] [[DocStringExtensions]]
deps = ["LibGit2", "Markdown", "Pkg", "Test"] deps = ["LibGit2", "Markdown", "Pkg", "Test"]
git-tree-sha1 = "50ddf44c53698f5e784bbebb3f4b21c5807401b1" git-tree-sha1 = "50ddf44c53698f5e784bbebb3f4b21c5807401b1"
...@@ -460,12 +454,6 @@ uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" ...@@ -460,12 +454,6 @@ uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
[[Logging]] [[Logging]]
uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"
[[LsqFit]]
deps = ["Distributions", "ForwardDiff", "LinearAlgebra", "NLSolversBase", "OptimBase", "Random", "StatsBase"]
git-tree-sha1 = "b32b5549461fcb93bce223e264d4a7ef0c9923fd"
uuid = "2fda8390-95c7-5789-9bda-21331edee243"
version = "0.11.0"
[[MKL_jll]] [[MKL_jll]]
deps = ["IntelOpenMP_jll", "Libdl", "Pkg"] deps = ["IntelOpenMP_jll", "Libdl", "Pkg"]
git-tree-sha1 = "eb540ede3aabb8284cb482aa41d00d6ca850b1f8" git-tree-sha1 = "eb540ede3aabb8284cb482aa41d00d6ca850b1f8"
...@@ -555,12 +543,6 @@ git-tree-sha1 = "8a8208d3a1b97994d15ebcdae235d332b4a69e78" ...@@ -555,12 +543,6 @@ git-tree-sha1 = "8a8208d3a1b97994d15ebcdae235d332b4a69e78"
uuid = "429524aa-4258-5aef-a3af-852621145aeb" uuid = "429524aa-4258-5aef-a3af-852621145aeb"
version = "1.2.3" version = "1.2.3"
[[OptimBase]]
deps = ["NLSolversBase", "Printf", "Reexport"]
git-tree-sha1 = "4c26a757fbb5b1893b7df19a44e21762d8f8e470"
uuid = "87e2bd06-a317-5318-96d9-3ecbac512eee"
version = "2.0.1"
[[Opus_jll]] [[Opus_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "f9d57f4126c39565e05a2b0264df99f497fc6f37" git-tree-sha1 = "f9d57f4126c39565e05a2b0264df99f497fc6f37"
...@@ -578,12 +560,6 @@ git-tree-sha1 = "1b556ad51dceefdbf30e86ffa8f528b73c7df2bb" ...@@ -578,12 +560,6 @@ git-tree-sha1 = "1b556ad51dceefdbf30e86ffa8f528b73c7df2bb"
uuid = "2f80f16e-611a-54ab-bc61-aa92de5b98fc" uuid = "2f80f16e-611a-54ab-bc61-aa92de5b98fc"
version = "8.42.0+4" version = "8.42.0+4"
[[PDMats]]
deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse", "Test"]
git-tree-sha1 = "95a4038d1011dfdbde7cecd2ad0ac411e53ab1bc"
uuid = "90014a1f-27ba-587c-ab20-58faa44d9150"
version = "0.10.1"
[[PGFPlotsX]] [[PGFPlotsX]]
deps = ["ArgCheck", "DataStructures", "Dates", "DefaultApplication", "DocStringExtensions", "MacroTools", "Parameters", "Requires", "Tables"] deps = ["ArgCheck", "DataStructures", "Dates", "DefaultApplication", "DocStringExtensions", "MacroTools", "Parameters", "Requires", "Tables"]
git-tree-sha1 = "1adde3d07cce96b6a3bb88572612db4bd9d6153b" git-tree-sha1 = "1adde3d07cce96b6a3bb88572612db4bd9d6153b"
...@@ -695,18 +671,6 @@ git-tree-sha1 = "cfbac6c1ed70c002ec6361e7fd334f02820d6419" ...@@ -695,18 +671,6 @@ git-tree-sha1 = "cfbac6c1ed70c002ec6361e7fd334f02820d6419"
uuid = "ae029012-a4dd-5104-9daa-d747884805df" uuid = "ae029012-a4dd-5104-9daa-d747884805df"
version = "1.1.2" version = "1.1.2"
[[Rmath]]
deps = ["Random", "Rmath_jll"]
git-tree-sha1 = "86c5647b565873641538d8f812c04e4c9dbeb370"
uuid = "79098fc4-a85e-5d69-aa6a-4863f24498fa"
version = "0.6.1"
[[Rmath_jll]]
deps = ["Libdl", "Pkg"]
git-tree-sha1 = "d76185aa1f421306dec73c057aa384bad74188f0"
uuid = "f50d1b31-88e8-58de-be2c-1cc44531875f"
version = "0.2.2+1"
[[Roots]] [[Roots]]
deps = ["Printf"] deps = ["Printf"]
git-tree-sha1 = "8f743e4f4368d1d753f3806bf635899dad6b4847" git-tree-sha1 = "8f743e4f4368d1d753f3806bf635899dad6b4847"
...@@ -770,22 +734,12 @@ git-tree-sha1 = "7bab7d4eb46b225b35179632852b595a3162cb61" ...@@ -770,22 +734,12 @@ git-tree-sha1 = "7bab7d4eb46b225b35179632852b595a3162cb61"
uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
version = "0.33.2" version = "0.33.2"
[[StatsFuns]]
deps = ["Rmath", "SpecialFunctions"]
git-tree-sha1 = "3b9f665c70712af3264b61c27a7e1d62055dafd1"
uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c"
version = "0.9.6"
[[StructArrays]] [[StructArrays]]
deps = ["DataAPI", "Tables"] deps = ["DataAPI", "Tables"]
git-tree-sha1 = "ad1f5fd155426dcc879ec6ede9f74eb3a2d582df" git-tree-sha1 = "ad1f5fd155426dcc879ec6ede9f74eb3a2d582df"
uuid = "09ab397b-f2b6-538f-b94a-2f83cf4a842a" uuid = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
version = "0.4.2" version = "0.4.2"
[[SuiteSparse]]
deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"]
uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9"
[[TableTraits]] [[TableTraits]]
deps = ["IteratorInterfaceExtensions"] deps = ["IteratorInterfaceExtensions"]
git-tree-sha1 = "b1ad568ba658d8cbb3b892ed5380a6f3e781a81e" git-tree-sha1 = "b1ad568ba658d8cbb3b892ed5380a6f3e781a81e"
......
...@@ -8,7 +8,6 @@ ADerrors = "5e92007d-7bf1-471c-8ceb-4591b8b567a9" ...@@ -8,7 +8,6 @@ ADerrors = "5e92007d-7bf1-471c-8ceb-4591b8b567a9"
BDIO = "375f315e-f2c4-11e9-2ef9-134f02f79e27" BDIO = "375f315e-f2c4-11e9-2ef9-134f02f79e27"
LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f" LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
LsqFit = "2fda8390-95c7-5789-9bda-21331edee243"
Optim = "429524aa-4258-5aef-a3af-852621145aeb" Optim = "429524aa-4258-5aef-a3af-852621145aeb"
PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee" PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee"
Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
module juobs module juobs
using ADerrors, PyPlot, LaTeXStrings, LinearAlgebra, LsqFit, Optim using ADerrors, PyPlot, LaTeXStrings, LinearAlgebra, Optim
import Statistics: mean import Statistics: mean
include("juobs_types.jl") include("juobs_types.jl")
......
...@@ -256,14 +256,17 @@ function fit_routine(model::Function, xdata::Array{<:Real}, ydata::Array{uwreal} ...@@ -256,14 +256,17 @@ function fit_routine(model::Function, xdata::Array{<:Real}, ydata::Array{uwreal}
yval = value.(ydata) yval = value.(ydata)
yer = err.(ydata) yer = err.(ydata)
chisq = gen_chisq(model, xdata, yer) chisq = gen_chisq(model, xdata, yer)
fit = curve_fit(model, xdata, yval, 1.0 ./ yer.^2, fill(0.5, param)) min_fun(t) = chisq(t, yval)
(upar, chi_exp) = isnothing(wpm) ? fit_error(chisq, coef(fit), ydata) : fit_error(chisq, coef(fit), ydata, wpm) sol = optimize(min_fun, fill(0.5, param), method=LBFGS())
par = Optim.minimizer(sol)
(upar, chi_exp) = isnothing(wpm) ? fit_error(chisq, par, ydata) : fit_error(chisq, par, ydata, wpm)
for i = 1:length(upar) for i = 1:length(upar)
isnothing(wpm) ? uwerr(upar[i]) : uwerr(upar[i], wpm) isnothing(wpm) ? uwerr(upar[i]) : uwerr(upar[i], wpm)
print("\n Fit parameter: ", i, ": ") print("\n Fit parameter: ", i, ": ")
details(upar[i]) details(upar[i])
end end
println("Chisq / chiexp: ", chisq(coef(fit), ydata), " / ", chi_exp, " (dof: ", dof(fit),")") println("Chisq / chiexp: ", sol.minimum, " / ", chi_exp, " (dof: ", length(yval) - param,")")
return upar return upar
end end
function gen_chisq(f::Function, x::Array{<:Real}, err::Vector{Float64}) function gen_chisq(f::Function, x::Array{<:Real}, err::Vector{Float64})
...@@ -293,19 +296,22 @@ function fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal} ...@@ -293,19 +296,22 @@ function fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}
data = vcat(xdata, ydata) data = vcat(xdata, ydata)
#guess #guess
fit = curve_fit(model, xval, yval, 1.0 ./ yer.^2, fill(0.5, param)) chisq = gen_chisq(model, xval, yer)
min_fun_cons(t) = chisq(t, yval)
sol_cons = optimize(min_fun_cons, fill(0.5, param), method=LBFGS())
par_cons = Optim.minimizer(sol_cons)
if covar if covar
C = [ADerrors.cov([xdata[k], ydata[k]]) for k = 1:length(ydata)] C = [ADerrors.cov([xdata[k], ydata[k]]) for k = 1:length(ydata)]
chisq_full_cov(p, d) = get_chi2_cov(model, d, C, p) chisq_full_cov(p, d) = get_chi2_cov(model, d, C, p)
min_fun_cov(t) = chisq_full_cov(t, dat) min_fun_cov(t) = chisq_full_cov(t, dat)
sol = optimize(min_fun_cov, vcat(fit.param, xval), method=LBFGS()) sol = optimize(min_fun_cov, vcat(par_cons, xval), method=LBFGS())
(upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full_cov, Optim.minimizer(sol), data) : fit_error(chisq_full_cov, Optim.minimizer(sol), data, wpm) (upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full_cov, Optim.minimizer(sol), data) : fit_error(chisq_full_cov, Optim.minimizer(sol), data, wpm)
else else
chisq_full(p, d) = get_chi2(model, d, ddat, p) chisq_full(p, d) = get_chi2(model, d, ddat, p)
min_fun(t) = chisq_full(t, dat) min_fun(t) = chisq_full(t, dat)
sol = optimize(min_fun, vcat(fit.param, xval), method=LBFGS()) sol = optimize(min_fun, vcat(par_cons, xval), method=LBFGS())
(upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full, Optim.minimizer(sol), data) : fit_error(chisq_full, Optim.minimizer(sol), data, wpm) (upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full, Optim.minimizer(sol), data) : fit_error(chisq_full, Optim.minimizer(sol), data, wpm)
end end
...@@ -317,7 +323,7 @@ function fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal} ...@@ -317,7 +323,7 @@ function fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}
print("\n Fit parameter: ", i, ": ") print("\n Fit parameter: ", i, ": ")
details(upar[i]) details(upar[i])
end end
println("Chisq / chiexp: ", sol.minimum, " / ", chi2_exp, " (dof: ", dof(fit),")") println("Chisq / chiexp: ", sol.minimum, " / ", chi2_exp, " (dof: ", length(ydata) - param,")")
return upar return upar
end end
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment