How do you determine the probabilities of outcomes on a die? Suppose you have a die, but you do not know the probabilities of the individual sides. How can you determine the individual probabilities by rolling the dice repeatedly and independently?
 A: Throw it many times and draw the histogram. The more times you throw it the more accurately (law of the  large numbers). The resulting shape will approximate the probability mass function describing the random variable dice-side.
This is at the core of Monte carlo methods which use randomness to approximate a stochastic phenomenon.
A: Consider a $k$-sided die and suppose we observe the counts $\mathbf{n} \equiv (n_1,...,n_k)$ with $n = \sum n_i$ total rolls.  Suppose that the sequence of die-rolls is exchangeable, so that the counts of the outcome come from a multinomial distribution with probability vector $\mathbf{p} \equiv (p_1,...,p_k)$.  Under these conditions we have $\mathbf{n}/n \rightarrow \mathbf{p}$ (from the law of large numbers) and so the way we determine the probabilities is to roll the die a large number of times and take the sample proportions as estimates of the probabilities.
If you would like to do this with some formal statistical modelling, this inference problem is handled quite well by using Bayesian analysis with a Dirichlet prior for the true probabilities for your die, which leads to the Multinomial-Dirichlet model.  Since there is no prior information suggesting bias in the die, we can stipulate a prior that is exchangeable with $n_0 > 0$ "pseudo data-points" in the prior.  The model is:
$$\begin{align}
\mathbf{n}|\mathbf{p}  &\sim \text{Multinomial}(n, \mathbf{p}), \\[12pt]
\mathbf{p} &\sim \text{Dirichlet} \Big( \frac{n_0}{k} \cdot \mathbf{1} \Big). \\[6pt]
\end{align}$$
Given an observed vector of counts $\mathbf{n}$ we then have the posterior distribution:
$$\mathbf{p}|\mathbf{n} \sim \text{Dirichlet} \Big( \frac{n_0}{k} \cdot \mathbf{1} + \mathbf{n} \Big).$$
For a reasonable choice of the prior strength $n_0$, and a substantial sample size, this posterior ought to give you a reasonable inference for the probability vector $\mathbf{p}$.  The Dirichlet distribution is programming in R in various packages (see e.g., here), so it is simple to obtain the density values for any input and generate the posterior density as a function.

Example: Suppose we have a six-sided die with true underlying probabilities:
$$\mathbf{p} = (0.12, \ 0.11, \ 0.13, \ 0.12, \ 0.15, \ 0.37).$$
Here we will conduct a simulation where we roll this die $n=400$ times and derive the resulting posterior distribution, which we program as the function posterior.  (Here we use a prior with $n_0 = 10$ pseudo-data points.)
#Generate simulated data
set.seed(1)
n      <- 400
PROB   <- c(0.12, 0.11, 0.13, 0.12, 0.15, 0.37)
k      <- length(PROB)
DATA   <- sample.int(k, size = n, replace = TRUE, prob = PROB)
COUNTS <- table(DATA)

#Determine posterior parameter
n0     <- 10
PAR    <- n0*rep(1, k)/k + COUNTS

#Determine the posterior distribution
posterior <- function(x) { DirichletReg::ddirichlet(x, alpha = PAR) }

A: Let's look at a formal test of the fairness of a die:
Suppose that the die is biased so that
probabilities of faces $1$ through $6$ are $(1,2,2,2,2,3)/12.$
If you roll the die $n = 600$ times, do you have a good chance
of detecting the bias?
set.seed(121)
faces = sample(1:6, 600, rep=T, p=c(1,2,2,2,2,3))
freq = tabulate(faces);  freq
[1]  33 111  87 109 109 151

Here is a histogram of the results of 600 rolls of this die.

We do have fewer $1$s and more $6$s than we would expect
from a fair die. Does a chi-squared 'goodness-of-fit' test detect the unevenness of the results.
We are testing the null hypothesis that all faces are
equally likely (probability $1/6$ each) against the alternative that they
are not. The expected counts for each face are $E_i = 100.$
The chi-squared test statistic is $Q = \sum_{i=1}^6 \frac{(X_i - E_i)^2}{E_i},$ where $X_i$ are the observed frequencies.
[This is called a chi-squared 'goodness-of-fit' statistic. Notice that large values of $Q$ indicate a bad fit to fairness.]
If the null hypothesis is true, then $Q \stackrel{aprx}{\sim}\mathsf{Chisq}(\nu=6-1=5).$ We would reject the null hypothesis at the 5% level if $Q \ge 11.07.$
qchisq(.95, 5)
[1] 11.0705

In R, the test is shown below: $Q = 75.42 > 11.07$ (X-squared in output).
The P-value is the tiny probability of getting a more extreme result if the die were fair.
chisq.test(freq)

        Chi-squared test for given probabilities

data:  freq
X-squared = 75.42, df = 5, p-value = 7.603e-15

This one simulated rolling of the die 600 times is not just a one-time
random accident. The following simulation uses 100,000 such tests to
show that we will almost always detect such an unfair die with data
from $n=600$ rolls.
set.seed(2021)
prob = c(1,2,2,2,2,3);  n = 600
pv = replicate(10^5, 
      chisq.test(tabulate(sample(1:6,n,rep=T,p=prob)) )$p.val)
mean(pv <= .05)
[1] 0.99999

However, $n = 60$ rolls would not be sufficient for vetting such a die.
With only 60 rolls of the die, we would detect unfairness only about 35% of the time. [But $n = 300$ roll would suffice for most purposes: 98% rejection. (Not shown.)]
set.seed(121)
prob = c(1,2,2,2,2,3);  n = 60
pv = replicate(10^5, 
      chisq.test(tabulate(sample(1:6,n,rep=T,p=prob)) )$p.val)
mean(pv <= .05)
[1] 0.34642

