---

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# ------------------------------------------------------------
# 0. Загрузка данных с правильными параметрами
# ------------------------------------------------------------
gipers = pd.read_csv(r"/Гиперспектр кукурузы.csv", 
                     sep=';',               # разделитель столбцов
                     decimal=',')           # десятичный разделитель

print("Названия столбцов:", gipers.columns.tolist())
print("Первые 5 строк:")
print(gipers.head())

# Извлекаем данные (обратите внимание на пробел в ' Spectr')
x = gipers['wavelength'].values
Y = gipers[' Spectr'].values   # или переименуйте столбец: gipers.rename(columns={' Spectr':'Spectr'}, inplace=True)

print(f"\nРазмер выборки: {len(x)} точек")
print(f"Диапазон x: [{x.min():.2f}, {x.max():.2f}]")
print(f"Диапазон Y: [{Y.min():.2f}, {Y.max():.2f}]")

# ------------------------------------------------------------
# 1. Аналитическое решение (МНК)
# ------------------------------------------------------------
a1 = ((x - x.mean()) * (Y - Y.mean())).mean() / ((x - x.mean())**2).mean()
a0 = Y.mean() - a1 * x.mean()
print(f"\nАналитическая модель: Y^ = {a0:.4f} + {a1:.4f} * x")

x_space = np.linspace(x.min(), x.max(), 100)
Y_pred_analyt = a0 + a1 * x_space

plt.figure(figsize=(10,6))
plt.scatter(x, Y, alpha=0.6, label='Обучающая выборка')
plt.plot(x_space, Y_pred_analyt, 'r', linewidth=2, label='МНК прогноз')
plt.title("Аналитическая парная линейная регрессия")
plt.xlabel("wavelength")
plt.ylabel("Spectr")
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# ------------------------------------------------------------
# 2. Класс для градиентного спуска
# ------------------------------------------------------------
class SimpleRegression:
    def __init__(self):
        self.a0 = 0
        self.a1 = 0
        self.history = []

    def predict(self, x):
        return self.a0 + self.a1 * x

    def MSE(self, x, Y):
        return ((Y - self.predict(x)) ** 2).mean()

    def MAE(self, x, Y):
        return np.abs(Y - self.predict(x)).mean()

    def fit(self, x, Y, alpha=0.001, epsylon=0.01, max_steps=5000):
        steps, errors = [], []
        step = 0
        self.history = []
        for _ in range(max_steps):
            dT_a0 = -2 * sum(Y - self.predict(x))
            dT_a1 = -2 * sum((Y - self.predict(x)) * x)
            self.a0 -= alpha * dT_a0
            self.a1 -= alpha * dT_a1
            new_error = self.MSE(x, Y)
            step += 1
            steps.append(step)
            errors.append(new_error)
            self.history.append(new_error)
            if new_error < epsylon:
                break
        return steps, errors

# ------------------------------------------------------------
# 3. Запуск градиентного спуска с разными alpha
# ------------------------------------------------------------
print("\n" + "="*60)
print("ГРАДИЕНТНЫЙ СПУСК (сравнение alpha)")
print("="*60)

# Подбираем alpha под масштаб данных (x ~ 400..1000)
alphas = [1e-9, 1e-8, 1e-7]   # пробуем разные порядки
models = []
mse_list = []

for alpha in alphas:
    regr = SimpleRegression()
    steps, errors = regr.fit(x, Y, alpha=alpha, epsylon=0.01, max_steps=5000)
    mse = regr.MSE(x, Y)
    models.append(regr)
    mse_list.append(mse)
    print(f"alpha = {alpha:.1e} -> MSE = {mse:.6f}, шагов: {len(steps)}")

best_idx = np.argmin(mse_list)
best_regr = models[best_idx]
best_alpha = alphas[best_idx]
print(f"\nЛучшая модель: alpha = {best_alpha:.1e}, MSE = {mse_list[best_idx]:.6f}")

Y_pred_gd = best_regr.predict(x_space)

plt.figure(figsize=(10,6))
plt.scatter(x, Y, alpha=0.6, label='Обучающая выборка')
plt.plot(x_space, Y_pred_gd, 'g', linewidth=2, label=f'GD прогноз (alpha={best_alpha:.1e})')
plt.title("Парная регрессия (градиентный спуск)")
plt.xlabel("wavelength")
plt.ylabel("Spectr")
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# ------------------------------------------------------------
# 4. Кривые обучения
# ------------------------------------------------------------
plt.figure(figsize=(10,6))
for i, (regr, alpha) in enumerate(zip(models, alphas)):
    plt.plot(range(1, len(regr.history)+1), regr.history, 
             label=f'alpha={alpha:.1e}', linewidth=2)
plt.title("Кривые обучения (MSE от шага)")
plt.xlabel("Шаг")
plt.ylabel("MSE")
plt.legend()
plt.grid(True, alpha=0.3)
plt.yscale('log')
plt.show()

# ------------------------------------------------------------
# 5. Сравнение моделей
# ------------------------------------------------------------
plt.figure(figsize=(12,6))
plt.scatter(x, Y, color='blue', alpha=0.4, s=30, label='Обучающая выборка')
plt.plot(x_space, Y_pred_analyt, 'r', linewidth=2.5, label='Аналитическое решение (МНК)')
plt.plot(x_space, Y_pred_gd, 'g--', linewidth=2.5, label=f'Градиентный спуск (alpha={best_alpha:.1e})')
plt.title("Сравнение аналитической модели и градиентного спуска")
plt.xlabel("wavelength")
plt.ylabel("Spectr")
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# ------------------------------------------------------------
# 6. Метрики качества
# ------------------------------------------------------------
mse_analyt = ((Y - (a0 + a1*x))**2).mean()
mae_analyt = np.abs(Y - (a0 + a1*x)).mean()
r2_analyt = 1 - mse_analyt / Y.var()

mse_gd = best_regr.MSE(x, Y)
mae_gd = np.abs(Y - best_regr.predict(x)).mean()
r2_gd = 1 - mse_gd / Y.var()

print("\n" + "="*60)
print("СРАВНЕНИЕ КАЧЕСТВА МОДЕЛЕЙ")
print("="*60)
print(f"{'Метрика':<15} {'Аналитическая':<20} {'Градиентный спуск':<20}")
print("-"*55)
print(f"{'MSE':<15} {mse_analyt:<20.6f} {mse_gd:<20.6f}")
print(f"{'MAE':<15} {mae_analyt:<20.6f} {mae_gd:<20.6f}")
print(f"{'R²':<15} {r2_analyt:<20.6f} {r2_gd:<20.6f}")

print(f"\nКоэффициенты:")
print(f"  Аналитически: a0 = {a0:.4f}, a1 = {a1:.4f}")
print(f"  GD (alpha={best_alpha:.1e}): a0 = {best_regr.a0:.4f}, a1 = {best_regr.a1:.4f}")