import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats
# Set random seed for reproducibility
np.random.seed(42)
# Create the figure and axis
plt.figure(figsize=(10, 6))
# Generate data points for x-axis
x = np.linspace(-4, 4, 1000)
# Calculate normal distribution values (PDF)
mean = 0
std_dev = 1
y = stats.norm.pdf(x, mean, std_dev)
# Plot the normal distribution curve
plt.plot(x, y, 'b-', linewidth=2, label=f'Normal Distribution\nμ={mean}, σ={std_dev}')
# Add a histogram of random samples from this distribution
samples = np.random.normal(mean, std_dev, 1000)
plt.hist(samples, bins=30, density=True, alpha=0.5, color='skyblue', label='Random Samples')
# Customize the plot
plt.title('Normal Distribution', fontsize=16)
plt.xlabel('Value', fontsize=12)
plt.ylabel('Probability Density', fontsize=12)
plt.grid(True, alpha=0.3, linestyle='--')
plt.legend()
# Mark the mean and standard deviations
plt.axvline(mean, color='red', linestyle='-', alpha=0.7, label='Mean')
plt.axvline(mean + std_dev, color='green', linestyle='--', alpha=0.7, label='μ ± σ (68%)')
plt.axvline(mean - std_dev, color='green', linestyle='--', alpha=0.7)
plt.axvline(mean + 2*std_dev, color='orange', linestyle='--', alpha=0.7, label='μ ± 2σ (95%)')
plt.axvline(mean - 2*std_dev, color='orange', linestyle='--', alpha=0.7)
# Add shaded regions for standard deviations
plt.fill_between(x, 0, y, where=((x >= mean-std_dev) & (x <= mean+std_dev)),
color='green', alpha=0.1)
plt.fill_between(x, 0, y, where=((x >= mean-2*std_dev) & (x <= mean+2*std_dev)),
color='orange', alpha=0.05)
plt.tight_layout()
plt.show()