Something not right with population sd estimate of gamma data?

[ProfessorPeregrine]

I am a frequentist dipping my toes into Bayes, trying to calculate credible intervals.

I am following along with your blog post and the notes in the code here, trying to put a general-purpose credible interval into an app.

I feel I am doing something wrong, because the estimates for sd are consistently off of what the actual sd is of a sample. I think the expected value of the posterior distribution should be the population standard deviation of the sample.

#' the data was generated from a random sample of the gamma distribution: samp<-rgamma(n=10,shape=5,scale=1)
samp<-c(5.353057075,
8.440596321,
2.253503205,
5.826479980,
5.711322996,
3.052636834,
4.008683698,
6.792975319,
3.553139157,
8.302714648)
samp<-data.frame(samp)
names(samp)<-"stress"
require(bayesboot)
#mean works as expected
bb_mean<-bayesboot(samp$stress,statistic=weighted.mean,use.weights=TRUE)
#keeping in mind we need the population sd:

pop.sd <- function(x) {
  n <- length(x)
  sd(x) * sqrt( (n - 1) / n)
}

bb_sd<-bayesboot(samp$stress,statistic = pop.sd)
summary(bb_sd)

I am consistently getting an estimate around 1.9 for the sd, when the theoretical distribution is 2.23 and the actual sd of the sample is:

sd(samp$stress)

is 2.123934876 and

pop.sd(samp$stress)

is 2.014941543

If I look at:

plot(bb_sd)

I see the posterior, but it is really not near 2.23 or 2.12 or 2.0. If I repeat the bootstrap, I never get near any of these, I am consistently around 1.9. I would have thought that multiple runs of bayesbootstrap() should give me numbers around 2.0 (although 2.12 has always been in the credible interval). I guess I would have expected that the posterior's mean would be above and below the actual sample pop.sd().

Am I doing something wrong?