So what about sparsity?
The Zig-Zag allows (that is what our paper is about) spike and slab priors for variable selection. Let’s demo this with Soss:
Let’s use MeasureTheory.jl 's SpikeMixture (a very spiky spike and slab with a Dirac measure as spike) to indicate that we have ground to believe the coefficients are actually zero:
m = @model x begin
α ~ SpikeMixture(Normal(), 0.2) # 0.2*Normal() + 0.8*Dirac(0)
β ~ SpikeMixture(Normal(), 0.2)
yhat = α .+ β .* x
y ~ For(eachindex(x)) do j
Normal(yhat[j], 2.0)
end
end
x = randn(20);
obs = -0.1 .+ 2x + 1randn(20);
T = 100.0
posterior = m(x=x) | (y=obs,)
sampling with the sticky ZigZag
include("sparsezigzag.jl")
trace, final, (num, acc) = @time sparse_zigzag(posterior, T, c=50)
ts, xs = ZigZagBoomerang.sep(discretize(trace, 0.1))
xs = xform(posterior).(xs) # make named tuples `NamedTuple{(:β, :α)`
now we can answer for example what is the posterior probability that the parameters are zero
julia> mean(getindex.(xs, :α) .== 0)
0.7211394302848576
julia> mean(getindex.(xs, :β) .== 0)
0.028985507246376812
The nice thing is that this scales very well…
