Statistics of dropping dice

[Sivarias]

Hello!

I'm working on a helpful program for me and my buddies that gives percentages for dice rolls. I can piece together simple stats, multiply fractions etc.
The issue is that I have no idea how to make a formula for dice rolls where the lowest dice is dropped. E.G. 4d6d1. How about 5d6d2? Or even worse. I roll 1d4+1d8+1d12 but I only keep the highest 2, what is the percentage breakdown of the dice rolls?
Right now I have a brute force solution where I simply roll the dice set 1000 times and use those results, but it's inelegant at best and inaccurate at worst.
Can you guys help out with the formula? Or link me to some place that can teach me what I need to know? I can't seem to find an online resource that is useful for this particular problem.

[Knaight]

You're repeating the work of Anydice to some extent, and the guy who made and released that is also contactable - I'd try reaching out.

Beyond that, brute forcing it is just fine. The 1000 roll method sounds like a Monte Carlo approximation, which I wouldn't use. Instead, just generate every result of absolute dice rolls and their probability, then convert that to a number. This can be all at the end or as you generate them, either way works.

Take 1d4+1d8+1d12, highest two. There are 4*8*12=384 different raw dice rolls, to be generated in three nested loops. Then those are converted to numbers, and you do data analysis on the outlet set. that 5d6 only generates 6^5=7776 individual die rolls, and that's without using combinatorics to shrink that ahead of time. For a computer, that's not really many calculations.

I might whip something up in Octave as example code here.

[Sivarias]

See the analysis of that number set is what I don't know how to do. What did you mean by that?

[JeenLeen]

I wrote some code to roll a d20 with advantage (from D&D 5e), which is basically 2d20b1. I think it could be expanded to what you want.

How I do it is
1) create a matrix with 1 column per die
2) in each column, put the total number of results of said die (e.g., 1-20), but have it so they line up... with each other... like 1-1, then 1-2, then 1-3... to 2-1, 2-2, 2-3... and so on. In other words, column 1 goes 1-20, 1-20, 1-20 ,and so on. Column 2 goes with 20 values of 1, then 20 values of 2, and so on. That way you get every combination.
NOTE: this gives a big dataset if you get a good handful of dice.
3) calculate the result for each row. For your case, you'd exclude the lowest out of each row.
This gives you a matrix, where the sum of each row is the true population of possible results. Create a frequency table of the results (PROC FREQ in SAS, or use some program to make a histogram but also output the counts and/or percentages.) If say, 8 out of 150 possible results were a 2, there's 8/150 % chance of getting a 2. (Those numbers were arbitrarily chosen.)

Here's a more technical explanation, based on the R programming language. I'm copying from a class project I did in grad school. In this description, I was trying to figure out the average roll, but you could easily create a histogram (or output the results of a histogram) to get the percentage from each result. That is, instead of taking the max or min (to emulate advantage or disadvantage), sum them minus the lowest.

It is more complex to calculate the true mean when rolling with advantage or disadvantage. One method, which is the method utilized by the plotrolladv() and plotrolldis() functions, is to first list every possible result of rolling two of the given die. This can be accomplished in R by creating a matrix containing every combination for a given die size x:
dmatrix=cbind(rep(1:x,x), rep(1:x, each=x))
Then, for advantage, take the maximum of each pair; for disadvantage, take the minimum. Then take the mean of the resultant vector. This can be accomplished in R by setting the variable adv to TRUE when we want to roll with advantage and taking the mean of the maximum of each row of dmatrix. For disadvantage, take the mean of the minimum of each row. For example:
if (adv) truemean=mean(pmax(dmatrix[,1],dmatrix[,2]))
else truemean=mean(pmin(dmatrix[,1],dmatrix[,2]))
This yields the true mean, as the variables truemean, of that given die with advantage or disadvantage.


If I have time, I'll try to brute force it in SAS and post something this week, but I got a couple deadlines so I might not remember.


Here's the R code. It first does a Monte Carlo simulation to give a histogram, but then calculates the true mean. (It would have made more sense just to create the actual distribution, but I think I didn't know how at the time and the professor liked Monte Carlo stuff.)

This doesn't do what you need, but if you know R, it might help.


Here's the R code to roll advantage.


#functions to be used by other functions--not intended for use by the user, although the user is welcome to use them
`%d%`<- function(x, y) rep(y, x) #creates a new operator %d%, so that 3%d%8 returns a vector of (8, 8, 8)
dicevector=function(dString)
{
#changes normal D&D nomenclature for writing dice rolls into a vector of dice sizes. No modifiers (no +/-x)
#example: turns "1d20" into 20. "1d20+3d6" to c(20,6,6,6). Also takes input using comma instead of plus, i.e., dicevector("1d20,3d6")
dString=gsub(pattern=",", x=dString, replacement="+")
`+` <- function(x, y) c(x,y) #overwrites + operator so that it concatenates vectors instead of adding numbers
dString=gsub(pattern="d", x=dString, replacement="%d%")

returnThis=eval(parse(text=dString))

rm("+") #This rm() call effectively restores the original meaning of +
return(returnThis)
}

