Le nombre prévu de lancers de dés nécessite de faire une somme supérieure ou égale à K?

9

Un dé à 6 faces est lancé de manière itérative. Quel est le nombre attendu de rouleaux requis pour faire une somme supérieure ou égale à K?

Avant de modifier

P(Sum>=1 in exactly 1 roll)=1
P(Sum>=2 in exactly 1 roll)=5/6
P(Sum>=2 in exactly 2 rolls)=1/6
P(Sum>=3 in exactly 1 roll)=5/6
P(Sum>=3 in exactly 2 rolls)=2/6
P(Sum>=3 in exactly 3 rolls)=1/36
P(Sum>=4 in exactly 1 roll)=3/6
P(Sum>=4 in exactly 2 rolls)=3/6
P(Sum>=4 in exactly 3 rolls)=2/36
P(Sum>=4 in exactly 4 rolls)=1/216

Après modification

P(Sum>=1 in atleast 1 roll)=1
P(Sum>=2 in atleast 1 roll)=5/6
P(Sum>=2 in atleast 2 rolls)=1
P(Sum>=3 in atleast 1 roll)=4/6
P(Sum>=3 in atleast 2 rolls)=35/36
P(Sum>=3 in atleast 3 rolls)=1
P(Sum>=4 in atleast 1 roll)=3/6
P(Sum>=4 in atleast 2 rolls)=33/36
P(Sum>=4 in atleast 3 rolls)=212/216
P(Sum>=4 in atleast 4 rolls)=1

Je ne suis pas sûr que ce soit correct tout d'abord et mais je pense que cette probabilité est liée au nombre de rouleaux attendu?

Mais je ne sais pas comment continuer. Suis-je dans la bonne direction?

Suspect habituel
la source
Comment avez-vous obtenu ? P(S2 in 2 rolls)
Glen_b -Reinstate Monica
@Glen_b Vous devez obtenir un nombre inférieur à 2 dans le premier lancer qui est 1. La probabilité d'obtenir 1 est donc 1/6 et le deuxième lancer peut être n'importe quel nombre. si vous obtenez un nombre supérieur ou égal à 2 lors du premier lancer, vous n'irez pas pour un deuxième lancer.
Suspect habituel du
1
Ah, je vois ce qui se passe. Vous ne décrivez pas cela comme "P (S \ geq 2 en 2 rouleaux)"; cette expression implique que le nombre de rouleaux est fixe. Ce que vous voulez, c'est "P (exactement 2 rouleaux requis pour obtenir )" ou "P (au moins 2 rouleaux requis pour obtenir )". S 2S2S2
Glen_b -Reinstate Monica
@Glen_b Ouais c'est la confusion. P (exactement 2 rouleaux nécessaires pour obtenir S> 2) je suppose. Tout ce que je veux finalement calculer, c'est le nombre de rouleaux prévu pour atteindre une somme supérieure à K?
Suspect habituel du
@Glen_b dois-je utiliser au moins ou exactement à cette fin? Et comment calculer le nombre attendu de rouleaux pour une somme plus importante comme 10000?
Suspect habituel

Réponses:

2

Il ne s'agit pour l'instant que de quelques idées pour une autre approche, plus exacte, basée sur le même constat que ma première réponse. Avec le temps, je prolongerai cela ...

Tout d'abord, une notation. Soit un entier donné (grand) positif. Nous voulons que la distribution de , qui est le nombre minimum de lancers d'un dé ordinaire pour obtenir la somme d' au moins . Donc, nous définissons d'abord comme le résultat du lancer de dés , et . Si nous pouvons trouver la distribution de pour tout alors nous pouvons trouver la distribution de en utilisant et nous sommes terminé.N K X i i X ( n ) = X 1 + + X n X ( n ) n NKNKXiiX(n)=X1++XnX(n)nN

P(Nn)=P(X1++XnK),

Maintenant, les valeurs possibles pour sont , et pour dans cette plage, pour trouver la probabilité , nous besoin de trouver le nombre total de façons d'écrire comme une somme d'exactement entiers, tous dans la plage . Mais cela s'appelle une composition entière restreinte, un problème bien étudié en combinatoire. Quelques questions connexes sur les mathématiques SE se trouve par https://math.stackexchange.com/search?q=integer+compositions n , n + 1 , n + 2 , , 6 n k P ( X 1 + + X n = k ) k n 1 , 2 , , 6X1++Xnn,n+1,n+2,,6nkP(X1++Xn=k)kn1,2,,6

