python - サンプル分布シミュレーションが正常にならない

Question

Pythonを使用して「サンプル比率のサンプリング分布」をシミュレートしようとしていました。ここの例のようにベルヌーイ変数で試しました

重要なのは、多数のガムボールのうち、真の比率が 0.6 の黄色のボールがあることです。サンプル (あるサイズ、たとえば 10) を取得し、それを平均してプロットすると、正規分布が得られるはずです。

私はpythonでやろうとしましたが、常に均一な分布しか得られません（または中央で平らになります）。何が欠けているのか理解できません。

プログラム：

from SDSP import create_bernoulli_population, get_frequency_df
from random import shuffle, choices
from bi_to_nor_demo import get_metrics, bare_minimal_plot
import matplotlib.pyplot as plt


N = 10000  # 10000 balls
p = 0.6    # probability of yellow ball is 0.6, and others (1-0.6)=>0.4
n_pickups = 1000       # sample size
n_experiments = 100  # I dont know what this is called 


# generate population
population = create_bernoulli_population(N,p)
theor_df = get_frequency_df(population)
theor_df

# choose sample, take mean and add to X_mean_list. Do this for n_experiments times
X_hat = []
X_mean_list = []
for each_experiment in range(n_experiments):
    X_hat = choices(population, k=n_pickups)  # this method is with replacement
    shuffle(population)
    X_mean = sum(X_hat)/len(X_hat)
    X_mean_list.append(X_mean)

# plot X_mean_list as bar graph
stats_df = get_frequency_df(X_mean_list)
fig, ax = plt.subplots(1,1, figsize=(5,5))
X = stats_df['x'].tolist()
P = stats_df['p(x)'].tolist()    
ax.bar(X, P, color="C0") 

plt.show()

依存関数:
bi_to_nor_demo
SDSP

出力：

更新: 以下のように均一な配布を試みましたが、同様の出力が得られました。通常に収束していません :(. (create_bernoulli_population の代わりに以下の関数を使用)

def create_uniform_population(N, Y=[]):
    """
    Given the total size of population N, 
    this function generates list of those outcomes uniformly distributed
    population list
    N - Population size, eg N=10000
    p - probability of interested outcome  
    Returns the outcomes spread out in population as a list
    """
    uniform_p = 1/len(Y)
    print(uniform_p)
    total_pops = []
    for i in range(0,len(Y)):
        each_o = [i]*(int(uniform_p*N))
        total_pops += each_o
    shuffle(total_pops)    
    return total_pops