Need a Programmer for Number-crunching

[Abisashi]

Alright, I fixed the code, and now my percentage of superheroes differs from Alex319 by only .01% of the total population. I updated my data above to be correct.

 

[genshou]

Minor note, there's no reason to store all those people. Instead, generate a person, run the math, increment the relevant counters, then generate the next person. Run the loop till you get bored.

This is why you're the programmer, and I'm a truck driver.

Thanks so much, guys! I should be able to have some fun figuring out facts about my D&D settings now. I'm trying to get a good feel for what "default D&D" really is like according to these sorts of statistics, compared to what people (myself included) think it's like based on gameplay experiences and the PCs' perspectives on the world at large.

 

[babomb]

Well, I went ahead and crunched the numbers for all possible combinations. I'm reporting these as probabilities, rather than # of people in a population.

*********Full population **********
*********Highest Stat Probabilities **********
3:	9.8464e-15
4:	4.0321e-11
5:	9.80607e-09
6:	6.20323e-07
7:	1.74701e-05
8:	0.000285572
9:	0.00247724
10:	0.0128441
11:	0.0439796
12:	0.105591
13:	0.181022
14:	0.212014
15:	0.194225
16:	0.141452
17:	0.0786341
18:	0.0274583
*********Point Buy Probabilities **********
-66:	9.8464e-15
-63:	1.77235e-13
-60:	1.68373e-12
-58:	5.90784e-13
-57:	1.06341e-11
-56:	8.86176e-13
-55:	1.01024e-11
-54:	5.06597e-11
-53:	1.48878e-11
-52:	9.10989e-11
-51:	1.97499e-10
-50:	1.46278e-10
-49:	5.45884e-10
-48:	7.24094e-10
-47:	1.00817e-09
-46:	2.50241e-09
-45:	2.85579e-09
-44:	5.42827e-09
-43:	9.67557e-09
-42:	1.23972e-08
-41:	2.38284e-08
-40:	3.5134e-08
-39:	5.2566e-08
-38:	9.03586e-08
-37:	1.26366e-07
-36:	2.03067e-07
-35:	3.11105e-07
-34:	4.493e-07
-33:	7.04594e-07
-32:	1.02198e-06
-31:	1.51338e-06
-30:	2.25833e-06
-29:	3.22787e-06
-28:	4.75332e-06
-27:	6.81073e-06
-26:	9.7205e-06
-25:	1.39073e-05
-24:	1.95538e-05
-23:	2.75096e-05
-22:	3.84192e-05
-21:	5.31432e-05
-20:	7.33872e-05
-19:	0.000100377
-18:	0.000136571
-17:	0.000184676
-16:	0.000248016
-15:	0.000330771
-14:	0.000438796
-13:	0.000577172
-12:	0.000755343
-11:	0.000981007
-10:	0.00126533
-9:	0.00162207
-8:	0.00206266
-7:	0.00260695
-6:	0.00326818
-5:	0.0040687
-4:	0.00502449
-3:	0.00615798
-2:	0.00748389
-1:	0.00902245
0:	0.0107839
1:	0.0127779
2:	0.0150095
3:	0.0174652
4:	0.0201393
5:	0.0229892
6:	0.0259826
7:	0.0290478
8:	0.032115
9:	0.0350865
10:	0.0378621
11:	0.0403307
12:	0.0423847
13:	0.0439331
14:	0.0448897
15:	0.0452203
16:	0.0448897
17:	0.0439331
18:	0.0423847
19:	0.0403307
20:	0.0378621
21:	0.0350865
22:	0.032115
23:	0.0290478
24:	0.0259826
25:	0.0229892
26:	0.0201393
27:	0.0174652
28:	0.0150095
29:	0.0127779
30:	0.0107839
31:	0.00902245
32:	0.00748389
33:	0.00615798
34:	0.00502449
35:	0.0040687
36:	0.00326818
37:	0.00260695
38:	0.00206266
39:	0.00162207
40:	0.00126533
41:	0.000981007
42:	0.000755343
43:	0.000577172
44:	0.000438796
45:	0.000330771
46:	0.000248016
47:	0.000184676
48:	0.000136571
49:	0.000100377
50:	7.33872e-05
51:	5.31432e-05
52:	3.84192e-05
53:	2.75096e-05
54:	1.95538e-05
55:	1.39073e-05
56:	9.7205e-06
57:	6.81073e-06
58:	4.75332e-06
59:	3.22787e-06
60:	2.25833e-06
61:	1.51338e-06
62:	1.02198e-06
63:	7.04594e-07
64:	4.493e-07
65:	3.11105e-07
66:	2.03067e-07
67:	1.26366e-07
68:	9.03586e-08
69:	5.2566e-08
70:	3.5134e-08
71:	2.38284e-08
72:	1.23972e-08
73:	9.67557e-09
74:	5.42827e-09
75:	2.85579e-09
76:	2.50241e-09
77:	1.00817e-09
78:	7.24094e-10
79:	5.45884e-10
80:	1.46278e-10
81:	1.97499e-10
82:	9.10989e-11
83:	1.48878e-11
84:	5.06597e-11
85:	1.01024e-11
86:	8.86176e-13
87:	1.06341e-11
88:	5.90784e-13
90:	1.68373e-12
93:	1.77235e-13
96:	9.8464e-15

