set.seed(1)
ploidy <- 8
## Simulate under the null
p <- stats::runif(ploidy / 2 + 1)
p <- p / sum(p)
q <- stats::convolve(p, rev(p), type = "open")
nvec <- c(stats::rmultinom(n = 1, size = 100, prob = q))
gl <- simgl(nvec)
gibbs_gl(gl = gl, alpha = rep(1, ploidy / 2 + 1), lg = TRUE)
Run the code above in your browser using DataLab