Нормальное распределение в Питоне: полное руководство для анализа

f(x) = (1 / (σ * sqrt(2π))) * e^(-(x-μ)²/(2σ²))

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

# Генерация 10000 случайных чисел из нормального распределения
# со средним 0 и стандартным отклонением 1
np_normal = np.random.normal(loc=0, scale=1, size=10000)

# Альтернативный способ с использованием NumPy Random Generator (для версий NumPy >= 1.17)
rng = np.random.default_rng()
np_normal_new = rng.normal(loc=0, scale=1, size=10000)

# Настройка стиля графиков
plt.figure(figsize=(12, 8))
sns.set_style('whitegrid')

# Гистограмма
plt.subplot(2, 2, 1)
plt.hist(np_normal, bins=50, alpha=0.7, color='skyblue', edgecolor='black')
plt.title('Гистограмма нормального распределения')
plt.xlabel('Значение')
plt.ylabel('Частота')

# Кривая плотности вероятности
plt.subplot(2, 2, 2)
sns.kdeplot(np_normal, fill=True, color='forestgreen')
plt.title('Кривая плотности вероятности')
plt.xlabel('Значение')
plt.ylabel('Плотность')

# QQ-график для проверки нормальности
plt.subplot(2, 2, 3)
from scipy import stats
stats.probplot(np_normal, plot=plt)
plt.title('QQ-график')

# Ящик с усами
plt.subplot(2, 2, 4)
sns.boxplot(x=np_normal, color='coral')
plt.title('Диаграмма размаха (Box plot)')
plt.xlabel('Значение')

plt.tight_layout()
plt.show()

# Визуализация влияния параметров на форму распределения
x = np.linspace(-10, 10, 1000)
plt.figure(figsize=(12, 6))

# Влияние среднего значения (μ)
plt.subplot(1, 2, 1)
for mu in [-3, 0, 3]:
y = stats.norm.pdf(x, mu, 1)
plt.plot(x, y, label=f'μ = {mu}, σ = 1')
plt.title('Влияние среднего значения (μ)')
plt.xlabel('x')
plt.ylabel('Плотность вероятности')
plt.legend()

# Влияние стандартного отклонения (σ)
plt.subplot(1, 2, 2)
for sigma in [0\.5, 1, 2]:
y = stats.norm.pdf(x, 0, sigma)
plt.plot(x, y, label=f'μ = 0, σ = {sigma}')
plt.title('Влияние стандартного отклонения (σ)')
plt.xlabel('x')
plt.ylabel('Плотность вероятности')
plt.legend()

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats

# Создадим три набора данных для сравнения
# Нормальное распределение
normal_data = np.random.normal(0, 1, 1000)
# Экспоненциальное распределение (асимметричное)
exp_data = np.random.exponential(1, 1000)
# Равномерное распределение
uniform_data = np.random.uniform(-3, 3, 1000)

# Визуальное сравнение распределений
plt.figure(figsize=(15, 10))

datasets = [normal_data, exp_data, uniform_data]
titles = ['Нормальное распределение', 'Экспоненциальное распределение', 'Равномерное распределение']

for i, (data, title) in enumerate(zip(datasets, titles)):
# Гистограмма с наложенной кривой плотности
plt.subplot(3, 3, i+1)
sns.histplot(data, kde=True, stat='density')
plt.title(f'{title} – Гистограмма')

# QQ-график
plt.subplot(3, 3, i+4)
stats.probplot(data, plot=plt)
plt.title(f'{title} – QQ-график')

# Ящик с усами
plt.subplot(3, 3, i+7)
sns.boxplot(x=data)
plt.title(f'{title} – Диаграмма размаха')

plt.tight_layout()
plt.show()

# Применение статистических тестов нормальности
def test_normality(data, name):
"""Проверка данных на нормальность с помощью различных тестов"""
print(f"\nРезультаты тестов нормальности для {name}:")

# Тест Шапиро-Уилка (лучший выбор для малых выборок, n < 50)
stat, p = stats.shapiro(data)
print(f"Тест Шапиро-Уилка: статистика={stat:.4f}, p-значение={p:.4e}")
if p > 0.05:
print("Данные выглядят нормально (не отвергаем нулевую гипотезу)")
else:
print("Данные не выглядят нормально (отвергаем нулевую гипотезу)")

