Número esperado en el que estaré después de robar cartas hasta que obtenga un as, 2, 3, etc.

12

Tengo algunos problemas para resolver lo siguiente.

Robas cartas de una baraja estándar de 52 cartas sin reemplazo hasta que obtengas un as. Sacas de lo que queda hasta que obtienes un 2. Continúas con 3. ¿Cuál es el número esperado en el que estarás después de que se acabe todo el mazo?

Era natural dejar

  • Ti=first position of card whose value is i
  • Ui=last position of card whose value is i

Entonces, el problema esencialmente equivale a calcular la probabilidad de que esté en k cuando se acabe el mazo, a saber:

Pr(T1<<TkUk+1<Tk)

Puedo ver eso

Pr(T1<<Tk)=1/k!andPr(Uk+1<Tk)=1/70

pero no pudo llegar más lejos ...

cuenta
fuente
1
¿Qué sucede si ya has sacado todos los s para cuando sacas tu primer as? 2
gung - Restablece a Monica
¿El número "esperado" realmente significa el número "más probable"?
Whuber
Este es un problema interesante, pero no estoy seguro acerca de las matemáticas que escribe después de "el problema esencialmente es". En la primera declaración, ¿ escribir ∩ en lugar de ? Incluso entonces, sin embargo, no estoy seguro de que la declaración sea correcta. Considere un comienzo de secuencia . Tenemos y, por lo tanto, , pero si entiendo correctamente la descripción del texto, ¿aún podemos elegir el As en la segunda posición y luego el 2 en la quinta posición? ¿Y por lo tanto no es una condición necesaria? 2AAA2T 1 > T 2 T 1 < T 2T1=2,T2=1T1>T2T1<T2
TooTone
@TooTone Oh, quise decir como dijiste, y tienes razón; no es una condición necesaria ...T 1 < T 2T1<T2
factura
@gung En ese caso, tu mazo se acabará y todavía estarás en 2.
factura

Respuestas:

0

Siguiendo la idea de @ Gung, ¿creo que el valor esperado sería 5.84? y según mi interpretación de los comentarios, supongo que "A" es un valor casi imposible (a menos que las últimas cuatro cartas del mazo sean todos ases). Aquí están los resultados de una simulación de Monte Carlo de 100,000 iteraciones

results
    2     3     4     5     6     7     8     9     J     K     Q     T 
 1406  7740 16309 21241 19998 15127  9393  4906   976   190   380  2334 

y aquí está el código R en caso de que quieras jugar con él ...

# monte carlo card-drawing functions from here
# http://streaming.stat.iastate.edu/workshops/r-intro/lectures/5-Rprogramming.pdf

# create a straightforward deck of cards
create_deck <-
    function( ){
        suit <- c( "H" , "C" , "D" , "S" )
        rank <- c( "A" , 2:9 , "T" , "J" , "Q" , "K" )
        deck <- NULL
        for ( r in rank ) deck <- c( deck , paste( r , suit ) )
        deck
    }

# construct a function to shuffle everything
shuffle <- function( deck ){ sample( deck , length( deck ) ) }

# draw one card at a time
draw_cards <-
    function( deck , start , n = 1 ){
        cards <- NULL

        for ( i in start:( start + n - 1 ) ){
            if ( i <= length( deck ) ){
                cards <- c( cards , deck[ i ] )
            }
        }

        return( cards )
    }

# create an empty vector for your results
results <- NULL