roll = function(dice="1d20", mod=0, fakeDice=FALSE, noprint=FALSE)
{
#VALID INPUT: dice can be accepted as "1d8+3d6" or "1d8,3d6". Also accepts numeric values like 10 to mean 1d10 or 'vectors' written as characters, i.e, "c(8,6,6,6)". Does not accept matrices. Accepts vectors, but likely will not work as intended.
#defaults to 1d20, with no modifier
#fakeDice option suppresses warnings about using non-standard dice. Defaults to FALSE (no suppression)
#noprint option suppresses printing the actual dice results, before summing/adding modifiers. Defaults to FALSE (no suppression)

#DECLARE INTERNAL FUNCTIONS--these functions are used by the roll function
internal.roll = function(d)
{
returnThis = sample(d, 1) #sample is a RNG build into R. Used instead of round() on runif() since round() rounds to even
return(returnThis)
}
internal.checkDie=function(d)
{ #checks if a standard die size was used. (could also be programmed via the | (or) operator)
returnThis=""

logical.obj=(any(c(d==2, d==4, d==6, d==8, d==10, d==12, d==20, d==100)))
if (logical.obj==FALSE) {returnThis="Warning: Non-standard die was rolled."}
return(returnThis)
}
#END INTERNAL FUNCTIONS. BEGIN ACTUAL ROLL FUNCTION CODE

d=dicevector(dice)
warnMsg=""
printVector=rep(0,length(d)) #create a vector to later be populated with the dice roll results

#the if-else below precludes most invalid input and determines if 1 die or multiple dice are being rolled
if (length(d)>1 & is.numeric(d) & !(is.matrix(dice))) #if more than 1 die
{
for (i in 1:length(d))
{
if (internal.checkDie(d[i]) != "") {warnMsg=internal.checkDie(d[i])}
printVector[i]=internal.roll(d[i])
}
}
else
{
if (is.numeric(d) & !(is.matrix(dice))) #if one die
{
if (internal.checkDie(d) != "") {warnMsg=internal.checkDie(d)}
printVector[1]=internal.roll(d)
}
else #if invalid dice value (ex. matrix, character string, etc.)
{
warning("invalid dice entry: returning NA")
return(NA_real_)
}
}

if (fakeDice==FALSE & warnMsg != "") {print(warnMsg)}
if (!noprint) {print(paste(c("Dice results:",printVector),collapse=" "), quote=FALSE)}
sum=sum(printVector)+mod
return(sum)
}



plotrolladv= function(mod=0, adv=TRUE, dice="1d20", details=FALSE, m=1000)
{
#as plotroll, but integrates advantage/disadvantage (that is, roll twice and use the best result for advantage and worst result for disadvantage)
#to simulate rolling with disadvantage, set adv to FALSE
#dice is not the first variable listed in the function-call because advantage/disadvantage almost if not always applies to a single d20 roll, so the user should not need to input the 'dice' variable. It is kept as an option
#some users may experience a lag when using large dice pools. Program ran successfully even at dice="100d100"

dicepool=dicevector(dice)
numDice=length(dicepool)

#simulation of m rolls of dice with advantage or disadvantage
resultvector=rep(0, m)
for (i in 1:m)
{
for (z in 1:numDice) #this loop is needed when the dice pool contains more than 1 die
{
roll1 = roll(dice=dicepool[z], mod=0, fakeDice=TRUE, noprint=TRUE) #roll works when 'dice' is equal to a single number instead of a string
roll2 = roll(dice=dicepool[z], mod=0, fakeDice=TRUE, noprint=TRUE)

if (adv) resultvector[i]=resultvector[i]+pmax(roll1, roll2)
else resultvector[i]=resultvector[i]+pmin(roll1,roll2)
}
}
resultvector=resultvector+mod

#calculate true mean
meanvector=rep(NA,numDice) #create a vector to house the true mean of each individual die in a dice pool
#loop through each die in the dice pool, determining the true mean
for (i in 1:numDice)
{
dsize = dicepool[i]
dmatrix = cbind(rep(1:dsize,dsize),rep(1:dsize,each=dsize))
if (adv) meanvector[i]=mean(pmax(dmatrix[,1],dmatrix[,2]))
else meanvector[i]=mean(pmin(dmatrix[,1],dmatrix[,2]))
}
truemean=sum(meanvector) + mod

#graph
if (adv) {advstring=" with Advantage"} else {advstring=" with Disadvantage"}
if (mod==0) {modstring=""} else {modstring=paste(c("+",mod), collapse="")}
titlestring=paste(c("Results for ",dice,modstring,advstring), collapse="")
subtitle=paste(c("red line is empirical mean (",mean(resultvector),"); blue line is the true mean (",truemean,")"), collapse="")
ylabel=paste(c("Frequency during",m,"Simulated Rolls"), collapse=" ")
plotvector=table(resultvector)
plot(plotvector, main=titlestring, sub=subtitle, xlab="Results of the Roll", ylab=ylabel)
abline(v=mean(resultvector), lwd=2, lty=2, col="red")
abline(v=truemean, lwd=2, lty=3, col="blue")

#output
if (!details) return(mean(resultvector)) #returns average result (mean)
else return(summary(resultvector)) #returns five-number summary
}



Also, here's R code to do something like 4d6b3, for rolling stats. It has code to sort the results of a vector, which could be used to sort the row of a column to drop the lowest value.


roll4d6b3 = function(n=6)
{ #simulates the practice of "roll 4d6 best 3" for rolling for ability scores during character creation
#n is how many times it rolls 4d6 and calculates the sum of the highest 3
#defaults to 6 for the six ability scores: Str, Dex, Con, Wis, Int, Cha
dice=sample(6, n*4, replace=TRUE) #creates vector of size n*4, each simulating a random roll of a d6
returnThis=rep(NA, n)

d=1
for (i in 1:n)
{ #look at each set of 4 dice (ex. 1-4) to determine which 3 are the best three, then sum them into returnThis[i]
tempD=c(dice[d],dice[d+1],dice[d+2],dice[d+3])
tempD=sort(tempD)
tempD[1]=0
returnThis[i]=sum(tempD)
d=d+4
}
return(returnThis)
}

[Knaight]

See the analysis of that number set is what I don't know how to do. What did you mean by that?

Say you're rolling 5d6 best three, and you have [1,3,4,2,3] as your individual numbers rolled. You rolled a 10.