***************heroes*************
Probability of getting a hero: 0.746229
*********Highest Stat Probabilities **********
12:	0.0654969
13:	0.154835
14:	0.237749
15:	0.227849
16:	0.178433
17:	0.0997931
18:	0.0358443
*********Point Buy Probabilities **********
-2:	2.67998e-07
-1:	4.21899e-06
0:	1.88056e-05
1:	4.27796e-05
2:	9.86008e-05
3:	0.000372585
4:	0.00115789
5:	0.00264444
6:	0.00484085
7:	0.00785077
8:	0.0134853
9:	0.0227038
10:	0.0336737
11:	0.0446837
12:	0.0525628
13:	0.0569206
14:	0.0592542
15:	0.0601603
16:	0.059987
17:	0.0588354
18:	0.0567932
19:	0.054046
20:	0.0507378
21:	0.0470183
22:	0.0430364
23:	0.0389261
24:	0.0348185
25:	0.0308071
26:	0.0269881
27:	0.0234046
28:	0.0201138
29:	0.0171233
30:	0.0144512
31:	0.0120907
32:	0.0100289
33:	0.00825213
34:	0.00673317
35:	0.00545235
36:	0.00437959
37:	0.0034935
38:	0.00276412
39:	0.00217368
40:	0.00169563
41:	0.00131462
42:	0.00101221
43:	0.000773452
44:	0.000588017
45:	0.000443256
46:	0.000332359
47:	0.000247479
48:	0.000183015
49:	0.000134513
50:	9.83441e-05
51:	7.12156e-05
52:	5.14845e-05
53:	3.68648e-05
54:	2.62035e-05
55:	1.86368e-05
56:	1.30262e-05
57:	9.12685e-06
58:	6.36979e-06
59:	4.32557e-06
60:	3.02632e-06
61:	2.02803e-06
62:	1.36952e-06
63:	9.44206e-07
64:	6.02093e-07
65:	4.16902e-07
66:	2.72124e-07
67:	1.6934e-07
68:	1.21087e-07
69:	7.04421e-08
70:	4.7082e-08
71:	3.19317e-08
72:	1.66132e-08
73:	1.29659e-08
74:	7.27427e-09
75:	3.82696e-09
76:	3.35341e-09
77:	1.35102e-09
78:	9.70338e-10
79:	7.31524e-10
80:	1.96023e-10
81:	2.64663e-10
82:	1.22079e-10
83:	1.99506e-11
84:	6.78876e-11
85:	1.35379e-11
86:	1.18754e-12
87:	1.42505e-11
88:	7.91692e-13
90:	2.25632e-12
93:	2.37508e-13
96:	1.31949e-14

