Relevancy of order statistics to the roll-and-keep dice mechanic?

Question

Throw $x$ dice, each of which has $z$ sides, keep the $y$ highest values rolled, and find their sum $s$. In roleplaying games, this is called a "roll and keep dice mechanic" and notated "$x$d$z$k$y$" or perhaps just "$x$k$y$" if the number of sides $z$ is understood. Let (capital) $S$ be the random variable for the sum of the $y$ highest dice, and let (lower case) $s$ be the observed sum of the $y$ highest dice.

A combinatorial formula for the probability mass funtion $f$ for this distribution was derived by a user called "techmologist" over at PhysicsForums: Puzzling "roll X dice, choose Y highest" problem.

The pmf is

$ f(s) = $ $$\sum_{r=1}^{z}\sum_{i=0}^{x-y}\sum_{j=x-y-i+1}^{x-i}{x \choose i,j,x-i-j}(r-1)^i N(s+i+j-x-ry; x-i-j, z-r-1) $$

where $N$ is the number of ways to distribute $B$ indistinguishable balls among $C$ cells so that each cell has between 0 and $D$ balls, namely,

$$ N(B;C,D) = \sum_{k=0}^{\lfloor B/(D+1) \rfloor} (-1)^k \binom{C}{k} \binom{B-k(D+1) + C-1}{C-1} $$

My intuition is that there might be a simpler way to express this mechanic using order statistics. Let $Z_1,Z_2,\ldots,Z_x$ be the unordered random variables for the $x$ iid dice throws, which are each discrete uniform distributions on the integer values $\{1, 2, \ldots, z \}$. And let $1 \leq Z_{1:x} \leq Z_{2:x} \leq \ldots \leq Z_{x:x} \leq z$ be the order statistics, so that for example, $Z_{x:x}$ is the r.v. representing maximum roll of the $x$ dice.

In support of my intuition, it can be shown that

$ E(S) = E(Z_{(x-y+1):x} + Z_{(x-y+2):x} + \ldots + Z_{(x-1):x} + Z_{x:x}) $

For instance, for the 4d6k3 distribution, the mean is 15869/1296 (or 12.2446). And this can be verified (quite easily in Mathematica) by using either techmologist's cominatorial formula or by using the expected value of the sum of the three highest order statistics shown above or by enumerating all outcomes or by simulating a large number of dice throws.

However,

$ S \not= Z_{(x-y+1):x} + Z_{(x-y+2):x} + \ldots + Z_{(x-1):x} + Z_{x:x} $ (WRONG)

EDIT: Thanks to Douglas Zare for explaining that the LHS and RHS are in fact equal... This helped me realize that my real confusion was that the pmf of S cannot be calculated easily by taking the convolution of the pmfs of the second through fourth order statistics, because the order statistics are not independent RVs. And that leads to a follow-up question: In a case like this, how would you find the convolution of the pmfs of the order statistics on the RHS, given that order statistics are not independent RVs? I will think about a better way to phrase this question and post it separately. Thanks again!

The question is whether my intuition is correct: Is there some way to express the pmf $f$ of $S$ in terms of the pmfs of the order statistics? And if not, what's wrong with my intuition that there should be some easy way to calculate the pmf $f$ from the pmf of the order statistics?

BTW, this isn't homework. I'm just a roleplayer, not a statistician, trying to satisfy my curiosity whether this dice mechanic can be analyzed and expressed more simply using order statistics.

Answer 1

The relationship is simple: $S = Z_{(x-y+1):x} + ... + Z_{x:x}$. This should make logical sense: The sum of the top $k$ dice is the highest die plue the second highest die etc. down to the $k$th highest die. If you still believe $S \ne Z_{(x-y+1):x} + ... + Z_{x:x}$, try to exhibit a roll of the dice so that the sum of the top $k$ dice is not equal to the highest die plus the second highest etc. plus the $k$th highest.

Order statistics are only independent for trivial constant distributions. They are not independent here. Since the summands are not independent, you can't use TransformedDistribution in Mathematica. See the documentation which says "TransformedDistribution represents a transformed distribution where , , ... are independent and follow the distributions , , ...." This is why the distribution you calculate for the right hand side is not correct.

Because of the dependence, you can't determine the distribution of the sum from the distributions of the summands. The same is true in much simpler cases. If $X_0$ and $X_1$ are both $1$ with probability $1/2$ and $0$ with probability $1/2$, then it is possible that $X_0 = 1-X_1$ so that $X_0 + X_1$ is the constant $1$. It is also possible that $X_0 = X_1$ so that $X_0+X_1$ is $0$ with probability $1/2$ and $2$ with probability $1/2$. It is also possible that $X_0$ and $X_1$ are independent so that $X_0+X_1$ takes the values $0,1,2$ with probabilities $1/4,1/2,1/4$, respectively.

Nevertheless, $S = Z_{(x-y+1):x} + ... + Z_{x:x}$. Expectation is linear regardless of whether the terms are independent, so $E(S) = E(Z_{(x-y+1):x}) + ... + E(Z_{x:x})$.

I've thought about expressing the distribution of $S$, and I keep getting the same expression that techmologist did (with the corrected upper bound I edited in). Order statistics for discrete distributions are messy, so I don't expect to find a big simplification by using them.