How do I simulate flip of biased coin in python?

问题:

In unbiased coin flip H or T occurs 50% of times.

But I want to simulate coin which gives H with probability 'p' and T with probability '(1-p)'.

something like this:

def flip(p):
   '''this function return H with probability p'''
   # do something
   return result

>> [flip(0.8) for i in xrange(10)]
[H,H,T,H,H,H,T,H,H,H]

回答1:

random.random() returns a uniformly distributed pseudo-random floating point number in the range [0, 1). This number is less than a given number p in the range [0,1) with probability p. Thus:

def flip(p):
    return 'H' if random.random() < p else 'T'

Some experiments:

>>> N = 100
>>> flips = [flip(0.2) for i in xrange(N)]
>>> float(flips.count('H'))/N
0.17999999999999999  # Approximately 20% of the coins are heads

>>> N = 10000
>>> flips = [flip(0.2) for i in xrange(N)]
>>> float(flips.count('H'))/N
0.20549999999999999  # Better approximation 

回答2:

Do you want the "bias" to be based in symmetric distribuition? Or maybe exponential distribution? Gaussian anyone?

Well, here are all the methods, extracted from random documentation itself.

First, an example of triangular distribution:

print random.triangular(0, 1, 0.7)

random.triangular(low, high, mode):

Return a random floating point number N such that low <= N < high and with the specified mode between those bounds. The low and high bounds default to zero and one. The mode argument defaults to the midpoint between the bounds, giving a symmetric distribution.

random.betavariate(alpha, beta):

Beta distribution. Conditions on the parameters are alpha > 0 and beta > 0. Returned values range between 0 and 1.

random.expovariate(lambd):

Exponential distribution. lambd is 1.0 divided by the desired mean. It should be nonzero. (The parameter would be called “lambda”, but that is a reserved word in Python.) Returned values range from 0 to positive infinity if lambd is positive, and from negative infinity to 0 if lambd is negative.

random.gammavariate(alpha, beta):

Gamma distribution. (Not the gamma function!) Conditions on the parameters are alpha > 0 and beta > 0.

random.gauss(mu, sigma):

Gaussian distribution. mu is the mean, and sigma is the standard deviation. This is slightly faster than the normalvariate() function defined below.

random.lognormvariate(mu, sigma):

Log normal distribution. If you take the natural logarithm of this distribution, you’ll get a normal distribution with mean mu and standard deviation sigma. mu can have any value, and sigma must be greater than zero.

random.normalvariate(mu, sigma):

Normal distribution. mu is the mean, and sigma is the standard deviation.

random.vonmisesvariate(mu, kappa):

mu is the mean angle, expressed in radians between 0 and 2*pi, and kappa is the concentration parameter, which must be greater than or equal to zero. If kappa is equal to zero, this distribution reduces to a uniform random angle over the range 0 to 2*pi.

random.paretovariate(alpha):

Pareto distribution. alpha is the shape parameter.

random.weibullvariate(alpha, beta)

Weibull distribution. alpha is the scale parameter and beta is the shape parameter.

回答3:

import random
def flip(p):
    return (random.random() < p)

That returns a boolean which you can then use to choose H or T (or choose between any two values) you want. You could also include the choice in the method:

def flip(p):
    if random.random() < p:
        return 'H'
    else:
        return 'T'

but it'd be less generally useful that way.

回答4:

How about:

import numpy as np
n, p = 1, .33  # n = coins flipped, p = prob of success
s = np.random.binomial(n, p, 100)

回答5:

• Import a random number between 0 - 1 (you can use randrange function)

• If the number is above (1-p), return tails.

• Else, return heads

回答6:

One can sample from the X ~ Bernoulli(p) distribution nsamples times using sympy too:

from sympy.stats import Bernoulli, sample_iter
list(sample_iter(Bernoulli('X', 0.8), numsamples=10)) # p = 0.8 and nsamples=10
# [1, 1, 0, 1, 1, 0, 1, 1, 1, 1]

Return 'H' or 'T' instead using

def flip(p, n):
    return list(map(lambda x: 'H' if x==1 else 'T', sample_iter(Bernoulli('X', p), numsamples=n)))

print(flip(0.8, 10)) # p = 0.8 and nsamples=10
# ['H', 'H', 'T', 'H', 'H', 'T', 'H', 'H', 'H', 'H']

回答7:

import random
def flip():
    return ["H" if random.randint(0,3) <= 2 else "T" for i in range(10)]

Right now probability of Head is 75% and tails is 25% (0,1,2 are all Heads and only 3 is Tails) . By using random.randint() you could have any probability of bias while still maintaining randomness.