Se lanza un dado justo 1.000 veces. ¿Cuál es la probabilidad de que salga el mismo número 5 veces seguidas? ¿Cómo se resuelve este tipo de pregunta para un número variable de lanzamientos y de repeticiones?
Respuestas
¿Demasiados anuncios?A continuación, calculamos la probabilidad de cuatro maneras:
Computation with Markov Chain 0.473981098314993
Computation with generating function 0.473981098314988
Estimation false method 0.536438013618686
Estimation correct method 0.473304632462677
Los dos primeros son métodos exactos y difieren sólo un poco (probablemente alguna ronda de errores), el tercer método es una estimación ingenua que no da el número correcto, el cuarto método es mejor y da un resultado muy cercano al método exacto.
Computacionalmente:
Cadena de Markov
Esto se puede modelar computacionalmente con un matriz de transición
Digamos que el vector columna $X_{k,j} = \lbrace x_1,x_2,x_3,x_4,x_5 \rbrace_{j}$ es la probabilidad de tener $k$ de los mismos números en una fila en el $j$ -a tirada de dados. Entonces (cuando se asume un dado de 6 caras)
$$X_{k,j} = M \cdot X_{k,j-1}$$ con
$$M = \begin{bmatrix} \frac{5}{6} & \frac{5}{6} & \frac{5}{6} & \frac{5}{6} & 0 \\ \frac{1}{6} & 0& 0 & 0 & 0 \\ 0& \frac{1}{6} & 0& 0 & 0 \\ 0 & 0& \frac{1}{6} & 0& 0 \\ 0&0 & 0& \frac{1}{6} & 1 \\ \end{bmatrix}$$
donde esta última entrada $M_{5,5} = 1$ se refiere a que 5 de los mismos seguidos son un estado de absorción en el que "paramos" el experimento.
Después de la primera tirada estará ciertamente en el estado 1 (ciertamente sólo hay 1 del mismo número en una fila).
$$X_{k,1} = \lbrace 1,0,0,0,0 \rbrace$$
Después de la $j$ -th roll esto se multiplicará con $M$ a $j-1$ veces
$$X_{k,j} = M^{j-1} \lbrace 1,0,0,0,0 \rbrace$$
Código R:
library(matrixcalc) ### allows us to use matrix.power
M <- matrix(c(5/6, 5/6, 5/6, 5/6, 0,
1/6, 0 , 0 , 0 , 0,
0, 1/6, 0 , 0 , 0,
0, 0 , 1/6, 0 , 0,
0, 0 , 0 , 1/6, 1),
5, byrow = TRUE)
start <- c(1,0,0,0,0)
matrix.power(M,999) %*% start
El resultado es $$X_{k,1000} = \begin{bmatrix} 0.438631855\\ 0.073152468\\ 0.012199943\\ 0.002034635\\ \color{red}{0.473981098}\end{bmatrix}$$
y esta última entrada 0,473981098 es la probabilidad de sacar el mismo número 5 veces seguidas en 1000 tiradas.
función generadora
Nuestra pregunta es:
- Cómo calcular la probabilidad de rodar cualquier número al menos $k$ veces seguidas, de $n$ ¿Intenta?
Esto equivale a la pregunta
- Cómo calcular la probabilidad de rodar el número 6 al menos $k-1$ veces seguidas, de $n-1$ ¿Intenta?
Puedes verlo como un seguimiento de si la tirada de dados $m$ es el mismo número que el de la tirada de dados $m-1$ (que tiene una probabilidad de 1/6). Y esto tiene que ocurrir $k-1$ veces seguidas (en nuestro caso 4 veces).
En este PREGUNTAS Y RESPUESTAS la pregunta alternativa se resuelve como un problema combinatorio: ¿De cuántas maneras podemos tirar los dados $n$ veces sin que aparezca el número "6 $k$ o más veces seguidas.
Esto se encuentra encontrando todas las combinaciones posibles de formas en que podemos combinar las cadenas "x", "x6", "x66", "x666" (donde "x" es cualquier número 1,2,3,4,5) en una cadena de longitud $n+1$ ( $n+1$ en lugar de $n$ porque en esta forma de construir cadenas la primera letra es siempre $x$ aquí). De esta manera contamos todas las posibilidades para hacer una cadena de longitud $n$ pero con sólo 1, 2 o 3 veces un 6 seguido (y no 4 o más veces).
Esas combinaciones se pueden encontrar utilizando un polinomio equivalente. Esto es muy similar al coeficientes binomiales que se relacionan con los coeficientes cuando expandimos la potencia $(x+y)^n$ pero también se relaciona con un combinación .
El polinomio es
$$\begin{array}{rcl} P(x) &=& \sum_{k=0}^\infty (5x+5x^2+5x^3+5x^4)^k\\ &=& \frac{1}{1-(5x+5x^2+5x^3+5x^4)} \\ &=& \frac{1}{1-5\frac{x-x^5}{1-x}}\\ &=& \frac{1-x}{1-6x+5x^5} \end{array}$$
El coeficiente del $x^n$ se refiere al número de maneras de ordenar los números 1,2,3,4,5,6 en una cadena de longitud $n-1$ sin 4 o más 6's seguidos. Este coeficiente se puede encontrar mediante una relación recursiva. $$P(x) (1-6x+5x^5) = 1-x$$ lo que implica que los coeficientes siguen la relación
$$a_n - 6a_{n-1} + 5 a_{n-5} = 0$$
y los primeros coeficientes se pueden calcular manualmente
$$a_1,a_2,a_3,a_4,a_5,a_6,a_7 = 5,30,180,1080,6475,38825,232800$$
Con esto, se puede calcular $a_{1000}$ y $1-a_{1000}/6^{999}$ será la probabilidad de sacar el mismo número 5 veces seguidas 5.
En el código R de abajo calculamos esto (e incluimos una división por 6 dentro de la recursión porque los números $a_{1000}$ y $6^{999}$ son demasiado grandes para calcularlas directamente). El resultado es $0.473981098314988$ , lo mismo que el cálculo con la cadena de Markov.
x <- 6/5*c(5/6,30/6^2,180/6^3,1080/6^4,6475/6^5,38825/6^6,232800/6^7)
for (i in 1:1000) {
t <- tail(x,5)
x <- c(x,(6/6*t[5]-5/6^5*t[1])) ### this adds a new number to the back of the vector x
}
1-x[1000]
Analítica/Estimación
Método 1: equivocado
Se podría pensar, que la probabilidad de tener en cualquier conjunto de 5 dados vecinos, 5 de los mismos números, es $\frac{1}{6^4} = \frac{1}{1296}$ y como hay 996 conjuntos de 5 dados vecinos la probabilidad de tener en al menos uno de estos conjuntos 5 dados iguales es:
$$ 1-(1-\frac{1}{6^4})^{996} \approx 0.536$$
Pero esto está mal. La razón es que los conjuntos 996 se superponen y no son independientes.
Método 2: correcto
Una forma mejor es aproximar la cadena de Markov que hemos calculado anteriormente. Al cabo de un tiempo conseguirás que la ocupación de los estados, con 1,2,3,4 del mismo número en una fila, sean más o menos estables y los cocientes serán aproximadamente $1/6,1/6^2,1/6^3,1/6^4$ (*). Por lo tanto, la fracción de tiempo que tenemos 4 en una fila es:
$$\text{frequency 4 in a row} = \frac{1/6^4}{1/6+1/6^2+1/6^3+1/6^4}$$
Si tenemos estos 4 seguidos entonces tenemos una probabilidad de 1/6 de terminar el juego. Así que la frecuencia de terminar el juego es
$$\text{finish-rate} = \frac{1}{6} \text{frequency 4 in a row} = \frac{1}{1554}$$
y la probabilidad de terminar después de $k$ pasos es aproximadamente
$$P_k \approx 1-(1-\frac{1}{1554})^{k-4} \underbrace{\approx 0.47330}_{\text{if $ k=1000 $}}$$
mucho más cerca del cálculo exacto.
(*) La ocupación en el estado $k$ durante el rodaje $j$ se relacionará con la ocupación en el estado $k-1$ durante el rodaje $j-1$ . Tendremos $x_{k,j} = \frac{1}{6} x_{k-1,j-1} \approx \frac{1}{6} x_{k-1,j}$ . Tenga en cuenta que esto requiere que usted tenga $x_{k-1,j} \approx x_{k-1,j-1}$ que se produce cuando la tasa de acabado es pequeña. Si no es el caso, se podría aplicar un factor para compensar, pero la suposición de una relación relativamente estable también será errónea.
Problemas relacionados
- Distribución límite asociada a los recuentos (problema combinatorio no trivial)
- Comprobar si una moneda es justa en función de la frecuencia de una subsecuencia
- ¿Cuál es la probabilidad de que salgan todas las caras de un dado después de un número n de tiradas?
- Probabilidad de una subsecuencia similar de longitud X en dos secuencias de longitud Y y Z
Este último problema relacionado da una aproximación diferente basada en los valores de las expectativas y estima la distribución como una distribución de Poisson sobredispersa. Dando una aproximación $1- \exp \left(-(1000-5+1)\left(\frac{1}{6^4}\right) /1.2 \right)\approx 0.4729354$ que tampoco está mal.
He obtenido un resultado diferente a la respuesta aceptada y me gustaría saber en qué me he equivocado.
Asumí un dado justo, de 6 caras, y simulé 1000 tiradas de 1000 tiros cada una. Cuando el resultado de una tirada coincide con los resultados de las 4 tiradas anteriores, se pone una bandera a TRUE. La media de esta columna de banderas y la media de las tiradas es entonces reportada. Obtengo ~0,07% como probabilidad de ver 5 tiradas seguidas del mismo número.
En R,
tibble(
run = rep(seq(1:1000), each = 1000),
roll = rep(seq(1:1000), 1000),
x = sample(1:6, 1000000, replace = T)
) %>%
group_by(run) %>%
mutate(
same_five = x == lag(x, 1) & x == lag(x, 2) & x == lag(x, 3) & x == lag(x, 4)
) %>%
summarize(
p_same_five = mean(same_five, na.rm = TRUE), .groups = "drop"
) %>%
summarize(mean(p_same_five)) * 100
mean(p_same_five)
1 0.07208702