#using Distributed
#@everywhere using MambaModels
using MambaModels

Data

globe_toss = Dict{Symbol, Any}(
  :w => [6, 7, 5, 6, 6],
  :n => [9, 9, 9, 9, 9]
)
globe_toss[:N] = length(globe_toss[:w]);

5

Model Specification

model = Model(
  w = Stochastic(1,
    (n, p, N) ->
      UnivariateDistribution[Binomial(n[i], p) for i in 1:N],
    false
  ),
  p = Stochastic(() -> Beta(1, 1))
);

Object of type "Model"
-------------------------------------------------------------------------------
w:
An unmonitored node of type "0-element ArrayStochastic{1}"
Float64[]
-------------------------------------------------------------------------------
p:
A monitored node of type "ScalarStochastic"
NaN

Initial Values

inits = [
  Dict(:w => globe_toss[:w], :n => globe_toss[:n], :p => 0.5),
  Dict(:w => globe_toss[:w], :n => globe_toss[:n], :p => rand(Beta(1, 1)))
];

2-element Array{Dict{Symbol,Any},1}:
 Dict(:w=>[6, 7, 5, 6, 6],:p=>0.5,:n=>[9, 9, 9, 9, 9])     
 Dict(:w=>[6, 7, 5, 6, 6],:p=>0.119402,:n=>[9, 9, 9, 9, 9])

Sampling Scheme

scheme = [NUTS(:p)]
setsamplers!(model, scheme);

Object of type "Model"
-------------------------------------------------------------------------------
w:
An unmonitored node of type "0-element ArrayStochastic{1}"
Float64[]
-------------------------------------------------------------------------------
p:
A monitored node of type "ScalarStochastic"
NaN

MCMC Simulations

chn = mcmc(model, globe_toss, inits, 10000, burnin=2500, thin=1, chains=2);

Object of type "ModelChains"

Iterations = 2501:10000
Thinning interval = 1
Chains = 1,2
Samples per chain = 7500

[0.571176; 0.686697; … ; 0.670241; 0.670241]

[0.625816; 0.695854; … ; 0.738604; 0.545944]

Describe draws

describe(chn)

Iterations = 2501:10000
Thinning interval = 1
Chains = 1,2
Samples per chain = 7500

Empirical Posterior Estimates:
    Mean       SD      Naive SE        MCSE      ESS
p 0.659196 0.0679238 0.0005545955 0.00067296178 7500

Quantiles:
   2.5%    25.0%     50.0%     75.0%      97.5%
p 0.5203 0.6154081 0.6616647 0.7062857 0.78600928

Convert to MCMCChains.Chains object

chn2 = MCMCChains.Chains(chn.value, String.(chn.names))

Object of type Chains, with data of type 7500×1×2 Array{Float64,3}

Iterations        = 1:7500
Thinning interval = 1
Chains            = 1, 2
Samples per chain = 7500
parameters        = p

2-element Array{ChainDataFrame,1}

Summary Statistics
. Omitted printing of 2 columns
│ Row │ parameters │ mean     │ std       │ naive_se    │ mcse        │
│     │ Symbol     │ Float64  │ Float64   │ Float64     │ Float64     │
├─────┼────────────┼──────────┼───────────┼─────────────┼─────────────┤
│ 1   │ p          │ 0.659196 │ 0.0679238 │ 0.000554596 │ 0.000672962 │

Quantiles

│ Row │ parameters │ 2.5%    │ 25.0%    │ 50.0%    │ 75.0%    │ 97.5%    │
│     │ Symbol     │ Float64 │ Float64  │ Float64  │ Float64  │ Float64  │
├─────┼────────────┼─────────┼──────────┼──────────┼──────────┼──────────┤
│ 1   │ p          │ 0.5203  │ 0.615408 │ 0.661665 │ 0.706286 │ 0.786009 │

Describe the MCMCChains

MCMCChains.describe(chn2)

2-element Array{ChainDataFrame,1}

Summary Statistics
. Omitted printing of 2 columns
│ Row │ parameters │ mean     │ std       │ naive_se    │ mcse        │
│     │ Symbol     │ Float64  │ Float64   │ Float64     │ Float64     │
├─────┼────────────┼──────────┼───────────┼─────────────┼─────────────┤
│ 1   │ p          │ 0.659196 │ 0.0679238 │ 0.000554596 │ 0.000672962 │

Quantiles

│ Row │ parameters │ 2.5%    │ 25.0%    │ 50.0%    │ 75.0%    │ 97.5%    │
│     │ Symbol     │ Float64 │ Float64  │ Float64  │ Float64  │ Float64  │
├─────┼────────────┼─────────┼──────────┼──────────┼──────────┼──────────┤
│ 1   │ p          │ 0.5203  │ 0.615408 │ 0.661665 │ 0.706286 │ 0.786009 │

Plot chn2

MCMCChains.plot(chn2)

End of 02/m2.1m.jl

This page was generated using Literate.jl.