---

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

# =======================================================
# 0. Загрузка данных
# =======================================================
# Загружаем датасет (файл должен лежать в той же папке, что и блокнот)
gipers = pd.read_csv("Гиперспектр кукурузы.csv")

# Извлекаем признаки в виде numpy-массивов
x = gipers['waveleng'].to_numpy(dtype=float)
Y = gipers['Spectr'].to_numpy(dtype=float)

# Сетка точек для гладкого построения линий регрессии
x_space = np.linspace(x.min(), x.max(), 300)


# =======================================================
# 1. Аналитическое вычисление коэффициентов (МНК)
# =======================================================
# Формулы коэффициентов парной линейной регрессии через ковариацию и дисперсию
a1 = ((x - x.mean()) * (Y - Y.mean())).mean() / (((x - x.mean()) ** 2).mean())
a0 = Y.mean() - a1 * x.mean()

print(f"1. Аналитическая модель: Y^ = {a0:.6f} + ({a1:.6f}) * x")

# Прогноз по аналитической модели
Y_pred_analyt = a0 + a1 * x_space

# График аналитической модели
plt.figure(figsize=(9, 5))
plt.scatter(x, Y, color='steelblue', alpha=0.4, label='Обучающая выборка')
plt.plot(x_space, Y_pred_analyt, color='crimson', linewidth=2, label=f'МНК: Y = {a0:.2f} + {a1:.4f}x')
plt.title('Парная линейная регрессия (Аналитический метод / МНК)')
plt.xlabel('waveleng')
plt.ylabel('Spectr')
plt.grid(True, linestyle='--', alpha=0.6)
plt.legend()
plt.show()


# =======================================================
# 2. Оценка параметров методом градиентного спуска (GD)
# =======================================================
class SimpleRegression(object):
    def __init__(self):
        self.a0 = 0.0
        self.a1 = 0.0

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

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

    def fit(self, x, Y, alpha=0.01, epsylon=1e-5, max_steps=4000):
        """
        Для численной устойчивости и быстрой сходимости градиентного спуска
        выполняется внутренняя стандартизация признака x.
        """
        x_mean = x.mean()
        x_std = x.std()
        x_norm = (x - x_mean) / x_std
        n = len(x)

        # Начальная инициализация параметров
        w0 = Y.mean()
        w1 = 0.0

        steps, errors = [], []

        for step in range(1, max_steps + 1):
            y_pred = w0 + w1 * x_norm
            diff = Y - y_pred

            # Градиенты MSE по параметрам
            grad_w0 = -2 * np.sum(diff) / n
            grad_w1 = -2 * np.sum(diff * x_norm) / n

            # Шаг градиентного спуска
            w0 -= alpha * grad_w0
            w1 -= alpha * grad_w1

            # Перевод весов обратно к исходному масштабу x
            self.a1 = w1 / x_std
            self.a0 = w0 - (w1 * x_mean) / x_std

            current_mse = self.MSE(x, Y)
            steps.append(step)
            errors.append(current_mse)

            # Проверка критерия остановки
            if len(errors) > 1 and abs(errors[-2] - errors[-1]) < epsylon:
                break

        return steps, errors

regr = SimpleRegression()
steps, errors = regr.fit(x, Y, alpha=0.01, epsylon=1e-6, max_steps=4000)

print(f"2. Модель градиентного спуска: Y^ = {regr.a0:.6f} + ({regr.a1:.6f}) * x")
print(f"   Количество итераций: {len(steps)}")
print(f"   Финальная MSE: {regr.MSE(x, Y):.6f}")

# График регрессии после градиентного спуска
Y_pred_gd = regr.predict(x_space)

plt.figure(figsize=(9, 5))
plt.scatter(x, Y, color='steelblue', alpha=0.4, label='Обучающая выборка')
plt.plot(x_space, Y_pred_gd, color='forestgreen', linewidth=2, label=f'GD: Y = {regr.a0:.2f} + {regr.a1:.4f}x')
plt.title('Парная линейная регрессия (Градиентный спуск)')
plt.xlabel('waveleng')
plt.ylabel('Spectr')
plt.grid(True, linestyle='--', alpha=0.6)
plt.legend()
plt.show()

# График кривой обучения (сходимость функции потерь)
plt.figure(figsize=(9, 4))
plt.plot(steps, errors, color='darkorange', linewidth=2)
plt.title('Кривая обучения (изменение MSE в процессе спуска)')
plt.xlabel('Номер шага')
plt.ylabel('MSE')
plt.grid(True, linestyle='--', alpha=0.6)
plt.show()


# =======================================================
# 3. Сравнение качества построенных моделей
# =======================================================
# Расчет метрик качества на выборке
mse_analyt = ((Y - (a0 + a1 * x)) ** 2).mean()
mse_gd = regr.MSE(x, Y)

mae_analyt = np.abs(Y - (a0 + a1 * x)).mean()
mae_gd = np.abs(Y - regr.predict(x)).mean()

# Итоговое сопоставление
print("--- Сравнение качества моделей ---")
print(f"MSE  (Аналитическая): {mse_analyt:.6f} | MAE: {mae_analyt:.6f}")
print(f"MSE  (Градиентный спуск): {mse_gd:.6f} | MAE: {mae_gd:.6f}")
print(f"Разница в MSE: {abs(mse_analyt - mse_gd):.8f}")

# Совмещенный график выборки и обоих прогнозов
plt.figure(figsize=(10, 6))
plt.scatter(x, Y, color='gray', alpha=0.35, s=20, label='Обучающая выборка')
plt.plot(x_space, Y_pred_analyt, color='red', linewidth=2.5, label='Аналитическое решение (МНК)')
plt.plot(x_space, Y_pred_gd, color='lime', linestyle='--', linewidth=2, label='Градиентный спуск')
plt.title('Сравнение прогнозов моделей на обучающей выборке')
plt.xlabel('waveleng')
plt.ylabel('Spectr')
plt.grid(True, linestyle='--', alpha=0.6)
plt.legend()
plt.show()