pythonprobabilitydice

Probabilities and rpg dice goes wrong for high values


I made this short code to calculate the chances of a success rolling dice, and it worked very well... but not in big numbers. Se the code, I'll explain better below.

def calc_dados(f_sucessos = 1, faces = 6, n_dados = 1):

  p_max = ((f_sucessos/faces)**n_dados) #chance de todos

  fator = 1
  p_meio = 0

  for i in range(n_dados-1):
    p_meio += (((f_sucessos/faces)**(n_dados-fator) * ((faces-f_sucessos)/faces)**(n_dados-(n_dados-fator))) * n_dados)
    fator += 1

  p = p_max + p_meio
  return p*100

So, ok, it works, why not see how my chances are better in function of adding dice? More the dice, better the chance. So I made this tiny table with pandas:

f_sucessos = 1    # how many faces are success
faces = 2         # faces of the dice
n_dados = 10      # n de dados lançados

suc_list = []
for i in range(0,n_dados): suc_list.append(f_sucessos)
fac_list = []
for i in range(0,n_dados): fac_list.append(faces)
cha_list = []
for i in range(0,n_dados): cha_list.append(calc_dados(f_sucessos, faces, i+1))

df = pd.DataFrame(
{
"n_dados" : range(1,n_dados+1),
"faces" : fac_list,
"sucessos" : suc_list,
"chance" : cha_list
 }
)
df

The results were very strange... So I wrote an coin probability table and tested as the coin was an 2 faced dice. The right table is this: table of right brute force tested results

But if you use my code to create this table the result will be this: table of the results of my code

Please, anybody can help me to understood why in a certain moment the probabilities just fall when they should be higher? For example:The chance of at least 1 'head' in 4 coins should be 93,75%, but my code says it is 81,25%...


Solution

  • Hello people, it's finally solved!

    First I would like to thank Raibek, who solved it using combinations, I didn't realize it was solved when he did it and below I'll tell you how and why.

    If you are not following the history of this code, you just need to know that it is used to calculate the probability of getting at least ns successes when rolling d amount of dice. Solution codes are at the end of this answer.

    I found out how to solve the problem by talking to a friend, Eber, who pointed me to an alternative to check the data, anydice.com. I quickly realized that my visual check, assembling tables in Excel/Calc was wrong, but why?

    Well, here comes my friend who, reading the table of large numbers with 7d6, where the error was already very evident, shows me that although at the beginning the account worked, my table did not have all the possible combinations. And the more possibilities there were, the more my accounts failed, with the odds getting smaller as more dice were added to the roll.

    This is the combinations I was considering, in this example on 7d6 case.

    table of 7d6 combinations

    In the first code the account was:

    successes**factor *failures**factor *d
    

    The mistake is in assuming that the number of possible combinations was equal to d (which is a coincidence up to 3 dice for the tests I did before thanks to factorials of 1 = 1 and factorial of 2 = 2).

    Now notice that, in 7d6 example, in the exact 3 block there are some missing possible combinations in yellow:

    table of 7d6 combinations where we can find exact 3 of 1 face

    The correct account for this term of the equation is:

    factorial(d) / factorial (failures) * factorial (successes)
    

    With this account we can find out what the chance of exactly n faces rolling is, and then if we want, for example, to know the chance of at least once getting the number 1 in 3d6, we just need to add the chances of getting exactly 1 time, 2 times and 3 times. What the code already did well.

    Finally, let's get to the code:

    Daniel-Eber solution:

    def dado(fs=1,ft=6,d=1,ns=1,exato=False):
        '''
        fs = faces sucesso
        ft = faces totais
        d  = n de dados
        ns - n de sucessos esperados modificados por exato
        exato = True: chance de exatamente ns ocorrerem, False: chance de pelo menos ns ocorrerem
        '''
        from math import factorial
        
        s =  fs/ft
        f =  (ft-fs)/ft
        d =  d
        ns =  ns
        
        p_max =  s**d
        
        falhas = 1
        po_soma = 0
        
        if exato == False:
            for i in range(d-1):
                po =  ( (s**(d-falhas)) * (f**(falhas))) * (factorial(d)/(factorial(falhas)*factorial((d-falhas))))
                po_soma += po
                falhas += 1
        else:
            p_max = 0
            falhas = d-ns
            po_soma =  ( (s**(d-falhas)) * (f**(falhas))) * (factorial(d)/(factorial(falhas)*factorial((d-falhas))))
    
        return f'{(p_max + po_soma)*100:.2f}%'
    
    print(dado(1,6,6,1))
    

    Raibek solution:

    from scipy.special import comb
    
    
    def calc_dados(f_sucessos=1, faces=6, n_dados=1):
    
        assert f_sucessos <= faces
    
        outcomes_total = faces ** n_dados
    
        outcomes_success = 0
    
        f_fail = faces - f_sucessos
    
        for i in range(1, n_dados + 1):
    
            one_permutation = (f_sucessos ** i) * (f_fail ** (n_dados - i))
            n_permutations = comb(n_dados, i)
    
            outcomes_success += one_permutation * n_permutations
    
        p = outcomes_success / outcomes_total
    
        return f'{(p)*100:.2f}%'