Is there a Julia package for sampling from a uniform distribution over a convex polytope?
I’m trying to sample a categorical distribution from an ambiguity set of categorical distributions subject to linear constraints. Here’s a minimum working example to sample a categorical distribution with k=3
using the hitandrun
R package (source and documentation):
# sample a categorical probability distribution subject to polytopic constraints
using RCall
ndimension = 3 # number of dimensions for the categorical probability distribution
plower = [0.2, 0.1, 0.4] # lower interval limits for p1, p2, p3
pupper = [0.8, 0.9, 0.6] # upper interval limits for p1, p2, p3
A = [-1.0 0.0 0.0;
1.0 0.0 0.0;
0.0 -1.0 0.0;
0.0 1.0 0.0;
0.0 0.0 -1.0;
0.0 0.0 1.0;
1.0 1.0 1.0] # constraint matrix LHS
b = [-0.2, 0.8, -0.1, 0.9, -0.4, 0.6, 1.0] # constraint matrix RHS
d = vcat(fill("<=", 2n),"=") # sense
nsample = 10 # number of samples
rseed = 9590 # random seed
rout = R"""
library(hitandrun)
set.seed($rseed)
constr <- list(constr = $A, rhs = $b, dir = $d)
samples <- hitandrun(constr, n.samples = $nsample, thin = ($ndimension) ^ 3)
"""
samples = rcopy(rout) # each row is a categorical probability distribution
julia> samples
10×3 Array{Float64,2}:
0.418657 0.132267 0.449077
0.228864 0.304185 0.466951
0.379767 0.190112 0.430121
0.293098 0.159821 0.547081
0.236806 0.245027 0.518167
0.36011 0.156859 0.483031
0.378486 0.162136 0.459377
0.332935 0.134134 0.532931
0.342103 0.129685 0.528212
0.252632 0.248207 0.499161
I found the Mamba package, which seems like a thorough, high quality package. This problem seems relatively straightforward, so I was expecting one of the Mamba methods to apply, but as far as I can tell none of the methods are applicable. The closest ones I found are the slice simplex and the shrinkage slice, but I couldn’t formulate either one to address my problem. Am I missing something? Can one of those (or any other method in Mamba
) be used for the same purpose as hitandrun
?
Using RCall is pain-free workaround (amazing how easy it worked), but I just want to make sure that I’m not missing something in the Mamba package that would let me do it natively. Maybe this is more appropriate for the Mamba github page as an issue/suggestion?
Thanks in advance for any advice!