Donc, en recherchant et en étudiant cette littérature combinatoire, nous pouvons obtenir des résultats précis et silencieux. J'y reviendrai, mais plus tard ...

kjetil b halvorsen
la source
2

Il existe une formule fermée simple en termes de racines d'un polynôme de degré 6.

Il est en fait un peu plus facile d'envisager un dé de foire général avec faces étiquetées avec les nombresd21,2,,d.

Soit le nombre attendu de rouleaux nécessaires pour égaler ou dépasser Pour Sinon, l'attente est une de plus que l'attente du nombre de rouleaux pour atteindre la valeur immédiatement précédente, qui serait parmi oùekk.k0, ek=0.kd,kd+1,,k1,

(1)ek=1+1d(ekd+ekd+1++ek1).

Cette relation de récurrence linéaire a une solution sous la forme

(2)ek=2kd+1+i=1daiλik

où les sont les racines complexes du polynômeλid

(3)Td1d(Td1+Td2++T+1).

Les constantes sont trouvées en appliquant la solution aux valeurs où dans tous les cas. Cela donne un ensemble d' équations linéaires dans les constantes et il a une solution unique. Que la solution fonctionne peut être démontrée en vérifiant la récurrence utilisant le fait que chaque racine satisfaitai(2)k=(d1),(d2),,1,0ek=0dd(1)(3):

1+1dj=1dekj=1+1dj=1d(2(kj)d+1+i=1daiλikj)=2kd+1+i=1daiλikd[1d(1+λi++λid1)]=2kd+1+i=1daiλikdλid=2kd+1+i=1daiλik=ek.

