二项分布模拟 python

Binomial distribution simulation python

假设 2 支球队 A 和 B 正在进行一系列比赛,第一场 赢得 4 场比赛的球队赢得系列赛。假设 A 队有 55% 赢得每场比赛的机会,每场比赛的结果是 独立。

(a) A 队赢得系列赛的概率是多少?给一个 准确结果并通过模拟确认。

(b) 预计玩了多少场比赛?给个确切的 结果并通过模拟确认。

(c) 给定 A 队的预期比赛次数是多少 赢得系列赛?给出一个准确的结果并通过模拟确认。

(d) 现在假设我们只知道A队更有可能获胜 每场比赛,但不知道确切的概率。如果最 可能玩的游戏数是 5,这意味着什么 A队每场比赛获胜的概率?

这就是我所做的,但没有得到它..需要一些输入。谢谢

import numpy as np

probs = np.array([.55 ,.45])
nsims = 500000

chance = np.random.uniform(size=(nsims, 7))

teamAWins = (chance > probs[None, :]).astype('i4')
teamBWins = 1 - teamAWins

teamAwincount = {}
teamBwincount = {}
for ngames in range(4, 8):
    afilt = teamAWins[:, :ngames].sum(axis=1) == 4
    bfilt = teamBWins[:, :ngames].sum(axis=1) == 4

    teamAwincount[ngames] = afilt.sum()
    teamBwincount[ngames] = bfilt.sum()

    teamAWins = teamAWins[~afilt]
    teamBWins = teamBWins[~bfilt]

teamAwinprops = {k : 1. * count/nsims for k, count in teamAwincount.iteritems()}
teamBwinprops = {k : 1. * count/nsims for k, count in teamBwincount.iteritems()}

好的,这是让你继续前进的想法和代码。

我相信这是 Negative Binomial Distribution,这很容易实现,并计算出热门和失败者的概率。

使用该代码,您可以定义整组事件,概率总和为 1。由此您可以:

  1. 获取准确答案
  2. 根据概率检查模拟

模拟代码为许多事件和单个事件模拟器添加了计数器。到目前为止,看起来它的概率与 负二项式

代码,Python 3.8 x64 Win10

import numpy as np
import scipy.special

# Negative Binomial as defined in
# https://mathworld.wolfram.com/NegativeBinomialDistribution.html
def P(x, r, p):
    return scipy.special.comb(x+r-1, r-1)*p**r*(1.0-p)**x

def single_event(p, rng):
    """
    Simulates single up-to-4-wins event,
    returns who won and how many opponent got
    """
    f = 0
    u = 0
    while True:
        if rng.random() < p:
            f += 1
            if f == 4:
                return (True,u) # favorite won
        else:
            u += 1
            if u == 4:
                return (False,f) # underdog won

def sample(p, rng, N):
    """
    Simulate N events and count all possible outcomes
    """

    f = np.array([0, 0, 0, 0], dtype=np.float64) # favorite counter
    u = np.array([0, 0, 0, 0], dtype=np.float64) # underdog counter

    for _ in range(0, N):
        w, i = single_event(p, rng)
        if w:
            f[i] += 1
        else:
            u[i] += 1

    return (f/float(N), u/float(N)) # normalization

def expected_nof_games(p, rng, N):
    """
    Simulate N events and computes expected number of games
    """

    Ng = 0
    for _ in range(0, N):
        w, i = single_event(p, rng)

        Ng += 4 + i # 4 games won by winner and i by loser

    return float(Ng)/float(N)


p = 0.55

# favorite
p04 = P(0, 4, p)
p14 = P(1, 4, p)
p24 = P(2, 4, p)
p34 = P(3, 4, p)

print(p04, p14, p24, p34, p04+p14+p24+p34)

# underdog
x04 = P(0, 4, 1.0-p)
x14 = P(1, 4, 1.0-p)
x24 = P(2, 4, 1.0-p)
x34 = P(3, 4, 1.0-p)
print(x04, x14, x24, x34, x04+x14+x24+x34)

# total probability
print(p04+p14+p24+p34+x04+x14+x24+x34)

# simulation of the games
rng = np.random.default_rng()
f, u = sample(p, rng, 200000)
print(f)
print(u)

# compute expected number of games

print("expected number of games")
E_ng = 4*p04 + 5*p14 + 6*p24 + 7*p34 + 4*x04 + 5*x14 + 6*x24 + 7*x34
print(E_ng)
# same result from Monte Carlo
print(expected_nof_games(p, rng, 200000))