Ever since I read the book Statistics with Confidence by Doug Altman and colleagues(Altman et al. 2013) I have been a very enthusiastic user of confidence intervals. What I specifically want from my students is for them to think about the process of estimation and inference and to realise that point estimates are rarely valid and an interval is a much better way to present the data.
I am also not part of the statistics wars between the frequentists and Bayesians. I come from a Biochemistry background and in Biology Bayes should be the dominant paradigm. Not because of some philosophical reason but because we know so little and there is so much complexity.
I have been reading Sewall Wright’s Evolution and the Genetics of Populations(Wright 1984) and he makes it very clear when deriving the formulae in the first volume that a Bayesian approach is required.
As a simple example I am going to use STAN to calculate the Bayesian Credible interval based on the height data that I collected from my students over many years. I am using STAN and not BUGS because of the improved sampling in STAN which reduces the computational cost.
I have split the dataset based on gender for this particular example.
library(haven)library(rstan)library(ggplot2)data <-read_sav("Expanded Student Data Female.sav")options(mc.cores=parallel::detectCores())Y <- data$Heightfit <-stan("female_heights.stan",iter=1000,chains=4, data=list(Y=Y))
Inference for Stan model: anon_model.
4 chains, each with iter=1000; warmup=500; thin=1;
post-warmup draws per chain=500, total post-warmup draws=2000.
mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
mu 1.63 0.00 0.00 1.62 1.63 1.64 2192 1.00
sigma 0.07 0.00 0.00 0.07 0.07 0.08 691 1.00
aMax_indicator 0.63 0.01 0.48 0.00 1.00 1.00 1863 1.00
aMin_indicator 0.35 0.01 0.48 0.00 0.00 1.00 1941 1.00
lp__ 1076.05 0.04 1.09 1073.02 1076.38 1077.11 775 1.01
Samples were drawn using NUTS(diag_e) at Mon Dec 15 13:30:27 2025.
For each parameter, n_eff is a crude measure of effective sample size,
and Rhat is the potential scale reduction factor on split chains (at
convergence, Rhat=1).
mu <-data.frame(extract(fit,"mu")[[1]])colnames(mu) <-c("Mean")sigma <-extract(fit,"sigma")[[1]]ggplot(mu, aes(x=Mean))+geom_histogram()+labs(title="Histogram of the STAN Simulated Mean Heights for Female Students", x="Mean Height (m)")
From the summary at the end of the Stan simulations the 95% credible interval is from 1.62 to 1.64 m.
This compares to the 95% confidence interval which is 1.62 and 1.64 m in SPSS (which also gave a Bayesian credible interval of 1.62 to 1.64 m).
This is a lot of computational work to get very little benefit which is why I am a pragmatic Bayesian. In this case it was a simple question and we had a large sample which did not contain any messy data. In this case the Bayesian Credible Interval and Confidence Interval are the same. The difference is in the validity of the interpretation. We cannot interpret the confidence interval in the same way as the Bayesian Credible Interval, but as they converge to the same result I consider that it is rather pedantic than practical.
data <-read_sav("Expanded Student Data Male.sav")options(mc.cores=parallel::detectCores())Y <- data$Heightfit <-stan("male_heights.stan",iter=1000,chains=4, data=list(Y=Y))
Inference for Stan model: anon_model.
4 chains, each with iter=1000; warmup=500; thin=1;
post-warmup draws per chain=500, total post-warmup draws=2000.
mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
mu 1.77 0.00 0.01 1.75 1.77 1.78 2038 1
sigma 0.09 0.00 0.01 0.08 0.09 0.11 1454 1
aMax_indicator 0.07 0.01 0.26 0.00 0.00 1.00 1899 1
aMin_indicator 0.30 0.01 0.46 0.00 0.00 1.00 1692 1
lp__ 250.23 0.03 1.02 247.44 250.57 251.23 918 1
Samples were drawn using NUTS(diag_e) at Mon Dec 15 13:31:07 2025.
For each parameter, n_eff is a crude measure of effective sample size,
and Rhat is the potential scale reduction factor on split chains (at
convergence, Rhat=1).
mu <-as.data.frame(extract(fit,"mu")[[1]])colnames(mu) <-c("Mean")sigma <-extract(fit,"sigma")[[1]]ggplot(mu, aes(x=Mean))+geom_histogram()+labs(title="Histogram of the STAN Simulated Mean Heights for Male Students", x="Mean Height (m)")
From the summary at the end of the Stan simulations the 95% credible interval is from 1.75 to 1.78 m.
This compares to the 95% confidence interval which is 1.75 and 1.79 m in SPSS (which also gave a Bayesian credible interval of 1.75 to 1.79 m).
The mean was 1.77 in SPSS and it is 1.76 from the Bayesian calculations.
Small differences are starting to appear once you get noisier data and smaller sample sizes (137 compared to 509).
