Roll an N sided die K times. Let S be the side that appeared most often. What is the expected number of times S appeared?

For example, consider a 6 sided die rolled 10 times. Based on the following monte-carlo simulation, I get that the side that appears most will appear 3.44 times on average.

n = 6
k = 10
samples = 10000
results = []

for _ in range(samples):
    counts = {s:0 for s in range(n)}
    for _ in range(k):
        s = randint(0, n-1)
        counts[s] += 1


print sum(results)/float(len(results))

But I can’t figure out how to get this in a closed form for any particular N and K.

Solutions Collecting From Web of "Roll an N sided die K times. Let S be the side that appeared most often. What is the expected number of times S appeared?"

Not an answer, but it’s worth trying some example with small $N$.

Let $M$ be your number.

For $n=2$, and any $k/2 < m \leq k$ you have $P(M=m)=\frac{\binom{k}{m}}{2^{k-1}}$. If $k$ even, then $P(M=k/2)=\frac{\binom{k}{k/2}}{2^k}$.

So, for $k=2j+1$ odd, you have:

$$\begin{align}E(M)&=\frac{1}{2^{2j}}\sum_{m=j+1}^{2j+1} m\binom{2j+1}{m}\\

For $N=2$ and $k=2j$, you get: $$

It’s gonna get messier when trying it for $N=3$.

Here is a closed form, we will need more sophisticated methods for the

Following the notation introduced at this MSE
link we suppose
that the die has $m$ faces and is rolled $n$ times. Rolling the die
with the most occured value being $q$ and instances of this size being
marked yields the species

+ \mathfrak{P}_{=1}(\mathcal{Z})
+ \cdots
+ \mathcal{V}\mathfrak{P}_{=q}(\mathcal{Z})).$$

This has generating function

$$G(z,v) =
\left(\sum_{r=0}^{q-1} \frac{z^r}{r!} + v\frac{z^q}{q!}\right)^m.$$

Subtracting the values where sets of size $q$ did not occur we obtain
the generating function

$$H_{q}(z) =
\left(\sum_{r=0}^{q} \frac{z^r}{r!}\right)^m
– \left(\sum_{r=0}^{q-1} \frac{z^r}{r!}\right)^m.$$

This also follows more or less by inspection.

We then obtain for the desired quantity the closed form

$$\bbox[5px,border:2px solid #00A000]{
[z^n] \sum_{q=1}^n q H_q(z).}$$


$$L_{q}(z) =
\left(\sum_{r=0}^{q} \frac{z^r}{r!}\right)^m$$

we thus have

$$\frac{n!}{m^n} [z^n] \sum_{q=1}^n q (L_{q}(z) – L_{q-1}(z)).$$

This is

$$\frac{n!}{m^n} [z^n]
\left(n L_n(z) – \sum_{q=0}^{n-1} L_q(z)\right).$$

We also have for

$$[z^n] L_q(z) =
\sum_{k=0}^{\min(q, n)} \frac{1}{k!} [z^{n-k}]
\left(\sum_{r=0}^{q} \frac{z^r}{r!}\right)^{m-1}$$

Furthermore we obtain for $m=1$

$$[z^n] L_q(z) =
[[n \le q]] \times \frac{1}{n!}.$$

With these we can implement a recursion, which in fact on being coded
proved inferior to Maple’s fast polynomial multiplication routines. It
is included here because it memoizes coefficients of $L_q(z)$, thereby
providing a dramatic speed-up of the plots at the cost of allocating
more memory.

All of this yields the following graph where we have scaled the plot
by a factor of $n/m.$ This is it for a six-sided die:

  3+  H                                                                                           
   |  H                                                                                           
   +   H                                                                                          
   +    H                                                                                         
   +     HH                                                                                       
   |      HH                                                                                      
  2+        HH                                                                                    
   |         HHH                                                                                  
   +            HHHH                                                                              
   +                HHHHHHH                                                                       
   |                       HHHHHHHHHHH                                                            
   +                                  HHHHHHHHHHHHHHHHHHHHHHHHH                                   
   |                                                           HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH   
   0              20             40             60             80             100            120  

And here is the plot for a twelve-sided die. (I consider it worth
observing that we have the exact value for the expectation in the
case of $120$ rolls of this die, a case count that has $130$ digits,
similar to what appeared in the companion post.)

   + H                                                                                            
   |  H                                                                                           
   +  H                                                                                           
  4+   HH                                                                                         
   +     H                                                                                        
   |      H                                                                                       
   +      HH                                                                                      
   +        HH                                                                                    
   |          HHHH                                                                                
   +              HHHHHH                                                                          
  2+                    HHHHHHHHHHH                                                               
   |                               HHHHHHHHHHHHHHHHHHHHHHHHH                                      
   +                                                        HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH   
   0              20             40             60             80             100            120  

This was the Maple code.


proc(n, m)
option remember;
local rolls, res, ind, counts, least, most;
    res := 0;

    for ind from m^n to 2*m^n-1 do
        rolls := convert(ind, base, m);

        counts := map(mel->op(2, mel),
                      convert(rolls[1..n], `multiset`));

        res := res + max(counts);


L := (m, rmax) -> add(z^r/r!, r=0..rmax)^m;

X :=
proc(n, m)
    option remember;
    local H;

    H := q -> expand(L(m,q)-L(m,q-1));

    coeff(add(q*H(q), q=1..n), z, n);

LCF :=
    option remember;

    if n < 0 then return 0 fi;

    if m = 1 then
       if n <= q then return 1/n! fi;
       return 0;


(m, q)  -> add(LCF(n, m, q)*z^n, n=0..q*m);

XX :=
proc(n, m)
option remember;
local res;

    res :=
    n*LCF(n,m,n) - add(LCF(n,m,q), q=0..n-1);


proc(nmx, m)
local pts;

    pts := [seq([n, XX(n,m)/(n/m)], n=1..nmx)];