Instantly share code, notes, and snippets.

# robbymeals/Y_samp.R Created Oct 31, 2012

What would you like to do?
Binomial Proportion Post Gists
 betaplot <- function(a,b){ theta = seq(0,1,0.005) p_theta = dbeta(theta, a, b) p <- qplot(theta, p_theta, geom='line') p <- p + theme_bw() p <- p + ylab(expression(paste('p(',theta,')', sep = ''))) p <- p + xlab(expression(theta)) return(p)}
 ## Generate population N = sample(seq(100000,400000),1) A = round(theta_true*N) B = N - A Zpop <- sample(c(rep(1,A),rep(0,B)))
 m = 0.95 n = 100 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)
 m = 0.05 n = 100 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)
 ### Pull random sample from population N_samp <- 500 Zsamp <- sample(Zpop,N_samp)
 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)
 ######### 1. Unknown Probability of Success ######## ## Using built-in R pseudo-random number generator theta_true <- runif(1,0,1)
 m = 0.5 n = 2 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)
 m = 0.5 n = 10 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)
 Y_samp <- sum(Zsamp)