Eräs noppapeli

Pelaaja saa aluksi yhden nopan. Joka kierros pelaaja heittää kaikki noppansa ja menettää ne joiden osoittamaa arvoa tuli useampi kuin yksi. Kierroksen päätteeksi pelaaja saa aina yhden uuden nopan. Tavoite on saada seitsemän noppaa (tai yleisemmin d+1 noppaa, kun nopassa on d tahkoa). Toisin sanoen peli päättyy, kun pelaajalla on kuusi noppaa ja niillä tulee jokaisella eri arvo. Mikä on yhden pelin kierrosten määrän (N) odotusarvo?

Ratkaisu

Peli on Markovin ketju, jossa tilat ovat noppien määrät 1, 2, \dots, d+1. Meidän täytyy selvittää siirtymätodennäköisyydet näiden välillä. Pienille arvoille, kuten d = 6, tämä voidaan tehdä tutkimalla kaikki mahdolliset heitot ja kirjaamalla kuinka moni niistä säilyttää kuinkakin monta noppaa. Seuraava Sage-koodi tekee tämän:

from itertools import product
from collections import Counter

def getTransitionsProbsBrute1(diceAmounts, dieFaces):
    n = dieFaces+1
    ret = [0]*n
    for dice in product(range(dieFaces), repeat=diceAmounts):
        keep = sum(x for x in Counter(dice).values() if x==1)
        #we're adding 1 die after round, but 1 is the
        #first state, so keep is the correct index
        ret[keep] += 1
    tot = dieFaces**diceAmounts
    return [x/tot for x in ret]

Laskentaa voidaan hieman nopeuttaa käyttämällä kombinaatioita (takaisin panolla):

from itertools import combinations_with_replacement as cmb_rep
from collections import Counter

def getType(dice):
    return list(Counter(dice).values())

def getTransitionsProbsBrute2(diceAmount, dieFaces):
    n = dieFaces+1
    ret = [0]*n
    for dice in cmb_rep(range(dieFaces), r=diceAmount):
        keep = sum(x for x in Counter(dice).values() if x==1)
        #multinomial of the type gives the amount
        #of orderings this corresponds to when account for all
        weight = multinomial(getType(dice))
        ret[keep] += weight
    tot = dieFaces**diceAmount
    return [x/tot for x in ret]

Tässä tarvittiin tietää kombinaatiosta sen tyyppi eli siinä esiintyvien lukujen esiintymien lukumäärät (koottuna vektoriksi).

Nyt voidaan muodostaa siirtymämatriisi ja käyttää tavanomaisia keinoja \mathbb{E}[N]:n selvittämiseen:

def getTransitionMatrixBrute(dieFaces):
    arr = [getTransitionsProbsBrute2(i, dieFaces)
           for i in range(1, dieFaces+2)]
    return Matrix(arr)

def getFundamentalMatrix(pMat):
    n = pMat.dimensions()[0]-1
    arr = [[(1 if i==j else 0)-pMat[i][j]
            for j in range(n)]
           for i in range(n)]
    return Matrix(arr)**-1

def getExpectedGameLen(dieFaces):
    P = getTransitionMatrixBrute(dieFaces)
    fundMat = getFundamentalMatrix(P)
    return sum(fundMat[0])


print (getExpectedGameLen(6))

Ja näin saadaan tulos

\begin{aligned} \mathbb{E}[N] = \frac{133191861}{12500} = 10655.34888 \end{aligned}

