1
$\begingroup$

I came across a nice answer by Sasha (I believe a Wolfram engineer) on mathematics stack exchange regarding the distribution of the sum of the top $k$ order statistics. (This was an answer to an RPG dice rolling question, where such distributions often arise).

His answer is copied below, edited to just the top $k$ case and for result finalization:

TotalHighestScoreDistribution[{x_, z_}, dist_] := 
  Block[{v}, 
   TransformedDistribution[Total[Array[v, z]], 
    Distributed[Array[v, z], 
     OrderDistribution[{dist, x}, Range[x - z + 1, x]]]]];

(d = TotalHighestScoreDistribution[{numdice, highest}, 
     DiscreteUniformDistribution[{1, numfaces}]];
   
   pgf = Expectation[z^t, t \[Distributed] d];
   
   cl = CoefficientList[
     Series[(1 - pgf)/(1 - z), {z, 0, highest numfaces}], z];
   
   sumdist = 
    cl // Differences[Reverse@#] & // Reverse // 
     Append[#, 1 - Tr@#][[highest ;;]] &); // AbsoluteTiming

I found the generating function approach interesting.

The same result can be had more compactly as:

shorter = 
   PDF[TransformedDistribution[Tr[Array[x, numdice][[-highest ;;]]], 
     Array[x, numdice][[-highest ;;]] \[Distributed] 
      OrderDistribution[{DiscreteUniformDistribution[{1, numfaces}], 
        numdice}, Range[numdice - highest + 1, numdice]]], 
    Range[highest, highest numfaces]]; // AbsoluteTiming

Depending on the parameters, one or the other is usually markedly faster.

In any case, on larger cases, both flounder. For example, with

{numdice, numfaces, highest} = {12, 6, 6}

They both take on the order of 2 to 4 minutes to complete.

Is there a more efficient means to derive the desired exact distribution using Mathematica?

$\endgroup$

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.