py · September 16, 2013 17:33 · Sep 16, 2013 · Sep 16, 2013 · Sep 16, 2013
diff --git a/bayesian_updating_distributions.R b/bayesian_updating_distributions.R
@@ -7,36 +7,46 @@ library('ggplot2')
 
 # Define functions 
 normalize <- function(pos) {
+  # Normalizes the possibilities so that sum = 1
   pos_sum <- sum(pos)
   pos <- pos/pos_sum
   return(pos)
 }
 
 # Update function - heads
 update_success <- function(pos) {
-  pos <- pos * ind/100
+  # Run when heads is flipped. Updates the possibilities by multiplying each 
+  # current possibility by its hypothesis value, then normalizing.
+  pos <- pos * hyp
   pos <- normalize(pos)
   return(pos)
 }
 
 # Update function - tails
 update_failure <- function(pos) {
-  pos <- pos * (1 - ind/100)
+  # Run when tails is flipped. Updates the possibilities by multiplying each
+  # current possibility by the probability the hypothesis is false (1-hyp).
+  pos <- pos * (1 - hyp)
   pos <- normalize(pos)
   return(pos)
 }
 
 # Plotting function - Plot Distribution
 pd <- function(pos) {
-  qplot(x=ind/100, y=pos, geom="bar", stat="identity") + 
+  # Function creates a qplot to display likelihood of hypothesis vs.
+  # hypothesis for chance of heads.
+  # y-axis limits set to .05 for consistency across plots, but may need
+  # to adjust if different prior or really unfair coin used.
+  qplot(x=hyp, y=pos, geom="bar", stat="identity", fill = I("blue")) +
     geom_vline(x=.50, color="red") +
     xlab("Hypothesis for chance of heads") + 
-    ylab("Likelihood of hypothesis")
+    ylab("Likelihood of hypothesis") + 
+    ylim(0,0.05)
 }
 
 # Initial state
-pos <- rep(1,101)
-ind <- 0:100
+pos <- rep(1,101) # Possibilities (starting with uniform distribution)
+hyp <- 0:100/100  # Hypotheses for chance of heads (0% to 100%)
 pos <- normalize(pos)
 pd(pos)
 

diff --git a/bayesian_updating_distributions → bayesian_updating_distributions.R b/bayesian_updating_distributions → bayesian_updating_distributions.R
diff --git a/bayesian_updating_distributions b/bayesian_updating_distributions
@@ -0,0 +1,62 @@
+#------------------------------------------------------------------------------#
+# Bayesian update with new information
+# Based on http://www.databozo.com/2013/09/15/Bayesian_updating_of_probability_distributions.html
+#------------------------------------------------------------------------------#
+# Load packages
+library('ggplot2')
+
+# Define functions 
+normalize <- function(pos) {
+  pos_sum <- sum(pos)
+  pos <- pos/pos_sum
+  return(pos)
+}
+
+# Update function - heads
+update_success <- function(pos) {
+  pos <- pos * ind/100
+  pos <- normalize(pos)
+  return(pos)
+}
+
+# Update function - tails
+update_failure <- function(pos) {
+  pos <- pos * (1 - ind/100)
+  pos <- normalize(pos)
+  return(pos)
+}
+
+# Plotting function - Plot Distribution
+pd <- function(pos) {
+  qplot(x=ind/100, y=pos, geom="bar", stat="identity") + 
+    geom_vline(x=.50, color="red") +
+    xlab("Hypothesis for chance of heads") + 
+    ylab("Likelihood of hypothesis")
+}
+
+# Initial state
+pos <- rep(1,101)
+ind <- 0:100
+pos <- normalize(pos)
+pd(pos)
+
+# Start flipping
+# Heads (h=1, t=0)
+pos <- update_success(pos)
+pd(pos)
+
+# Tails (h=1, t=1)
+pos <- update_failure(pos)
+pd(pos)
+
+# Heads (h=2, t=1)
+pos <- update_success(pos)
+pd(pos)
+
+# Head (h=3, t=1)
+pos <- update_success(pos)
+pd(pos)
+
+# Tails (h=3, t=2)
+pos <- update_failure(pos)
+pd(pos)
No results found