一類骰子游戲中的概率計算

埠默笙聲聲聲脈發表於2022-01-20

一個骰子,一個跑道,停在某個格子上有獎勵。包含這種玩法的遊戲不要太多,拿“大富翁”作個圖示:

在玩的時候時常在問自己:

  • 我停在前方第n格的概率是多少?
  • 我停在前方第n格的期望擲骰子數是多少?
  • 感性上,我停在前方第100格的概率,應該和我停在前方第1000格的概率是一樣的,那麼這個概率是多少?

不妨就來程式設計解決這些疑問!

Q1:停在前方第n格的概率是多少?

不妨先考慮簡單的情形:

n = 1,那麼我們至多能擲1次骰子,且僅值為1時滿足條件,那麼概率是 1 / 6 = 0.166667。

n = 2,那麼我們至多能擲2次骰子,點數為 (1, 1) 或 (2) 時滿足條件。所有擲骰子的情況是:{ (1, 1) (1, 2) (1, 3) (1, 4) (1, 5) (1, 6) (2) (3) (4) (5) (6) } 一共11種情況,那麼概率是 2 / 11 = 0.181818。

我們可以總結出這樣的公式:

  • 停在前方第n格的概率 = 正好停在第n格的情況數 ÷ 給定最多投擲數的情況下,停留的位置大於等於第n格的情況數 = F(n) / G(n)

 Step1:計算F(n) 

方法一:動態規劃

顯然,停在n格的情況數 = 停在n-1格的情況數 + 停在n-2格的情況數 + … + 停在n-6格的情況數,所以有:

  • 轉移方程:F(n) = F(n-1) + F(n-2) + F(n-3) + F(n-4) + F(n-5) + F(n-6)
  • 邊界條件:F(0) = 1,F(n<0) = 0

程式碼如下:

def calF(n):
    F_list = [0 for i in range(n+1)]
    for i in range(n+1):
        for j in range(1,7):
            if i - j == 0: F_list[i] += 1 # 即F(0)=1
            elif i - j < 0: F_list[i] += 0 # 即F(n<0)=0
            else: F_list[i] += F_list[i-j]
    return F_list[-1]

方法二:遞迴演算法

動態規劃與遞迴都有著“分而治之”的思想,在某些情況下是能相互轉換的。遞迴演算法是這麼看這個問題的:於當前位置,重複地執行擲骰子這一動作。若正好停在第n格,則計數器+1,函式返回;若大於第n格,則函式返回;若小於第n格,則繼續擲骰子。

程式碼如下:

def calF_with_RA(n, counter): 
    # n:到目標的距離
    # counter:達成條件的情況計數器
    for i in range(1,7):
        if i == n: 
            counter += 1
            return counter
        elif i < n:
            counter = calF_with_RA(n-i, counter)
    return counter

Step2:計算G(n)

如果在Step1中採用了遞迴演算法,那麼只要改寫一下函式中計數器的達成條件,和函式退出的位置即可,程式碼如下:

def calG_with_RA(n, counter): 
    # n:到目標的距離
    # counter:達成條件的情況計數器
    for i in range(1,7):
        if i >= n: 
            counter += 1
        else: counter = calG_with_RA(n-i, counter)
    return counter

Step3:計算並畫圖

於是計算概率的函式有:

def calF_divide_G(n):
    F = calF_with_DP(n)
    G = calG_with_RA(n, 0)
    return F/G

把n=1到20的概率畫出來:

一類骰子游戲中的概率計算
if __name__ == '__main__': 
    prob = []
    for i in range(1,21):
        prob.append(calF_divide_G(i))
    n = np.arange(1, 21).astype(dtype=np.str_) # 轉換資料型別 否則plot中會顯示浮點數
    plt.plot(n,prob,linewidth=2,color='r',marker='o', markerfacecolor='blue',markersize=8)
    plt.xlabel('n')
    plt.ylabel('P(n)')
    plt.show()
View Code

可以看出確實有收斂的趨勢(n取大於25時,程式會變得異常難算,難以展示),收斂到的概率是0.196717(保留6位小數),這也正好回答了Q3。另外,n=6時概率最大,值為0.198758(保留6位小數)。

Q2:停在前方第n格的期望擲骰子數是多少?

一個思路是,記錄所有到達第n格的擲骰子“軌跡”,然後對所有軌跡進行統計。由於遞迴演算法本身就有列印軌跡的能力,順便再統計點相關的資料自然不在話下,於是改寫遞迴函式:

def throw_dice(n, track, record):
    # n:到目標的距離
    # track:骰子的軌跡
    # record:記錄所有軌跡的長度(即擲骰子次數)
    for i in range(1,7):
        if i == n:
            track.append(i)
            print(track)
            record.append(len(track))
            return record
        elif i < n:
            track_ = track.copy()
            track_.append(i)
            throw_dice(n-i, track_, record)
    return record

遍歷所有擲骰子軌跡,求出擲骰子數的概率分佈和期望:

