diff --git a/examples/anasol/runtests.jl b/examples/anasol/runtests.jl index 50788f3e..ee9dca38 100644 --- a/examples/anasol/runtests.jl +++ b/examples/anasol/runtests.jl @@ -41,7 +41,7 @@ if !haskey(ENV, "MADS_NO_GADFLY") Mads.spaghettiplots(md, 2) Mads.rmfile(joinpath(workdir, "w01short-ax-2-spaghetti.svg")) Mads.rmfile(joinpath(workdir, "w01short-vx-2-spaghetti.svg")) - Mads.plotmass([1.,2.],[10.,20.],[5.,10.], joinpath(workdir, "plotmass-test.svg")) + Mads.plotmass([1.0, 2.0], [10.0, 20.0], [5.0, 10.0], joinpath(workdir, "plotmass-test.svg")) Mads.rmfile(joinpath(workdir, "plotmass-test.svg")) Mads.graphon() end @@ -99,7 +99,7 @@ Mads.wellon!(md, "w1a") Mads.savemadsfile(md, joinpath(workdir, "test.mads")) Mads.savemadsfile(md, Mads.getparamdict(md)) Mads.savemadsfile(md, Mads.getparamdict(md), joinpath(workdir, "test.mads")) -Mads.savemadsfile(md, Mads.getparamdict(md), joinpath(workdir, "test.mads"), explicit=true) +Mads.savemadsfile(md, Mads.getparamdict(md), joinpath(workdir, "test.mads"); explicit=true) Mads.rmfile(joinpath(workdir, "w01shortexp-rerun.mads")) Mads.rmfile(joinpath(workdir, "test.mads")) Mads.setmadsinputfile("test.mads") @@ -128,7 +128,7 @@ mdbad = deepcopy(md) Mads.setparamsdistnormal!(mdbad, fill(1, length(m1)), fill(2, length(m1))) Mads.setparamsdistuniform!(mdbad, fill(1, length(m1)), fill(2, length(m1))) -Mads.computemass("w01lambda", time=50, path=workdir) +Mads.computemass("w01lambda"; time=50, path=workdir) Mads.rmdir(joinpath(workdir, "mass_reduced.svg")) d = joinpath(workdir, "test_results") @@ -163,35 +163,35 @@ good_targetkeys = JLD2.load(joinpath(d, "targetkeys.jld2"), "tk") Test.@testset "Anasol" begin # Test Mads.forward(md; all=true) - ssr = 0. + ssr = 0.0 for obskey in union(Set(keys(forward_results)), Set(keys(good_forward_results))) ssr += (forward_results[obskey] - good_forward_results[obskey])^2 end - Test.@test isapprox(ssr, 0., atol=1e-8) + Test.@test isapprox(ssr, 0.0, atol=1e-8) - # Test Mads.getparamrandom(md, 5, init_dist=true) + # Test Mads.getparamrandom(md, 5; init_dist=true) for obskey in union(Set(keys(paramrandom_true)), Set(keys(good_paramrandom_true))) Test.@test isapprox(paramrandom_true[obskey], good_paramrandom_true[obskey], atol=1e-8) end - # Test Mads.getparamrandom(md, 5, init_dist=false) + # Test Mads.getparamrandom(md, 5; init_dist=false) for obskey in union(Set(keys(paramrandom_false)), Set(keys(good_paramrandom_false))) Test.@test isapprox(paramrandom_false[obskey], good_paramrandom_false[obskey], atol=1e-8) end # Test Mads.getparamdict(md) - ssr = 0. + ssr = 0.0 for obskey in union(Set(keys(pd)), Set(keys(good_pd))) ssr += (pd[obskey] - good_pd[obskey])^2 end - Test.@test isapprox(ssr, 0., atol=1e-8) + Test.@test isapprox(ssr, 0.0, atol=1e-8) # Test computeconcentrations(paramdict) - ssr = 0. + ssr = 0.0 for obskey in union(Set(keys(forward_preds)), Set(keys(good_forward_preds))) ssr += (forward_preds[obskey] - good_forward_preds[obskey])^2 end - Test.@test isapprox(ssr, 0., atol=1e-8) + Test.@test isapprox(ssr, 0.0, atol=1e-8) # Test Mads.of(md) Test.@test isapprox(madsOf, 628963.6972820368, atol=1e-8) #test diff --git a/examples/model_analysis/bayes_weight_analsis.jl b/examples/model_analysis/bayes_weight_analsis.jl index edba1a2e..4cfebb4a 100644 --- a/examples/model_analysis/bayes_weight_analsis.jl +++ b/examples/model_analysis/bayes_weight_analsis.jl @@ -7,29 +7,29 @@ Mads.mkdir("bayes_results") @info("Bayesian analysis with different initial parameter sets and observation weights (standard deviation errors)") @info("Bayesian analysis for the initial parameter guesses:") -for w = (1000000, 1000, 1) +for w in (1000000, 1000, 1) Mads.setobsweights!(md, w) mcmcchain = Mads.bayessampling(md; nsteps=1000000, burnin=100000, thinning=100, seed=2016) Mads.scatterplotsamples(md, permutedims(mcmcchain.value), "bayes_results/bayes_init_w$w.png") o = Mads.forward(md, mcmcchain.value) - Mads.spaghettiplot(md, o, filename="bayes_results/bayes_init_w$(w)_spaghetti.png") - Printf.@printf "Init: Observation Weight %d StdDev %f ->`o5` prediction: min = %f max = %f\n" w 1/w min(o[:,5]...) max(o[:,5]...) - f = Gadfly.plot(x=o[:,5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) + Mads.spaghettiplot(md, o; filename="bayes_results/bayes_init_w$(w)_spaghetti.png") + Printf.@printf "Init: Observation Weight %d StdDev %f ->`o5` prediction: min = %f max = %f\n" w 1 / w min(o[:, 5]...) max(o[:, 5]...) + f = Gadfly.plot(; x=o[:, 5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) Gadfly.draw(Gadfly.PNG("bayes_results/bayes_init_w$(w)_o5.png", 6Gadfly.inch, 4Gadfly.inch), f) end pinit = Dict(zip(Mads.getparamkeys(md), Mads.getparamsinit(md))) n = 100 @info("Calibration using $n random initial guesses for model parameters") -r = Mads.calibraterandom(md, n, all=true, seed=2016, save_results=false) -pnames = collect(keys(r[1,3])) -p = hcat(map(i->collect(values(r[i,3])), 1:n)...)' +r = Mads.calibraterandom(md, n; all=true, seed=2016, save_results=false) +pnames = collect(keys(r[1, 3])) +p = hcat(map(i -> collect(values(r[i, 3])), 1:n)...)' np = length(pnames) @info("Identify the 3 different global optima with different model parameter estimates") -ind_n0 = abs(p[:,4]) .< 0.1 +ind_n0 = abs(p[:, 4]) .< 0.1 in0 = findall(ind_n0 .== true)[1] -ind_n1 = abs(p[:,4]-1) .< 0.1 +ind_n1 = abs(p[:, 4] - 1) .< 0.1 in1 = findall(ind_n1 .== true)[1] ind_n01 = !(ind_n0 | ind_n1) in01 = findall(ind_n01 .== true)[1] @@ -39,15 +39,15 @@ v = [in0, in1, in01] @info("Bayesian analysis for the 3 different global optima") for i = 1:3 - Mads.setparamsinit!(md, r[v[i],3]) - for w = (1000000, 1000, 1) + Mads.setparamsinit!(md, r[v[i], 3]) + for w in (1000000, 1000, 1) Mads.setobsweights!(md, w) mcmcchain = Mads.bayessampling(md; nsteps=1000000, burnin=100000, thinning=100, seed=2016) Mads.scatterplotsamples(md, permutedims(mcmcchain.value), "bayes_results/bayes_opt_$(optnames[i])_w$w.png") o = Mads.forward(md, mcmcchain.value) - Mads.spaghettiplot(md, o, filename="bayes_results/bayes_opt_$(optnames[i])_w$(w)_spaghetti.png") - Printf.@printf "O%-3s: Observation Weight %d StdDev %f -> `o5` prediction: min = %f max = %f\n" optnames[i] w 1/w min(o[:,5]...) max(o[:,5]...) - f = Gadfly.plot(x=o[:,5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) + Mads.spaghettiplot(md, o; filename="bayes_results/bayes_opt_$(optnames[i])_w$(w)_spaghetti.png") + Printf.@printf "O%-3s: Observation Weight %d StdDev %f -> `o5` prediction: min = %f max = %f\n" optnames[i] w 1 / w min(o[:, 5]...) max(o[:, 5]...) + f = Gadfly.plot(; x=o[:, 5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) Gadfly.draw(Gadfly.PNG("bayes_results/bayes_opt_$(optnames[i])_w$(w)_o5.png", 6Gadfly.inch, 4Gadfly.inch), f) end end diff --git a/examples/model_analysis/emcee_walkers_analsis.jl b/examples/model_analysis/emcee_walkers_analsis.jl index 946b45b8..7abba4c9 100644 --- a/examples/model_analysis/emcee_walkers_analsis.jl +++ b/examples/model_analysis/emcee_walkers_analsis.jl @@ -6,13 +6,13 @@ md = Mads.loadmadsfile(joinpath("models", "internal-polynomial.mads")) Mads.mkdir("emcee_results") @info("AffineInvariantMCMC (EMCEE) Bayesian analysis with different number of walkers:") -for nw = (100, 50, 20, 10, 5, 1) +for nw in (100, 50, 20, 10, 5, 1) Mads.setobsweights!(md, 1000000) @time chain, llhoods = Mads.emceesampling(md; numwalkers=nw, nsteps=1000000, burnin=100000, thinning=100, seed=2016) Mads.scatterplotsamples(md, chain', "emcee_results/emcee_init_nw$nw.png") o = Mads.forward(md, chain') - Mads.spaghettiplot(md, o, filename="emcee_results/emcee_init_nw$(nw)_spaghetti.png") - Printf.@printf "Init: Number of walkers %d ->`o5` prediction: min = %f max = %f\n" nw min(o[:,5]...) max(o[:,5]...) - f = Gadfly.plot(x=o[:,5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) + Mads.spaghettiplot(md, o; filename="emcee_results/emcee_init_nw$(nw)_spaghetti.png") + Printf.@printf "Init: Number of walkers %d ->`o5` prediction: min = %f max = %f\n" nw min(o[:, 5]...) max(o[:, 5]...) + f = Gadfly.plot(; x=o[:, 5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) Gadfly.draw(Gadfly.PNG("emcee_results/emcee_init_nw$(nw)_o5.png", 6Gadfly.inch, 4Gadfly.inch), f) -end \ No newline at end of file +end diff --git a/examples/model_analysis/emcee_weight_analsis.jl b/examples/model_analysis/emcee_weight_analsis.jl index 6ac25c42..ab59e4dd 100644 --- a/examples/model_analysis/emcee_weight_analsis.jl +++ b/examples/model_analysis/emcee_weight_analsis.jl @@ -7,29 +7,29 @@ Mads.mkdir("emcee_results") @info("AffineInvariantMCMC (EMCEE) Bayesian analysis with different initial parameter sets and observation weights (standard deviation errors)") @info("AffineInvariantMCMC (EMCEE) Bayesian analysis for the initial parameter guesses:") -for w = (1000000, 1000, 1) +for w in (1000000, 1000, 1) Mads.setobsweights!(md, w) chain, llhoods = Mads.emceesampling(md; numwalkers=100, nsteps=1000000, burnin=100000, thinning=100, seed=2016) Mads.scatterplotsamples(md, chain', "emcee_results/emcee_init_w$w.png") o = Mads.forward(md, chain') - Mads.spaghettiplot(md, o, filename="emcee_results/emcee_init_w$(w)_spaghetti.png") - Printf.@printf "Init: Observation Weight %d StdDev %f ->`o5` prediction: min = %f max = %f\n" w 1/w min(o[:,5]...) max(o[:,5]...) - f = Gadfly.plot(x=o[:,5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) + Mads.spaghettiplot(md, o; filename="emcee_results/emcee_init_w$(w)_spaghetti.png") + Printf.@printf "Init: Observation Weight %d StdDev %f ->`o5` prediction: min = %f max = %f\n" w 1 / w min(o[:, 5]...) max(o[:, 5]...) + f = Gadfly.plot(; x=o[:, 5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) Gadfly.draw(Gadfly.PNG("emcee_results/emcee_init_w$(w)_o5.png", 6Gadfly.inch, 4Gadfly.inch), f) end pinit = Dict(zip(Mads.getparamkeys(md), Mads.getparamsinit(md))) n = 100 @info("Calibration using $n random initial guesses for model parameters") -r = Mads.calibraterandom(md, n, all=true, seed=2016, save_results=false) -pnames = collect(keys(r[1,3])) -p = hcat(map(i->collect(values(r[i,3])), 1:n)...)' +r = Mads.calibraterandom(md, n; all=true, seed=2016, save_results=false) +pnames = collect(keys(r[1, 3])) +p = hcat(map(i -> collect(values(r[i, 3])), 1:n)...)' np = length(pnames) @info("Identify the 3 different global optima with different model parameter estimates") -ind_n0 = abs(p[:,4]) .< 0.1 +ind_n0 = abs(p[:, 4]) .< 0.1 in0 = findall(ind_n0 .== true)[1] -ind_n1 = abs(p[:,4]-1) .< 0.1 +ind_n1 = abs(p[:, 4] - 1) .< 0.1 in1 = findall(ind_n1 .== true)[1] ind_n01 = !(ind_n0 | ind_n1) in01 = findall(ind_n01 .== true)[1] @@ -39,15 +39,15 @@ v = [in0, in1, in01] @info("AffineInvariantMCMC (EMCEE) Bayesian analysis for the 3 different global optima") for i = 1:3 - Mads.setparamsinit!(md, r[v[i],3]) - for w = (1000000, 1000, 1) + Mads.setparamsinit!(md, r[v[i], 3]) + for w in (1000000, 1000, 1) Mads.setobsweights!(md, w) chain, llhoods = Mads.emceesampling(md; numwalkers=100, nsteps=1000000, burnin=100000, thinning=100, seed=2016) Mads.scatterplotsamples(md, chain', "emcee_results/emcee_opt_$(optnames[i])_w$w.png") o = Mads.forward(md, chain') - Mads.spaghettiplot(md, o, filename="emcee_results/emcee_opt_$(optnames[i])_w$(w)_spaghetti.png") - Printf.@printf "O%-3s: Observation Weight %d StdDev %f -> `o5` prediction: min = %f max = %f\n" optnames[i] w 1/w min(o[:,5]...) max(o[:,5]...) - f = Gadfly.plot(x=o[:,5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) + Mads.spaghettiplot(md, o; filename="emcee_results/emcee_opt_$(optnames[i])_w$(w)_spaghetti.png") + Printf.@printf "O%-3s: Observation Weight %d StdDev %f -> `o5` prediction: min = %f max = %f\n" optnames[i] w 1 / w min(o[:, 5]...) max(o[:, 5]...) + f = Gadfly.plot(; x=o[:, 5], Gadfly.Guide.xlabel("o5"), Gadfly.Geom.histogram()) Gadfly.draw(Gadfly.PNG("emcee_results/emcee_opt_$(optnames[i])_w$(w)_o5.png", 6Gadfly.inch, 4Gadfly.inch), f) end end