# Тест Д'Агостино-Пирсона (комбинирует тесты на асимметрию и эксцесс)
stat, p = stats.normaltest(data)
print(f"Тест Д'Агостино-Пирсона: статистика={stat:.4f}, p-значение={p:.4e}")
if p > 0.05:
print("Данные выглядят нормально (не отвергаем нулевую гипотезу)")
else:
print("Данные не выглядят нормально (отвергаем нулевую гипотезу)")

# Тест Колмогорова-Смирнова
stat, p = stats.kstest(data, 'norm', args=(np.mean(data), np.std(data, ddof=1)))
print(f"Тест Колмогорова-Смирнова: статистика={stat:.4f}, p-значение={p:.4e}")
if p > 0.05:
print("Данные выглядят нормально (не отвергаем нулевую гипотезу)")
else:
print("Данные не выглядят нормально (отвергаем нулевую гипотезу)")

# Применяем тесты к нашим наборам данных
for data, name in zip(datasets, titles):
test_normality(data, name)

# Попробуем трансформировать экспоненциально распределенные данные
plt.figure(figsize=(15, 5))

# Исходные данные
plt.subplot(1, 3, 1)
sns.histplot(exp_data, kde=True, stat='density')
plt.title('Исходные данные (экспоненциальные)')

# Логарифмическая трансформация
log_data = np.log1p(exp_data) # log(1+x) предотвращает проблемы с нулевыми значениями
plt.subplot(1, 3, 2)
sns.histplot(log_data, kde=True, stat='density')
plt.title('Логарифмическая трансформация')

# Трансформация Box-Cox
from scipy import special
# Находим оптимальный параметр лямбда
lambda_opt, _ = stats.boxcox_normplot(exp_data + 1, -10, 10)
boxcox_data = special.boxcox(exp_data + 1, lambda_opt)
plt.subplot(1, 3, 3)
sns.histplot(boxcox_data, kde=True, stat='density')
plt.title(f'Box-Cox трансформация (λ={lambda_opt:.2f})')

plt.tight_layout()
plt.show()

# Проверяем нормальность трансформированных данных
test_normality(log_data, "Экспоненциальные данные (лог-трансформация)")
test_normality(boxcox_data, "Экспоненциальные данные (Box-Cox трансформация)")

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

# Параметры нормального распределения
mu = 100 # среднее значение, например IQ
sigma = 15 # стандартное отклонение

# Вероятность, что случайный человек имеет IQ между 85 и 115
# (в пределах одного стандартного отклонения от среднего)
prob_within_one_sigma = stats.norm.cdf(115, mu, sigma) – stats.norm.cdf(85, mu, sigma)
print(f"Вероятность IQ между 85 и 115: {prob_within_one_sigma:.4f} ({prob_within_one_sigma*100:.1f}%)")

# Вероятность, что IQ выше 130 (высокий интеллект)
prob_above_130 = 1 – stats.norm.cdf(130, mu, sigma)
print(f"Вероятность IQ выше 130: {prob_above_130:.4f} ({prob_above_130*100:.2f}%)")

# Построение 95% доверительного интервала
confidence_level = 0.95
margin_of_error = stats.norm.ppf((1 + confidence_level) / 2) * sigma / np.sqrt(100) # для выборки размером 100
ci_lower = mu – margin_of_error
ci_upper = mu + margin_of_error
print(f"95% доверительный интервал для среднего IQ: [{ci_lower:.2f}, {ci_upper:.2f}]")

# Визуализация
x = np.linspace(mu – 3*sigma, mu + 3*sigma, 1000)
y = stats.norm.pdf(x, mu, sigma)

plt.figure(figsize=(12, 6))
plt.plot(x, y, 'r-', lw=2)
plt.fill_between(x, y, where=(x >= 85) & (x <= 115), color='blue', alpha=0.3, label='IQ 85-115')
plt.fill_between(x, y, where=x >= 130, color='green', alpha=0.3, label='IQ > 130')
plt.axvline(ci_lower, color='purple', linestyle='--', label=f'95% CI: {ci_lower:.2f}')
plt.axvline(ci_upper, color='purple', linestyle='--', label=f'95% CI: {ci_upper:.2f}')
plt.title('Нормальное распределение IQ в популяции')
plt.xlabel('Значение IQ')
plt.ylabel('Плотность вероятности')
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# Имитация данных о ежедневных продажах с выбросами
np.random.seed(42)
sales_data = np.random.normal(5000, 800, 100) # 100 дней продаж
# Добавление аномальных значений
sales_data[20] = 9000 # Праздничный день
sales_data[50] = 8500 # Распродажа
sales_data[75] = 2000 # Технические проблемы

