library(MCMCpack) p = rdirichlet(100,rep(1,4))# create some artificial data ###### Shouldn't have to change anything from this point on # N is the number of stomachs / scat samples you have N = dim(p)[1] # prey is the number of items prey = dim(p)[2] findAlpha = function(a) { NLL = -sum(log(ddirichlet(p,a))) return(NLL) } o = optim(runif(prey),findAlpha) # these are the estimates of alpha for the prior o$par