Also here's the Octave/Matlab code. I wrote it as a function so you can just substitute in the dice and which you keep.

function [retval] = DicePoolDiscard (input1, input2)
%This figures out the probability distribution of a roll and keep system.
%Input 1 is a string of dice, expressed in number of sides. [2,4,6,8,10] would be 1d2, 1d4, 1d8, 1d8, 1d10.
%Input 2 expresses which dice are kept in an ordered list, using 1s and 0s. [0,1,1,1,0] would keep the middle three dice.

%CHECK SIZE LENGTH
if abs(length(input1)-length(input2))>0.5;
error('Input lengths are different')
end
%VERIFY INTEGERS
input1=int32(round(input1));

%RECURSIVE COUNTER
Counter = 1+zeros(1,length(input1)); %Counts up by 1, where each "digit" is a seperated out and the base is per digit by die.
Startdigit = length(input1); %Looks at last digit first. Essentially a pointer.

%TOTAL COUNT
Count = 1;
for i = 1:length(input1);
Count = Count*input1(i);
end
Tracker = zeros (1,Count);

%COUNTING
for j = 1:Count
Incremented = 0; %Tracks whether or not the counter has been updated for the roll.
Tracker(j) = sum(input2.*sort(Counter)); %Tracks individual roll.
if j == Count
Incremented = 1; %Forced break.
end
%Updates Counter
UseDigit = Startdigit;
while Incremented == 0;
if input1(UseDigit)-Counter(UseDigit)>0.5;
Counter(UseDigit)=Counter(UseDigit)+1;
Incremented = 1;
else
Counter(UseDigit)=1;
UseDigit=UseDigit-1;
end
end
end

%TRACKER TO FREQUENCY DATA
Tracker = int32(sort(Tracker)); %The int32 is to allow the next loop to work.
min=Tracker(1);
max=Tracker(Count);
ValDist = linspace (min, max, (max-min+1));
Freq=zeros(1,(max-min+1));

for k = 1:Count
ValCurrent = Tracker(k);
Freq(ValCurrent-min+1) = Freq(ValCurrent-min+1)+1;
end

Freq=Freq./sum(Freq);
PercentageFreq=Freq*100;
retval = [ValDist;PercentageFreq];
endfunction


EDIT: This doesn't calculate the mean, just the probability of every possible result. That output set can then be used to easily calculate means, standard deviations, whatever.

EDIT 2: I think there's a slight bug that crops up when using long values here. If it were more than thrown together demonstration code I'd fix it, but that's what it is.

[jayem]

For high numbers of dice rolls you can pretty much assume that the lowest number will be a 1.

There's 0.9^10 of the numbers all being greater than 1 (35%)
There's 0.8^10 of all the numbers being greater than 2 (10.7%)
There's 0.7^10 of all the numbers being greater than 2 (2.8%)
etc...
So the chances that the lowest number is:
1=1.0^10*(1-0.9^10) (65%)
2=0.9^10*(1-0.8^10) (31.3%)
3=0.8^10*(1-0.7^10) (10.4%)
etc... [though from here you can kind of guess the 2nd half is pretty close to 1]
And for the average 'lowest roll' you can just add them up.
(1*0.65+2*0.313+3*0.104+4*0.7^10+5*0.6^10+6*0.5^10 +7*0.4^10+8*0.3^10+9*0.2^10+0.1^10)

So I get the actual average as 1.1.

So at that point for totals you can basically take 1 off the 'average'.

This of course is calculated independently of the total die roll (5.5*10). So although you can take the difference for an average thing, when dealing with individual cases or distributions you have to be a bit more cautious.

However in general it will be almost certainly to be 1 at the lower end (exactly one where total is less than 20), and will be higher at the far end (ending up at exactly 10 when the total is 100). But for much of the range the lowest is going to be between 1&2.