*************superheroes**********
Probability of getting a superhero: 0.136984
Probability of a hero being a superhero: 0.183568
*********Highest Stat Probabilities **********
15:	0.225784
16:	0.334358
17:	0.301997
18:	0.137861
*********Point Buy Probabilities **********
25:	0.14258
26:	0.130526
27:	0.117279
28:	0.103604
29:	0.0900222
30:	0.0770776
31:	0.065106
32:	0.054321
33:	0.0448431
34:	0.0366472
35:	0.0296952
36:	0.0238574
37:	0.0190311
38:	0.0150578
39:	0.0118413
40:	0.0092371
41:	0.00716149
42:	0.00551411
43:	0.00421344
44:	0.00320327
45:	0.00241467
46:	0.00181055
47:	0.00134816
48:	0.000996988
49:	0.000732768
50:	0.000535738
51:	0.000387953
52:	0.000280466
53:	0.000200824
54:	0.000142746
55:	0.000101526
56:	7.09611e-05
57:	4.97193e-05
58:	3.47e-05
59:	2.35639e-05
60:	1.64861e-05
61:	1.10479e-05
62:	7.46059e-06
63:	5.14364e-06
64:	3.27995e-06
65:	2.27111e-06
66:	1.48242e-06
67:	9.22493e-07
68:	6.59631e-07
69:	3.83739e-07
70:	2.56483e-07
71:	1.73951e-07
72:	9.05017e-08
73:	7.0633e-08
74:	3.96272e-08
75:	2.08477e-08
76:	1.8268e-08
77:	7.35981e-09
78:	5.286e-09
79:	3.98504e-09
80:	1.06785e-09
81:	1.44177e-09
82:	6.65035e-10
83:	1.08683e-10
84:	3.69823e-10
85:	7.3749e-11
86:	6.46921e-12
87:	7.76306e-11
88:	4.31281e-12
90:	1.22915e-11
93:	1.29384e-12
96:	7.18802e-14

 

[genshou]

Wow, that must have taken a while to calculate! Thanks!

I guess computers have come a long way; I calculated the probability of each score with 4d6dl back in high school, and it took ten minutes just to figure out 1,296 possible rolls.

 

[genshou]

babomb, I just noticed in your results that "superhero" point-buy probabilities do not contain point-buys of 25. Did you run it for results greater than 25--instead of greater than or equal to 25--to generate the other "superhero" statistics (such as percentage of "superheroes" among the common rabble) as well? If it's not too much trouble, could you include the 25 point-buys for factoring the "superhero" statistics?

I only noticed because I looked really hard to figure out why your "superhero" percentage was a few points lower than the people who did it by randomly generating a few million results.

 

[babomb]

Ha, yeah I used > 25. I updated the stats in the above post to include >= 25.

As for how long it takes to compute, ~20s for the whole thing, and I could probably make it run a lot faster if I had a mind, using knowledge of permutations and such. But the straightforward way is easier to implement, and fast enough for me.

The key that makes this a lot faster than computing all 6^18 possibilities is that I don't care if you got 10 by rolling 4, 4, 2 or 8, 1, 1. So I can pre-compute the probability of rolling a 10 on 3d6 (27/216).

Then I can say, "What's the probability of an all 10s character?" Well that's (27/216)^6. What about 3,3,3,3,3,4? That's (1/216)^5 * (3/216), etc.
In this way, I only have to count 16^6 or 16777216 possibilities. That's still a lot, but it's manageable for a computer. I could substantially reduce this by taking into account the fact that, e.g., 3,3,3,3,3,4 and 3,3,3,3,4,3 are the same thing for our purposes. I don't currently. For reference, my C++ code is below.

#include <algorithm>
#include <iostream>
#include <cstdlib>
#include <cmath>
#include <map>

using namespace std;

