Sa 16 Juli 2016 — under bayes, statistics
There is a lot of ink spent on the topic by most prominent scholars of all time so I really dont want to add anything to the debate.
The 1 paper which influenced me in getting a practical overview of the topic and showing tangiblly was by Jake Vanderplass , where he explains with 2 examples
It establishes from a practical prespective that Bayesian answer to the problem is near to the real world. But it made me thinking how do we know the real world answer? Either we need to play the table example for 1000 times or toss the coin 1000 times. The answer to the question is YES again which was provided by Jake himself in another (talk)[]. So in this work I try to experimentally validate both the experiments by running the examples 1000 times and finding the real values.
Where frequentist and Bayesian differ. [http://www.rpgroup.caltech.edu/courses/PBoC_CSHL_2015/files_2015/articles/eddy%20bayesian%20methods.pdf]
Problem is to find the conditional probability of Bob wining given that Alice has won 5 and Bob 3. The game is won whoever reaches 6 points first!
import scipy.stats as stats
import numpy as np
prob_bob_wining = np.array([])
for j in range(100):
count_alice_win = 0
count_bob_win = 0
count_cond_hit = 0
for i in range(10000):
# Select random p
p = np.random.uniform(0,1)
alice = 0
bob = 0
flag_condition = False
# Loop till anyone wins. The game ends when anyone reaches 6 points
while (alice < 6 and bob < 6):
#The outcome of a single trial game
outcome = np.random.binomial(1, p)
#Alice wins if outcome is 1
alice += outcome
#Bob winsif outcome is 0
bob += (1 - outcome)
#The condition which is provided in the problem
if (alice == 5 and bob == 3):
flag_condition = True
count_cond_hit += 1
#Someone has won the game
#checking if the condition was hit
#in such condition how many times did bob win
if flag_condition:
if 6 == alice:
count_alice_win += 1
else:
count_bob_win += 1
#The divison of condition_hit to the time bob win is the
#condition probability of Bob wining given P(Bob| A=5, B=3)
prob_bob_wining = np.append(prob_bob_wining, count_bob_win/ count_cond_hit)
print ("P : ",p ," condition hit : ", count_cond_hit)
print ("Bob win after condition : ", count_bob_win)
print ("Alice win after condition : ", count_alice_win)
print ("Bob chance of wining :: ", count_bob_win/count_cond_hit)
%matplotlib inline
import matplotlib.pyplot as plt
#Plotting the conditional probability of p(Bob| A=5, B=3)
plt.hist(prob_bob_wining)
{See the document for the proof}
This is a special case of Bayesian analysis in which we can can have a analytic solution. (Not true for all other cases of bayesian analysis)
$$ P( Bob | A=5, B=3 ) = \frac {B(5+a,3+b+3) }{B(5+a,3+b)} $$import math
a_1 = 1 #prior
a_2 = 1 #prior for uniform beta distribution. Non informative prior
def beta_1(a,b):
'''uses gamma function or inbuilt math.gamma() to compute values of beta function'''
beta = math.gamma(a)*math.gamma(b)/math.gamma(a+b)
return beta
print (beta_1(5 + a_1, 3+a_2+3) / beta_1(5 + a_1, 3 + a_2 ))
http://www.behind-the-enemy-lines.com/2008/01/are-you-bayesian-or-frequentist-or.html
Need to find
count_2_heads = 0
count_condition_hit = 0
percentage_of_2_heads = np.array([])
#for p in np.arange(0.1, 1, 0.1):
for j in range(100000):
#for i in range(10000):
#randomly find a p value
p = np.random.uniform(0,1)
#Toss the coint 14 times
outcome = np.random.binomial(1,p, 14)
#Check if the coin toss had 10 heads
if (10 == outcome.sum() ):
count_condition_hit += 1
#Find the probability of 2 heads
heads_2 = np.random.binomial(1,p, 2)
# Check if there are 2 heads
if (2 == heads_2.sum()):
count_2_heads += 1
#print ("For P : ", p ,"Condition Hit : ", count_condition_hit, "Count 2 heads : ", count_2_heads)
#if (count_condition_hit):
# print ("Percentage of 2 Heads : ", (count_2_heads/ count_condition_hit))
#else:
# print ("Precentage of 2 Heads : ", 0)
if (count_condition_hit):
percentage_of_2_heads = np.append(percentage_of_2_heads, (count_2_heads/ count_condition_hit))
print ("Percentage of 2 Heads : ", (count_2_heads/ count_condition_hit))
$$ P ( HH | data ) = \frac {B(10+a+2,4+b) } {B(10+a,4+b)} $$
import math
def beta_1(a,b):
'''uses gamma function or inbuilt math.gamma() to compute values of beta function'''
beta = math.gamma(a)*math.gamma(b)/math.gamma(a+b)
return beta
print (beta_1(10 + 1 + 2, 4+1) / beta_1(10 + 1 , 4 + 1 ))
All content copyright 2014-2016 Deebul Nair unless otherwise noted. Licensed under Creative Commons.