____________________________-For two die
For each dice, the effect of rolling with advantage is the same as having a 36 sided die with 11 (6's), 9 (5's) 7 (4's), 5 (3's), 3(2's),1 (1)
You can draw the table and count them [as described above], and see how the pattern generalises.

This looks linear (so I suspect 3 will be cubic or quartic), will have a look for a few small numbers but am busy.

[Errata]

See the analysis of that number set is what I don't know how to do. What did you mean by that?

The brute force method that is exact, not an approximation, is in your example of 1d4+1d8+1d12 drop 1, instead of calculating 1000 possible random rolls, to systematically go through all 4 * 8 * 12 possible rolls. In this case, it works out to 384 possible rolls, so it's actually cheaper than your method. If you get to extreme cases like rolling hundreds of dice, then this method will break down as too computationally expensive. One possibility, would be to figure out how much you can afford to calculate, and have two methods. One systematic brute force method if it's below, say 1 million possibilities, and the monte carlo random approximation if it's over that many possible rolls.

To systematically go through all possible rolls, you'd just be doing:
1, 1, 1. 1, 1, 2. 1, 1, 3... etc all the way up to 4, 8, 12. In each of those it's easy to determine which one you'd drop and what the value is. Each permutation has odds of exactly 1/384, so it's easy to add them up for whatever statistics you're trying to generate. To do the permutations, you may want to initialize one array with the base of all the dice (e.g. 2d4 + 3d8 would be 2, 2, 8, 8, 8), and one array with the current permutation (e.g. 1, 1, 1, 1, 1), then you increment the lowest position until it overflows, then reset it and increment the next position.

[Errata]

I thought more about more efficient ways to calculate this. There are probably better references out there that have delved into exactly this and have a more efficient method, but it's an interesting problem to work through, so here is my attempt.

For 20d6 drop 8:

Odds of the lowest roll being 6 is (1/6)^20.
Odds of the lowest roll being 5 is (2/6)^20 - (1/6)^20.
Odds of the lowest roll being 4 is (3/6)^20 - (2/6)^20 etc.

So take each case and stipulate temporarily that that is the lowest roll. Then multiply whatever result you get with that assumption by the odds of that being true.
If you assume the lowest roll is 4, you know you can't roll lower than 4, so they're not really d6 anymore. You rolled at least one 4. So in this case, it's the odds of the lowest being 4 odds of the result being the same as 19d3 drop 7 + 12*3
The degenerate case of the lowest being 6 is just 12*6.

Now you need to calculate 19d2 drop 7, 19d3 drop 7, 19d4 drop 7, 19d5 drop 7, and 19d6 drop 7. At first glance, it seems like it would grow exponentially for a large number of dropped dice, but I think not, since these will end up needing to calculate the same few expressions. You just need to make sure that they are re-used rather than recalculated, and it should be efficient.

It gets a little more complex with heterogeneous dice, but not that much more.

[jayem]

Brute force has lots of pluses,
you can simply change the 'score' function for the die set, and get any date you want
it will be perfectly accurate
but it can be expensive

Monte Carlo has the flexibility, and can deal with crazy amounts of die but will get errors and can't quite be sure how much (though if you know roughly what it looks like you can get a pretty good idea)

If you can keep the appropriate distribution tabled then that saves a lot of time over brute force (rather than doing the calculation for 1&6, 2&5, 3&4, 4&3, 5&2, 6&1 you can calculate it once). However you lose some information. For this case we also need to keep the minimum.

If you can deal with an approximate distribution you can really take it to silly numbers that monte-carlo struggles with. This works best if you can split it into 'independent' functions, which will be hard to do in this case. My "for large numbers of die the lowest number is almost certainly 1" is a very crude method in this family.


We need to keep the probabilities for each combination of total and the mimimum
For 1d6 the probabilities of each (total, minimum) options are trivially P(1,1)=P(2,2)=P(3,3)=P(4,4)=P(5,5)=P(6,6)=1/6 . Anything else=0

Given two distributions, P&Q the probabilities for the combined distribution R can be found as follows
R(t,m)=sum of all P(t1,m1)*Q(t2*m2) where t=t1+t2 and (m=m1 and m2>m1) or (m=m2 and m1>m2)



first={(1,1):1.0/6,(2,2):1.0/6,(3,3):1.0/6,(4,4):1.0/6,(5,5):1.0/6,(6,6):1.0/6}
count=1

def Combine(P,Q):
output={} #use a dictionary to store combinations (not an efficient way of doing it)
for (t1,m1) in P.keys():
for (t2,m2) in Q.keys():
res=(t1+t2, min(m1,m2))
global count
count=count+1 #keep record of loop calls
#Add to relevant total
if res in output:
output[res]=output[res]+P[(t1,m1)]*Q[(t2,m2)]
else:
output[res]=P[(t1,m1)]*Q[(t2,m2)]
return output

#When complete we can find the boxes that have the same Advantaged toal
def conclude(P):
output={}
for (t,m) in P.keys():
if (t-m) in output:
output[t-m]=output[t-m]+P[(t,m)]
else:
output[t-m]=P[(t,m)]
total=0.0
for i in range(1,100000):
if i in output:
print (str(i)+":"+str(output[i]),)
total=total+output[i]
print total #should be 1

dists=[first]
for t in range(0,10):
print t
dists.append(Combine(dists[t],dists[t]))
conclude(dists[2])
print count



For 4 Die I get
('3:0.000771604938272',)
('4:0.00308641975309',)
('5:0.00771604938272',)
('6:0.0162037037037',)
('7:0.0293209876543',)
('8:0.0478395061728',)
('9:0.0702160493827',)
('10:0.0941358024691',) (3Die Average)
('11:0.114197530864',)
('12:0.128858024691',)
('13:0.132716049383',) Peak
('14:0.123456790123',) (4 die Average)
('15:0.101080246914',)
('16:0.0725308641975',)
('17:0.0416666666667',)
('18:0.0162037037037',)
1.0

Note of course that totals over the penultimate die are impossible, and it's slightly right shifted compared to 3die. This probably is quite an easy correction factor to apply.

And for 32 die, comparing with loss to without loss
('100:0.021707013425--0.0192154853948',)
('101:0.0242586764818--0.0217067732967',)
('102:0.0268211110944--0.0242583239806',)
('103:0.0293391593033--0.0268206039946',)
('104:0.0317538474--0.0293384442643',)
('105:0.0340045759349--0.0317528589834',)
('106:0.0360315747804--0.0340032363284',)
('107:0.0377784890122--0.0360297945553',)
('108:0.03919494263--0.0377761692339',)
('109:0.0402389197576--0.0391919785996',)
('110:0.0408788082362--0.0402352065533',)
('111:0.0410949685986--0.0408742480647',)
('112:0.040880721116--0.0410894797497',)
('113:0.0402426826143--0.0408742480647',)
('114:0.0392004297442--0.0402352065533',)
('115:0.0377855124168--0.0391919785996',)
('116:0.0360398860242--0.0377761692339',)
('117:0.0340138699218--0.0360297945553',)
('118:0.0317637691675--0.0340032363284',)
('119:0.02934931438--0.0317528589834',)
('120:0.0268310796842--0.0293384442643',)
('121:0.0242680312008--0.0268206039946',)
('122:0.0217153397571--0.0242583239806',)
('123:0.0192225637806--0.0217067732967',)
0.784415286462


Again note the peak value is one less for the 'with subtraction', [but it is clearly a different shape with a flatter local surface. ETA actually it's pretty similar if I don't mix up the non subtraction version]
Note also that it is 1 less than the total die, you probably also want to compare it with the case where you never rolled the lost die, where the average will be 2.5 higher (I think).

It also takes minutes to do 1024 die (there are 36000 possible batches outcomes), and the results are effectively 0 on the tail (which I think will mean the centre will be very very similar, but I don't know).

[Knaight]

While the brute force situation can get expensive (e.g. hundreds of dice) in terms of systems that actually exist it doesn't seem particularly likely to - the most I've seen is roll and keep with ten dice, which is a hundred million results. That's somewhere like 30 seconds of calculation on a modest computer with a program that isn't total garbage but isn't necessarily well optimized either.

[jayem]

While the brute force situation can get expensive (e.g. hundreds of dice) in terms of systems that actually exist it doesn't seem particularly likely to - the most I've seen is roll and keep with ten dice, which is a hundred million results. That's somewhere like 30 seconds of calculation on a modest computer with a program that isn't total garbage but isn't necessarily well optimized either.
I'm pretty sure it would be pretty pointless much above ten. It would be interesting to find where the boundary is, and confirm this. But brute force is clearly viable in most practical cases.

I think about 65% of NdN's outcomes contain at least a 1. Which seems like it ought to be pretty near boring.
For 10D10 the lowest total that need not is of course 20, and the highest total that might is of course 91 (but this range contains almost all the outcomes, so we don't really get anywhere).

[Errata]

I think about 65% of NdN's outcomes contain at least a 1. Which seems like it ought to be pretty near boring.

That's assuming there is only 1 dropped value. If there are a large number of dice rolled, there might be a larger number of dice dropped as well. The algorithm I outlined above should give a way to calculate it efficiently and precisely even for large numbers of dice rolled and large numbers dropped, without relying on approximations.

[jayem]

That's assuming there is only 1 dropped value. If there are a large number of dice rolled, there might be a larger number of dice dropped as well. The algorithm I outlined above should give a way to calculate it efficiently and precisely even for large numbers of dice rolled and large numbers dropped, without relying on approximations.

Good point. I still think it will tend to a very narrow distribution.
And algorithm looks promising. I thought to do multiple drops you'd pretty much have to keep that number of die in memory (though actually just recording the number for each result would also work?! *) which would get costly.

So you could build up your table. It's probably easier to build it from the bottom?
I think you'd need (total, die count, die size (implicit minimum), die dropped).
So it's not going to be easy or quick, but it's probably the best we can do.

*yes but it's also very costly. For aDb you need a^(b-1) boxes to hold the final distribution, and each requires calculation from distributions below.

[factotum]

*yes but it's also very costly. For aDb you need a^(b-1) boxes to hold the final distribution, and each requires calculation from distributions below.

You don't need to store the result of every die roll? You only need to store the actual total. e.g. for 4d6 drop 1, you don't store all 1,296 possible sets of die rolls, you only store the 16 possible results (3 to 18) with a number in each indicating how many times they came up. It's pretty easy to calculate the maximum possible roll and allocate an array appropriately.

[jayem]

You don't need to store the result of every die roll? You only need to store the actual total. e.g. for 4d6 drop 1, you don't store all 1,296 possible sets of die rolls, you only store the 16 possible results (3 to 18) with a number in each indicating how many times they came up. It's pretty easy to calculate the maximum possible roll and allocate an array appropriately.
Exactly, and you can use that to save time as well, rather than going through each combinations. if you store previous data and cross multiplying them all. So for 4D6 rather than doing 1296 individual throws. You take the 12 results for 2D6, and do 144 calculations (noting e.g. that the 6 ways of getting 7 mean that there are 36 ways of getting 14 as a result of getting a 7 with each pair, and then adding 8*6 etc into the box to get the total number of 14s.)

But for dropping 1, you need to keep (6*) more data, as it now matters what the lowest is as well.

And for drop 2,3,4 you need to start keeping ever more information as they interact at the end (I mean you get to 10ish dropped dice in memory, but computationally 1000-1D6 took minutes, so it would take days for 1000-4D6, which is still miles better than naive brute force where 20D6 takes days). In the limit of dropping all the die, this would of course be keeping every roll.

Errata's method, however, I think gives a lease of life beyond that.

I'd hoped just skipping straight to storing how many of each number you get, you'd have a perfectly flexible compromise. And it kind of does, you get perfect data, and can do about 15 D6 in the time it takes to do 10 by brute force. And it's exponential wrt die size not die count, but that still won't get very far.

[Knaight]

You don't need to store the result of every die roll? You only need to store the actual total. e.g. for 4d6 drop 1, you don't store all 1,296 possible sets of die rolls, you only store the 16 possible results (3 to 18) with a number in each indicating how many times they came up. It's pretty easy to calculate the maximum possible roll and allocate an array appropriately.

It's also easy to know how many unmodified rolls you need and preallocate an array for that, updating as necessary. I made sure to do it that way in my example code, mostly because even small test cases potentially get out of hand if I'm repeatedly copying and deleting arrays that long when it's totally unnecessary.

[factotum]

But for dropping 1, you need to keep (6*) more data, as it now matters what the lowest is as well.


But it doesn't? It only matters insofar as it alters the total you'll get from that dice roll, so you just calculate the total and increment the appropriate array element. So, if you roll 4, 4, 3 and 2, the total is 11 and you increment array element 11. You neither know nor care that this result had a lowest roll of 2, because it's irrelevant in the long run.

[jayem]

But it doesn't? It only matters insofar as it alters the total you'll get from that dice roll, so you just calculate the total and increment the appropriate array element. So, if you roll 4, 4, 3 and 2, the total is 11 and you increment array element 11. You neither know nor care that this result had a lowest roll of 2, because it's irrelevant in the long run.

At the very end. But if you want to use your mid result you need to keep it (otherwise you have to redo it every time).
Without minimums You can merge (8&5) not caring about whether it came from (2&6,2&3) (3&5,2&3) etc... which saves time (to add an extra die you have to do 6 steps per element in the range(I.E 6*N), rather than 6 times the outcomes (6^N) ) and doesn't cost much space.

With minimums you also need to keep that with you, you can keep (5,2,1)(4,3,1)(3,4,1)(4,4,1)(1,2,5) etc... in the same box (Total 8, Minimum 1) and use that to calculate 4 die with.
And so put (5,2,1,1)(4,3,1,1)(3,4,1,1)(4,4,1,1)(1,2,5,1) etc... straight into the (T9,M1),(T10,M2) etc.. boxes in six steps.
Then add all the (T8,M2) into the (T10,M2)
But if (2,2,2,2) was in the same box then that is different (unless it gets merged with a 1 later).
It's only when you get the final total (and then take the minimum off) that you can lose the minimum. Which if storage were the critical issue could easily be done by doing them each individually.

[Knaight]

At the very end. But if you want to use your mid result you need to keep it (otherwise you have to redo it every time).
Without minimums You can merge (8&5) not caring about whether it came from (2&6,2&3) (3&5,2&3) etc... which saves time (to add an extra die you have to do 6 steps per element in the range(I.E 6*N), rather than 6 times the outcomes (6^N) ) and doesn't cost much space.

With minimums you also need to keep that with you, you can keep (5,2,1)(4,3,1)(3,4,1)(4,4,1)(1,2,5) etc... in the same box (Total 8, Minimum 1) and use that to calculate 4 die with.
And so put (5,2,1,1)(4,3,1,1)(3,4,1,1)(4,4,1,1)(1,2,5,1) etc... straight into the (T9,M1),(T10,M2) etc.. boxes in six steps.
Then add all the (T8,M2) into the (T10,M2)
But if (2,2,2,2) was in the same box then that is different (unless it gets merged with a 1 later).
It's only when you get the final total (and then take the minimum off) that you can lose the minimum. Which if storage were the critical issue could easily be done by doing them each individually.

You'd want to vary the minimum to an array as determined by the initial drop conditions.

That said, this involves some unnecessary storage - you don't ever need to store more than two rolls at one time, in any state of completion - one as a counter, and one as raw data. Sort it, add the ones you're using, call it a day. Your method appears to involve building up all of them in parallel, which is faster but demands a lot more active memory.

Though I still maintain that brute force should work fine for the sort of rolls people actually make.

[jayem]

Your method appears to involve building up all of them in parallel, which is faster but demands a lot more active memory.

Though I still maintain that brute force should work fine for the sort of rolls people actually make.

Exactly, though by the time memory is a problem the processor cant cope.

[crayzz]

So brute force is fun and always works, but I wanted to see if I could come up with something more elegant.

So, simplest case: 2d4 drop lowest. only 16 possibilities.




1
2
3
4


1
1
2
3
4


2
2
2
3
4


3
3
3
3
4


4
4
4
4
4



You end up with seven ways to get 4, five ways to get 3, three ways to get 2, and only one way to get 1. If you continue on to larger dice, you'll find that the pattern holds: there are always 2n-1 ways to get n from a 2dXd1 roll. Since there are X2 possibilities, the probability of rolling a given number is (2n-1)/X2.

Complicating this a bit: 1d4 + 1d6 drop lowest. Now there's 24 possibilities.




1
2
3
4
5
6


1
1
2
3
4
5
6


2
2
2
3
4
5
6


3
3
3
3
4
5
6


4
4
4
4
4
5
6



So interestingly, 2n-1 ways to roll a value holds for all values of 4 and under. For 5 and 6, there are only four ways to roll. I'm having a hard time summarizing that distribution in a neat, concise way, but extrapolating to two dice of size X1 and X2, where X1<=X2:



P(n) =
(2n-1)/(X1X2)
n<=X1



1/X2
n>X1



Where X1 = X2, the above just simplifies to the 1st case.

[Knaight]

Looking at that data visualization, that suggests that there is a recursive geometric perspective that can be taken to an arbitrarily high number of dimensions (one per die), where the Nth dimension has a (N-1)th dimensional shell. Defining a hypervolume V(a,b,c,d,e,f,...) where the lowercase letters are individual dice you can essentially create shells of a lower dimension V/a, V/b, V/c, etc. then cascade them inwards, while dropping yet a further dimension down (V/ab, V/ac, V/ad) etc. to eliminate intersections in lower dimensions. It's a top down structure that cuts down on the calculations a bit.

That said, this is a first look at a high level problem solving strategy which needs a lot of detail work, and could potentially break down if weird phenomena start showing up at higher dimensions. Still, three's a potential strategy there.

It is however a lot less robust in terms of cases handled - it's a short cut for the single drop case, and while there are probably symmetries for higher drop cases as well that's an area where it can easily get out of hand quickly.

[crayzz]

Looking at that data visualization, that suggests that there is a recursive geometric perspective that can be taken to an arbitrarily high number of dimensions (one per die), where the Nth dimension has a (N-1)th dimensional shell. Defining a hypervolume V(a,b,c,d,e,f,...) where the lowercase letters are individual dice you can essentially create shells of a lower dimension V/a, V/b, V/c, etc. then cascade them inwards, while dropping yet a further dimension down (V/ab, V/ac, V/ad) etc. to eliminate intersections in lower dimensions. It's a top down structure that cuts down on the calculations a bit.

That said, this is a first look at a high level problem solving strategy which needs a lot of detail work, and could potentially break down if weird phenomena start showing up at higher dimensions. Still, three's a potential strategy there.

It is however a lot less robust in terms of cases handled - it's a short cut for the single drop case, and while there are probably symmetries for higher drop cases as well that's an area where it can easily get out of hand quickly.

"Recursive Geometry" was basically was I was going for before I realized things were gonna get out of hand. Rolling an arbitrary number of dice of arbitrary size and keeping only the highest roll (e.g. 1d4 + 1d8 + 1d20 drop 2) ends up being relatively easy to calculate. Making it truly arbitrary gets complicated, fast.

[Errata]

So you could build up your table. It's probably easier to build it from the bottom?
I think you'd need (total, die count, die size (implicit minimum), die dropped).
So it's not going to be easy or quick, but it's probably the best we can do.

I think if you were doing say XdY drop Z, you would start from all the drop 0 result, then use that to calculate drop 1 results, then you would no longer need the drop 0 results and could free that memory. So you would only need 2 different drop amounts in memory at any one time.

You would array with buckets for the frequency of all X*Y possible total roll values, times all Y different die size (approximately), times 2 drop values at a a time. So for example, if you wanted to calculate 200d100 drop 100, with odds of each roll to a precision of 64-bit double precision floating point, that would be around 32MB of RAM. Which is not bad at all for such an unreasonably large about of dice. For 20d20 drop 10, it would be only 128k. A more reasonable roll like 6d6 drop 3 would be around 2k RAM, though at that point brute force is perfectly fine.

And that's assuming that you want to know the odds of every possible dice roll (which would be hard to even display the full result properly for very large dice values). If you have a more limited interest like the average or median, then you need to track a lot less data for that.

ETA: I overestimated those by a factor of around 2, because it's not X*Y possible roll values, it's (X-Z)*Y.

[Knaight]

"Recursive Geometry" was basically was I was going for before I realized things were gonna get out of hand. Rolling an arbitrary number of dice of arbitrary size and keeping only the highest roll (e.g. 1d4 + 1d8 + 1d20 drop 2) ends up being relatively easy to calculate. Making it truly arbitrary gets complicated, fast.

It definitely looked like something that was less an approach to a die app and more the beginning steps of a serious research paper in discrete math.

[Sivarias]

Hey guys! So I've been trying to follow along. And it seems my best bet that I understand is to not bother with the approximation since the values we'll be dealing with are small enough to actually calculate the true value.

To help out a bit more I can't fathom a situation where more that 6 dice will be rolled ever, and you only ever keep the highest or two highest ever. Dice will be from 4 to 12 so the bigger dice (20 and 100) won't be there to eat up the memory.

So in case something comes up that I'm not foreseeing as long as the program can handle up to 8 dice without taking a long time it'll be fine.

The table solution seems the easiest if I'm understanding correctly. Create a 2-D array (Or parallel arrays in this case because of C++) check for the two highest values and sum them. Shave off a few seconds if I wanted to get complicated and have it check to see if it's got the two highest possible values.

[AMFV]

One of the things that you have to realize is that probability generators won't work to give you exact percentages for the die. Because they aren't perfectly milled, sanded, and are possibly worn. Not to mention that the numbers written on them involve the removal of different amounts of material on them that would affect the weighting. Now this probably doesn't matter too much in a single roll, but once you start getting into large scale probabilities it would.

[Sivarias]

Of course!

But I need a general idea of how different abilities interact to make sure thier all on par.

[Errata]

The table solution seems the easiest if I'm understanding correctly. Create a 2-D array (Or parallel arrays in this case because of C++) check for the two highest values and sum them. Shave off a few seconds if I wanted to get complicated and have it check to see if it's got the two highest possible values.

Not really. If you're referring to crayzz's table, that is for the very special case of rolling 2 dice and dropping 1. It gets a lot more complicated to generalize.

Since you've decided that you only care about fewer than 10 dice being rolled at a time, then by far the simplest solution is the brute force method of systematically going through every permutation. That solution gets slow for very large numbers of dice, but works perfectly fine for small numbers of dice and is very simple to implement.

[Sivarias]

I was actually referring to jeenlean's idea because that seems to be the way to calculate it manually. You create an array with every possible combination and go through checking for the highest values.

[JeenLeen]

I was actually referring to jeenlean's idea because that seems to be the way to calculate it manually. You create an array with every possible combination and go through checking for the highest values.

Yeah, that's what my idea does (or should do.)

It made the array as big as you needed, so if you have 10d6, it'd be ten columns (one for each die) and as many rows as you needed to do every possible combination.
In my idea, I tested the code using 100d100, and it took a noticeable time to run*, but it still ran fine. This was using R on an about-5 or 6 year old computer (though the computer was decent when new. Not top of the line, but pretty good).
*I forget the details, but about a minute to a few minutes. It wasn't something like a "leave it running for 20 minutes" job.

For this purpose (getting a distribution of dropping lowest die, as opposed to mine of getting a mean), you'd need to add another loop to take out the lowest value and then store the sums somewhere else to get your distribution, but I reckon that wouldn't be computationally expensive.

[Errata]

I guess it depends on your environment what implementation makes sense, but in general there is no reason to keep that much data around at once in a big table like that.

For the brute force method, you only need a single 1 dimensional array with a number of elements equal to the maximum possible roll (obtained by adding up the sides of all the dice if you rolled the highest value on each). Initialize the array to 0 in every element. Go through every possible permutation of dice rolls, and for each permutation, calculate what the sum of the roll would be after applying rules like drop N, then increment the array element for that value by 1. In the end you have an array that tells you the relative frequency for every possible roll. To convert them into a fraction, you just divide the integer in each element by the total number of permutations you went through.

As for how you go through every possible permutation, I discussed one easy way to implement that that requires two small arrays equal to the number of dice rolled.

[Calthropstu]

In my experience, dropping dice results in lost dice.

Hope this helps!

[jayem]

See above for how to do with normal amounts of die.


I can't really remember. I think I'd been planning to build up from even powers of the total Die count in which case you needed the lower drops kept (but you could in theory discard the lower die counts once used).



Taking the initial algorithm Errata gave.


For 20d6 drop 8:

Odds of the lowest roll being 6 is (1/6)^20.
Odds of the lowest roll being 5 is (2/6)^20 - (1/6)^20.
Odds of the lowest roll being 4 is (3/6)^20 - (2/6)^20 etc.

So take each case and stipulate temporarily that that is the lowest roll. Then multiply whatever result you get with that assumption by the odds of that being true.

If you assume the lowest roll is 4, you know you can't roll lower than 4, so they're not really d6 anymore. You rolled at least one 4. So in this case, it's the odds of the lowest being 4 odds of the result being the same as 19d3 drop 7 + 12*3
The degenerate case of the lowest being 6 is just 12*6.

Using N continuous (0-1) dice temporarily, so there is a M'th least die'. We can skip straight to that.
The odds of the m die being less than an X is x^m*(1-x)^(n-m)*(n!)/(m!)/(n-m!)
[note that n&m are constant for any one 'game']

If the die are all greater than X, then we can consider the situation where they can roll between X&1 which is easy to calculate the average = (x+1)/2

The average overall is thus just the integral of
x^m*(1-x)^(n-m)*(n!)/(m!)/(n-m!) * (x+1)/2

Which has to be pretty near doable, though you may need to break it into sections and approximate each.
For discrete die you'd be able to sum exactly, but then you'd need to deal with die on the 'border'

I'll have to put some more thought in it

[DavidSh]

See above for how to do with normal amounts of die.


Taking the initial algorithm Errata gave.

Using N continuous (0-1) dice temporarily, so there is a M'th least die'. We can skip straight to that.
The odds of the m die being less than an X is x^m*(1-x)^(n-m)*(n!)/(m!)/(n-m!)
[note that n&m are constant for any one 'game']

If the die are all greater than X, then we can consider the situation where they can roll between X&1 which is easy to calculate the average = (x+1)/2

The average overall is thus just the integral of
x^m*(1-x)^(n-m)*(n!)/(m!)/(n-m!) * (x+1)/2

Which has to be pretty near doable, though you may need to break it into sections and approximate each.
For discrete die you'd be able to sum exactly, but then you'd need to deal with die on the 'border'

I'll have to put some more thought in it


I don't think those spoilerified formulas can be correct. If I plug in n=2, m=1, the integral you assert becomes the integral (from 0 to 1) of x(1-x)*2 * (x+1)/2. This is the integral of x - x^3, which is 1/4.

I'm pretty sure that taking two independent uniform random variables over [0,1], and dropping the smaller, leaves you with a value that has expectation greater than 1/2, and I think the expectation should be 2/3.

I can't think of a general formula that doesn't have summations in it.

[Errata]

My formula was a lot simpler, and I'm pretty sure it's valid.

If you roll XdY dice, and the lowest value is Y, then every single one of those dice has to be a Y, so the odds are (1/Y)^X. Very simple.
For the lowest value to be Y-1, then every single dice has to be a Y-1 or a Y, so the odds are (2/Y)^X, minus the odds of the lowest being Y, which is a subset of the possible conditions in which they are all Y-1 or Y. So it's (2/Y)^X - (1/Y)^X, and so on.

I don't think you need to worry about factorials for this specific situation.

ETA: I see, you're trying not to calculate the odds of the lowest dice, but the Nth lowest dice. I don't follow the math on how to do that.

[jayem]

My formula was a lot simpler, and I'm pretty sure it's valid.

If you roll XdY dice, and the lowest value is Y, then every single one of those dice has to be a Y, so the odds are (1/Y)^X. Very simple.
For the lowest value to be Y-1, then every single dice has to be a Y-1 or a Y, so the odds are (2/Y)^X, minus the odds of the lowest being Y, which is a subset of the possible conditions in which they are all Y-1 or Y. So it's (2/Y)^X - (1/Y)^X, and so on.

I don't think you need to worry about factorials for this specific situation.

ETA: I see, you're trying not to calculate the odds of the lowest dice, but the Mth lowest dice. I don't follow the math on how to do that.

It's not so bad, based on your work (and I think I see where I went firstly wrong)
Each die has a probability X of being less than X (the die goes from 0-1), Our case is when M-1 die are less than X, so multiply them together. NB the Mth die is of course exactly X (which is where the dX comes from), this is my error.
Each die has a probability (1-X) of being more than X. Again There are N-M die in the situation.
We can then rearrange all the dice, in both groups, as in every other binomial case (hence the factorials), except when corrected it's actually a trinomial distribution with a single element in the middle block [for the discrete die case the middle block can also vary in the same 'game', so you need the double sum you had].

x^(m-1)*(1-x)^(n-m)*(n!)/(m-1!)/(n-m!) * (x+1)/2
Which seems much closer to the predicted 2/3rds for the 2dice drop 1 situation and giving about 3/4 for 3dice drop 2 [ETA which seems believable with Monte Carlo].

[DavidSh]

Well, I went and looked up a few things. The term for values like the k-th smallest of a set of independent identically distributed random variables is order statistic. The Wikipedia article on Order Statistics gives the distribution for the k-th smallest of n random variables uniform over the unit interval [0,1] as Beta(k,n+1,k), and says that is has mean k/(n+1). Not too surprising. Then the expected value of taking n random uniform [0,1] variables, dropping the m smallest, and summing the rest, is the sum{k=m+1 to n} of k/(n+1), which I calculate to be (n+m+1)(n-m)/(2(n+1)).

(Even though the order statistics are not independent from one another, I can still calculate the mean of a sum as the sum of means. Means are nice that way.)

At least, this evaluates to 2/3 when n=2 and m=1, so it passes one sanity check.