伯努利样品和伯努利pmf的图密度直方图

问题描述

问题摘要：

为什么我的样本密度与pmf如此不同，我如何执行此模拟以使pmf和样本估计值相似。

问题：

我使用scipy模拟了独立的伯努利试验的样本。我现在尝试获取我创建的样本的密度直方图，并将其与pmf（概率质量函数）进行比较。我希望密度直方图显示两个bin，每个bin都徘徊在pmf附近，但我在5的pmf值上方有2个bin。请问有人可以告诉我如何创建对Bernoulli不这样做的密度直方图吗？我尝试了一些其他发行版的类似仿真，但似乎效果很好。我在这里想念的是什么？您能告诉我如何操纵我的代码来完成这项工作吗？

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats

trials = 10**3
p = 0.5


sample_bernoulli = stats.bernoulli.rvs(p,size=trials) # Generate benoulli RV
plt.plot((0,1),stats.bernoulli.pmf((0,p),'bo',ms=8,label='bernoulli pmf')

# Density histogram of generated values
plt.hist(sample_bernoulli,density=True,alpha=0.5,color='steelblue',edgecolor='none')
plt.show()

如果这是一个简单或琐碎的问题，我必须道歉，但我找不到在线解决方案，却发现这个问题很有趣。任何帮助将不胜感激。

解决方法

原因是plt.hist主要用于连续分布。如果不提供明确的bin边界，则plt.hist只会在最小值和最大值之间创建10个等距的bin。这些垃圾箱大多数将是空的。在只有两个可能的数据值的情况下，应该只有两个bin，所以有3个边界：

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats

trials = 10**3
p = 0.5

sample_bernoulli = stats.bernoulli.rvs(p,size=trials) # Generate benoulli RV
plt.plot((0,1),stats.bernoulli.pmf((0,p),'bo',ms=8,label='bernoulli pmf')

# Density histogram of generated values
plt.hist(sample_bernoulli,density=True,alpha=0.5,color='steelblue',edgecolor='none',bins=np.linspace(-0.5,1.5,3))
plt.show()

这是默认箱边界以及样本如何放入箱的可视化。请注意，使用density=True，对直方图进行了归一化处理，以使所有条形的面积总和为1。在这种情况下，两个条形的宽度为0.1，高约5.0，而其他8条具有高度为零。因此，总面积为2*0.1*5 + 8*0.0 = 1。

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats

trials = 10 ** 3
p = 0.5

sample_bernoulli = stats.bernoulli.rvs(p,size=trials)  # Generate benoulli RV

# Density histogram of generated values with default bins
values,binbounds,bars = plt.hist(sample_bernoulli,alpha=0.2,edgecolor='none')
# show the bin boundaries
plt.vlines(binbounds,max(values) * 1.05,color='crimson',ls=':')
# show the sample values with a random displacement
plt.scatter(sample_bernoulli * 0.9 + np.random.uniform(0,0.1,trials),np.random.uniform(0,max(values),color='lime')
# show the index of each bin
for i in range(len(binbounds) - 1):
    plt.text((binbounds[i] + binbounds[i + 1]) / 2,max(values) / 2,i,ha='center',va='center',fontsize=20,color='crimson')
plt.show()