Testing the collection of Bernoulli variables Suppose that I am throwing a triple dice (i.e, three dices simultaneously) $n$ times. Then number of times that the sum of these three dices $k=3,\ldots,18$ appears is my sample $X_1,\ldots,X_n$ is the random variable $S_k$. Now I have derived that $S_k$ follows binomial distribution $B(n,p_k)$ where $p_k$ are the probabilities derived from $p_k := P(\sum_{i=1}^{3}X^{(i)} = k)$.
Given a sample $x_1,\ldots,x_n$, how do I test that $S_k$ follows the theoretical distribution? Note that for each $k=3,\ldots,18$, I have one value $S_k$ for the whole sample.
 A: Here is a worked example in R:
First simulate the dice, say 10,000 times:
set.seed(2022)
numsides <- 6
numdice <- 3 
cases <- 10000
simsumdice <- function(sides,dice){
  sum(sample(sides, dice, replace=TRUE))
  }
sims <- replicate(cases, simsumdice(sides=numsides, dice=numdice))
results <- table(factor(sims, levels = numdice:(numsides*numdice))
results

to get
   3    4    5    6    7    8    9   10   11   12   13   14   15   16   17   18 
  43  136  287  458  668  957 1191 1269 1205 1163 1005  708  467  255  132   56

You might have expected about $\frac{10000}{216} \approx 46.3$ values of $3$ and similarly of $18$, but this is a simulation.
Then compare these simulation results to the theoretical probabilities, i.e. to the expected numbers of each sum, for example using a chi-squared test
theoryprob <- c(1,3,6,10,15,21,25,27,27,25,21,15,10,6,3,1)/216
chisq.test(results, p=theoryprob)

and see
        Chi-squared test for given probabilities

data:  results
X-squared = 10.459, df = 15, p-value = 0.7899

and this $p$-value fails to reject the hypothesis that the simulation comes from this theoretical distribution.
