二项分布模拟python

问题描述

假设A和B两支球队进行了一系列比赛，第一支球队 球队赢得4连胜。假设A队有55％ 赢得每场比赛的机会，并且每场比赛的结果都是 独立的。

（a）A队赢得系列赛的概率是多少？给一个 准确的结果，并通过仿真进行确认。

（b）预期打多少场比赛？给出确切的 结果并通过仿真进行确认。

（c）鉴于甲队，预计将打多少场比赛 赢得系列赛？给出准确的结果并通过仿真确认。

（d）现在假设我们只知道A队更有可能获胜 每个游戏，但不知道确切的概率。如果最 可能玩的游戏数是5，这意味着什么 甲队赢得每场比赛的几率？

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

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],dtype=np.float64) # favorite counter
    u = np.array([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,N):
    """
    Simulate N events and computes expected number of games
    """

    Ng = 0
    for _ in range(0,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,p)
p24 = P(2,p)
p34 = P(3,p)

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

# underdog
x04 = P(0,1.0-p)
x14 = P(1,1.0-p)
x24 = P(2,1.0-p)
x34 = P(3,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,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,200000))