CODEWITHSUNDEEP
python
hardlyphil
hardlyphil

Reputation: 21

Simulation for conditional probabilty problem in python

I am trying to simulate a simple conditional probability problem. You hae two boxes. If you open A you have a 50% change of wining the prize, If you open B you have a 75% chance of winning. With some simple (bad) python I have tired But the appending doesn't work. Any thoughts on a neater way of doing this?

import random
import numpy as np
def liveORdie(prob):
    #Takes an argument of the probability of survival
    live = 0
    for i in range(100):
        if random.random() <= prob*1.0:
            live =1
    return live

def simulate(n):
    trials = np.array([0])
    for i in range(n):
        if random.random() <= 0.5:
            np.append(trials,liveORdie(0.5))
            print(trials)
        else:
            np.append(trials,liveORdie(0.75))
    return(sum(trials)/n)


simulate(10)

Upvotes: 2

Views: 480

Answers (4)

Sam Mason
Sam Mason

Reputation: 16184

I'm pretty sure that just reduces to:

from numpy import mean
from numpy.random import choice
from scipy.stats import bernoulli

def simulate(n):
  probs = choice([0.5, 0.75], n)
  return 1-mean(bernoulli.rvs((1-probs)**100))

which as others have pointed out will basically always return 1 — 0.5**100 is ~1e-30.

Upvotes: 0

digitalarbeiter
digitalarbeiter

Reputation: 2335

The loop in liveORdie() (please consider PEP8 for naming conventions) will cause the probability of winning to increase: each pass of the loop has a prop chance of winning, and you give it 100 tries, so with 50% resp. 75% you are extremely likely to win.

Unless I really misunderstood the problem, you probably just want

def live_or_die(prob):
    return random.random() < prob

Upvotes: 0

Richard
Richard

Reputation: 61289

You could make the code tighter by using list comprehensions and numpy's array operations, like so:

import random
import numpy as np

def LiveOrDie():
    prob = 0.5 if random.random()<=0.5 else 0.75
    return np.sum(np.random.random(100)<=prob)

def simulate(n):
    trials = [LiveOrDie() for x in range(n)]
    return(sum(trials)/n)

Simulate(10)

Upvotes: 1

Prune
Prune

Reputation: 77857

append is a list operation; you're forcing it onto a numpy array, which is not the same thing. Stick with a regular list, since you're not using any extensions specific to an array.

def simulate(n):
    trials = []
    for i in range(n):
        if random.random() <= 0.5:
            trials.append(liveORdie(0.5))

Now look at your liveORdie routine. I don't think this is what you want: you loop 100 times to produce a single integer ... and if any one of your trials comes up successful, you return a 1. Since you haven't provided documentation for your algorithms, I'm not sure what you want, but I suspect that it's a list of 100 trials, rather than the conjunction of all 100. You need to append here, as well.

Better yet, run through a tutorial on list comprehension, and use those.

Upvotes: 0

Related Questions