Bernoulli样本和Bernoulli pmf的图密度直方图

Question

问题摘要：

为什么我的密度从我的样本如此不同的pmf和我如何进行这个模拟，使pmf和样本的估计是相似的。

问题：

我用scipy模拟了一个独立的伯努利试验的样本。现在，我试图获取我创建的样本的密度直方图，并将其与pmf (概率质量函数)进行比较。我希望密度直方图能显示每个在pmf附近徘徊的两个桶，但相反，我在pmf值之上有2个桶，在5。有人能告诉我如何创建一个密度直方图，而不是针对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,1), 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()

​

​

如果这是一个简单或琐碎的问题，我必须道歉，但我在网上找不到解决办法，而且觉得这个问题很有趣。任何帮助都将不胜感激。

Answer 1

原因是plt.hist主要用于连续发行版。如果不提供明确的bin边界，plt.hist只在最小值和最大值之间创建10个间隔相等的回收箱。这些垃圾箱大部分都是空的。如果只有两个可能的数据值，就应该只有两个回收箱，所以有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,1), 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, density=True, alpha=0.2, color='steelblue', edgecolor='none')
# show the bin boundaries
plt.vlines(binbounds, 0, 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), trials), 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()

​

​