# run your simulation this many times..
for ( i in seq( 100000 ) ){
    # create a new deck
    sdeck <- shuffle( create_deck() )

    d <- sdeck[ grep('A|2' , sdeck ) ]
    e <- identical( grep( "2" , d ) , 1:4 )

    # loop through ranks in this order
    rank <- c( "A" , 2:9 , "T" , "J" , "Q" , "K" )

    # start at this position
    card.position <- 0

    # start with a blank current.draw
    current.draw <- ""

    # start with a blank current rank
    this.rank <- NULL

    # start with the first rank
    rank.position <- 1

    # keep drawing until you find the rank you wanted
    while( card.position < 52 ){

        # increase the position by one every time
        card.position <- card.position + 1

        # store the current draw for testing next time
        current.draw <- draw_cards( sdeck , card.position )

        # if you draw the current rank, move to the next.
        if ( grepl( rank[ rank.position ] , current.draw ) ) rank.position <- rank.position + 1

        # if you have gone through every rank and are still not out of cards,
        # should it still be a king?  this assumes yes.
        if ( rank.position == length( rank ) ) break        

    }

    # store the rank for this iteration.
    this.rank <- rank[ rank.position ]

    # at the end of the iteration, store the result
    results <- c( results , this.rank )

}

# print the final results
table( results )

# make A, T, J, Q, K numerics
results[ results == 'A' ] <- 1
results[ results == 'T' ] <- 10
results[ results == 'J' ] <- 11
results[ results == 'Q' ] <- 12
results[ results == 'K' ] <- 13
results <- as.numeric( results )

# and here's your expected value after 100,000 simulations.
mean( results )
Anthony Damico
fuente
¿Por qué es Aimposible? Considere la secuencia de 48 cartas seguidas por, AAAApor ejemplo.
TooTone
tienes razón ... es uno de 270725 - o con código R1/prod( 48:1 / 52:5 )
Anthony Damico
1
Esta respuesta es incorrecta. Considere la cuenta para "2": como esto puede resultar solo cuando se encuentran todos los 2 antes que cualquiera de los 1, su probabilidad es uno en cada y, por lo tanto, su expectativa en su simulación es con un error estándar de . Su rendimiento de es más de seis errores estándar demasiado altos, por lo que es casi seguro que es erróneo. Un valor exacto para la media (basado en una simulación diferente con iteraciones) es . 105/ ( 8(84)=7037.516601065.833±0.004105/(84)1428.637.516601065.833±0.004
Whuber
1
Desafortunadamente, su código muy documentado es varias veces más largo y más lento de lo necesario. Demostré que su salida es incorrecta; aunque desearía tener tiempo para depurar su código, no lo hago y no es mi tarea hacer eso. Mi argumento es el siguiente: seguirá trabajando en "2" al final si y solo si todos los "2" preceden a todos los "A". Entre las formas igualmente probables de organizar las cuatro "2" sy las cuatro "A", exactamente una de ellas satisface este criterio. Por lo tanto, su valor bajo el encabezado "2" debe estar cerca de , pero no lo es. 105/70=1,429(4+44)=70results105/70=1429
whuber
1
Incluso los moderadores no pueden eliminar los votos de otras personas :-). Una prueba de chi-cuadrado ahora sugiere que sus resultados concuerdan con los míos, pero sería bueno saber cómo probó su simulación, porque eso mejoraría la confianza en su respuesta. De hecho, según una edición que hiciste en el primer párrafo de tu respuesta, ahora ambos resultados están equivocados: como he interpretado tu pregunta, nunca es posible seguir trabajando en un as cuando todas las cartas están agotadas.
Whuber
7

Para una simulación es crucial ser correcto y rápido. Ambos objetivos sugieren escribir código que apunte a las capacidades centrales del entorno de programación, así como un código que sea lo más breve y simple posible, porque la simplicidad brinda claridad y la claridad promueve la corrección. Aquí está mi intento de lograr ambos en R:

#
# Simulate one play with a deck of `n` distinct cards in `k` suits.
#
sim <- function(n=13, k=4) {
  deck <- sample(rep(1:n, k)) # Shuffle the deck
  deck <- c(deck, 1:n)        # Add sentinels to terminate the loop
  k <- 0                      # Count the cards searched for
  for (j in 1:n) {
    k <- k+1                          # Count this card
    deck <- deck[-(1:match(j, deck))] # Deal cards until `j` is found
    if (length(deck) < n) break       # Stop when sentinels are reached
  }
  return(k)                   # Return the number of cards searched
}