# Определение выбросов с использованием Z-оценки
z_scores = (sales_data – np.mean(sales_data)) / np.std(sales_data)
outliers = np.abs(z_scores) > 2 # выбросы: более 2 стандартных отклонений от среднего

print(f"Количество выявленных выбросов: {np.sum(outliers)}")
print(f"Индексы выбросов: {np.where(outliers)[0]}")
print(f"Значения выбросов: {sales_data[outliers]}")

# Визуализация данных и выбросов
plt.figure(figsize=(12, 6))
plt.scatter(range(len(sales_data)), sales_data, c=['red' if x else 'blue' for x in outliers])
plt.axhline(np.mean(sales_data), color='green', linestyle='-', label='Среднее')
plt.axhline(np.mean(sales_data) + 2*np.std(sales_data), color='orange', linestyle='--', 
label='Верхняя граница (μ+2σ)')
plt.axhline(np.mean(sales_data) – 2*np.std(sales_data), color='orange', linestyle='--', 
label='Нижняя граница (μ-2σ)')
plt.title('Выявление выбросов с помощью Z-оценки')
plt.xlabel('День')
plt.ylabel('Продажи')
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# Имитация результатов A/B-теста для конверсии
np.random.seed(42)
# Группа A: 5% конверсия, 1000 пользователей
conversions_A = np.random.binomial(1, 0.05, 1000).sum()
# Группа B: 6.5% конверсия, 1000 пользователей
conversions_B = np.random.binomial(1, 0.065, 1000).sum()

# Расчет конверсии
conversion_rate_A = conversions_A / 1000
conversion_rate_B = conversions_B / 1000
print(f"Конверсия в группе A: {conversion_rate_A:.1%}")
print(f"Конверсия в группе B: {conversion_rate_B:.1%}")
print(f"Абсолютное улучшение: {(conversion_rate_B – conversion_rate_A):.1%}")
print(f"Относительное улучшение: {((conversion_rate_B – conversion_rate_A) / conversion_rate_A):.1%}")

# Проведение Z-теста для сравнения пропорций
from statsmodels.stats.proportion import proportions_ztest
count = np.array([conversions_A, conversions_B])
nobs = np.array([1000, 1000])
stat, p_value = proportions_ztest(count, nobs)

print(f"Z-статистика: {stat:.4f}")
print(f"P-значение: {p_value:.4f}")
if p_value < 0.05:
print("Результат статистически значим (p < 0.05)")
else:
print("Результат статистически не значим (p >= 0.05)")

# Визуализация распределения разницы конверсий
def simulate_ab_test(p_a, p_b, n_a, n_b, n_simulations=10000):
"""Симуляция множества A/B-тестов для построения распределения разницы"""
results_a = np.random.binomial(n_a, p_a, n_simulations) / n_a
results_b = np.random.binomial(n_b, p_b, n_simulations) / n_b
differences = results_b – results_a
return differences

differences = simulate_ab_test(0.05, 0.065, 1000, 1000)

plt.figure(figsize=(10, 6))
plt.hist(differences, bins=50, alpha=0.7, color='skyblue')
plt.axvline(0, color='red', linestyle='--', label='Нулевая гипотеза (нет разницы)')
plt.axvline(conversion_rate_B – conversion_rate_A, color='green', linestyle='-', 
label=f'Наблюдаемая разница: {(conversion_rate_B – conversion_rate_A):.1%}')
plt.title('Распределение разницы конверсий между группами A и B')
plt.xlabel('Разница конверсий (B – A)')
plt.ylabel('Частота')
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# Оценка риска инвестиционного портфеля
# Предположим, годовая доходность имеет нормальное распределение
annual_return_mean = 0.08 # 8% средняя доходность
annual_return_std = 0.15 # 15% волатильность