Kombinaatioiden (eritoteen takaisinpanon kanssa) määrä kuitenkin kasvaa melko nopeasti, kun d kasvaa. Edellä jo mainittua tyyppi-ideaa voidaan jatkojalostaa: se kuinka monta noppaa heittotuloksesta menetetään riippuu vain tyypistä (kaikki, joiden arvoa yli yksi kappale menetetään). Tyypit ovat noppien lukumäärän ositukset (paitsi viimeisen tilan ositus [1,1,\dots,1] ei kelpaa, sillä se tarkoittaa d+1 eri arvoa, mutta nopassa on vain d eri arvoa. (Tästä vielä huomio-osioissa.) Nyt täytyy enää laskea kuinka moni noppien heitto vastaa kutakin ositusta. Tätä on selvitetty koodin kommenteissa:

def getFactorOfType(t, dieFaces):
    #first choose a value for each number in the type
    ret = factorial(dieFaces)/factorial(dieFaces-len(t))
    #divide away the orderings in the type
    #(we need the type's type)
    ret /= factorialProd(getType(t))
    #now we can multiply by the total number of orderings
    ret *= multinomial(t)
    return ret

def getTransitionMatrix(dieFaces):
    n = dieFaces+1
    arr = [[0]*n for i in range(n)]
    for i in range(1, n+1):
        for t in Partitions(i):
            if (len(t)>dieFaces): continue
            toState = i-sum(x for x in t if x>1)
            weight = getFactorOfType(t, dieFaces)
            arr[i-1][toState] += weight
    for i in range(n):
        tot = dieFaces**(i+1)
        for j in range(n):
            arr[i][j] /= tot
    return Matrix(arr)

Tässä vielä lopullinen tyypeillä laskeva koodi kokonaisuudessaan ja taulukko tuloksesta muutamilla d:n arvoilla:

from collections import Counter

def getType(dice):
    return list(Counter(dice).values())

def factorialProd(arr):
    return reduce(lambda acc,x: acc*factorial(x), arr, 1)

def getFactorOfType(t, dieFaces):
    ret = factorial(dieFaces)/factorial(dieFaces-len(t))
    ret /= factorialProd(getType(t))
    ret *= multinomial(t)
    return ret

def getTransitionMatrix(dieFaces):
    n = dieFaces+1
    arr = [[0]*n for i in range(n)]
    for i in range(1, n+1):
        for t in Partitions(i):
            if (len(t)>dieFaces): continue
            toState = i-sum(x for x in t if x>1)
            weight = getFactorOfType(t, dieFaces)
            arr[i-1][toState] += weight
    for i in range(n):
        tot = dieFaces**(i+1)
        for j in range(n):
            arr[i][j] /= tot
    return Matrix(arr)

def getFundamentalMatrix(pMat):
    n = pMat.dimensions()[0]-1
    arr = [[(1 if i==j else 0)-pMat[i][j]
            for j in range(n)]
           for i in range(n)]
    return Matrix(arr)**-1

def getExpectedGameLen(dieFaces):
    P = getTransitionMatrix(dieFaces)
    fundMat = getFundamentalMatrix(P)
    return sum(fundMat[0])


res = [[d, float(getExpectedGameLen(Integer(d)))]
      for d in range(1, 21)]
print (res)
dE[N]
11
24
315
480.96296
5718.61364
610655.34888
7261903.04770
810636087.48690
9713183947.85983
107.89467e+10
111.44259e+13
124.35102e+15
132.16594e+18
141.77942e+21
152.41246e+24
165.39710e+27
171.99232e+31
181.21348e+35
191.21944e+39
202.02174e+43

Huomioita

  • Viimeisen tilan siirtymiä ei tarvitsisi laskea (siinä mielessä, että peli jatkuisi) vaan se asetetaan absorboivaksi eikä sitä tarvita odotusarvon laskussa ollenkaan. Itse asiassa sen laskeminen on juuri kaikista matriisin riveistä vaativin, sillä tarvittavia objekteja (permutaatioita takaisinpanolla, kombinaatioita takaisin panolla tai osituksia) on eniten.
  • Myös ensimmäinen tila on käytännössä turha: yhtä noppaa on turha edes heittää, koska sitä ei voi menettää; ottaa vaan suoraan toisen nopan. Mutta pidetään se mukana yksinkertaisuuden vuoksi ja lasketaan pelin pituuteen mukaan myös nämä deterministiset askeleet.

Tehtäviä

  • Muokkaa koodia, niin että turhia viimeisen tilan siirtymiä ei lasketa. Aseta viimeiseen tilaan itsesiirtymä todennäköisyydellä 1 tai unohda koko viimeinen tila, sillä sitä ei laskuissa tarvita.
  • Laske N:n varianssi. Vinkki: katso Wikipedian artikkeli.
  • Entä jos sovitaan, että pelaaja häviää pelin, mikäli hän jollain kierroksella menettää kaikki noppansa. Mikä on häviämisen todennäköisyys?
  • Voidaan sopia mitä erilaisimpia sääntöjä siitä miten pelaaja saa lisänopan /-noppia tai menettää niitä riippuen hänen heitostaan. Tutki jotain tällaista muokattua peliä.

Vastaa

Täytä tietosi alle tai klikkaa kuvaketta kirjautuaksesi sisään:

WordPress.com-logo

Olet kommentoimassa WordPress.com -tilin nimissä. Log Out /  Muuta )

Google photo

Olet kommentoimassa Google -tilin nimissä. Log Out /  Muuta )

Twitter-kuva

Olet kommentoimassa Twitter -tilin nimissä. Log Out /  Muuta )

Facebook-kuva

Olet kommentoimassa Facebook -tilin nimissä. Log Out /  Muuta )

Muodostetaan yhteyttä palveluun %s