m = 0.5 n = 500 pr <- prior(m,n) lk <- likelihood(N_samp,Y_samp) po <- posterior(m,n,N_samp,Y_samp) model_plot <- data.frame('Dist'=c(rep('Prior',nrow(pr)), rep('Likelihood',nrow(lk)), rep('Posterior',nrow(po))), rbind(pr,lk,po)) with(model_plot, Dist <- factor(Dist, levels = c('Prior', 'Likelihood', 'Posterior'), ordered = TRUE)) mean_po <- mean_of_posterior(m,n,N_samp,Y_samp) mode_po <- mode_of_posterior(m,n,N_samp,Y_samp) sd_po <- sd_of_posterior(m,n,N_samp,Y_samp)</code>