# Вероятность отрицательной доходности за год
prob_negative_return = stats.norm.cdf(0, annual_return_mean, annual_return_std)
print(f"Вероятность отрицательной годовой доходности: {prob_negative_return:.2%}")

# Вероятность потери более 10% капитала
prob_loss_10percent = stats.norm.cdf(-0.10, annual_return_mean, annual_return_std)
print(f"Вероятность потери более 10% капитала: {prob_loss_10percent:.2%}")

# Расчет Value at Risk (VaR) с уровнем доверия 95%
VaR_95 = -stats.norm.ppf(0.05, annual_return_mean, annual_return_std)
print(f"Value at Risk (95% CI): {VaR_95:.2%}")

# Расчет Expected Shortfall (ES) / Conditional VaR с уровнем доверия 95%
def expected_shortfall(mean, std, confidence_level):
"""Расчет ожидаемых потерь при превышении VaR"""
alpha = 1 – confidence_level
var = -stats.norm.ppf(alpha, mean, std)
pdf_value = stats.norm.pdf(stats.norm.ppf(alpha))
es = -(mean – std * pdf_value / alpha)
return es

ES_95 = expected_shortfall(annual_return_mean, annual_return_std, 0.95)
print(f"Expected Shortfall (95% CI): {ES_95:.2%}")

# Визуализация доходностей и риск-метрик
x = np.linspace(-0.5, 0.5, 1000)
y = stats.norm.pdf(x, annual_return_mean, annual_return_std)

plt.figure(figsize=(12, 6))
plt.plot(x, y, 'b-', lw=2)
plt.fill_between(x, y, where=x < 0, color='red', alpha=0.3, label='Отрицательная доходность')
plt.fill_between(x, y, where=x < -0.10, color='darkred', alpha=0.5, label='Потери > 10%')
plt.axvline(-VaR_95, color='purple', linestyle='--', label=f'VaR 95%: {VaR_95:.2%}')
plt.axvline(-ES_95, color='magenta', linestyle='-.', label=f'ES 95%: {ES_95:.2%}')
plt.title('Распределение годовой доходности инвестиционного портфеля')
plt.xlabel('Годовая доходность')
plt.ylabel('Плотность вероятности')
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

import numpy as np
import time

# Пример: расчет Z-оценок для большого массива данных

# Неоптимизированный подход (циклы)
def z_scores_loop(data):
z_scores = []
mean = sum(data) / len(data)
std = (sum((x – mean) ** 2 for x in data) / len(data)) ** 0.5
for x in data:
z_scores.append((x – mean) / std)
return z_scores

# Оптимизированный подход (векторизация)
def z_scores_numpy(data):
return (data – np.mean(data)) / np.std(data)

# Генерация большого массива данных
large_data = np.random.normal(100, 15, 1_000_000)

# Сравнение производительности
start_time = time.time()
z_scores_loop_result = z_scores_loop(large_data[:10000]) # Используем меньшую выборку для цикла
loop_time = time.time() – start_time
print(f"Время выполнения с циклами (для 10,000 элементов): {loop_time:.4f} секунд")

start_time = time.time()
z_scores_numpy_result = z_scores_numpy(large_data)
numpy_time = time.time() – start_time
print(f"Время выполнения с NumPy (для 1,000,000 элементов): {numpy_time:.4f} секунд")
print(f"Ускорение (с учетом размера данных): ~{(loop_time / numpy_time) * (1_000_000 / 10_000):.1f}x")

from joblib import Parallel, delayed
import multiprocessing
import numpy as np
from scipy import stats
import time

# Функция для симуляции множества тестов нормальности
def simulate_normality_test(n_samples, sample_size):
"""Симулирует тесты Шапиро-Уилка для n_samples выборок размера sample_size"""
results = []
for _ in range(n_samples):
# Генерация случайной выборки из нормального распределения
sample = np.random.normal(0, 1, sample_size)
# Выполнение теста нормальности
_, p_value = stats.shapiro(sample)
results.append(p_value)
return results

# Последовательный подход
start_time = time.time()
sequential_results = simulate_normality_test(100, 500)
sequential_time = time.time() – start_time
print(f"Последовательное выполнение: {sequential_time:.4f} секунд")

# Параллельный подход
def single_test(i, sample_size):
sample = np.random.normal(0, 1, sample_size)
_, p_value = stats.shapiro(sample)
return p_value

