Pkg.build("Mamba")
INFO: Building SpecialFunctions INFO: Building Cairo INFO: Building Rmath
using Mamba,GraphViz
# there is an error on GraphViz.jl at next.juliabox.org. this error is occured by missing download url of graphviz.tar.gz. Jun.03.2018
INFO: Precompiling module Distributions. INFO: Recompiling stale cache file /root/.julia/lib/v0.6/Calculus.ji for module Calculus. INFO: Recompiling stale cache file /root/.julia/lib/v0.6/Measures.ji for module Measures. INFO: Recompiling stale cache file /root/.julia/lib/v0.6/Compose.ji for module Compose. INFO: Precompiling module Gadfly. WARNING: StatsBase.WeightVec is deprecated, use StatsBase.Weights{#1#S, #2#T, #3#V} where #3#V<:AbstractArray{#2#T, 1} where #2#T<:Real where #1#S<:Real instead. likely near /root/.julia/v0.6/KernelDensity/src/univariate.jl:80 INFO: Precompiling module LightGraphs.
## Model and User-Defined Sampler Specifications
model = Model(
y = Stochastic(1,
(mu, s2) -> MvNormal(mu, sqrt(s2)),
false
),
mu = Logical(1,
(xmat, beta) -> xmat * beta,
false
),
beta = Stochastic(1,
() -> MvNormal(2, sqrt(1000))
),
s2 = Stochastic(
() -> InverseGamma(0.001, 0.001)
)
)
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "0-element Mamba.ArrayStochastic{1}" Float64[] ------------------------------------------------------------------------------- y: An unmonitored node of type "0-element Mamba.ArrayStochastic{1}" Float64[] ------------------------------------------------------------------------------- mu: An unmonitored node of type "0-element Mamba.ArrayLogical{1}" Float64[] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" NaN
Gibbs_beta = Sampler([:beta],
(beta, s2, xmat, y) ->
begin
beta_mean = mean(beta.distr)
beta_invcov = invcov(beta.distr)
Sigma = inv(Symmetric(xmat' * xmat / s2 + beta_invcov))
mu = Sigma * (xmat' * y / s2 + beta_invcov * beta_mean)
rand(MvNormal(mu, Sigma))
end
)
An object of type "Mamba.Sampler{Dict{Any,Any}}" Sampling Block Nodes: Symbol[:beta] CodeInfo(:(begin SSAValue(3) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:beta))) SSAValue(2) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:s2))) SSAValue(1) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:xmat))) SSAValue(0) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:y))) $(Expr(:inbounds, false)) # meta: location In[4] #17 4 beta_mean = (Main.mean)((Core.getfield)(SSAValue(3), :distr)) # line 5: beta_invcov = (Main.invcov)((Core.getfield)(SSAValue(3), :distr)) # line 6: Sigma = (Main.inv)((Main.Symmetric)((Main.Ac_mul_B)(SSAValue(1), SSAValue(1)) / SSAValue(2) + beta_invcov)::Symmetric{_,_} where _ where _) # line 7: mu = Sigma * ((Main.Ac_mul_B)(SSAValue(1), SSAValue(0)) / SSAValue(2) + beta_invcov * beta_mean) # meta: pop location $(Expr(:inbounds, :pop)) return (Main.rand)((Main.MvNormal)(mu, Sigma)) end))=>Any
Gibbs_s2 = Sampler([:s2],
(mu, s2, y) ->
begin
a = length(y) / 2.0 + shape(s2.distr)
b = sumabs2(y - mu) / 2.0 + scale(s2.distr)
rand(InverseGamma(a, b))
end
)
An object of type "Mamba.Sampler{Dict{Any,Any}}" Sampling Block Nodes: Symbol[:s2] CodeInfo(:(begin SSAValue(2) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:mu))) SSAValue(1) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:s2))) SSAValue(0) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:y))) $(Expr(:inbounds, false)) # meta: location In[5] #21 4 a = (Main.length)(SSAValue(0)) / 2.0 + (Main.shape)((Core.getfield)(SSAValue(1), :distr)) # line 5: b = (Main.sumabs2)(SSAValue(0) - SSAValue(2)) / 2.0 + (Main.scale)((Core.getfield)(SSAValue(1), :distr)) # meta: pop location $(Expr(:inbounds, :pop)) return (Main.rand)((Main.InverseGamma)(a, b)::Distributions.InverseGamma{_} where _)::Float64 end))=>Float64
## Hybrid No-U-Turn and Slice Sampling Scheme
scheme1 = [NUTS(:beta),
Slice(:s2, 3.0)]
2-element Array{Mamba.Sampler,1}: An object of type "Mamba.Sampler{Mamba.NUTSTune}" Sampling Block Nodes: Symbol[:beta] CodeInfo(:(begin f = Mamba.#319 return (f)(model, block) end))=>Any An object of type "Mamba.Sampler{Mamba.SliceTune{Distributions.Multivariate}}" Sampling Block Nodes: Symbol[:s2] CodeInfo(:(begin f = Mamba.#337 return (f)(model, block) end))=>Any
## No-U-Turn Sampling Scheme
scheme2 = [NUTS([:beta, :s2])]
1-element Array{Mamba.Sampler{Mamba.NUTSTune},1}: An object of type "Mamba.Sampler{Mamba.NUTSTune}" Sampling Block Nodes: Symbol[:beta, :s2] CodeInfo(:(begin f = Mamba.#319 return (f)(model, block) end))=>Any
## User-Defined Sampling Scheme
scheme3 = [Gibbs_beta, Gibbs_s2]
2-element Array{Mamba.Sampler{Dict{Any,Any}},1}: An object of type "Mamba.Sampler{Dict{Any,Any}}" Sampling Block Nodes: Symbol[:beta] CodeInfo(:(begin SSAValue(3) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:beta))) SSAValue(2) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:s2))) SSAValue(1) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:xmat))) SSAValue(0) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:y))) $(Expr(:inbounds, false)) # meta: location In[4] #17 4 beta_mean = (Main.mean)((Core.getfield)(SSAValue(3), :distr)) # line 5: beta_invcov = (Main.invcov)((Core.getfield)(SSAValue(3), :distr)) # line 6: Sigma = (Main.inv)((Main.Symmetric)((Main.Ac_mul_B)(SSAValue(1), SSAValue(1)) / SSAValue(2) + beta_invcov)::Symmetric{_,_} where _ where _) # line 7: mu = Sigma * ((Main.Ac_mul_B)(SSAValue(1), SSAValue(0)) / SSAValue(2) + beta_invcov * beta_mean) # meta: pop location $(Expr(:inbounds, :pop)) return (Main.rand)((Main.MvNormal)(mu, Sigma)) end))=>Any An object of type "Mamba.Sampler{Dict{Any,Any}}" Sampling Block Nodes: Symbol[:s2] CodeInfo(:(begin SSAValue(2) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:mu))) SSAValue(1) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:s2))) SSAValue(0) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:y))) $(Expr(:inbounds, false)) # meta: location In[5] #21 4 a = (Main.length)(SSAValue(0)) / 2.0 + (Main.shape)((Core.getfield)(SSAValue(1), :distr)) # line 5: b = (Main.sumabs2)(SSAValue(0) - SSAValue(2)) / 2.0 + (Main.scale)((Core.getfield)(SSAValue(1), :distr)) # meta: pop location $(Expr(:inbounds, :pop)) return (Main.rand)((Main.InverseGamma)(a, b)::Distributions.InverseGamma{_} where _)::Float64 end))=>Float64
## Sampling Scheme Assignment
setsamplers!(model, scheme1)
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "0-element Mamba.ArrayStochastic{1}" Float64[] ------------------------------------------------------------------------------- y: An unmonitored node of type "0-element Mamba.ArrayStochastic{1}" Float64[] ------------------------------------------------------------------------------- mu: An unmonitored node of type "0-element Mamba.ArrayLogical{1}" Float64[] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" NaN
## Graph Representation of Nodes
draw(model)
# draw(model, filename="lineDAG.dot")
digraph MambaModel { "beta" [shape="ellipse"]; "beta" -> "mu"; "y" [shape="ellipse", style="filled", fillcolor="gray85"]; "mu" [shape="diamond", style="filled", fillcolor="gray85"]; "mu" -> "y"; "s2" [shape="ellipse"]; "s2" -> "y"; "xmat" [shape="box", style="filled", fillcolor="gray85"]; "xmat" -> "mu"; }
display(Graph(graph2dot(model)))
## Data
line = Dict{Symbol, Any}(
:x => [1, 2, 3, 4, 5],
:y => [1, 3, 3, 3, 5]
)
Dict{Symbol,Any} with 2 entries: :y => [1, 3, 3, 3, 5] :x => [1, 2, 3, 4, 5]
line[:xmat] = [ones(5) line[:x]]
5×2 Array{Float64,2}: 1.0 1.0 1.0 2.0 1.0 3.0 1.0 4.0 1.0 5.0
## Set Random Number Generator Seed
srand(123)
MersenneTwister(UInt32[0x0000007b], Base.dSFMT.DSFMT_state(Int32[1464307935, 1073116007, 222134151, 1073120226, -290652630, 1072956456, -580276323, 1073476387, 1332671753, 1073438661 … 138346874, 1073030449, 1049893279, 1073166535, -1999907543, 1597138926, -775229811, 32947490, 382, 0]), [1.75595, 1.62839, 1.58743, 1.23614, 1.10934, 1.63961, 1.69764, 1.24159, 1.94223, 1.21538 … 1.54803, 1.39513, 1.10156, 1.0517, 1.96602, 1.02523, 1.14843, 1.92493, 1.43046, 1.71541], 382)
## Initial Values
inits = [
Dict{Symbol, Any}(
:y => line[:y],
:beta => rand(Normal(0, 1), 2),
:s2 => rand(Gamma(1, 1))
)
for i in 1:3
]
3-element Array{Dict{Symbol,Any},1}: Dict{Symbol,Any}(Pair{Symbol,Any}(:beta, [1.19027, 2.04818]),Pair{Symbol,Any}(:y, [1, 3, 3, 3, 5]),Pair{Symbol,Any}(:s2, 1.63439)) Dict{Symbol,Any}(Pair{Symbol,Any}(:beta, [0.459416, -0.396679]),Pair{Symbol,Any}(:y, [1, 3, 3, 3, 5]),Pair{Symbol,Any}(:s2, 0.140438)) Dict{Symbol,Any}(Pair{Symbol,Any}(:beta, [-0.0754831, 0.273815]),Pair{Symbol,Any}(:y, [1, 3, 3, 3, 5]),Pair{Symbol,Any}(:s2, 0.372091))
## MCMC Simulations
setsamplers!(model, scheme1)
sim1 = mcmc(model, line, inits, 10000, burnin=250, thin=2, chains=3)
MCMC Simulation of 10000 Iterations x 3 Chains... Chain 1: 0% [0:30:33 of 0:30:35 remaining] Chain 1: 10% [0:00:29 of 0:00:32 remaining] Chain 1: 20% [0:00:18 of 0:00:23 remaining] Chain 1: 30% [0:00:14 of 0:00:20 remaining] Chain 1: 40% [0:00:11 of 0:00:18 remaining] Chain 1: 50% [0:00:08 of 0:00:16 remaining] Chain 1: 60% [0:00:06 of 0:00:16 remaining] Chain 1: 70% [0:00:05 of 0:00:15 remaining] Chain 1: 80% [0:00:03 of 0:00:15 remaining] Chain 1: 90% [0:00:01 of 0:00:15 remaining] Chain 1: 100% [0:00:00 of 0:00:15 remaining] Chain 2: 0% [0:00:12 of 0:00:12 remaining] Chain 2: 10% [0:00:13 of 0:00:14 remaining] Chain 2: 20% [0:00:11 of 0:00:14 remaining] Chain 2: 30% [0:00:10 of 0:00:14 remaining] Chain 2: 40% [0:00:08 of 0:00:14 remaining] Chain 2: 50% [0:00:07 of 0:00:13 remaining] Chain 2: 60% [0:00:05 of 0:00:14 remaining] Chain 2: 70% [0:00:04 of 0:00:14 remaining] Chain 2: 80% [0:00:03 of 0:00:14 remaining] Chain 2: 90% [0:00:01 of 0:00:14 remaining] Chain 2: 100% [0:00:00 of 0:00:14 remaining] Chain 3: 0% [0:00:10 of 0:00:10 remaining] Chain 3: 10% [0:00:11 of 0:00:12 remaining] Chain 3: 20% [0:00:10 of 0:00:12 remaining] Chain 3: 30% [0:00:08 of 0:00:12 remaining] Chain 3: 40% [0:00:07 of 0:00:12 remaining] Chain 3: 50% [0:00:06 of 0:00:12 remaining] Chain 3: 60% [0:00:05 of 0:00:12 remaining] Chain 3: 70% [0:00:04 of 0:00:12 remaining] Chain 3: 80% [0:00:02 of 0:00:12 remaining] Chain 3: 90% [0:00:01 of 0:00:12 remaining] Chain 3: 100% [0:00:00 of 0:00:12 remaining]
Object of type "Mamba.ModelChains" Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 [0.918604 1.08849 0.719733; 0.291973 1.51878 0.629041; … ; 0.480511 1.04643 0.783623; 0.308753 0.843837 0.772927] [0.31982 0.35836 0.997317; 0.175972 0.331401 0.769161; … ; 0.475919 0.647507 0.865755; 0.308566 0.647507 0.865755] [0.539945 1.2289 0.722182; 0.543182 -0.0625791 1.07568; … ; 1.21675 2.40923 0.38377; 3.39815 1.94268 0.472419]
setsamplers!(model, scheme2)
sim2 = mcmc(model, line, inits, 10000, burnin=250, thin=2, chains=3)
MCMC Simulation of 10000 Iterations x 3 Chains... Chain 1: 0% [0:01:19 of 0:01:19 remaining] Chain 1: 10% [0:00:19 of 0:00:21 remaining] Chain 1: 20% [0:00:16 of 0:00:20 remaining] Chain 1: 30% [0:00:14 of 0:00:20 remaining] Chain 1: 40% [0:00:12 of 0:00:20 remaining] Chain 1: 50% [0:00:10 of 0:00:20 remaining] Chain 1: 60% [0:00:08 of 0:00:20 remaining] Chain 1: 70% [0:00:06 of 0:00:20 remaining] Chain 1: 80% [0:00:04 of 0:00:20 remaining] Chain 1: 90% [0:00:02 of 0:00:20 remaining] Chain 1: 100% [0:00:00 of 0:00:20 remaining] Chain 2: 0% [0:00:19 of 0:00:19 remaining] Chain 2: 10% [0:00:16 of 0:00:18 remaining] Chain 2: 20% [0:00:13 of 0:00:16 remaining] Chain 2: 30% [0:00:11 of 0:00:15 remaining] Chain 2: 40% [0:00:09 of 0:00:15 remaining] Chain 2: 50% [0:00:07 of 0:00:14 remaining] Chain 2: 60% [0:00:06 of 0:00:14 remaining] Chain 2: 70% [0:00:04 of 0:00:14 remaining] Chain 2: 80% [0:00:03 of 0:00:14 remaining] Chain 2: 90% [0:00:01 of 0:00:14 remaining] Chain 2: 100% [0:00:00 of 0:00:14 remaining] Chain 3: 0% [0:00:48 of 0:00:48 remaining] Chain 3: 10% [0:00:18 of 0:00:20 remaining] Chain 3: 20% [0:00:15 of 0:00:19 remaining] Chain 3: 30% [0:00:14 of 0:00:20 remaining] Chain 3: 40% [0:00:12 of 0:00:20 remaining] Chain 3: 50% [0:00:09 of 0:00:19 remaining] Chain 3: 60% [0:00:08 of 0:00:19 remaining] Chain 3: 70% [0:00:06 of 0:00:19 remaining] Chain 3: 80% [0:00:04 of 0:00:19 remaining] Chain 3: 90% [0:00:02 of 0:00:19 remaining] Chain 3: 100% [0:00:00 of 0:00:19 remaining]
Object of type "Mamba.ModelChains" Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 [0.373215 0.144735 0.795453; 0.360428 0.375819 0.790722; … ; 0.5952 1.18019 0.532691; 0.603331 1.39923 0.485838] [0.574327 1.76252 0.42196; 1.84602 0.24232 1.00588; … ; 0.414817 1.11963 0.86145; 0.414817 1.11963 0.86145] [2.79335 -0.16188 1.02612; 0.654844 -1.01357 1.23578; … ; 1.00966 1.97517 0.329111; 0.828978 1.9125 0.562846]
setsamplers!(model, scheme3)
sim3 = mcmc(model, line, inits, 10000, burnin=250, thin=2, chains=3)
MCMC Simulation of 10000 Iterations x 3 Chains...
WARNING: sumabs2(x) is deprecated, use sum(abs2, x) instead. Stacktrace: [1] depwarn(::String, ::Symbol) at ./deprecated.jl:70 [2] sumabs2(::Array{Float64,1}) at ./deprecated.jl:57 [3] #21 at ./In[5]:5 [inlined] [4] (::##23#24)(::Mamba.Model, ::Int64) at ./<missing>:0 [5] sample!(::Mamba.Model, ::Int64) at /root/.julia/v0.6/Mamba/src/model/simulation.jl:99 [6] sample!(::Mamba.Model) at /root/.julia/v0.6/Mamba/src/model/simulation.jl:94 [7] mcmc_worker!(::Array{Any,1}) at /root/.julia/v0.6/Mamba/src/model/mcmc.jl:75 [8] _collect(::Array{Array{Any,1},1}, ::Base.Generator{Array{Array{Any,1},1},Mamba.#mcmc_worker!}, ::Base.EltypeUnknown, ::Base.HasShape) at ./array.jl:454 [9] map(::Function, ::Array{Array{Any,1},1}) at ./abstractarray.jl:1865 [10] pmap2(::Function, ::Array{Array{Any,1},1}) at /root/.julia/v0.6/Mamba/src/utils.jl:99 [11] mcmc_master!(::Mamba.Model, ::UnitRange{Int64}, ::Int64, ::Int64, ::UnitRange{Int64}, ::Bool) at /root/.julia/v0.6/Mamba/src/model/mcmc.jl:52 [12] #mcmc#29(::Int64, ::Int64, ::Int64, ::Bool, ::Function, ::Mamba.Model, ::Dict{Symbol,Any}, ::Array{Dict{Symbol,Any},1}, ::Int64) at /root/.julia/v0.6/Mamba/src/model/mcmc.jl:32 [13] (::Mamba.#kw##mcmc)(::Array{Any,1}, ::Mamba.#mcmc, ::Mamba.Model, ::Dict{Symbol,Any}, ::Array{Dict{Symbol,Any},1}, ::Int64) at ./<missing>:0 [14] include_string(::String, ::String) at ./loading.jl:515 [15] include_string(::Module, ::String, ::String) at /root/.julia/v0.6/Compat/src/Compat.jl:174 [16] execute_request(::ZMQ.Socket, ::IJulia.Msg) at /root/.julia/v0.6/IJulia/src/execute_request.jl:154 [17] (::Compat.#inner#16{Array{Any,1},IJulia.#execute_request,Tuple{ZMQ.Socket,IJulia.Msg}})() at /root/.julia/v0.6/Compat/src/Compat.jl:496 [18] eventloop(::ZMQ.Socket) at /root/.julia/v0.6/IJulia/src/eventloop.jl:8 [19] (::IJulia.##14#17)() at ./task.jl:335 while loading In[18], in expression starting on line 2
Chain 1: 0% [0:38:06 of 0:38:08 remaining] Chain 1: 10% [0:00:24 of 0:00:27 remaining] Chain 1: 20% [0:00:12 of 0:00:15 remaining] Chain 1: 30% [0:00:08 of 0:00:11 remaining] Chain 1: 40% [0:00:06 of 0:00:10 remaining] Chain 1: 50% [0:00:04 of 0:00:08 remaining] Chain 1: 60% [0:00:03 of 0:00:08 remaining] Chain 1: 70% [0:00:02 of 0:00:07 remaining] Chain 1: 80% [0:00:01 of 0:00:07 remaining] Chain 1: 90% [0:00:01 of 0:00:06 remaining] Chain 1: 100% [0:00:00 of 0:00:06 remaining] Chain 2: 0% [0:00:06 of 0:00:06 remaining] Chain 2: 10% [0:00:03 of 0:00:04 remaining] Chain 2: 20% [0:00:03 of 0:00:04 remaining] Chain 2: 30% [0:00:03 of 0:00:04 remaining] Chain 2: 40% [0:00:02 of 0:00:04 remaining] Chain 2: 50% [0:00:02 of 0:00:04 remaining] Chain 2: 60% [0:00:02 of 0:00:04 remaining] Chain 2: 70% [0:00:01 of 0:00:04 remaining] Chain 2: 80% [0:00:01 of 0:00:04 remaining] Chain 2: 90% [0:00:00 of 0:00:04 remaining] Chain 2: 100% [0:00:00 of 0:00:04 remaining] Chain 3: 0% [0:00:04 of 0:00:04 remaining] Chain 3: 10% [0:00:03 of 0:00:04 remaining] Chain 3: 20% [0:00:03 of 0:00:04 remaining] Chain 3: 30% [0:00:03 of 0:00:04 remaining] Chain 3: 40% [0:00:02 of 0:00:04 remaining] Chain 3: 50% [0:00:02 of 0:00:04 remaining] Chain 3: 60% [0:00:02 of 0:00:04 remaining] Chain 3: 70% [0:00:01 of 0:00:04 remaining] Chain 3: 80% [0:00:01 of 0:00:04 remaining] Chain 3: 90% [0:00:00 of 0:00:04 remaining] Chain 3: 100% [0:00:00 of 0:00:04 remaining]
Object of type "Mamba.ModelChains" Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 [0.558365 0.754007 0.721287; 0.504818 1.94074 0.395813; … ; 0.325436 0.299898 0.855316; 1.02423 0.520895 0.897936] [2.42999 0.820647 1.09475; 4.26356 -0.0745875 0.535085; … ; 2.26708 2.51055 0.127477; 2.89595 -1.59997 1.22045] [58.1332 -0.421533 1.89501; 2.59648 4.18922 0.0347426; … ; 0.247849 0.427886 0.7765; 0.270629 0.396874 0.779213]
## Gelman, Rubin, and Brooks Convergence Diagnostic
gelmandiag(sim1, mpsrf=true, transform=true) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Gelman, Rubin, and Brooks Diagnostic: PSRF 97.5% s2 1.002 1.003 beta[1] 1.002 1.003 beta[2] 1.003 1.003 Multivariate 1.000 NaN
gelmandiag(sim2, mpsrf=true, transform=true) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Gelman, Rubin, and Brooks Diagnostic: PSRF 97.5% s2 1.005 1.015 beta[1] 1.000 1.000 beta[2] 1.000 1.000 Multivariate 1.007 NaN
gelmandiag(sim3, mpsrf=true, transform=true) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Gelman, Rubin, and Brooks Diagnostic: PSRF 97.5% s2 1.000 1.000 beta[1] 1.002 1.002 beta[2] 1.001 1.002 Multivariate 1.000 NaN
## Geweke Convergence Diagnostic
gewekediag(sim1) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Geweke Diagnostic: First Window Fraction = 0.1 Second Window Fraction = 0.5 Z-score p-value s2 -2.448 0.0144 beta[1] -0.450 0.6528 beta[2] 0.657 0.5112 Z-score p-value s2 -1.361 0.1736 beta[1] 0.266 0.7904 beta[2] -0.291 0.7714 Z-score p-value s2 -1.677 0.0936 beta[1] -0.022 0.9821 beta[2] 0.283 0.7771
## Geweke Convergence Diagnostic
gewekediag(sim2) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Geweke Diagnostic: First Window Fraction = 0.1 Second Window Fraction = 0.5 Z-score p-value s2 -0.901 0.3675 beta[1] 0.646 0.5184 beta[2] -0.476 0.6339 Z-score p-value s2 1.985 0.0472 beta[1] 0.228 0.8195 beta[2] 0.073 0.9421 Z-score p-value s2 -1.742 0.0816 beta[1] 0.636 0.5245 beta[2] -0.323 0.7468
## Geweke Convergence Diagnostic
gewekediag(sim2) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Geweke Diagnostic: First Window Fraction = 0.1 Second Window Fraction = 0.5 Z-score p-value s2 -0.901 0.3675 beta[1] 0.646 0.5184 beta[2] -0.476 0.6339 Z-score p-value s2 1.985 0.0472 beta[1] 0.228 0.8195 beta[2] 0.073 0.9421 Z-score p-value s2 -1.742 0.0816 beta[1] 0.636 0.5245 beta[2] -0.323 0.7468
## Heidelberger-Welch Diagnostic
heideldiag(sim1) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Heidelberger and Welch Diagnostic: Target Halfwidth Ratio = 0.1 Alpha = 0.05 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.5998 1.15916075 0.176484028 0 beta[1] 251 1 0.6081 0.56689198 0.062442773 0 beta[2] 251 1 0.6189 0.80665808 0.016980313 1 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.2885 1.13947629 0.155614707 0 beta[1] 251 1 0.5622 0.57205135 0.057180269 1 beta[2] 251 1 0.4376 0.80577830 0.015814972 1 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.4505 1.29747773 0.265610771 0 beta[1] 251 1 0.1861 0.54408853 0.064555202 0 beta[2] 251 1 0.1812 0.81327110 0.018693595 1
heideldiag(sim2) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Heidelberger and Welch Diagnostic: Target Halfwidth Ratio = 0.1 Alpha = 0.05 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.8491 1.80264129 0.613183420 0 beta[1] 251 1 0.6533 0.57201198 0.100010924 0 beta[2] 251 1 0.6978 0.80564502 0.029361227 1 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.1238 1.50022282 0.196508215 0 beta[1] 251 1 0.7523 0.58188054 0.104408366 0 beta[2] 251 1 0.7694 0.80354294 0.029188774 1 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.8393 1.67274868 0.308027488 0 beta[1] 251 1 0.9750 0.57719994 0.090231219 0 beta[2] 251 1 0.9738 0.80692656 0.025309048 1
heideldiag(sim3) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Heidelberger and Welch Diagnostic: Target Halfwidth Ratio = 0.1 Alpha = 0.05 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.3739 1.4855501 0.163476834 0 beta[1] 251 1 0.3771 0.6299113 0.034193149 1 beta[2] 251 1 0.6818 0.7873776 0.010856256 1 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.1608 1.60896183 0.395504511 0 beta[1] 251 1 0.8160 0.62212856 0.040120777 1 beta[2] 251 1 0.4684 0.79213694 0.012274618 1 Burn-in Stationarity p-value Mean Halfwidth Test s2 251 1 0.3840 1.4719120 0.156612791 0 beta[1] 251 1 0.4325 0.6224132 0.033551099 1 beta[2] 251 1 0.2163 0.7923469 0.010592564 1
## Raftery-Lewis Convergence Diagnostic
rafterydiag(sim1) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Raftery and Lewis Diagnostic: Quantile (q) = 0.025 Accuracy (r) = 0.005 Probability (s) = 0.95 Thinning Burn-in Total Nmin Dependence Factor s2 2 257 8.6890×10³ 3746 2.3195408 beta[1] 2 273 2.3015×10⁴ 3746 6.1438868 beta[2] 2 265 1.6555×10⁴ 3746 4.4193807 Thinning Burn-in Total Nmin Dependence Factor s2 2 257 8.4090×10³ 3746 2.2447944 beta[1] 2 271 2.2509×10⁴ 3746 6.0088094 beta[2] 2 265 1.4793×10⁴ 3746 3.9490123 Thinning Burn-in Total Nmin Dependence Factor s2 2 255 7.8770×10³ 3746 2.1027763 beta[1] 2 275 2.5215×10⁴ 3746 6.7311799 beta[2] 2 265 1.5069×10⁴ 3746 4.0226909
rafterydiag(sim2) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Raftery and Lewis Diagnostic: Quantile (q) = 0.025 Accuracy (r) = 0.005 Probability (s) = 0.95 Thinning Burn-in Total Nmin Dependence Factor s2 2 361 1.20301×10⁵ 3746 32.1145222 beta[1] 2 269 1.90070×10⁴ 3746 5.0739455 beta[2] 4 271 2.63790×10⁴ 3746 7.0419114 Thinning Burn-in Total Nmin Dependence Factor s2 2 337 9.5029×10⁴ 3746 25.368126 beta[1] 2 271 2.1093×10⁴ 3746 5.630806 beta[2] 4 275 2.3971×10⁴ 3746 6.399092 Thinning Burn-in Total Nmin Dependence Factor s2 2 371 1.30115×10⁵ 3746 34.7343833 beta[1] 2 273 2.46350×10⁴ 3746 6.5763481 beta[2] 4 279 3.11390×10⁴ 3746 8.3126001
rafterydiag(sim3) |> showall
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Raftery and Lewis Diagnostic: Quantile (q) = 0.025 Accuracy (r) = 0.005 Probability (s) = 0.95 Thinning Burn-in Total Nmin Dependence Factor s2 2 255 7625 3746 2.0355045 beta[1] 2 255 8137 3746 2.1721837 beta[2] 2 257 9127 3746 2.4364656 Thinning Burn-in Total Nmin Dependence Factor s2 2 255 7877 3746 2.1027763 beta[1] 2 257 8547 3746 2.2816337 beta[2] 2 257 8547 3746 2.2816337 Thinning Burn-in Total Nmin Dependence Factor s2 2 255 7625 3746 2.0355045 beta[1] 2 257 8273 3746 2.2084891 beta[2] 2 257 8409 3746 2.2447944
## Summary Statistics
describe(sim1)
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Empirical Posterior Estimates: Mean SD Naive SE MCSE ESS s2 1.19870492 1.56864002 0.012971057 0.0600219481 683.009 beta[1] 0.56101062 1.15536100 0.009553660 0.0156473377 4875.000 beta[2] 0.80856916 0.34726763 0.002871550 0.0043911549 4875.000 Quantiles: 2.5% 25.0% 50.0% 75.0% 97.5% s2 0.17591020 0.391720912 0.67219297 1.29812099 5.8893323 beta[1] -1.80132114 -0.055811678 0.56674250 1.19351175 2.8785273 beta[2] 0.13335967 0.623366952 0.80857307 0.99002252 1.5021781
describe(sim2)
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Empirical Posterior Estimates: Mean SD Naive SE MCSE ESS s2 1.65853760 7.67749431 0.0634850646 0.122786322 3909.6603 beta[1] 0.57703082 1.32603422 0.0109649535 0.030772632 1856.8644 beta[2] 0.80537151 0.39995784 0.0033072443 0.008675861 2125.2168 Quantiles: 2.5% 25.0% 50.0% 75.0% 97.5% s2 0.194999834 0.412008359 0.70378644 1.34486155 7.8372754 beta[1] -1.869579355 -0.008477116 0.55332126 1.17379551 3.1352361 beta[2] 0.056019003 0.628740742 0.81490848 0.98495395 1.5351725
describe(sim3)
Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 Empirical Posterior Estimates: Mean SD Naive SE MCSE ESS s2 1.52214133 7.58640005 0.0627318078 0.0778954259 4875 beta[1] 0.62481770 1.26263455 0.0104407027 0.0099656668 4875 beta[2] 0.79062047 0.38313403 0.0031681285 0.0031775716 4875 Quantiles: 2.5% 25.0% 50.0% 75.0% 97.5% s2 0.17038427 0.386716322 0.66829529 1.2984785 7.1679905 beta[1] -1.76718183 0.029597673 0.62215365 1.2074545 3.0674552 beta[2] 0.04950653 0.616954117 0.79309708 0.9704145 1.5102962
## Highest Posterior Density Intervals
hpd(sim1) |> show
95% Lower 95% Upper s2 0.08105333 4.0471892 beta[1] -1.82356032 2.8269553 beta[2] 0.13379775 1.5023687
hpd(sim2) |> show
95% Lower 95% Upper s2 0.120542196 4.5846500 beta[1] -1.993136241 3.0001487 beta[2] 0.013429415 1.4853466
hpd(sim3) |> show
95% Lower 95% Upper s2 0.071760184 4.3348775 beta[1] -1.862824151 2.9422668 beta[2] 0.070223727 1.5260982
## Cross-Correlations
cor(sim1) |> show
s2 beta[1] beta[2] s2 1.000000000 -0.05849085 0.031009024 beta[1] -0.058490849 1.00000000 -0.902530376 beta[2] 0.031009024 -0.90253038 1.000000000
cor(sim2) |> show
s2 beta[1] beta[2] s2 1.000000000 -0.03673615 -0.004026096 beta[1] -0.036736149 1.00000000 -0.905842663 beta[2] -0.004026096 -0.90584266 1.000000000
cor(sim3) |> show
s2 beta[1] beta[2] s2 1.00000000 0.08159652 -0.13687408 beta[1] 0.08159652 1.00000000 -0.90489107 beta[2] -0.13687408 -0.90489107 1.00000000
## Lag-Autocorrelations
autocor(sim1) |> show
Lag 2 Lag 10 Lag 20 Lag 100 s2 0.76778414 0.35282562 0.131592764 0.044746958 beta[1] 0.32363233 -0.06878183 -0.010020497 -0.004079708 beta[2] 0.25761327 -0.06478366 -0.025873051 0.012010870 Lag 2 Lag 10 Lag 20 Lag 100 s2 0.72262408 0.2915667603 0.132101886 -0.037133178 beta[1] 0.28503104 0.0148294471 0.005551656 0.019832008 beta[2] 0.22748197 -0.0040102710 0.008956909 0.028501515 Lag 2 Lag 10 Lag 20 Lag 100 s2 0.85240357 0.545525422 0.2687135047 0.0165931570 beta[1] 0.33115909 -0.007933394 0.0195832385 -0.0064098384 beta[2] 0.28378709 0.004359849 -0.0022212270 -0.0148778720
autocor(sim2) |> show
Lag 2 Lag 10 Lag 20 Lag 100 s2 0.18488454 0.035202057 -0.003940754 -0.0054302732 beta[1] 0.44605699 0.146730076 0.016926444 -0.0408891861 beta[2] 0.36545891 0.139713177 0.016543257 -0.0262060509 Lag 2 Lag 10 Lag 20 Lag 100 s2 0.25756416 0.044745078 0.017576072 -0.004635488 beta[1] 0.54635868 0.109015883 0.044875654 0.044184269 beta[2] 0.46683183 0.088258890 0.029745194 0.029725365 Lag 2 Lag 10 Lag 20 Lag 100 s2 0.28006727 0.027014925 0.0017233265 -0.006291950 beta[1] 0.47993193 0.036792864 0.0197812905 0.022203871 beta[2] 0.41089419 0.044369596 0.0254449988 0.005718589
autocor(sim3) |> show
Lag 2 Lag 10 Lag 20 Lag 100 s2 0.084059883 -0.001807798 0.0087606658 -0.008382782 beta[1] 0.009371122 -0.033369868 0.0145553415 0.015623240 beta[2] 0.011278466 -0.027321121 -0.0017182889 0.010721298 Lag 2 Lag 10 Lag 20 Lag 100 s2 0.02625439178 -0.0027369637 0.0022959722 -0.0013664788 beta[1] -0.00017040833 0.0054179462 -0.0069941660 0.0169695433 beta[2] 0.01452332938 0.0185093311 -0.0042258046 0.0164201209 Lag 2 Lag 10 Lag 20 Lag 100 s2 0.0506426039 -0.0023500907 -0.0029839163 0.006845880 beta[1] -0.0049671116 -0.0061117412 -0.0039590427 0.024093818 beta[2] 0.0062727269 -0.0061675828 -0.0006290705 0.031324733
## State Space Change Rate (per Iteration)
changerate(sim1) |> show
Change Rate s2 1.000 beta[1] 0.876 beta[2] 0.876 Multivariate 1.000
changerate(sim2) |> show
Change Rate s2 0.828 beta[1] 0.828 beta[2] 0.828 Multivariate 0.828
changerate(sim3) |> show
Change Rate s2 1 beta[1] 1 beta[2] 1 Multivariate 1
## Deviance Information Criterion
dic(sim1) |> show
MCMC Processing of 4875 Iterations x 3 Chains... Chain 1: 0% [0:01:43 of 0:01:43 remaining] Chain 1: 10% [0:00:02 of 0:00:02 remaining] Chain 1: 20% [0:00:01 of 0:00:01 remaining] Chain 1: 30% [0:00:01 of 0:00:01 remaining] Chain 1: 40% [0:00:00 of 0:00:01 remaining] Chain 1: 50% [0:00:00 of 0:00:01 remaining] Chain 1: 60% [0:00:00 of 0:00:00 remaining] Chain 1: 70% [0:00:00 of 0:00:00 remaining] Chain 1: 80% [0:00:00 of 0:00:00 remaining] Chain 1: 90% [0:00:00 of 0:00:00 remaining] Chain 1: 100% [0:00:00 of 0:00:00 remaining] Chain 2: 0% [0:03:02 of 0:03:03 remaining] Chain 2: 10% [0:00:04 of 0:00:04 remaining] Chain 2: 20% [0:00:02 of 0:00:02 remaining] Chain 2: 30% [0:00:01 of 0:00:01 remaining] Chain 2: 40% [0:00:01 of 0:00:01 remaining] Chain 2: 50% [0:00:00 of 0:00:01 remaining] Chain 2: 60% [0:00:00 of 0:00:01 remaining] Chain 2: 70% [0:00:00 of 0:00:01 remaining] Chain 2: 80% [0:00:00 of 0:00:01 remaining] Chain 2: 90% [0:00:00 of 0:00:01 remaining] Chain 2: 100% [0:00:00 of 0:00:01 remaining] Chain 3: 0% [0:04:16 of 0:04:16 remaining] Chain 3: 10% [0:00:05 of 0:00:05 remaining] Chain 3: 20% [0:00:02 of 0:00:03 remaining] Chain 3: 30% [0:00:01 of 0:00:02 remaining] Chain 3: 40% [0:00:01 of 0:00:01 remaining] Chain 3: 50% [0:00:01 of 0:00:01 remaining] Chain 3: 60% [0:00:00 of 0:00:01 remaining] Chain 3: 70% [0:00:00 of 0:00:01 remaining] Chain 3: 80% [0:00:00 of 0:00:01 remaining] Chain 3: 90% [0:00:00 of 0:00:01 remaining] Chain 3: 100% [0:00:00 of 0:00:01 remaining] DIC Effective Parameters pD 14.194649 1.3814665 pV 22.608880 5.5885820
dic(sim2) |> show
MCMC Processing of 4875 Iterations x 3 Chains... Chain 1: 0% [0:00:00 of 0:00:00 remaining] Chain 1: 10% [0:00:00 of 0:00:00 remaining] Chain 1: 20% [0:00:00 of 0:00:00 remaining] Chain 1: 30% [0:00:00 of 0:00:00 remaining] Chain 1: 40% [0:00:00 of 0:00:00 remaining] Chain 1: 50% [0:00:00 of 0:00:00 remaining] Chain 1: 60% [0:00:00 of 0:00:00 remaining] Chain 1: 70% [0:00:00 of 0:00:00 remaining] Chain 1: 80% [0:00:00 of 0:00:00 remaining] Chain 1: 90% [0:00:00 of 0:00:00 remaining] Chain 1: 100% [0:00:00 of 0:00:00 remaining] Chain 2: 0% [0:01:08 of 0:01:08 remaining] Chain 2: 10% [0:00:01 of 0:00:02 remaining] Chain 2: 20% [0:00:01 of 0:00:01 remaining] Chain 2: 30% [0:00:00 of 0:00:01 remaining] Chain 2: 40% [0:00:00 of 0:00:01 remaining] Chain 2: 50% [0:00:00 of 0:00:00 remaining] Chain 2: 60% [0:00:00 of 0:00:00 remaining] Chain 2: 70% [0:00:00 of 0:00:00 remaining] Chain 2: 80% [0:00:00 of 0:00:00 remaining] Chain 2: 90% [0:00:00 of 0:00:00 remaining] Chain 2: 100% [0:00:00 of 0:00:00 remaining] Chain 3: 0% [0:02:23 of 0:02:24 remaining] Chain 3: 10% [0:00:03 of 0:00:03 remaining] Chain 3: 20% [0:00:01 of 0:00:02 remaining] Chain 3: 30% [0:00:01 of 0:00:01 remaining] Chain 3: 40% [0:00:01 of 0:00:01 remaining] Chain 3: 50% [0:00:00 of 0:00:01 remaining] Chain 3: 60% [0:00:00 of 0:00:01 remaining] Chain 3: 70% [0:00:00 of 0:00:01 remaining] Chain 3: 80% [0:00:00 of 0:00:01 remaining] Chain 3: 90% [0:00:00 of 0:00:00 remaining] Chain 3: 100% [0:00:00 of 0:00:00 remaining] DIC Effective Parameters pD 13.114703 0.21530766 pV 26.573505 6.94470864
dic(sim3) |> show
MCMC Processing of 4875 Iterations x 3 Chains... Chain 1: 0% [0:00:00 of 0:00:00 remaining] Chain 1: 10% [0:00:00 of 0:00:00 remaining] Chain 1: 20% [0:00:00 of 0:00:00 remaining] Chain 1: 30% [0:00:00 of 0:00:00 remaining] Chain 1: 40% [0:00:00 of 0:00:00 remaining] Chain 1: 50% [0:00:00 of 0:00:00 remaining] Chain 1: 60% [0:00:00 of 0:00:00 remaining] Chain 1: 70% [0:00:00 of 0:00:00 remaining] Chain 1: 80% [0:00:00 of 0:00:00 remaining] Chain 1: 90% [0:00:00 of 0:00:00 remaining] Chain 1: 100% [0:00:00 of 0:00:00 remaining] Chain 2: 0% [0:01:13 of 0:01:13 remaining] Chain 2: 10% [0:00:01 of 0:00:02 remaining] Chain 2: 20% [0:00:01 of 0:00:01 remaining] Chain 2: 30% [0:00:00 of 0:00:01 remaining] Chain 2: 40% [0:00:00 of 0:00:01 remaining] Chain 2: 50% [0:00:00 of 0:00:00 remaining] Chain 2: 60% [0:00:00 of 0:00:00 remaining] Chain 2: 70% [0:00:00 of 0:00:00 remaining] Chain 2: 80% [0:00:00 of 0:00:00 remaining] Chain 2: 90% [0:00:00 of 0:00:00 remaining] Chain 2: 100% [0:00:00 of 0:00:00 remaining] Chain 3: 0% [0:02:28 of 0:02:29 remaining] Chain 3: 10% [0:00:03 of 0:00:03 remaining] Chain 3: 20% [0:00:01 of 0:00:02 remaining] Chain 3: 30% [0:00:01 of 0:00:01 remaining] Chain 3: 40% [0:00:01 of 0:00:01 remaining] Chain 3: 50% [0:00:00 of 0:00:01 remaining] Chain 3: 60% [0:00:00 of 0:00:01 remaining] Chain 3: 70% [0:00:00 of 0:00:01 remaining] Chain 3: 80% [0:00:00 of 0:00:01 remaining] Chain 3: 90% [0:00:00 of 0:00:00 remaining] Chain 3: 100% [0:00:00 of 0:00:00 remaining] DIC Effective Parameters pD 13.326638 0.49244858 pV 25.339951 6.49910490
## Subset Sampler Output
sim = sim1[1000:5000, ["beta[1]", "beta[2]"], :]
describe(sim)
Iterations = 1000:5000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 2001 Empirical Posterior Estimates: Mean SD Naive SE MCSE ESS beta[1] 0.53682913 1.10492601 0.0142609687 0.025219002 1919.5994 beta[2] 0.81491586 0.33253903 0.0042919876 0.007129130 2001.0000 Quantiles: 2.5% 25.0% 50.0% 75.0% 97.5% beta[1] -1.64673860 -0.059486655 0.55755745 1.15954184 2.6330071 beta[2] 0.19097135 0.634450846 0.81168506 0.98640910 1.4749587
## Write to and Read from an External File
write("sim1.jls", sim1)
sim1 = read("sim1.jls", ModelChains)
Object of type "Mamba.ModelChains" Iterations = 252:10000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 4875 [0.918604 1.08849 0.719733; 0.291973 1.51878 0.629041; … ; 0.480511 1.04643 0.783623; 0.308753 0.843837 0.772927] [0.31982 0.35836 0.997317; 0.175972 0.331401 0.769161; … ; 0.475919 0.647507 0.865755; 0.308566 0.647507 0.865755] [0.539945 1.2289 0.722182; 0.543182 -0.0625791 1.07568; … ; 1.21675 2.40923 0.38377; 3.39815 1.94268 0.472419]
## Restart the Sampler
sim = mcmc(sim1, 5000)
describe(sim)
MCMC Simulation of 5000 Iterations x 3 Chains... Chain 1: 0% [0:00:10 of 0:00:10 remaining] Chain 1: 10% [0:00:05 of 0:00:06 remaining] Chain 1: 20% [0:00:06 of 0:00:07 remaining] Chain 1: 30% [0:00:05 of 0:00:07 remaining] Chain 1: 40% [0:00:04 of 0:00:06 remaining] Chain 1: 50% [0:00:03 of 0:00:06 remaining] Chain 1: 60% [0:00:02 of 0:00:06 remaining] Chain 1: 70% [0:00:02 of 0:00:06 remaining] Chain 1: 80% [0:00:01 of 0:00:06 remaining] Chain 1: 90% [0:00:01 of 0:00:07 remaining] Chain 1: 100% [0:00:00 of 0:00:07 remaining] Chain 2: 0% [0:00:06 of 0:00:06 remaining] Chain 2: 10% [0:00:05 of 0:00:06 remaining] Chain 2: 20% [0:00:05 of 0:00:06 remaining] Chain 2: 30% [0:00:05 of 0:00:07 remaining] Chain 2: 40% [0:00:04 of 0:00:07 remaining] Chain 2: 50% [0:00:03 of 0:00:07 remaining] Chain 2: 60% [0:00:03 of 0:00:07 remaining] Chain 2: 70% [0:00:02 of 0:00:07 remaining] Chain 2: 80% [0:00:01 of 0:00:07 remaining] Chain 2: 90% [0:00:01 of 0:00:07 remaining] Chain 2: 100% [0:00:00 of 0:00:07 remaining] Chain 3: 0% [0:00:12 of 0:00:12 remaining] Chain 3: 10% [0:00:05 of 0:00:06 remaining] Chain 3: 20% [0:00:05 of 0:00:06 remaining] Chain 3: 30% [0:00:04 of 0:00:06 remaining] Chain 3: 40% [0:00:04 of 0:00:06 remaining] Chain 3: 50% [0:00:03 of 0:00:06 remaining] Chain 3: 60% [0:00:03 of 0:00:06 remaining] Chain 3: 70% [0:00:02 of 0:00:07 remaining] Chain 3: 80% [0:00:01 of 0:00:07 remaining] Chain 3: 90% [0:00:01 of 0:00:07 remaining] Chain 3: 100% [0:00:00 of 0:00:07 remaining] Iterations = 252:15000 Thinning interval = 2 Chains = 1,2,3 Samples per chain = 7375 Empirical Posterior Estimates: Mean SD Naive SE MCSE ESS s2 1.32570551 2.1403772 0.0143895985 0.097462523 482.28675 beta[1] 0.58688944 1.1955036 0.0080372825 0.013073860 7375.00000 beta[2] 0.80096933 0.3593708 0.0024160235 0.003736499 7375.00000 Quantiles: 2.5% 25.0% 50.0% 75.0% 97.5% s2 0.17447190 0.39188898 0.67500933 1.32734455 7.1426628 beta[1] -1.77547058 -0.03033497 0.58313055 1.19879591 2.9924841 beta[2] 0.09708961 0.61962761 0.80333058 0.98158896 1.4994692
## Plotting
## Default summary plot (trace and density plots)
p1 = plot(sim1)
2×3 Array{Gadfly.Plot,2}: Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.Identity(), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000] y=[0.918604, 0.291973, 0.499291, 1.42713, 1.18368, 0.682828, 0.837245, 0.498076, 0.74477, 0.397509 … 1.2082, 1.28321, 0.634384, 1.23317, 1.1752, 1.11254, 1.2893, 1.31298, 1.21675, 3.39815] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:y, "y"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Iteration", :horizontal), Gadfly.Guide.YLabel("Value", :vertical), Gadfly.Guide.Title("s2")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:y, [0.918604, 0.291973, 0.499291, 1.42713, 1.18368, 0.682828, 0.837245, 0.498076, 0.74477, 0.397509 … 1.2082, 1.28321, 0.634384, 1.23317, 1.1752, 1.11254, 1.2893, 1.31298, 1.21675, 3.39815]),Pair{Symbol,Any}(:x, [252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000]))) … Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.Identity(), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000] y=[0.719733, 0.629041, 0.629041, 1.07921, 0.435366, 0.570586, 0.8951, 0.916515, 0.916515, 1.00009 … 0.686276, 0.729036, 0.966384, 0.709008, 0.958812, 1.4681, 0.279791, 0.161171, 0.38377, 0.472419] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:y, "y"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Iteration", :horizontal), Gadfly.Guide.YLabel("Value", :vertical), Gadfly.Guide.Title("beta[2]")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:y, [0.719733, 0.629041, 0.629041, 1.07921, 0.435366, 0.570586, 0.8951, 0.916515, 0.916515, 1.00009 … 0.686276, 0.729036, 0.966384, 0.709008, 0.958812, 1.4681, 0.279791, 0.161171, 0.38377, 0.472419]),Pair{Symbol,Any}(:x, [252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000]))) Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.DensityStatistic(256, -Inf), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[0.918604, 0.291973, 0.499291, 1.42713, 1.18368, 0.682828, 0.837245, 0.498076, 0.74477, 0.397509 … 1.2082, 1.28321, 0.634384, 1.23317, 1.1752, 1.11254, 1.2893, 1.31298, 1.21675, 3.39815] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Value", :horizontal), Gadfly.Guide.YLabel("Density", :vertical), Gadfly.Guide.Title("s2")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:x, [0.918604, 0.291973, 0.499291, 1.42713, 1.18368, 0.682828, 0.837245, 0.498076, 0.74477, 0.397509 … 1.2082, 1.28321, 0.634384, 1.23317, 1.1752, 1.11254, 1.2893, 1.31298, 1.21675, 3.39815]))) Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.DensityStatistic(256, -Inf), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[0.719733, 0.629041, 0.629041, 1.07921, 0.435366, 0.570586, 0.8951, 0.916515, 0.916515, 1.00009 … 0.686276, 0.729036, 0.966384, 0.709008, 0.958812, 1.4681, 0.279791, 0.161171, 0.38377, 0.472419] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Value", :horizontal), Gadfly.Guide.YLabel("Density", :vertical), Gadfly.Guide.Title("beta[2]")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:x, [0.719733, 0.629041, 0.629041, 1.07921, 0.435366, 0.570586, 0.8951, 0.916515, 0.916515, 1.00009 … 0.686276, 0.729036, 0.966384, 0.709008, 0.958812, 1.4681, 0.279791, 0.161171, 0.38377, 0.472419])))
p2 = plot(sim2)
2×3 Array{Gadfly.Plot,2}: Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.Identity(), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000] y=[0.373215, 0.360428, 0.360428, 0.385009, 0.385009, 1.49714, 0.514227, 0.471849, 0.471849, 0.623983 … 1.02593, 1.54352, 0.656545, 0.721253, 0.80178, 1.19915, 0.636863, 0.628456, 1.00966, 0.828978] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:y, "y"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Iteration", :horizontal), Gadfly.Guide.YLabel("Value", :vertical), Gadfly.Guide.Title("s2")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:y, [0.373215, 0.360428, 0.360428, 0.385009, 0.385009, 1.49714, 0.514227, 0.471849, 0.471849, 0.623983 … 1.02593, 1.54352, 0.656545, 0.721253, 0.80178, 1.19915, 0.636863, 0.628456, 1.00966, 0.828978]),Pair{Symbol,Any}(:x, [252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000]))) … Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.Identity(), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000] y=[0.795453, 0.790722, 0.790722, 0.861822, 0.861822, 1.2972, 0.995092, 0.873184, 0.873184, 0.563108 … 0.823057, 0.55214, 0.476842, 1.0473, 0.710604, 0.799717, 0.613455, 0.833058, 0.329111, 0.562846] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:y, "y"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Iteration", :horizontal), Gadfly.Guide.YLabel("Value", :vertical), Gadfly.Guide.Title("beta[2]")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:y, [0.795453, 0.790722, 0.790722, 0.861822, 0.861822, 1.2972, 0.995092, 0.873184, 0.873184, 0.563108 … 0.823057, 0.55214, 0.476842, 1.0473, 0.710604, 0.799717, 0.613455, 0.833058, 0.329111, 0.562846]),Pair{Symbol,Any}(:x, [252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000]))) Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.DensityStatistic(256, -Inf), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[0.373215, 0.360428, 0.360428, 0.385009, 0.385009, 1.49714, 0.514227, 0.471849, 0.471849, 0.623983 … 1.02593, 1.54352, 0.656545, 0.721253, 0.80178, 1.19915, 0.636863, 0.628456, 1.00966, 0.828978] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Value", :horizontal), Gadfly.Guide.YLabel("Density", :vertical), Gadfly.Guide.Title("s2")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:x, [0.373215, 0.360428, 0.360428, 0.385009, 0.385009, 1.49714, 0.514227, 0.471849, 0.471849, 0.623983 … 1.02593, 1.54352, 0.656545, 0.721253, 0.80178, 1.19915, 0.636863, 0.628456, 1.00966, 0.828978]))) Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.DensityStatistic(256, -Inf), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[0.795453, 0.790722, 0.790722, 0.861822, 0.861822, 1.2972, 0.995092, 0.873184, 0.873184, 0.563108 … 0.823057, 0.55214, 0.476842, 1.0473, 0.710604, 0.799717, 0.613455, 0.833058, 0.329111, 0.562846] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Value", :horizontal), Gadfly.Guide.YLabel("Density", :vertical), Gadfly.Guide.Title("beta[2]")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:x, [0.795453, 0.790722, 0.790722, 0.861822, 0.861822, 1.2972, 0.995092, 0.873184, 0.873184, 0.563108 … 0.823057, 0.55214, 0.476842, 1.0473, 0.710604, 0.799717, 0.613455, 0.833058, 0.329111, 0.562846])))
p3 = plot(sim3)
2×3 Array{Gadfly.Plot,2}: Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.Identity(), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000] y=[0.558365, 0.504818, 1.68034, 0.449903, 0.434914, 0.676003, 0.420836, 0.79273, 1.15149, 0.512384 … 0.886454, 0.246371, 0.765125, 0.140835, 1.71716, 0.213482, 0.206701, 4.09362, 0.247849, 0.270629] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:y, "y"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Iteration", :horizontal), Gadfly.Guide.YLabel("Value", :vertical), Gadfly.Guide.Title("s2")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:y, [0.558365, 0.504818, 1.68034, 0.449903, 0.434914, 0.676003, 0.420836, 0.79273, 1.15149, 0.512384 … 0.886454, 0.246371, 0.765125, 0.140835, 1.71716, 0.213482, 0.206701, 4.09362, 0.247849, 0.270629]),Pair{Symbol,Any}(:x, [252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000]))) … Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.Identity(), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000] y=[0.721287, 0.395813, 0.83418, 1.01035, 0.810929, 0.69633, 0.706221, 1.0258, 0.56517, 0.56345 … 0.526292, 0.64138, 0.956293, 0.689051, 1.45857, 0.831149, 0.759107, 1.38993, 0.7765, 0.779213] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:y, "y"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Iteration", :horizontal), Gadfly.Guide.YLabel("Value", :vertical), Gadfly.Guide.Title("beta[2]")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:y, [0.721287, 0.395813, 0.83418, 1.01035, 0.810929, 0.69633, 0.706221, 1.0258, 0.56517, 0.56345 … 0.526292, 0.64138, 0.956293, 0.689051, 1.45857, 0.831149, 0.759107, 1.38993, 0.7765, 0.779213]),Pair{Symbol,Any}(:x, [252, 254, 256, 258, 260, 262, 264, 266, 268, 270 … 9982, 9984, 9986, 9988, 9990, 9992, 9994, 9996, 9998, 10000]))) Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.DensityStatistic(256, -Inf), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[0.558365, 0.504818, 1.68034, 0.449903, 0.434914, 0.676003, 0.420836, 0.79273, 1.15149, 0.512384 … 0.33802, 0.886454, 0.246371, 0.765125, 1.71716, 0.213482, 0.206701, 4.09362, 0.247849, 0.270629] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Value", :horizontal), Gadfly.Guide.YLabel("Density", :vertical), Gadfly.Guide.Title("s2")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:x, [0.558365, 0.504818, 1.68034, 0.449903, 0.434914, 0.676003, 0.420836, 0.79273, 1.15149, 0.512384 … 0.33802, 0.886454, 0.246371, 0.765125, 1.71716, 0.213482, 0.206701, 4.09362, 0.247849, 0.270629]))) Gadfly.Plot(Gadfly.Layer[Gadfly.Layer(nothing, Dict{Any,Any}(), Gadfly.StatisticElement[], Gadfly.Geom.LineGeometry(Gadfly.Stat.DensityStatistic(256, -Inf), false, 2, Symbol("")), nothing, 0)], nothing, Data( x=[0.721287, 0.395813, 0.83418, 1.01035, 0.810929, 0.69633, 0.706221, 1.0258, 0.56517, 0.56345 … 0.526292, 0.64138, 0.956293, 0.689051, 1.45857, 0.831149, 0.759107, 1.38993, 0.7765, 0.779213] color=[1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3] titles=Dict{Symbol,AbstractString}(Pair{Symbol,AbstractString}(:color, "color"),Pair{Symbol,AbstractString}(:x, "x")) ) , Gadfly.ScaleElement[Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true)], Gadfly.StatisticElement[], nothing, Gadfly.GuideElement[Gadfly.Guide.ColorKey("Chain"), Gadfly.Guide.XLabel("Value", :horizontal), Gadfly.Guide.YLabel("Density", :vertical), Gadfly.Guide.Title("beta[2]")], Gadfly.Theme(LCHab{Float32}(70.0,60.0,240.0), 0.9mm, 0.45mm, 1.8mm, Function[Compose.circle, Gadfly.square, Gadfly.diamond, Gadfly.cross, Gadfly.xcross, Gadfly.utriangle, Gadfly.dtriangle, Gadfly.star1, Gadfly.star2, Gadfly.hexagon, Gadfly.octagon, Gadfly.hline, Gadfly.vline], 0.3mm, :solid, nothing, nothing, 1.0, nothing, 5.0mm, RGB{N0f8}(0.816,0.816,0.878), Measures.Length{:mm,Float64}[0.5mm, 0.5mm], RGB{N0f8}(0.627,0.627,0.627), 0.2mm, "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.424,0.376,0.42), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.337,0.29,0.333), "'PT Sans Caption','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 3.88056mm, RGB{N0f8}(0.212,0.165,0.208), "'PT Sans','Helvetica Neue','Helvetica',sans-serif", 2.82222mm, RGB{N0f8}(0.298,0.251,0.294), 40, -0.05mm, 1.0mm, 3.0mm, Gadfly.default_stroke_color, 0.3mm, Gadfly.default_discrete_highlight_color, Gadfly.default_continuous_highlight_color, Gadfly.default_lowlight_color, 0.6, Gadfly.default_middle_color, 0.6mm, :left, :square, :none, nothing, 2.0mm, 1000, 10.0, 0.5, 0.2, 4, Gadfly.Scale.DiscreteColorScale(Gadfly.Scale.default_discrete_colors, nothing, nothing, true), Gadfly.Scale.ContinuousColorScale(Gadfly.Scale.#71, Gadfly.Scale.ContinuousScaleTransform(identity, identity, Gadfly.Scale.identity_formatter), nothing, nothing)), Dict{Symbol,Any}(Pair{Symbol,Any}(:color, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1 … 3, 3, 3, 3, 3, 3, 3, 3, 3, 3]),Pair{Symbol,Any}(:x, [0.721287, 0.395813, 0.83418, 1.01035, 0.810929, 0.69633, 0.706221, 1.0258, 0.56517, 0.56345 … 0.526292, 0.64138, 0.956293, 0.689051, 1.45857, 0.831149, 0.759107, 1.38993, 0.7765, 0.779213])))
draw(p1)
draw(p2)
draw(p3)
#draw(plot(sim1))
## Write plot to file
draw(p1, filename="summaryplot.svg")
#draw(p, filename="summaryplot.pdf", fmt=:pdf)
## Autocorrelation and running mean plots
p = plot(sim1, [:autocor, :mean], legend=true)
draw(p, nrow=3, ncol=2)
#draw(p, nrow=3, ncol=2, filename="autocormeanplot.svg")
#draw(p, nrow=3, ncol=2, filename="autocormeanplot.pdf", fmt=:pdf)
p2 = plot(sim2, [:autocor, :mean], legend=true)
draw(p, nrow=3, ncol=2)
p3 = plot(sim3, [:autocor, :mean], legend=true)
draw(p, nrow=3, ncol=2)
## Pairwise contour plots
p1 = plot(sim1, :contour)
draw(p1, nrow=2, ncol=2)
#draw(p, nrow=2, ncol=2, filename="contourplot.svg")
#draw(p, nrow=2, ncol=2, filename="contourplot.pdf", fmt=:pdf)
p2 = plot(sim2, :contour)
draw(p2, nrow=2, ncol=2)
p3 = plot(sim3, :contour)
draw(p3, nrow=2, ncol=2)
## Development and Testing
setinputs!(model, line) # Set input node values
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "0-element Mamba.ArrayStochastic{1}" Float64[] ------------------------------------------------------------------------------- y: An unmonitored node of type "0-element Mamba.ArrayStochastic{1}" Float64[] ------------------------------------------------------------------------------- mu: An unmonitored node of type "0-element Mamba.ArrayLogical{1}" Float64[] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; … ; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" NaN
setinits!(model, inits[1]) # Set initial values
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "2-element Mamba.ArrayStochastic{1}" [1.19027, 2.04818] ------------------------------------------------------------------------------- y: An unmonitored node of type "5-element Mamba.ArrayStochastic{1}" [1.0, 3.0, 3.0, 3.0, 5.0] ------------------------------------------------------------------------------- mu: An unmonitored node of type "5-element Mamba.ArrayLogical{1}" [3.23845, 5.28663, 7.33481, 9.38299, 11.4312] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; … ; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" 1.634388973402134
setsamplers!(model, scheme1) # Set sampling scheme
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "2-element Mamba.ArrayStochastic{1}" [1.19027, 2.04818] ------------------------------------------------------------------------------- y: An unmonitored node of type "5-element Mamba.ArrayStochastic{1}" [1.0, 3.0, 3.0, 3.0, 5.0] ------------------------------------------------------------------------------- mu: An unmonitored node of type "5-element Mamba.ArrayLogical{1}" [3.23845, 5.28663, 7.33481, 9.38299, 11.4312] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; … ; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" 1.634388973402134
showall(model) # Show detailed node information
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "2-element Mamba.ArrayStochastic{1}" [1.19027, 2.04818] Distribution: ZeroMeanIsoNormal( dim: 2 μ: [0.0, 0.0] Σ: [1000.0 0.0; 0.0 1000.0] ) Function: CodeInfo(:(begin $(Expr(:inbounds, false)) # meta: location In[3] #3 16 SSAValue(0) = (Base.Math.sqrt_llvm)((Base.sitofp)(Float64, 1000)::Float64)::Float64 # meta: location /root/.julia/v0.6/Distributions/src/multivariate/mvnormal.jl Type 224 SSAValue(1) = (Base.mul_float)(SSAValue(0), SSAValue(0))::Float64 SSAValue(4) = (Base.div_float)((Base.sitofp)(Float64, 1)::Float64, SSAValue(1))::Float64 # meta: pop location # meta: pop location $(Expr(:inbounds, :pop)) return $(Expr(:new, Distributions.MvNormal{Float64,PDMats.ScalMat{Float64},Distributions.ZeroVector{Float64}}, :($(Expr(:new, Distributions.ZeroVector{Float64}, 2))), :($(Expr(:new, PDMats.ScalMat{Float64}, 2, SSAValue(1), SSAValue(4)))))) end))=>Distributions.MvNormal{Float64,PDMats.ScalMat{Float64},Distributions.ZeroVector{Float64}} Source Nodes: Symbol[] Target Nodes: Symbol[:mu, :y] ------------------------------------------------------------------------------- y: An unmonitored node of type "5-element Mamba.ArrayStochastic{1}" [1.0, 3.0, 3.0, 3.0, 5.0] Distribution: IsoNormal( dim: 5 μ: [3.23845, 5.28663, 7.33481, 9.38299, 11.4312] Σ: [1.63439 0.0 0.0 0.0 0.0; 0.0 1.63439 0.0 0.0 0.0; 0.0 0.0 1.63439 0.0 0.0; 0.0 0.0 0.0 1.63439 0.0; 0.0 0.0 0.0 0.0 1.63439] ) Function: CodeInfo(:(begin SSAValue(1) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:mu))) SSAValue(0) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:s2))) return (Main.MvNormal)(SSAValue(1), (Main.sqrt)(SSAValue(0))) end))=>Any Source Nodes: Symbol[:mu, :s2] Target Nodes: Symbol[] ------------------------------------------------------------------------------- mu: An unmonitored node of type "5-element Mamba.ArrayLogical{1}" [3.23845, 5.28663, 7.33481, 9.38299, 11.4312] Function: CodeInfo(:(begin SSAValue(1) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:xmat))) SSAValue(0) = $(Expr(:invoke, MethodInstance for getindex(::Mamba.Model, ::Symbol), :(Main.getindex), :(model), :(:beta))) return SSAValue(1) * SSAValue(0) end))=>Any Source Nodes: Symbol[:xmat, :beta] Target Nodes: Symbol[:y] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; 1.0 3.0; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" 1.634388973402134 Distribution: Distributions.InverseGamma{Float64}( invd: Distributions.Gamma{Float64}(α=0.001, θ=1000.0) θ: 0.001 ) Function: CodeInfo(:(begin return $(Expr(:invoke, MethodInstance for Distributions.InverseGamma{Float64}(::Float64, ::Float64), Distributions.InverseGamma{Float64}, 0.001, 0.001)) end))=>Distributions.InverseGamma{Float64} Source Nodes: Symbol[] Target Nodes: Symbol[:y]
logpdf(model, 1) # Log-density sum for block 1
-48.56941725488224
logpdf(model, 2) # Block 2
-47.22743779767911
logpdf(model) # All blocks
-55.97587603194255
sample!(model, 1) # Sample values for block 1
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "2-element Mamba.ArrayStochastic{1}" [0.517366, 1.47502] ------------------------------------------------------------------------------- y: An unmonitored node of type "5-element Mamba.ArrayStochastic{1}" [1.0, 3.0, 3.0, 3.0, 5.0] ------------------------------------------------------------------------------- mu: An unmonitored node of type "5-element Mamba.ArrayLogical{1}" [1.99239, 3.46741, 4.94244, 6.41746, 7.89248] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; … ; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" 1.634388973402134
sample!(model, 2) # Block 2
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "2-element Mamba.ArrayStochastic{1}" [0.517366, 1.47502] ------------------------------------------------------------------------------- y: An unmonitored node of type "5-element Mamba.ArrayStochastic{1}" [1.0, 3.0, 3.0, 3.0, 5.0] ------------------------------------------------------------------------------- mu: An unmonitored node of type "5-element Mamba.ArrayLogical{1}" [1.99239, 3.46741, 4.94244, 6.41746, 7.89248] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; … ; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" 2.3113776949290785
sample!(model) # All blocks
Object of type "Mamba.Model" ------------------------------------------------------------------------------- beta: A monitored node of type "2-element Mamba.ArrayStochastic{1}" [0.988215, 0.187305] ------------------------------------------------------------------------------- y: An unmonitored node of type "5-element Mamba.ArrayStochastic{1}" [1.0, 3.0, 3.0, 3.0, 5.0] ------------------------------------------------------------------------------- mu: An unmonitored node of type "5-element Mamba.ArrayLogical{1}" [1.17552, 1.36282, 1.55013, 1.73743, 1.92474] ------------------------------------------------------------------------------- xmat: [1.0 1.0; 1.0 2.0; … ; 1.0 4.0; 1.0 5.0] ------------------------------------------------------------------------------- s2: A monitored node of type "Mamba.ScalarStochastic" 3.8164143777752413