import numpy as np
#import statsmodels.api as sm
import scipy.stats as stats
import matplotlib.pyplot as plt
import seaborn as sns

# set the random seed:
np.random.seed(123456)

# set sample size:
n = 11

# initialize r to an array of length r=1000 to later store results:
r = 1000

sample_mean = np.empty(r)
sample_var = np.empty(r)

# repeat r times:
for j in range(r):
    # draw a sample and store the sample mean in pos. j=0,1,... of ybar:
    sample = stats.norm.rvs(10, 2, size=n)
    sample_mean[j] = np.mean(sample)
    sample_var[j] = np.var(sample, ddof=1)

mean = np.mean(sample_var)
variance = np.var(sample_var, ddof=1)
print("Mean of sampling distribution from Normal is :", mean)
print("Variance of sampling distribution from Normal is :", variance)

# simulated density:
#kde = sm.nonparametric.KDEUnivariate(sample_var)
#kde.fit()

fig, ax = plt.subplots()

sns.histplot(data=sample_var, x=None, kde=True).set(title='Histogram of Sampling_dist w/ n=11')
ax.set_xlim(0,15)
#plt.savefig('C:/BOOK/PyBasics/PyStat/code/b1-ch5-8.png')