Se puede aplicar esto de forma reproducible con la replicatefunción después de establecer la semilla de número aleatorio, como en

> set.seed(17);  system.time(d <- replicate(10^5, sim(13, 4)))
   user  system elapsed 
   5.46    0.00    5.46

Eso es lento, pero lo suficientemente rápido como para realizar simulaciones bastante largas (y por lo tanto precisas) repetidamente sin esperar. Hay varias formas en que podemos exhibir el resultado. Comencemos con su significado:

> n <- length(d)
> mean(d)
[1] 5.83488

> sd(d) / sqrt(n)
[1] 0.005978956

El último es el error estándar: esperamos que la media simulada esté dentro de dos o tres SE del valor verdadero. Eso coloca la verdadera expectativa en algún lugar entre y5.8535.8175.853 .

También podríamos querer ver una tabulación de las frecuencias (y sus errores estándar). El siguiente código embellece un poco la tabulación:

u <- table(d)
u.se <- sqrt(u/n * (1-u/n)) / sqrt(n)
cards <- c("A", "2", "3", "4", "5", "6", "7", "8", "9", "T", "J", "Q", "K")
dimnames(u) <- list(sapply(dimnames(u), function(x) cards[as.integer(x)]))
print(rbind(frequency=u/n, SE=u.se), digits=2)

Aquí está la salida:

                2       3      4      5      6      7       8       9       T       J       Q       K
frequency 0.01453 0.07795 0.1637 0.2104 0.1995 0.1509 0.09534 0.04995 0.02249 0.01009 0.00345 0.00173
SE        0.00038 0.00085 0.0012 0.0013 0.0013 0.0011 0.00093 0.00069 0.00047 0.00032 0.00019 0.00013

¿Cómo podemos saber que la simulación es correcta? Una forma es probarlo exhaustivamente para detectar problemas más pequeños. Por esa razón, este código fue escrito para atacar una pequeña generalización del problema, reemplazando cartas distintas con y palos con . Sin embargo, para la prueba es importante poder alimentar el código de una plataforma en un orden predeterminado. Escribamos una interfaz ligeramente diferente para el mismo algoritmo:413n4k

draw <- function(deck) {
  n <- length(sentinels <- sort(unique(deck)))
  deck <- c(deck, sentinels)
  k <- 0
  for (j in sentinels) {
    k <- k+1
    deck <- deck[-(1:match(j, deck))]
    if (length(deck) < n) break
  }
  return(k)
}

(Es posible usarlo drawen lugar de en simtodas partes, pero el trabajo adicional realizado al principio drawhace que sea el doble de lento que sim).

Podemos usar esto aplicándolo a cada barajado distinto de un mazo dado. Dado que el propósito aquí es solo unas pocas pruebas puntuales, la eficiencia en la generación de esas mezclas no es importante. Aquí hay una forma rápida de fuerza bruta:

n <- 4 # Distinct cards
k <- 2 # Number of suits
d <- expand.grid(lapply(1:(n*k), function(i) 1:n))
e <- apply(d, 1, function(x) var(tabulate(x))==0)
g <- apply(d, 1, function(x) length(unique(x))==n)
d <- d[e & g,]

Ahora des un marco de datos cuyas filas contienen todas las barajas. Aplicar drawa cada fila y contar los resultados:

d$result <- apply(as.matrix(d), 1, draw)
    (counts <- table(d$result))

La salida (que usaremos en una prueba formal momentáneamente) es

   2    3    4 
 420  784 1316 