def get_dice_distribution(n):
    record = throw_dice(n, [], [])
    Min = min(record)
    Max = max(record)
    distribution = [0 for i in range(Max - Min + 1)]
    for num in record:
        distribution[num - Min] += 1
    Sum = sum(distribution)
    for i, num in enumerate(distribution): # 求擲骰子數的概率分佈
        distribution[i] = num / Sum
    
    expectation = 0 # 求期望擲骰子數
    for i in range(Max - Min +1):
        expectation += (Min + i) * distribution[i]
    print(expectation)
    
    return Min, Max, distribution, expectation

畫出擲骰子數的概率分佈(以n=10為例):

一類骰子游戲中的概率計算
def plotDistribution(n):
    Min, Max, distribution, expectation = get_dice_distribution(n)
    N = np.arange(Min, Max+1).astype(dtype=np.str_)
    plt.bar(N, distribution)
    plt.xlabel('n')
    plt.ylabel('throwing times')
    plt.show()
View Code

畫出n=1到20的期望擲骰子數:

一類骰子游戲中的概率計算
def plotExpectation():
    expectation_list = []
    for i in range(1,21):
        Min, Max, distribution, expectation = get_dice_distribution(i)
        expectation_list.append(expectation)
    n = np.arange(1, 21).astype(dtype=np.str_)
    plt.plot(n, expectation_list, linewidth=2, color='r', marker='o', markerfacecolor='blue', markersize=8)
    plt.xlabel('n')
    plt.ylabel('E(n)')
    plt.show()
    return expectation_list
View Code

可以看出,大致上是一個線性的關係,列印出具體的資料:

可以看出,當n小於等於6時,是線性的;大於6時,大致上是線性的。

整體程式碼:

一類骰子游戲中的概率計算
import numpy as np
from matplotlib import pyplot as plt 


def calF_with_DP(n):
    F_list = [0 for i in range(n+1)]
    for i in range(n+1):
        for j in range(1,7):
            if i - j == 0: F_list[i] += 1 # 即F(0)=1
            elif i - j < 0: F_list[i] += 0 # 即F(n<0)=0
            else: F_list[i] += F_list[i-j]
    return F_list[-1]


def calF_with_RA(n, counter): 
    # n:到目標的距離
    # counter:達成條件的情況計數器
    for i in range(1,7):
        if i == n: 
            counter += 1
            return counter
        elif i < n:
            counter = calF_with_RA(n-i, counter)
    return counter


def calG_with_RA(n, counter): 
    # n:到目標的距離
    # counter:達成條件的情況計數器
    for i in range(1,7):
        if i >= n: 
            counter += 1
        else: counter = calG_with_RA(n-i, counter)
    return counter


def calF_divide_G(n):
    F = calF_with_RA(n, 0)
    G = calG_with_RA(n, 0)
    return F/G


def plotPn():
    prob = []
    for i in range(1,21):
        prob.append(calF_divide_G(i))
    n = np.arange(1, 21).astype(dtype=np.str_) # 轉換資料型別 否則plot中會顯示浮點數
    plt.plot(n, prob, linewidth=2, color='r', marker='o', markerfacecolor='blue', markersize=8)
    plt.xlabel('n')
    plt.ylabel('P(n)')
    plt.show()


def throw_dice(n, track, record):
    # n:到目標的距離
    # track:骰子的軌跡
    # record:記錄所有軌跡的長度(即擲骰子次數)
    for i in range(1,7):
        if i == n:
            track.append(i)
            record.append(len(track))
            return record
        elif i < n:
            track_ = track.copy()
            track_.append(i)
            throw_dice(n-i, track_, record)
    return record


def get_dice_distribution(n):
    record = throw_dice(n, [], [])
    Min = min(record)
    Max = max(record)
    distribution = [0 for i in range(Max - Min + 1)]
    for num in record:
        distribution[num - Min] += 1
    Sum = sum(distribution)
    for i, num in enumerate(distribution): # 求擲骰子數的概率分佈
        distribution[i] = num / Sum
    
    expectation = 0 # 求期望擲骰子數
    for i in range(Max - Min +1):
        expectation += (Min + i) * distribution[i]
    print(expectation)
    
    return Min, Max, distribution, expectation

    
def plotDistribution(n):
    Min, Max, distribution, expectation = get_dice_distribution(n)
    N = np.arange(Min, Max+1).astype(dtype=np.str_)
    plt.bar(N, distribution)
    plt.xlabel('n')
    plt.ylabel('throwing times')
    plt.show()


def plotExpectation():
    expectation_list = []
    for i in range(1,21):
        Min, Max, distribution, expectation = get_dice_distribution(i)
        expectation_list.append(expectation)
    n = np.arange(1, 21).astype(dtype=np.str_)
    plt.plot(n, expectation_list, linewidth=2, color='r', marker='o', markerfacecolor='blue', markersize=8)
    plt.xlabel('n')
    plt.ylabel('E(n)')
    plt.show()
    return expectation_list
    
if __name__ == '__main__':
    # plotPn()
    # plotDistribution(10)
    expectation_list = plotExpectation()
View Code

小結

使用動態規劃和遞迴演算法解決了問題(實際上可以完全依靠遞迴)。遞迴函式由於具有列印軌跡的能力,攜帶資訊多,可擴充性強,稍微改寫就能滿足需要。

 

相關文章