Jamol
Jamol

Reputation: 3908

Python probability

We have a six sided die, with sides numbered 1 through 6. The probability of first seeing a 1 on the n-th roll decreases as n increases. I want to find the smallest number of rolls such that this probability is less than some given limit.

def probTest(limit):
    prob = 1.0
    n = 1
    while prob > limit:
        prob = (1/6)**n
        n += 1        
    return n-1

What is wrong with my code?

Upvotes: 3

Views: 4903

Answers (6)

Sandipan Dey
Sandipan Dey

Reputation: 23099

The probability distribution of the number X of Bernoulli trials needed to get one success, follows a geometric distribution (https://en.wikipedia.org/wiki/Geometric_distribution), we can compute the corresponding probability X=n given p directly with the following formula:

p = 1./6 # prob of successes
geom_prob = []
for n in range(1,25):
    geom_prob.append((1-p)**(n-1)*p) 
print geom_prob 
# [0.16666666666666666, 0.1388888888888889, 0.11574074074074076, 0.09645061728395063, 0.08037551440329219, 0.06697959533607684, 0.05581632944673069, 0.04651360787227558, 0.038761339893562986, 0.032301116577969156, 0.02691759714830763, 0.022431330956923026, 0.018692775797435855, 0.015577313164529882, 0.012981094303774901, 0.010817578586479085, 0.009014648822065905, 0.0075122073517215875, 0.006260172793101323, 0.005216810660917769, 0.0043473422174314744, 0.003622785181192896, 0.0030189876509940797, 0.002515823042495067]
import matplotlib.pyplot as plt
plt.plot(geom_prob)
plt.xlabel('n')
plt.ylabel('probability')
plt.title('Probability of getting first head on nth roll of a die')
plt.show()

enter image description here

We can use simulation to find the probability as well, like the following:

import numpy as np
sim_geom = np.random.geometric(p=p, size=1000)
import seaborn as sns   
sns.distplot(sim_geom)
plt.show()

enter image description here

Upvotes: 0

james
james

Reputation: 11

Thanks. Please take in consideration the other factors.

probTest(25/216.0) return 4 rather then the correct result n=3

Upvotes: 1

jammy
jammy

Reputation: 21

the correct will be: prob = (5.0/6)**(n-1)*1/6.0

Upvotes: 2

chunky
chunky

Reputation: 1

def probTest(limit):
    n=1
    prob = 1.0
    while prob > limit:
        prob = prob * (1/6.0)*((5/6.0)**n-1)
        n +=1
    return n-1

Upvotes: -1

user1887513
user1887513

Reputation: 7

def probTest(limit):
    prob = 1.0
    n = 1
    while prob > limit:
        prob =  5/6^(n-1)*1/6.0
        n += 1        
    return n

Upvotes: -1

Jeff
Jeff

Reputation: 12418

The probably of rolling a one on the nth roll is 5/6^(n-1)*1/6, not 1/6^n.
1/6^n is the probability of rolling one on all n rolls.

The first n-1 rolls each have a 5/6 chance of not being one.
The nth roll has a 1/6th chance of being one.

Upvotes: 7

Related Questions