How do I simulate biased die in python?

Question

I want to simulate N-sided biased die?

def roll(N,bias):
     '''this function rolls N dimensional die with biasing provided'''
     # do something
     return result

>> N=6
>> bias=( 0.20,0.20,0.15,0.15,0.14,0.16,)
>> roll(N,bias)
   2

Answer 1

A little bit of math here.

A regular die will give each number 1-6 with equal probability, namely 1/6. This is referred to as uniform distribution (the discrete version of it, as opposed to the continuous version). Meaning that if X is a random variable describing the result of a single role then X~U[1,6] - meaning X is distributed equally against all possible results of the die roll, 1 through 6.

This is equal to choosing a number in [0,1) while dividing it into 6 sections: [0,1/6), [1/6,2/6), [2/6,3/6), [3/6,4/6), [4/6,5/6), [5/6,1).

You are requesting a different distribution, which is biased. The easiest way to achieve this is to divide the section [0,1) to 6 parts depending on the bias you want. So in your case you would want to divide it into the following: [0,0.2), [0.2,0.4), [0.4,0.55), 0.55,0.7), [0.7,0.84), [0.84,1).

If you take a look at the wikipedia entry, you will see that in this case, the cumulative probability function will not be composed of 6 equal-length parts but rather of 6 parts which differ in length according to the bias you gave them. Same goes for the mass distribution.

Back to the question, depending on the language you are using, translate this back to your die roll. In Python, here is a very sketchy, albeit working, example:

import random
sampleMassDist = (0.2, 0.1, 0.15, 0.15, 0.25, 0.15)
# assume sum of bias is 1
def roll(massDist):
    randRoll = random.random() # in [0,1]
    sum = 0
    result = 1
    for mass in massDist:
        sum += mass
        if randRoll < sum:
            return result
        result+=1

print(roll(sampleMassDist))

Answer 2

More language agnostic, but you could use a lookup table.

Use a random number in the range 0-1 and lookup the value in a table:

0.00 - 0.20   1
0.20 - 0.40   2
0.40 - 0.55   3
0.55 - 0.70   4
0.70 - 0.84   5
0.84 - 1.00   6

Answer 3

import random

def roll(sides, bias_list):
    assert len(bias_list) == sides
    number = random.uniform(0, sum(bias_list))
    current = 0
    for i, bias in enumerate(bias_list):
        current += bias
        if number <= current:
            return i + 1

The bias will be proportional.

>>> print roll(6, (0.20, 0.20, 0.15, 0.15, 0.14, 0.16))
6
>>> print roll(6, (0.20, 0.20, 0.15, 0.15, 0.14, 0.16))
2

Could use integers too (better):

>>> print roll(6, (10, 1, 1, 1, 1, 1))
5
>>> print roll(6, (10, 1, 1, 1, 1, 1))
1
>>> print roll(6, (10, 1, 1, 1, 1, 1))
1
>>> print roll(6, (10, 5, 5, 10, 4, 8))
2
>>> print roll(6, (1,) * 6)
4

Answer 4

It is a little surprising that the np.random.choice answer hasn't been given here.

from numpy import random 
def roll(N,bias):
    '''this function rolls N dimensional die with biasing provided'''
    return random.choice(np.range(N),p=bias)

The p option gives "the probabilities associated with each entry in a", where a is np.range(N) for us. "If not given the sample assumes a uniform distribution over all entries in a".

Answer 5

See the recipe for Walker's alias method for random objects with different probablities.
An example, strings A B C or D with probabilities .1 .2 .3 .4 --

abcd = dict( A=1, D=4, C=3, B=2 )
  # keys can be any immutables: 2d points, colors, atoms ...
wrand = Walkerrandom( abcd.values(), abcd.keys() )
wrand.random()  # each call -> "A" "B" "C" or "D"
                # fast: 1 randint(), 1 uniform(), table lookup

cheers
-- denis

Answer 6

Just to suggest a more efficient (and pythonic3) solution, one can use bisect to search in the vector of accumulated values — that can moreover be precomputed and stored in the hope that subsequent calls to the function will refer to the same "bias" (to follow the question parlance).

from bisect import bisect
from itertools import accumulate
from random import uniform

def pick( amplitudes ):
    if pick.amplitudes != amplitudes:
        pick.dist = list( accumulate( amplitudes ) )
        pick.amplitudes = amplitudes
    return bisect( pick.dist, uniform( 0, pick.dist[ -1 ] ) )
pick.amplitudes = None

In absence of Python 3 accumulate, one can just write a simple loop to compute the cumulative sum.

Answer 7

from random import random
biases = [0.0,0.3,0.5,0.99]
coins = [1 if random()<bias else 0 for bias in biases]

Answer 8

i have created a code for a dictionary giving a event and corresponding probability, it gives back the corresponding key ie the event of that probability.

import random


def WeightedDie(Probabilities):   


    high_p = 0   
    rand = random.uniform(0,1)

    for j,i in Probabilities.items():
        high_p = high_p + i
        if rand< high_p:
            return j

Answer 9

We could use numpy's multinomial distribution too

import numpy as np

bias = [0.10,0.10,0.15,0.15,0.14,0.16,0.05,0.06,0.04,0.05]   # a 10-sided biased die
np.where(np.random.multinomial(1, bias, size=1)[0]==1)[0][0]+1 # just 1 roll
# 4

If you want to roll the biased die (with the given bias probabilities) for n times, use the following function

def roll(probs, ntimes):   # roll a len(probs) sided biased die with bias probs for ntimes
   return np.apply_along_axis(lambda x: x.tolist().index(1)+1, 1, 
                       np.random.multinomial(1, bias,  size=10)).tolist()

roll(probs=bias, ntimes=10)    # 10 rolls
# [5, 6, 8, 4, 8, 3, 1, 5, 8, 6]

Answer 10

For python 3.6 and above, you can make use of random's choices() method already part of stdlib. To simulate the die in your example, the equivalent code would be:-

import random

def roll(N, bias_list):
   return random.choices(list(range(N)), weights=bias_list, k=1)[-1]