Using Bayesian methods, you could start with a conjugate Dirichlet prior for the probabilities of the six sides, update it with your observations, and then find the probability from the Dirichlet posterior that side five has the highest underlying probability of the six sides.
This will be affected slightly by the prior you choose and substantially by the actual observations. It may produce slightly counter-intuitive results for small numbers of observations. To take a simpler example with a biased coin,
- if you start with a uniform prior for the probability of it being heads then toss it once and see heads, the posterior probability of it being biased towards heads would be $0.75$ and towards tails $0.25$;
- if instead you tossed it $200$ times and see heads $101$ times then the posterior probability of it being biased towards heads would be about $0.556$;
- if you tossed it $200$ times and see heads $115$ times then the posterior probability of it being biased towards heads would be about $0.983$.
I do not see a simple way of doing the integration with six-sided dice to find the probability a given face is most probable, but simulation will get close enough. The following uses R and a so-called uniform Dirichlet prior for the biases, supposing you observed $21$ dice throws of $2$ ones, $3$ twos, $4$ threes, $5$ fours, $6$ fives and $1$ six:
probmostlikely <- function(obs, prior=rep(1, length(obs)),
posterior <- prior + obs
sims <- rdirichlet(cases, posterior)
table(apply(sims, 1, function(x) which(x == max(x)))) / cases
probmostlikely(c(2, 3, 4, 5, 6, 1))
# 1 2 3 4 5 6
# 0.027885 0.072102 0.152415 0.279509 0.460498 0.007591
so suggesting that the die is biased most towards five with posterior probability about $0.46$ (and most towards six with posterior probability just under $0.008$).
Seeing that pattern of observations ten times as often would increase the posterior probability that the die is biased most towards five to just under $0.82$ (and reduces those for one and six to something so small that they never appeared as most likely in a million simulations).
probmostlikely(c(20, 30, 40, 50, 60, 10))
# 2 3 4 5
# 0.000222 0.013702 0.167825 0.818251