我想在
python上计算二项式概率.我试着应用公式:
probability = scipy.misc.comb(n,k)*(p**k)*((1-p)**(n-k))
我得到的一些概率是无限的.我检查了一些p = inf的值.对于其中一个,n = 450,000,k = 17.该值必须大于1e302,这是浮点数处理的最大值.
然后我尝试使用sum(np.random.binomial(n,p,numberOfTrials)== valueOfInterest)/ numberOfTrials
这将绘制numberOfTrials样本并计算绘制valueOfInterest值的平均次数.
这不会带来任何无限的价值.但是,这是一种有效的方法吗?为什么这种方式不会提高任何无限值,而计算概率呢?
解决方法
在日志域中工作以计算组合和取幂函数,然后将它们提升为指数.
像这样的东西:
combination_num = range(k+1,n+1) combination_den = range(1,n-k+1) combination_log = np.log(combination_num).sum() - np.log(combination_den).sum() p_k_log = k * np.log(p) neg_p_K_log = (n - k) * np.log(1 - p) p_log = combination_log + p_k_log + neg_p_K_log probability = np.exp(p_log)
由于数字较大,删除了数字下溢/溢出.在n = 450000且p = 0.5,k = 17的示例中,它返回p_log = -311728.4,i.例如,最终概率的对数非常小,因此在获取np.exp时会发生下溢.但是,您仍然可以使用日志概率.