Cette solution sous forme fermée nous donne de bons moyens d'approximer la réponse ainsi que de l'évaluer avec précision. (Pour les valeurs petites à modestes de l'application directe de la récurrence est une technique de calcul efficace.) Par exemple, avec nous pouvons facilement calculerk,d=6

e1000000=285714.761905

Pour les approximations, il y aura une racine unique plus grande donc finalement (pour un suffisamment grand ) le terme dominera les termes dansL'erreur diminuera de façon exponentielle selon la deuxième plus petite norme des racines. Poursuivant l'exemple avec le coefficient de est et la norme la plus petite suivante est ( passant , les autres ont tendance à être très proches de en taille.) Ainsi, nous pouvons approximer la valeur précédente commeλ+=1kλ+kd(2).k = 6 , λ + a + = 0,4761905 0,7302500. a i 1k=6,λ+a+=0.47619050.7302500.ai1

e10000002×1066+1+0.4761905=285714.761905

avec une erreur de l'ordre de0.730250010610314368.


Pour démontrer à quel point cette solution est pratique, voici du Rcode qui retourne une fonction pour évaluer pour tout (dans le cadre des calculs en virgule flottante double précision) et pas trop grand (il s'embourbera une fois ):ekkdd100

die <- function(d, mult=1, cnst=1, start=rep(0,d)) {
  # Create the companion matrix (its eigenvalues are the lambdas).
  X <- matrix(c(0,1,rep(0,d-1)),d,d+1)
  X[, d] <- mult/d
  lambda <- eigen(X[, 1:d], symmetric=FALSE, only.values=TRUE)$values

  # Find the coefficients that agree with the starting values.
  u <- 2*cnst/(d+1)
  a <- solve(t(outer(lambda, 1:d, `^`)), start - u*((1-d):0))

  # This function assumes the starting values are all real numbers.
  f <- Vectorize(function(i) Re(sum(a * lambda ^ (i+d))) + u*i)

  list(f=f, lambda=lambda, a=a, multiplier=mult, offset=cnst)
}

À titre d'exemple de son utilisation, il calcule ici les attentes pourk=1,2,,16:

round(die(6)$f(1:10), 3)

1.000 1.167 1.361 1.588 1.853 2.161 2.522 2.775 3.043 3.324 3.613 3.906 4.197 4.476 4.760 5.046

L'objet qu'il renvoie inclut les racines et leurs multiplicateurs pour une analyse plus approfondie. Le premier composant du tableau des multiplicateurs est le coefficient utileλiaia+.

(Si vous êtes curieux de savoir à quoi dieservent les autres paramètres , exécutez die(2, 2, 0, c(1,0))$f(1:10)et voyez si vous reconnaissez la sortie ;-). Cette généralisation a aidé à développer et à tester la fonction.)

whuber
la source
+1. La fonction diedonne une erreur pour moi: object 'phi' not found.
COOLSerdash
1
@COOL Merci d'avoir vérifié. Un changement de dernière minute du nom de variable (de phià a) pour correspondre au texte était le coupable. Je l'ai corrigé (et vérifié).
whuber
1

il n'y a aucun moyen d'obtenir le nombre exact de rouleaux attendus en général, mais pour un K.

Soit N l'événement du roulement attendu pour obtenir somme => K.

pour K = 1, E (N) = 1

pour K = 2,E(N)=(56+21)/(56+1)=1711

etc.

Il sera difficile d'obtenir E (N) pour les gros K. Par exemple, pour K = 20, vous devrez vous attendre à (4 rouleaux, 20 rouleaux)

Le théorème de limite centrale sera plus bénéfique avec un certain% de confiance. comme nous le savons, l'occurrence est uniformément distribuée, pour les grandes valeurs de K. (Distribution normale)

K(Sum) follows N(3.5N,35N12)

Maintenant, vous avez besoin de "N" pour obtenir Sum au moins K .... nous le convertissons en distribution normale standard. où % Vous pouvez obtenir des valeurs Z à partir des "tables normales standard" ou d'ici par exemple

K3.5N35N12=Zα
α=1confidenceZ0.01=2.31,Z0.001=2.98

Vous connaissez K, Z (à toute erreur) ........ alors vous pouvez obtenir N = E (N) à un certain% de confiance en résolvant l'équation.

Hemant Rupani
la source
2
Comment avez-vous calculé ces probabilités? Comment en êtes-vous arrivé à cette équation E (N)?
Suspect habituel
@UsualSuspect P (Sum> = 2 in 1 roll) = 5/6 (vous savez) P (Sum> = 2 in 2 rolls) = 1 (car vous devez obtenir la somme d'au moins 2 de 2 rollings) et pour E (N ) ......... c'est juste une moyenne attendue
Hemant Rupani
Désolé je ne peux pas le mentionner. Ce n'est pas au moins, exactement 2 rouleaux. J'ai compris l'équation E (N) maintenant.
Suspect habituel du
@UsualSuspect ohh! au fait si vous avez besoin de E (N) pour un K particulier, alors je peux le faire :).
Hemant Rupani
j'ai besoin de k = 20 et k = 10000. C'est mieux si vous m'expliquez plutôt que de donner directement des réponses.
Suspect habituel
0

Je vais donner une méthode pour trouver une solution approximative. Soit d'abord la variable aléatoire, "résultat du lancer avec les dés" et soit le nombre de lancers nécessaires pour atteindre une somme d'au moins . Nous avons alors que donc pour trouver la distribution de nous devons trouver les convolutions des distributions de pour , pour tout . Ces convolutions peuvent être trouvées numériquement, mais pour les grandsXiiNk

P(Nn)=P(X1+X2++Xnk)
NXii=1,2,,nnncela pourrait être beaucoup de travail, alors nous essayons plutôt d'approximer la fonction de distribution cumulative pour les convolutions, en utilisant des méthodes de point de selle. Pour un autre exemple de méthodes de point de selle, voir ma réponse à la somme générique de variables aléatoires gamma

Nous utiliserons l'approximation de Lugannini-Rice pour le cas discret, et suivons R Butler: "Saddlepoint Approximations with Applications", page 18 (deuxième correction de continuité). Tout d'abord, nous avons besoin de la fonction de génération de moment de , qui est Ensuite, la fonction génératrice de cumul pour la somme de dés indépendants devient et nous avons également besoin des premières dérivées de , mais nous trouverons celles qui utilisent symboliquement R. Le code est le suivant:Xi

M(T)=EetXi=16(et+e2t+e3t+e4t+e5t+e6t)
n
Kn(t)=nlog(16i=16eit)
K

 DD <- function(expr, name, order = 1) {
        if(order < 1) stop("'order' must be >= 1")
        if(order == 1) D(expr, name)
        else DD(D(expr, name), name, order - 1)
     }

make_cumgenfun  <-  function() {
    fun0  <-  function(n, t) n*log(mean(exp((1:6)*t)))
    fun1  <-  function(n, t) {}
    fun2  <-  function(n, t) {}
    fun3  <-  function(n, t) {}
    d1  <-  DD(expression(n*log((1/6)*(exp(t)+exp(2*t)+exp(3*t)+exp(4*t)+exp(5*t)+exp(6*t)))),  "t", 1)
    d2  <-  DD(expression(n*log((1/6)*(exp(t)+exp(2*t)+exp(3*t)+exp(4*t)+exp(5*t)+exp(6*t)))),  "t", 2)
    d3  <-  DD(expression(n*log((1/6)*(exp(t)+exp(2*t)+exp(3*t)+exp(4*t)+exp(5*t)+exp(6*t)))),  "t", 3)
    body(fun1)  <-  d1
    body(fun2)  <-  d2
    body(fun3)  <-  d3
    return(list(fun0,  fun1,  fun2,  fun3))
}

Ensuite, nous devons résoudre l'équation du point de selle.

Cela se fait par le code suivant:

funlist  <-  make_cumgenfun()

# To solve the saddlepoint equation for n,  k:
solve_speq  <-   function(n, k)  {# note that n+1 <= k <= 6n is needed
    Kd  <-  function(t) funlist[[2]](n, t)
    k  <-  k-0.5
    uniroot(function(s) Kd(s)-k,  lower=-100,  upper=1,  extendInt="upX")$root
}

Notez que le code ci-dessus n'est pas très robuste, pour les valeurs de loin dans l'une ou l'autre queue de la distribution, il ne fonctionnera pas. Ensuite, un code pour calculer réellement la fonction de probabilité de queue, approximativement, par l'approximation de Luganini-Rice, après Butler, page 18, (deuxième correction de continuité):k

Fonction de retour de la probabilité de queue:

#

Ghelp  <-  function(n, k) {
    stilde  <-  solve_speq(n, k)
    K  <-  function(t) funlist[[1]](n, t)
    Kd <-  function(t) funlist[[2]](n, t)
    Kdd <- function(t) funlist[[3]](n, t)
    Kddd <- function(t) funlist[[4]](n, t)
    w2tilde  <-  sign(stilde)*sqrt(2*(stilde*(k-0.5)-K(stilde)))  
    u2tilde  <-  2*sinh(stilde/2)*sqrt(Kdd(stilde))
    mu  <-  Kd(0)
    result  <- if (abs(mu-(k-0.5)) <= 0.001) 0.5-Kddd(0)/(6*sqrt(2*pi)*Kdd(0)^(3/2))  else
    1-pnorm(w2tilde)-dnorm(w2tilde)*(1/w2tilde - 1/u2tilde)
    return(result)
}
G  <- function(n, k) {
      fun  <- function(k) Ghelp(n, k)
      Vectorize(fun)(k)
  }

Essayons ensuite de l'utiliser pour calculer un tableau de la distribution, basé sur la formule où est la fonction du code R ci-dessus.

P(Nn)=P(X1+X2++Xnk)=1P(X1++Xnk+1)=1G(n,k+1)
G

Maintenant, laissez-nous répondre à la question d'origine avec . Ensuite, le nombre minimum de rouleaux est de 4 et le nombre maximal de rouleaux est de 20. La probabilité que 20 rouleaux soient nécessaires est très petite et peut être calculée exactement à partir de la formule binomiale, je laisse cela au lecteur. (l'approximation ci-dessus ne fonctionnera pas pour ).K=20n=20

Donc, la probabilité que soit approximée parN19

> 1-G(20, 21)
[1] 2.220446e-16

La probabilité que soit approximée par:N10

> 1-G(10, 21)
[1] 0.002880649

Etc. En utilisant tout cela, vous pouvez obtenir vous-même une approximation de l'attente. Cela devrait être bien meilleur que les approximations basées sur le théorème de la limite centrale.

kjetil b halvorsen
la source