|
6 | 6 |
|
7 | 7 | import FigureHelpers as CM |
8 | 8 | using FigureHelpers |
| 9 | +using MCMCChains |
| 10 | + |
| 11 | + |
| 12 | +function CP.density_params(chns, pars=names(chns, :parameters); |
| 13 | + makie_config::MakieConfig=MakieConfig(), |
| 14 | + fig = figure_conf(cm2inch.((8.3,8.3/1.618)); makie_config), |
| 15 | + column = 1, xlims=nothing, |
| 16 | + labels=nothing, colors = nothing, ylabels = nothing, normalize = false, |
| 17 | + kwargs_axis = repeat([()],length(pars)), |
| 18 | + prange = (0.025, 0.975), # do not extend x-scale to outliers |
| 19 | + kwargs... |
| 20 | + ) |
| 21 | + n_chains = size(chns,3) |
| 22 | + n_samples = length(chns) |
| 23 | + labels_ch = isnothing(labels) ? string.(1:n_chains) : string.(labels) |
| 24 | + ylabels = isnothing(ylabels) ? string.(pars) : ylabels |
| 25 | + !isnothing(xlims) && (length(xlims) != length(pars)) && error( |
| 26 | + "Expected length(xlims)=$(length(xlims)) (each a Tuple or nothing) to be length(pars)=$(length(pars))") |
| 27 | + for (i, param) in enumerate(pars) |
| 28 | + ax = Axis(fig[i, column]; ylabel=ylabels[i], yaxisposition = :right, kwargs_axis[i]...) |
| 29 | + if isnothing(colors) |
| 30 | + pal = fig.scene.theme.palette |
| 31 | + #pal = ax.palette |
| 32 | + colors = pal.color[] |
| 33 | + end |
| 34 | + for i_chain in 1:n_chains |
| 35 | + _values = chns[:, param, i_chain] |
| 36 | + if !isnothing(prange) |
| 37 | + qmin,qmax = StatsBase.quantile(_values, prange) |
| 38 | + _values = _values[qmin .<= _values .<= qmax] |
| 39 | + end |
| 40 | + col = colors[i_chain] |
| 41 | + if normalize |
| 42 | + k = KernelDensity.kde(_values) |
| 43 | + md = maximum(k.density) |
| 44 | + lines!(ax, k.x, k.density ./ md; label=labels_ch[i_chain], color = col, kwargs...) |
| 45 | + else |
| 46 | + density!(ax, _values; label=labels_ch[i_chain], color = (col, 0.3), strokecolor = col, strokewidth = 1, |
| 47 | + #strokearound = true, |
| 48 | + kwargs...) |
| 49 | + end |
| 50 | + end |
| 51 | + xlim = CP.passnothing(getindex)(xlims, i) |
| 52 | + !isnothing(xlim) && xlims!(ax, xlim) |
| 53 | + #hideydecorations!(ax, ticklabels=false, ticks=false, grid=false) |
| 54 | + hideydecorations!(ax, label=false, ticklabels=true) |
| 55 | + # if i < length(params) |
| 56 | + # hidexdecorations!(ax; grid=false) |
| 57 | + # else |
| 58 | + # ax.xlabel = "Parameter estimate" |
| 59 | + # end |
| 60 | + end |
| 61 | + # axes = [only(contents(fig[i, 2])) for i in 1:length(params)] |
| 62 | + # linkxaxes!(axes...) |
| 63 | + #axislegend(only(contents(fig[2, column]))) |
| 64 | + fig |
| 65 | +end |
| 66 | + |
| 67 | +function histogram_params(chns, pars=names(chns, :parameters); |
| 68 | + makie_config::MakieConfig=MakieConfig(), |
| 69 | + fig = figure_conf(cm2inch.((8.3,8.3/1.618)); makie_config), |
| 70 | + column = 1, xlims=nothing, |
| 71 | + labels=nothing, colors = nothing, ylabels = nothing, normalization = :pdf, |
| 72 | + kwargs_axis = repeat([()],length(pars)), |
| 73 | + ylimits = (0,1), |
| 74 | + kwargs... |
| 75 | + ) |
| 76 | + n_chains = size(chns,3) |
| 77 | + n_samples = length(chns) |
| 78 | + labels_ch = isnothing(labels) ? string.(1:n_chains) : string.(labels) |
| 79 | + ylabels = isnothing(ylabels) ? string.(pars) : ylabels |
| 80 | + !isnothing(xlims) && (length(xlims) != length(pars)) && error( |
| 81 | + "Expected length(xlims)=$(length(xlims)) (each a Tuple or nothing) to be length(pars)=$(length(pars))") |
| 82 | + for (i, param) in enumerate(pars) |
| 83 | + ax = Axis(fig[i, column]; ylabel=ylabels[i], kwargs_axis[i]...) |
| 84 | + ylims!(ax, ylimits) |
| 85 | + if isnothing(colors) |
| 86 | + colors = ax.palette.color[] |
| 87 | + end |
| 88 | + for i_chain in 1:n_chains |
| 89 | + _values = chns[:, param, i_chain] |
| 90 | + col = colors[i_chain] |
| 91 | + hist!(ax, _values; label=labels_ch[i_chain], color = (col, 0.3), strokecolor = col, strokewidth = 1, |
| 92 | + normalization, |
| 93 | + #strokearound = true, |
| 94 | + kwargs...) |
| 95 | + end |
| 96 | + xlim = passnothing(getindex)(xlims, i) |
| 97 | + !isnothing(xlim) && xlims!(ax, xlim) |
| 98 | + #hideydecorations!(ax, ticklabels=false, ticks=false, grid=false) |
| 99 | + hideydecorations!(ax, label=false, ticklabels=true) |
| 100 | + # if i < length(params) |
| 101 | + # hidexdecorations!(ax; grid=false) |
| 102 | + # else |
| 103 | + # ax.xlabel = "Parameter estimate" |
| 104 | + # end |
| 105 | + end |
| 106 | + # axes = [only(contents(fig[i, 2])) for i in 1:length(params)] |
| 107 | + # linkxaxes!(axes...) |
| 108 | + #axislegend(only(contents(fig[2, column]))) |
| 109 | + fig |
| 110 | +end |
| 111 | + |
9 | 112 |
|
10 | 113 | end # module |
11 | 114 |
|
0 commit comments