const double g_prob[] = { 1, 3, 6, 10, 15, 21, 25, 27 };
const int   g_pointBuyValue[16] = { -11, -8, -5, -3, -1, 0, 1, 2, 3, 4, 5, 6, 8, 10, 13, 16 };

double prob3d6( int n )
{
    if( n < 3 || n > 18 )
    {
        return 0;
    }

    if( n > 10 )
    {
        n = 21 - n;
    }

    return g_prob[n - 3] / 216;
}

int getModifier( int score )
{
    return score / 2 - 5;
}

int main()
{
    map<int, double> pointBuyProbs;
    map<int, double> heroPointBuyProbs;
    map<int, double> superheroPointBuyProbs;

    map<int, double> maxStatProbs;
    map<int, double> heroMaxStatProbs;
    map<int, double> superheroMaxStatProbs;

    double heroProb = 0;
    double superheroProb = 0;

    for( int i = 3; i <= 18; ++i )
    {
        int iCost = g_pointBuyValue[ i -3 ];
        int iModifier = getModifier( i );
        double iProb = prob3d6( i );

        for( int j = 3; j <= 18; ++j )
        {
            int jCost = g_pointBuyValue[ j -3 ];
            int jModifier = getModifier( j );
            double jProb = prob3d6( j );

            for( int k = 3; k <= 18; ++k )
            {
                int kCost = g_pointBuyValue[ k -3 ];
                int kModifier = getModifier( k );
                double kProb = prob3d6( k );

                for( int l = 3; l <= 18; ++l )
                {
                    int lCost = g_pointBuyValue[ l -3 ];
                    int lModifier = getModifier( l );
                    double lProb = prob3d6( l );

                    for( int m = 3; m <= 18; ++m )
                    {
                        int mCost = g_pointBuyValue[ m -3 ];
                        int mModifier = getModifier( m );
                        double mProb = prob3d6( m );

                        for( int n = 3; n <= 18; ++n )
                        {
                            int nCost = g_pointBuyValue[ n - 3 ];
                            int nModifier = getModifier( n );
                            double nProb = prob3d6( n );
                            int cost = iCost + jCost + kCost + lCost + mCost + nCost;
                            int totalModifier = iModifier + jModifier + kModifier + lModifier + mModifier + nModifier;
                            double prob = iProb * jProb * kProb * lProb * mProb * nProb;
                            int maxStat = max( max( max(i, j), max(k, l)), max(m, n));

                            maxStatProbs[maxStat] += prob;
                            pointBuyProbs[cost] += prob;
                            if( totalModifier > -3 &&
                                maxStat >= 12 )
                            {
                                heroProb += prob;
                                heroMaxStatProbs[maxStat] += prob;
                                heroPointBuyProbs[cost] += prob;
                                if(cost >= 25 &&
                                    maxStat >= 15 )
                                {
                                    superheroProb += prob;
                                    superheroMaxStatProbs[maxStat] += prob;
                                    superheroPointBuyProbs[cost] += prob;
                                }
                            }
                        }
                    }
                }
            }
        }
    }

    // Done computing. Output stuff.
    cout << "*********Full population **********n";
    cout << "*********Highest Stat Probabilities **********n";
    for( map<int, double>::iterator i = maxStatProbs.begin(); i!= maxStatProbs.end(); ++i )
    {
        cout << i->first << ":t" << i->second << "n";
    }

    cout << "*********Point Buy Probabilities **********" << "n";
    for( map<int, double>::iterator i = pointBuyProbs.begin(); i!= pointBuyProbs.end(); ++i )
    {
        cout << i->first << ":t" << i->second << "n";
    }
    cout << "***************heroes*************" << endl;
    cout << "Probability of getting a hero: " << heroProb << "n";
    cout << "*********Highest Stat Probabilities **********n";
    for( map<int, double>::iterator i = heroMaxStatProbs.begin(); i!= heroMaxStatProbs.end(); ++i )
    {
        cout << i->first << ":t" << i->second / heroProb << "n";
    }
    cout << "*********Point Buy Probabilities **********" << "n";
    for( map<int, double>::iterator i = heroPointBuyProbs.begin(); i!= heroPointBuyProbs.end(); ++i )
    {
        cout << i->first << ":t" << i->second / heroProb << "n";
    }
    cout << "*************superheroes**********" << endl;
    cout << "Probability of getting a superhero: " << superheroProb << "n";
    cout << "Probability of a hero being a superhero: " << superheroProb / heroProb << "n";
    cout << "*********Highest Stat Probabilities **********n";
    for( map<int, double>::iterator i = superheroMaxStatProbs.begin(); i!= superheroMaxStatProbs.end(); ++i )
    {
        cout << i->first << ":t" << i->second / superheroProb << "n";
    }
    cout << "*********Point Buy Probabilities **********" << "n";
    for( map<int, double>::iterator i = superheroPointBuyProbs.begin(); i!= superheroPointBuyProbs.end(); ++i )
    {
        cout << i->first << ":t" << i->second / superheroProb << "n";
    }

    return 0;
}

 

