Hello everyone,
I am using Stan and Mamba to write an API to find the distribution of the parameters of the differential equations.
The API is
using OrdinaryDiffEq, ParameterizedFunctions, DiffEqBayes, Mamba
f1 = @ode_def_nohes LotkaVolterraTest1 begin
dx = a*x - b*x*y
dy = -c*y + d*x*y
end a=>1.5 b=1.0 c=3.0 d=1.0
u0 = [1.0,1.0]
tspan = (0.0,10.0)
prob1 = ODEProblem(f1,u0,tspan)
sol = solve(prob1,Tsit5())
t = collect(linspace(1,10,10))
using RecursiveArrayTools # for VectorOfArray
randomized = VectorOfArray([(sol(t[i]) + .01randn(2)) for i in 1:length(t)])
data = convert(Array,randomized)
bayesian_result = bayesian_inference(prob1,t,data;num_samples=1,num_warmup=1)
The source code of the API is here
bayesian_result is a struct with one of the members as chain_results which is of type Mamba.Chains
julia> bayesian_result.chain_results
Warmup took (0.00090, 0.0043, 0.0016, 0.00078) seconds, 0.0076 seconds total
Sampling took (0.13, 6.7, 0.025, 0.019) seconds, 6.8 seconds total
Mean MCSE StdDev 5% 50% 95% N_Eff N_Eff/s R_hat
lp__ -179 45 89 -248 -200 -49 4.0 0.59 nan
accept_stat__ 0.00 0.00 0.00 0.00 0.00 0.00 4.0 0.59 nan
stepsize__ 14 0.00 0.00 14 14 14 4.0 0.59 nan
treedepth__ 0.00 0.00 0.00 0.00 0.00 0.00 4.0 0.59 nan
n_leapfrog__ 1.0 0.00 0.00 1.0 1.0 1.0 4.0 0.59 nan
divergent__ 1.0 0.00 0.00 1.0 1.0 1.0 4.0 0.59 nan
energy__ 181 45 90 50 222 253 4.0 0.59 nan
sigma[1] 1.6 0.62 1.2 0.55 2.6 2.9 4.0 0.59 nan
sigma[2] 1.1 0.59 1.2 0.32 0.87 2.9 4.0 0.59 nan
theta[1] -0.55 0.69 1.4 -1.8 0.23 0.98 4.0 0.59 nan
Samples were drawn using hmc with nuts.
For each parameter, N_Eff is a crude measure of effective sample size,
and R_hat is the potential scale reduction factor on split chains (at
convergence, R_hat=1).
DiffEqBayes.StanModel{Int64,Mamba.Chains}(0, Object of type "Mamba.Chains"
Iterations = 1:1
Thinning interval = 1
Chains = 1,2,3,4
Samples per chain = 1
[-48.6759 0.0 … 0.869348 0.228172]
[-248.306 0.0 … 0.505688 0.979491]
[-199.57 0.0 … 2.88785 -1.79169]
[-218.627 0.0 … 0.32407 -1.6235])
I want to extract some useful information from this chain.
When I am doing:
bayesian_result.chain_results[1:100,["lp__"],:]
Object of type "Mamba.Chains"
Iterations = 1:1
Thinning interval = 1
Chains = 1,2,3,4
Samples per chain = 1
[-48.6759]
[-248.306]
[-199.57]
[-218.627]
All well. but when I am trying to extract sigma or theta[1] by
bayesian_result.chain_results[1:100,["theta[1]"],:]
I get the following error:
ERROR: BoundsError: attempt to access 10-element Array{AbstractString,1} at index [[0]]
Stacktrace:
[1] throw_boundserror(::Array{AbstractString,1}, ::Tuple{Array{Int64,1}}) at ./abstractarray.jl:433
[2] checkbounds at ./abstractarray.jl:362 [inlined]
[3] macro expansion at ./multidimensional.jl:441 [inlined]
[4] _getindex at ./multidimensional.jl:438 [inlined]
[5] getindex at ./abstractarray.jl:882 [inlined]
[6] getindex(::Mamba.Chains, ::UnitRange{Int64}, ::Array{String,1}, ::Int64) at /Users/ayush/.julia/v0.6/Mamba
How do I access the sigma and theta variables from the Mamba chain?
TIA 