(El valor de es fácil de entender, por cierto: todavía estaríamos trabajando en la tarjeta si y solo si todos los dos precedieron a todos los ases. La posibilidad de que esto suceda (con dos palos) es . De los shuffles distintos, tienen esta propiedad).2 1 / ( 2 + 2420225202520/6=4201/(2+22)=1/625202520/6=420

Podemos probar la salida con una prueba de chi-cuadrado. Con este fin, aplico veces a este caso de cartas distintas en palos:sim n = 4 k = 210,000n=4k=2

>set.seed(17)
>d.sim <- replicate(10^4, sim(n, k))
>print((rbind(table(d.sim) / length(d.sim), counts / dim(d)[1])), digits=3)

         2     3     4
[1,] 0.168 0.312 0.520
[2,] 0.167 0.311 0.522

> chisq.test(table(d.sim), p=counts / dim(d)[1])

    Chi-squared test for given probabilities

data:  table(d.sim) 
X-squared = 0.2129, df = 2, p-value = 0.899

Debido a que es tan alto, no encontramos diferencias significativas entre lo que dice y los valores calculados por enumeración exhaustiva. La repetición de este ejercicio para otros valores (pequeños) de y produce resultados comparables, lo que nos da una amplia razón para confiar cuando se aplica a y .n k n = 13psimnksimn=13k=4

Finalmente, una prueba de chi-cuadrado de dos muestras comparará la salida de simcon la salida informada en otra respuesta:

>y <- c(1660,8414,16973,21495,20021,14549,8957,4546,2087,828,313,109)
>chisq.test(cbind(u, y))

data:  cbind(u, y) 
X-squared = 142.2489, df = 11, p-value < 2.2e-16

La enorme estadística de chi cuadrado produce un valor p que es esencialmente cero: sin duda, simno está de acuerdo con la otra respuesta. Hay dos posibles resoluciones del desacuerdo: una (¡o ambas!) De estas respuestas es incorrecta o implementan diferentes interpretaciones de la pregunta. Por ejemplo, he interpretado que "después de que se acaba el mazo" significa que después de observar la última carta y, si es posible, actualizar el "número en el que estará" antes de finalizar el procedimiento. Es concebible que el último paso no esté destinado a ser dado. Quizás alguna diferencia de interpretación tan sutil explicará el desacuerdo, en cuyo punto podemos modificar la pregunta para aclarar lo que se está preguntando.

whuber
fuente
4

Hay una respuesta exacta (en forma de un producto matricial, presentado en el punto 4 a continuación). Existe un algoritmo razonablemente eficiente para calcularlo, derivado de estas observaciones:

  1. Se puede generar una combinación aleatoria de cartas barajando aleatoriamente tarjetas y luego intercalando aleatoriamente las cartas restantes dentro de ellas.N kN+kNk

  2. Al mezclar solo los ases, y luego (aplicando la primera observación) intercalando los dos, luego los tres, etc., este problema puede verse como una cadena de trece pasos.

  3. Necesitamos hacer un seguimiento de más del valor de la tarjeta que estamos buscando. Sin embargo, al hacer esto, no necesitamos tener en cuenta la posición de la marca en relación con todas las cartas, sino solo su posición en relación con las cartas de igual o menor valor.

    Imagina colocar una marca en el primer as, y luego marcar los dos primeros encontrados después, y así sucesivamente. (Si en algún momento el mazo se agota sin mostrar la carta que estamos buscando actualmente, dejaremos todas las cartas sin marcar). Deje que el "lugar" de cada marca (cuando exista) sea el número de cartas de igual o menor valor que se repartieron cuando se realizó la marca (incluida la propia tarjeta marcada). Los lugares contienen toda la información esencial.

  4. El lugar después de la marca es un número aleatorio. Para un mazo dado, la secuencia de estos lugares forma un proceso estocástico. De hecho, es un proceso de Markov (con matriz de transición variable). Por lo tanto, se puede calcular una respuesta exacta a partir de doce multiplicaciones matriciales.ith

Usando estas ideas, esta máquina obtiene un valor de (computación en coma flotante de doble precisión) en segundo. Esta aproximación del valor exacto es precisa para todos los dígitos mostrados.1 / 9 19826005792658947850269453319689390235225425695.83258855290199651/9

1982600579265894785026945331968939023522542569339917784579447928182134345929899510000000000

El resto de esta publicación proporciona detalles, presenta una implementación funcional (en R) y concluye con algunos comentarios sobre la pregunta y la eficiencia de la solución.


Generación aleatoria de barajas

En realidad, es más claro conceptualmente y matemáticamente no más complicado considerar un "mazo" (también conocido como multiset ) de cartas de las cuales hay de la denominación más baja, de la siguiente más baja, y así sucesivamente . (La pregunta formulada se refiere al mazo determinado por el vector .N=k1+k2++kmk1k213(4,4,,4)

¡Una "combinación aleatoria" de cartas es una permutación tomada uniforme y aleatoriamente de la permutaciones de las cartas. Estas mezclas se agrupan en grupos de configuraciones equivalentes porque permutar los "ases" entre ellos no cambia nada, permutar los "dos" entre ellos tampoco cambia nada, y así sucesivamente. Por lo tanto, cada grupo de permutaciones que se ven idénticas cuando se ignoran los palos de las cartas contienepermutaciones Estos grupos, cuyo número viene dado por el coeficiente multinomialNN!=N×(N1)××2×1Nk1k2k1!×k2!××km!

(Nk1,k2,,km)=N!k1!k2!km!,

se llaman "combinaciones" de la baraja.

Hay otra forma de contar las combinaciones. ¡Las primeras cartas solo pueden formar combinación. Dejan "espacios" entre y alrededor de ellos en los que se pueden colocar las siguientes cartas. Podríamos indicar esto con un diagrama donde " " designa una de las tarjetas y " " designa una ranura que puede contener entre y tarjetas adicionales:k1k1!/k1!=1k1+1k2k1_0k2

_____k1 stars

Cuando se tarjetas adicionales, el patrón de estrellas y nuevas tarjetas divide las tarjetas en dos subconjuntos. El número de subconjuntos distintos es .k2k1+k2(k1+k2k1,k2)=(k1+k2)!k1!k2!

Repitiendo este procedimiento con "tres", encontramos que hay formas de intercalarlos entre las primeras cartas. Por lo tanto, el número total de formas distintas de organizar las primeras tarjetas de esta manera es igualk3((k1+k2)+k3k1+k2,k3)=(k1+k2+k3)!(k1+k2)!k3!k1+k2k1+k2+k3

1×(k1+k2)!k1!k2!×(k1+k2+k3)!(k1+k2)!k3!=(k1+k2+k3)!k1!k2!k3!.

Después de terminar las últimas cartas y continuar multiplicando estas fracciones telescópicas, encontramos que el número de combinaciones distintas obtenidas es igual al número total de combinaciones contadas previamente, . Por lo tanto, no hemos pasado por alto ninguna combinación. Eso significa que este proceso secuencial de barajar las cartas captura correctamente las probabilidades de cada combinación, suponiendo que en cada etapa cada forma distinta posible de intercalar las nuevas cartas entre las viejas se tome con una probabilidad uniformemente igual.kn(Nk1,k2,,km)

El proceso del lugar

Inicialmente, hay ases y obviamente el primero está marcado. En etapas posteriores hay tarjetas, el lugar (si existe una tarjeta marcada) es igual a (algún valor de a ), y estamos a punto de intercalar cartas a su alrededor. Podemos visualizar esto con un diagrama comok1n=k1+k2++kj1p1nk=kj

_____p1 stars____np stars

donde " " designa el símbolo marcado actualmente. Condicional a este valor del lugar , deseamos encontrar la probabilidad de que el próximo lugar sea igual a (algún valor de a ; según las reglas del juego, el siguiente lugar debe venir después de , de donde ). Si podemos encontrar cuántas formas hay de intercalar las nuevas tarjetas en los espacios en blanco para que el siguiente lugar sea igual a , entonces podemos dividir por el número total de formas de intercalar estas cartas (igual a , como hemos visto) para obtener elpq1n+kpqp+1kq(n+kk)probabilidad de transición de que el lugar cambie de a . (También habrá una probabilidad de transición para que el lugar desaparezca por completo cuando ninguna de las nuevas cartas siga a la carta marcada, pero no hay necesidad de calcular esto explícitamente).pq

Actualicemos el diagrama para reflejar esta situación:

_____p1 starss stars | ____nps stars

La barra vertical " " muestra dónde aparece la primera carta nueva después de la carta marcada: por lo tanto, no pueden aparecer cartas nuevas entre el y el (y, por lo tanto, no se muestran espacios en ese intervalo). No sabemos cuántas estrellas hay en este intervalo, así que lo acabo de llamar (que puede ser cero). Lo desconocido desaparecerá una vez que encontremos la relación entre él y .||ssq

Supongamos, entonces, que intercalamos nuevas cartas alrededor de las estrellas antes del y luego, independientemente de eso, intercalamos las nuevas cartas restantes alrededor de las estrellas después del . Existenjkj1|

τn,k(s,p)=((p1)+jj)((nps)+(kj)1kj1)

maneras de hacer esto. Tenga en cuenta, sin embargo, esta es la parte más complicada del análisis, que el lugar de es igual a porque|p+s+j+1

  • Hay "viejas" tarjetas en o antes de la marca.p
  • Hay tarjetas antiguas después de la marca pero antes de .s|
  • Hay cartas nuevas antes de la marca.j
  • Existe la nueva tarjeta representada por misma.|

Por lo tanto, nos da información sobre la transición del lugar al lugar . Cuando rastreamos esta información cuidadosamente para todos los valores posibles de , y sumamos todas estas posibilidades (disjuntas), obtenemos la probabilidad condicional del lugar después del lugar ,τn,k(s,p)pq=p+s+j+1sqp

Prn,k(q|p)=(j(p1+jj)(n+kqkj1))/(n+kk)

donde la suma comienza en y termina en . (La longitud variable de esta suma sugiere que hay es poco probable que sea una fórmula cerrada en función de y , excepto en casos especiales).j=max(0,q(n+1))j=min(k1,q(p+1)n,k,q,p

El algoritmo

Inicialmente, existe la probabilidad que el lugar sea y la probabilidad que tenga cualquier otro valor posible en . Esto puede ser representado por un vector .1102,3,,k1p1=(1,0,,0)

Después de intercalar las siguientes tarjetas , el vector se actualiza a multiplicándolo (a la izquierda) por la matriz de transición . Esto se repite hasta que se hayan colocado tarjetas . En cada etapa , la suma de las entradas en el vector de probabilidad es la posibilidad de que se haya marcado alguna tarjeta. Lo que quede para hacer que el valor sea igual a por lo tanto, existe la posibilidad de que no quede ninguna tarjeta marcada después del pasok2p1p2(Prk1,k2(q|p),1pk1,1qk2)k1+k2++kmjpj1j. Por lo tanto, las diferencias sucesivas en estos valores nos dan la probabilidad de que no podamos encontrar una carta de tipo para marcar: esa es la distribución de probabilidad del valor de la carta que estábamos buscando cuando el mazo se agota al final del juego .j


Implementación

El siguiente Rcódigo implementa el algoritmo. Es paralela a la discusión anterior. Primero, el cálculo de las probabilidades de transición se realiza mediante t.matrix(sin normalización con la división por , lo que facilita el seguimiento de los cálculos al probar el código):(n+kk)

t.matrix <- function(q, p, n, k) {
  j <- max(0, q-(n+1)):min(k-1, q-(p+1))
  return (sum(choose(p-1+j,j) * choose(n+k-q, k-1-j))
}

Esto lo utiliza transitionpara actualizar a . Calcula la matriz de transición y realiza la multiplicación. También se encarga de calcular el vector inicial si el argumento es un vector vacío:pj1pjp1p

#
# `p` is the place distribution: p[i] is the chance the place is `i`.
#
transition <- function(p, k) {
  n <- length(p)
  if (n==0) {
    q <- c(1, rep(0, k-1))
  } else {
    #
    # Construct the transition matrix.
    #
    t.mat <- matrix(0, nrow=n, ncol=(n+k))
    #dimnames(t.mat) <- list(p=1:n, q=1:(n+k))
    for (i in 1:n) {
      t.mat[i, ] <- c(rep(0, i), sapply((i+1):(n+k), 
                                        function(q) t.matrix(q, i, n, k)))
    }
    #
    # Normalize and apply the transition matrix.
    #
    q <- as.vector(p %*% t.mat / choose(n+k, k))
  }
  names(q) <- 1:(n+k)
  return (q)
}

Ahora podemos calcular fácilmente las probabilidades sin marca en cada etapa para cualquier mazo:

#
# `k` is an array giving the numbers of each card in order;
# e.g., k = rep(4, 13) for a standard deck.
#
# NB: the *complements* of the p-vectors are output.
#
game <- function(k) {
  p <- numeric(0)
  q <- sapply(k, function(i) 1 - sum(p <<- transition(p, i)))
  names(q) <- names(k)
  return (q)
}

Aquí están para el mazo estándar:

k <- rep(4, 13)
names(k) <- c("A", 2:9, "T", "J", "Q", "K")
(g <- game(k))

La salida es

         A          2          3          4          5          6          7          8          9          T          J          Q          K 
0.00000000 0.01428571 0.09232323 0.25595013 0.46786622 0.66819134 0.81821790 0.91160622 0.96146102 0.98479430 0.99452614 0.99818922 0.99944610

Según las reglas, si se marcara un rey, no buscaríamos más cartas: esto significa que el valor de debe aumentarse a . Al hacerlo, las diferencias dan la distribución del "número en el que se encontrará cuando se acabe el mazo":0.99944611

> g[13] <- 1; diff(g)
          2           3           4           5           6           7           8           9           T           J           Q           K 
0.014285714 0.078037518 0.163626897 0.211916093 0.200325120 0.150026562 0.093388313 0.049854807 0.023333275 0.009731843 0.003663077 0.001810781

(Compare esto con el resultado que informo en una respuesta separada que describe una simulación de Monte-Carlo: parecen ser los mismos, hasta las cantidades esperadas de variación aleatoria).

El valor esperado es inmediato:

> sum(diff(g) * 2:13)
[1] 5.832589

En total, esto requirió solo una docena de líneas de código ejecutable. Lo he comparado con cálculos manuales para valores pequeños de (hasta ). Por lo tanto, si se observa alguna discrepancia entre el código y el análisis anterior del problema, confíe en el código (porque el análisis puede tener errores tipográficos).k3


Observaciones

Relaciones con otras secuencias.

Cuando hay una de cada tarjeta, la distribución es una secuencia de recíprocos de números enteros:

> 1/diff(game(rep(1,10)))
[1]      2      3      8     30    144    840   5760  45360 403200

El valor en el lugar es(comenzando en el lugar ). Esta es la secuencia A001048 en la Enciclopedia en línea de secuencias enteras . En consecuencia, podríamos esperar una fórmula cerrada para los mazos con constante (los mazos "adecuados") que generalizaría esta secuencia, que tiene algunos significados profundos. (Por ejemplo, cuenta los tamaños de las clases de conjugación más grandes en los grupos de permutación y también está relacionado con los coeficientes trinomiales ). (Desafortunadamente, los recíprocos en la generalización para no suelen ser enteros).ii!+(i1)!i=1kik>1

El juego como proceso estocástico

Nuestro análisis deja en claro que los coeficientes iniciales de los vectores , , son constantes. Por ejemplo, rastreemos la salida de mientras procesa cada grupo de tarjetas:ipjjigame

> sapply(1:13, function(i) game(rep(4,i)))

[[1]]
[1] 0

[[2]]
[1] 0.00000000 0.01428571

[[3]]
[1] 0.00000000 0.01428571 0.09232323

[[4]]
[1] 0.00000000 0.01428571 0.09232323 0.25595013

...

[[13]]
 [1] 0.00000000 0.01428571 0.09232323 0.25595013 0.46786622 0.66819134 0.81821790 0.91160622 0.96146102 0.98479430 0.99452614 0.99818922 0.99944610

Por ejemplo, el segundo valor del vector final (que describe los resultados con un mazo completo de 52 cartas) ya apareció después de que se procesó el segundo grupo (y es igual a ). Por lo tanto, si desea información solo sobre las marcas a través del valor de la tarjeta , solo tiene que realizar el cálculo para un mazo de cartas .1/(84)=1/70jthk1+k2++kj

Debido a que la posibilidad de no marcar una carta de valor se acerca rápidamente a medida que aumenta, después de tipos de cartas en cuatro palos casi hemos alcanzado un valor límite para la expectativa. De hecho, el valor límite es de aproximadamente (calculado para un mazo de cartas, en cuyo punto el error de redondeo de doble precisión evita ir más allá).j1j135.8333554×32

Sincronización

Al observar el algoritmo aplicado al vector , vemos que su sincronización debe ser proporcional a y, usando un límite superior bruto, no es peor que proporcional a . Por temporización todos los cálculos para a través de y a través de , y el análisis de solamente aquellos que toman tiempos relativamente largos ( segundo o más), calculo el tiempo de cálculo es de aproximadamente , apoyando esta evaluación de límite superior.( k , k , ... , k ) k 2 m 3 k = 1 7 n = 10 30 1 / 2 O ( k 2 n 2.9 )m(k,k,,k)k2m3k=17n=10301/2O(k2n2.9)

Un uso de estas asintóticas es proyectar tiempos de cálculo para problemas mayores. Por ejemplo, viendo que el caso toma aproximadamente segundos, estimaríamos que el caso (muy interesante) tomaría aproximadamente segundos. (En realidad, toma segundos).1,31 k = 1 , n = 100 1,31 ( 1 / 4 ) 2 ( 100 / 30 ) 2,92,7 2,87k=4,n=301.31k=1,n=1001.31(1/4)2(100/30)2.92.72.87

whuber
fuente
0

Hackeó un simple Monte Carlo en Perl y encontró aproximadamente .5.8329

#!/usr/bin/perl

use strict;

my @deck = (1..13) x 4;

my $N = 100000; # Monte Carlo iterations.

my $mean = 0;

for (my $i = 1; $i <= $N; $i++) {
    my @d = @deck;
    fisher_yates_shuffle(\@d);
    my $last = 0;
        foreach my $c (@d) {
        if ($c == $last + 1) { $last = $c }
    }
    $mean += ($last + 1) / $N;
}

print $mean, "\n";

sub fisher_yates_shuffle {
    my $array = shift;
        my $i = @$array;
        while (--$i) {
        my $j = int rand($i + 1);
        @$array[$i, $j] = @$array[$j, $i];
    }
}
zen
fuente
Dada la aguda discrepancia entre esta y todas las respuestas anteriores, incluidas dos simulaciones y una teórica (exacta), sospecho que está interpretando la pregunta de una manera diferente. En ausencia de cualquier explicación de su parte, solo tenemos que tomarlo como incorrecto. (Sospecho que puede estar contando uno menos, en cuyo caso su 4.8 debería compararse con 5.83258 ...; pero aun así, sus dos dígitos significativos de precisión no proporcionan información adicional sobre este problema.)
whuber
1
¡Sí! Hubo un error off-by-one.
Zen