From 156f788f3df8be557bbf6883dfefb6f437a19f08 Mon Sep 17 00:00:00 2001
From: baggepinnen <cont-frb@ulund.org>
Date: Sun, 7 Jan 2018 14:12:43 +0100
Subject: [PATCH] update

---
 flux/nn_prior/nn_prior.jl  | 34 +++++++++++++++++++---------------
 flux/nn_prior/utilities.jl |  2 +-
 2 files changed, 20 insertions(+), 16 deletions(-)

diff --git a/flux/nn_prior/nn_prior.jl b/flux/nn_prior/nn_prior.jl
index 08e71f8..2cd7bd1 100644
--- a/flux/nn_prior/nn_prior.jl
+++ b/flux/nn_prior/nn_prior.jl
@@ -9,7 +9,7 @@ using OrdinaryDiffEq, LTVModels
     num_params     = 20
     stepsize       = 0.05
     n_bootstrap    = 4
-    num_montecarlo = 30
+    num_montecarlo = 20
 end
 
 inoffice() = gethostname() ∈ ["billman", "Battostation"]
@@ -21,7 +21,7 @@ if !isdefined(:AbstractEnsembleSystem)
     @everywhere include("linear_sys.jl")
 end
 
-sys = TwoLinkSys(N=1000, h=0.02, σ0 = 0)#; warn("Added noise")
+sys = TwoLinkSys(N=1000, h=0.02, σ0 = 0.1)#; warn("Added noise")
 # sys = LinearSys(1,N=100, h=0.2, σ0 = 1, n = 3, ns = 3)
 
 @everywhere function fit_model(opt, loss, m, x, y, u, xv, yv, uv, sys, modeltype;
@@ -62,8 +62,9 @@ end
 
 @everywhere function fit_system(seed, sys, modeltype, jacprop = 0; doplot=false)
     x,y,u    = generate_data(sys, modeltype, seed)
-
+    iters = modeltype==System ? 3000 : 1500
     if jacprop > 0
+        iters = iters ÷ (jacprop+1)
         n,m = sys.ns,sys.n
         R1,R2 = 0.01eye(n*(n+m)), 10eye(n) # Fast
         P0 = 10000R1
@@ -90,7 +91,7 @@ end
 
         m          = modeltype(sys, num_params, activation)
         opt        = [ADAM(params(m), stepsize, decay=0.005); [expdecay(Param(p), wdecay) for p in params(m) if p isa AbstractMatrix]]
-        results    = fit_model(opt, loss(m), m, xt, yt, ut, xv, yv, uv, sys,modeltype, iters=modeltype==System ? 3000 : 1500, doplot=doplot, batch_size = sys.N)
+        results    = fit_model(opt, loss(m), m, xt, yt, ut, xv, yv, uv, sys,modeltype, iters=iters, doplot=doplot, batch_size = size(yt,2))
         @pack results = x, u, y, xv, uv, yv
         println("Done: ", it)
         results
@@ -177,30 +178,33 @@ open(file->serialize(file, res), "res","w") # Save
 error()
 
 ##
+sleep(60*25)
+using StatPlots
 res = open(file->deserialize(file), "res") # Load
+sys = res[1][1][1][:sys]
 nr = length(res[1]) ÷ 2
 labelvec = ["f" "g" "h"]
+infostring = @sprintf("Num hidden: %d, sigma: %2.2f, Montecarlo: %d", num_params, sys.σ0, num_montecarlo)
 
-pgfplots(size=(600,300), show=true)
-using StatPlots
+# pgfplots(size=(600,300), show=true)
 
 pred  = hcat([eval_pred.(get_res(res,i), true) for i ∈ 1:nr]...)
 sim   = hcat([eval_sim.(get_res(res,i), true) for i ∈ 1:nr]...)
 predp = hcat([eval_pred.(get_res(res,i), true) for i ∈ nr+1:2nr]...)
 simp  = hcat([eval_sim.(get_res(res,i), true) for i ∈ nr+1:2nr]...)
-vio1 = violin(log10.(pred), side=:left, lab=["Standard" "" ""],xticks=(1:nr,labelvec), ylabel="log(Prediction RMS)", reuse=false, c=:blue)
-violin!((1:nr)',log10.(predp), side=:right, lab=["Jacobian Propagation" "" ""],xticks=(1:nr,labelvec), c=:red)
-vio2 = violin(min.(log10.(sim),2), side=:left, lab=["Standard" "" ""],xticks=(1:nr,labelvec), ylabel="log(Simulation RMS)", c=:blue)
-violin!((1:nr)',min.(log10.(simp),2), side=:right, lab=["Jacobian Propagation" "" ""],xticks=(1:nr,labelvec), c=:red)
-plot(vio1,vio2); gui()
-savefig2("/local/home/fredrikb/papers/nn_prior/figs/valerr.tex")
+vio1 = violin(log10.(pred), side=:left, lab=["Standard" "" ""],xticks=(1:nr,labelvec), ylabel="log(Prediction RMS)", reuse=false, c=:red)
+violin!((1:nr)',log10.(predp), side=:right, lab=["Jacobian Propagation" "" ""],xticks=(1:nr,labelvec), c=:blue)
+vio2 = violin(min.(log10.(sim),2), side=:left, lab=["Standard" "" ""],xticks=(1:nr,labelvec), ylabel="log(Simulation RMS)", c=:red)
+violin!((1:nr)',min.(log10.(simp),2), side=:right, lab=["Jacobian Propagation" "" ""],xticks=(1:nr,labelvec), c=:blue, legend=false)
+plot(vio1,vio2,title=infostring); gui()
+# savefig2("/local/home/fredrikb/papers/nn_prior/figs/valerr.tex")
 
 # jacerrtrain = hcat([eval_jac.(get_res(res,i), false) for i ∈ 1:nr]...)
 jacerr = hcat([eval_jac.(get_res(res,i), true) for i ∈ 1:nr]...)
 jacerrp = hcat([eval_jac.(get_res(res,i), true) for i ∈ nr+1:2nr]...)
-violin(log10.(jacerr), side=:left, xticks=(1:nr,labelvec), title="Jacobian error (validation data)", ylabel="log(norm of Jacobian error)", reuse=true, lab=["Standard" "" ""], c=:red)
-violin!((1:nr)',log10.(jacerrp), side=:right, xticks=(1:nr,labelvec), ylabel="log(norm of Jacobian error)", reuse=true, lab=["Jacobian Propagation" "" ""], c=:blue); gui()
-savefig2("/local/home/fredrikb/papers/nn_prior/figs/jacerr.tex")
+violin(identity.(jacerr), side=:left, xticks=(1:nr,labelvec), title="Jacobian error (validation data)"*infostring, ylabel="norm of Jacobian error", reuse=true, lab=["Standard" "" ""], c=:red)
+violin!((1:nr)',identity.(jacerrp), side=:right, xticks=(1:nr,labelvec), ylabel="Jacobian RMS", reuse=true, lab=["Jacobian Propagation" "" ""], c=:blue); gui()
+# savefig2("/local/home/fredrikb/papers/nn_prior/figs/jacerr.tex")
 
 ##
 # TODO: ====================================================================================
diff --git a/flux/nn_prior/utilities.jl b/flux/nn_prior/utilities.jl
index 3f4848a..4744688 100644
--- a/flux/nn_prior/utilities.jl
+++ b/flux/nn_prior/utilities.jl
@@ -231,7 +231,7 @@ end
 
 function eval_jac(results, eval=false)
     Jm, Js, Jtrue = all_jacobians(results, eval)
-    vecnorm(Jm.-Jtrue)
+    sqrt(mean(abs2, Jm .- Jtrue))
 end
 
 function plotresults(results, eval=false)
-- 
GitLab