start_time = time.time()
n_cores = multiprocessing.cpu_count()
parallel_results = Parallel(n_jobs=n_cores)(
delayed(single_test)(i, 500) for i in range(100)
)
parallel_time = time.time() – start_time
print(f"Параллельное выполнение ({n_cores} ядер): {parallel_time:.4f} секунд")
print(f"Ускорение: {sequential_time / parallel_time:.2f}x")

import numpy as np
import pandas as pd
import time

# Создаем большой DataFrame с нормально распределенными данными
n_rows = 1_000_000
n_cols = 10
df = pd.DataFrame(
np.random.normal(0, 1, (n_rows, n_cols)),
columns=[f'feature_{i}' for i in range(n_cols)]
)

# Расчет статистик с помощью NumPy (требует извлечения данных из DataFrame)
start_time = time.time()
means_np = np.mean(df.values, axis=0)
stds_np = np.std(df.values, axis=0)
numpy_time = time.time() – start_time
print(f"NumPy время: {numpy_time:.4f} секунд")

# Расчет статистик с помощью оптимизированных методов pandas
start_time = time.time()
means_pd = df.mean()
stds_pd = df.std()
pandas_time = time.time() – start_time
print(f"Pandas время: {pandas_time:.4f} секунд")

# Использование описательных статистик pandas
start_time = time.time()
describe_pd = df.describe()
describe_time = time.time() – start_time
print(f"Pandas describe() время: {describe_time:.4f} секунд")

import numpy as np
import time

# Сравнение различных методов генерации случайных чисел из нормального распределения
n_samples = 10_000_000

# Метод 1: numpy.random.normal
start_time = time.time()
samples_1 = np.random.normal(0, 1, n_samples)
time_1 = time.time() – start_time
print(f"numpy.random.normal: {time_1:.4f} секунд")

# Метод 2: Random Generator API (быстрее для больших массивов в NumPy >= 1.17)
start_time = time.time()
rng = np.random.default_rng()
samples_2 = rng.normal(0, 1, n_samples)
time_2 = time.time() – start_time
print(f"numpy.random.default_rng().normal: {time_2:.4f} секунд")

# Метод 3: преобразование равномерного распределения через Box-Muller
start_time = time.time()
u1 = np.random.random(n_samples)
u2 = np.random.random(n_samples)
samples_3 = np.sqrt(-2 * np.log(u1)) * np.cos(2 * np.pi * u2)
time_3 = time.time() – start_time
print(f"Box-Muller transform: {time_3:.4f} секунд")

# Метод 4: использование Cupy для GPU-ускорения (если доступно)
try:
import cupy as cp
start_time = time.time()
samples_4 = cp.random.normal(0, 1, n_samples)
# Перенос результатов обратно в CPU для корректного сравнения
samples_4_cpu = cp.asnumpy(samples_4)
time_4 = time.time() – start_time
print(f"cupy.random.normal (GPU): {time_4:.4f} секунд")
except ImportError:
print("CuPy не установлен. Для GPU-ускорения установите: pip install cupy-cuda11x")

from functools import lru_cache
import numpy as np
from scipy import stats
import time

# Функция без кэширования
def normal_probability(a, b, mu, sigma):
"""Вычисляет вероятность P(a <= X <= b) для нормального распределения"""
return stats.norm.cdf(b, mu, sigma) – stats.norm.cdf(a, mu, sigma)

# Функция с кэшированием
@lru_cache(maxsize=128)
def normal_probability_cached(a, b, mu, sigma):
"""Кэшированная версия функции вычисления вероятности"""
return stats.norm.cdf(b, mu, sigma) – stats.norm.cdf(a, mu, sigma)

# Тестирование на повторяющихся запросах
queries = [(i/10, (i+1)/10, 0, 1) for i in range(-30, 30)] * 10

# Без кэширования
start_time = time.time()
results_no_cache = [normal_probability(*q) for q in queries]
no_cache_time = time.time() – start_time
print(f"Время без кэширования: {no_cache_time:.6f} секунд")

# С кэшированием
start_time = time.time()
results_cached = [normal_probability_cached(*q) for q in queries]
cache_time = time.time() – start_time
print(f"Время с кэшированием: {cache_time:.6f} секунд")
print(f"Ускорение: {no_cache_time / cache_time:.2f}x")