import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import kstest, norm

n_simulations = 10000  
time_steps = 1000  
T = 1 
dt = T / time_steps  
time_array = np.linspace(0, T, time_steps) 

# Генерация траекторий винеровского процесса W_t
dW = np.sqrt(dt) * np.random.randn(n_simulations, time_steps)
W_t = np.cumsum(dW, axis=1) 

# Вычисление интегралов Ито для каждого из трех процессов

# Процесс a
X_T_a = np.cumsum((time_array**2) * dW, axis=1)[:, -1]

# Процесс b
X_T_b = (1 - T) * np.cumsum((1 / (1 - time_array)) * dW, axis=1)[:, -1]

# Процесс c
X_T_c = np.cumsum((W_t**2) * dW, axis=1)[:, -1]

# Построение графиков
plt.figure(figsize=(15, 10))

for i, (X_T, title) in enumerate(zip([X_T_a, X_T_b, X_T_c], 
                                     ['Процесс a: X_T = ∫ t^2 dW_t', 
                                      'Процесс b: X_T = (1 - T) ∫ 1/(1 - t) dW_t', 
                                      'Процесс c: X_T = ∫ W_t^2 dW_t']), 
                                 start=1):
    plt.subplot(3, 1, i)
    for j in range(5):  
        plt.plot(time_array, np.cumsum((time_array**2) * dW[j, :], axis=0) if i == 1 else 
                 (1 - time_array) * np.cumsum((1 / (1 - time_array)) * dW[j, :], axis=0) if i == 2 else 
                 np.cumsum((W_t[j, :]**2) * dW[j, :], axis=0), 
                 label=f'Траектория {j+1}')
    plt.xlabel("Time")
    plt.ylabel("X_t")
    plt.title(f"5 произвольных траекторий для {title}")
    plt.legend()

plt.tight_layout()
plt.show()

# Вычисление X_0.5 для каждого процесса
T_half_index = int(0.5 / dt) 
X_half_a = np.cumsum((time_array[:T_half_index]**2) * dW[:, :T_half_index], axis=1)[:, -1]
X_half_b = (1 - 0.5) * np.cumsum((1 / (1 - time_array[:T_half_index])) * dW[:, :T_half_index], axis=1)[:, -1]
X_half_c = np.cumsum((W_t[:, :T_half_index]**2) * dW[:, :T_half_index], axis=1)[:, -1]

# Построение гистограмм и тест Колмогорова-Смирнова для нормальности
# распределения
plt.figure(figsize=(18, 6))

for i, (X_half, title) in enumerate(zip([X_half_a, X_half_b, X_half_c], 
                                        ['Процесс a: X_T = ∫ t^2 dW_t', 
                                         'Процесс b: X_T = (1 - T) ∫ 1/(1 - t) dW_t', 
                                         'Процесс c: X_T = ∫ W_t^2 dW_t']), 
                                    start=1):
    plt.subplot(1, 3, i)
    plt.hist(X_half, bins=50, density=True, alpha=0.6, color='blue')
    mu, std = np.mean(X_half), np.std(X_half)
    xmin, xmax = plt.xlim()
    x = np.linspace(xmin, xmax, 100)
    p = norm.pdf(x, mu, std)
    plt.plot(x, p, 'k', linewidth=2)
    plt.title(f"Гистограмма для {title} в T=0.5")

# Тест Колмогорова-Смирнова
    ks_stat, p_value = kstest(X_half, 'norm', args=(mu, std))
    print(f"{title} в T=0.5: K-S test statistic = {ks_stat}, p-value = {p_value}")
    if p_value > 0.05:
        print(f"{title} в T=0.5: Распределение близко к нормальному (не отвергаем гипотезу нормальности)")
    else:
        print(f"{title} в T=0.5: Распределение отличается от нормального (отвергаем гипотезу нормальности)")

plt.tight_layout()
plt.show()