[resistor]

Some of these aren't particularly hard to solve by hand (as in, on a sheet of paper), just using counting formulae.

For instance, here's the probabilities for highest stat being X:

3 -> 6 * (1/216) * (1/216)^5 ~= 0.0000000000000590784
4 -> 6 * (3/216) * (4/126)^5 ~= 0.00000000181489
5 -> 6 * (6/216) * (10/216)^5 ~= 0.000000035447
6 -> 6 * (10/216) * (20/216)^5 ~= 0.00000189051
7 -> 6 * (15/216) * (35/216)^5 ~= 0.0000465436
8 -> 6 * (21/216) * (56/216)^5 ~= 0.000683263
9 -> 6 * (25/216) * (81/216)^5 ~= 0.00514984
10 -> 6 * (27/216) * (108/216)^5 ~= 0.0234375
11 -> 6 * (27/216) * (135/216)^5 ~= 0.0715256
12 -> 6 * (25/216) * (160/216)^5 ~= 0.15487
13 -> 6 * (21/216) * (181/216)^5 ~= 0.241013
14 -> 6 * (15/216) * (196/216)^5 ~= 0.256331
15 -> 6 * (10/216) * (206/216)^5 ~= 0.219162
16 -> 6 * (6/216) * (212/216)^5 ~= 0.151796
17 -> 6 * (3/216) * (215/216)^5 ~= 0.0814221
18-> 6 * (1/216) * (216/216)^5 ~= 0.0277778

 

[Elric]

18-> 6 * (1/216) * (216/216)^5 ~= 0.0277778

This isn't right. The chance that your highest stat is an 18 is 1- (the chance that your highest stat is 17 or below).

The chance that one stat is 17 or below is 215/216. So the chance that all 6 stats are 17 or below is (215/216)^6. This means that the answer is 1- (215/216)^6= 0.027458.

What you've got here is the answer for "what is the average number of 18s a character has." The problem is that some characters have multiple 18s, so the expected number of 18s is higher than the percentage of characters with 18s.

 

[resistor]

What you've got here is the answer for "what is the average number of 18s a character has." The problem is that some characters have multiple 18s, so the expected number of 18s is higher than the percentage of characters with 18s.

Hrm? No.

What I computed is the odds of having exactly one X and five X-or-below, multiplied by six because the X could be any of the six stats.

To put it another way, your odds of having a highest stat of 18 is the odds of rolling at least one 18 in six rolls of 3d6. Each roll is independent, with a chance on 1/216 on each. So the probability of rolling at least one 18 is 6/216.

 

[arscott]

For instance, here's the probabilities for highest stat being X:

3 -> 6 * (1/216) * (1/216)^5 ~= 0.0000000000000590784
<snip>

That's pretty obviously wrong. The probability of the highest stat being 3 is the same thing as the probability of all stats being 3s. And that's just (1/216)^6. Your value is 6 times as high.

I'm gonna work on the math now and come up with better values.