Практика 08 · Градієнтний спуск

📖 Теорія: відкрий lecture.html у цій же теці. 📝 Домашнє завдання: homework.md · 🧪 Тест: quiz.html

Наскрізний приклад той самий, що в лекції: дошка оголошень про продаж вживаних телефонів. Підбираємо два числа моделі ціна = w · (рік − 2018) + b.

Що зробимо:

1. Зберемо дошку оголошень і подивимось на функцію втрат
2. Порахуємо похідну двома способами — чисельно й формулою — і звіримо
3. Напишемо градієнтний спуск з нуля, десятком рядків
4. Порівняємо знайдені w і b з аналітичною формулою і з LinearRegression
5. Проженемо три швидкості навчання й побачимо три долі
6. Додамо стохастичний спуск і порахуємо, скільки кроків він встигає зробити

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

from sklearn.linear_model import LinearRegression

# один і той самий генератор — щоб числа в тебе збіглися з числами тут
rng = np.random.default_rng(42)

print("numpy      ", np.__version__)
print("генератор випадкових чисел зафіксовано: default_rng(42)")

numpy       2.5.1
генератор випадкових чисел зафіксовано: default_rng(42)

1 · Дошка оголошень

Спершу «правда», якої в реальних даних не видно: справедлива ціна телефона. Вона залежить від року випуску, обсягу памʼяті й стану. Далі частина продавців — шахраї: вони ставлять ціну помітно нижчу за справедливу. Це і є шум, з яким доведеться жити нашій моделі.

КІЛЬКІСТЬ_ОГОЛОШЕНЬ = 700
ЧАСТКА_ШАХРАЇВ = 0.25

моделі_телефонів = np.array(["Galaxy A54", "iPhone 12", "Redmi Note 12", "Pixel 7"])
назви_станів = np.array(["потертий", "добрий", "новий"])
надбавка_за_стан = {"новий": 1500, "добрий": 0, "потертий": -1200}

рік = rng.integers(2018, 2024, КІЛЬКІСТЬ_ОГОЛОШЕНЬ)
памʼять_гб = rng.choice([64, 128, 256], КІЛЬКІСТЬ_ОГОЛОШЕНЬ)
стан = rng.choice(назви_станів, КІЛЬКІСТЬ_ОГОЛОШЕНЬ)
модель = rng.choice(моделі_телефонів, КІЛЬКІСТЬ_ОГОЛОШЕНЬ)

справедлива_ціна = (2900
                    + 1750 * (рік - 2018)
                    + 21 * (памʼять_гб - 64)
                    + np.array([надбавка_за_стан[c] for c in стан])
                    + rng.normal(0, 500, КІЛЬКІСТЬ_ОГОЛОШЕНЬ))

шахрай = (rng.random(КІЛЬКІСТЬ_ОГОЛОШЕНЬ) < ЧАСТКА_ШАХРАЇВ).astype(int)

# шахрайська знижка глибока, чесний торг — дрібний
знижка = np.where(шахрай == 1,
                  rng.uniform(0.35, 0.95, КІЛЬКІСТЬ_ОГОЛОШЕНЬ),
                  rng.uniform(0.85, 1.15, КІЛЬКІСТЬ_ОГОЛОШЕНЬ))
ціна_грн = np.round(справедлива_ціна * знижка, -1)

оголошення = pd.DataFrame({
    "модель": модель,
    "рік": рік,
    "стан": стан,
    "памʼять_гб": памʼять_гб,
    "ціна_грн": ціна_грн.astype(int),
    "шахрай": шахрай,
})

print(оголошення.head(6).to_string(index=False))
print()
print("усього оголошень:", len(оголошення))
print("шахрайських серед них:", int(шахрай.sum()))

       модель  рік     стан  памʼять_гб  ціна_грн  шахрай
      Pixel 7 2018    новий         128      6070       0
Redmi Note 12 2022 потертий         128     11070       0
   Galaxy A54 2021 потертий          64      2730       1
   Galaxy A54 2020 потертий          64      6410       0
      Pixel 7 2020 потертий          64      5320       0
    iPhone 12 2023 потертий         128      9180       1

усього оголошень: 700
шахрайських серед них: 173

2 · Дві змінні, з якими працюватимемо

Модель бачить лише одну ознаку — вік телефона. Ціни переводимо в тисячі гривень: з ними зручніше рахувати руками й дивитись на числа.

# x — скільки років телефону від 2018-го; y — ціна в тисячах гривень
x = (рік - 2018).astype(float)
y = ціна_грн / 1000.0

print("x: від", x.min(), "до", x.max(), "· середнє", round(x.mean(), 3))
print("y: від", round(y.min(), 2), "до", round(y.max(), 2),
      "тис. грн · середнє", round(y.mean(), 3))

x: від 0.0 до 5.0 · середнє 2.497
y: від 0.94 до 18.49 тис. грн · середнє 8.278

3 · Функція втрат

Функція втрат — це одне число, яке каже, наскільки погана конкретна пара (w, b). Беремо середній квадрат помилки:

$$L(w, b) = \frac{1}{n} \sum_i \left( y_i - (w x_i + b) \right)^2$$

Порахуємо її для кількох навмання взятих пар — просто щоб побачити, що це звичайне число, яке міняється разом із параметрами.

def втрата(w, b):
    "Середній квадрат помилки для пари параметрів (w, b)."
    помилки = y - (w * x + b)
    return np.mean(помилки ** 2)


for пара in [(0.0, 0.0), (1.0, 4.0), (1.6, 4.2), (3.0, 2.0)]:
    print(f"w = {пара[0]:>4}   b = {пара[1]:>4}   L = {втрата(*пара):8.3f}")

w =  0.0   b =  0.0   L =   83.036
w =  1.0   b =  4.0   L =   11.099
w =  1.6   b =  4.2   L =    6.811
w =  3.0   b =  2.0   L =   13.897

4 · Похідна двома способами

Це найважливіша клітинка всієї практики. Похідна відповідає на питання: якщо трохи зсунути w праворуч, втрата зросте чи впаде і наскільки різко?

Спосіб перший, чисельний. Буквально зсуваємо w на крихітну величину й дивимось, як змінилась втрата.

Спосіб другий, за формулою. Виведення є в лекції, результат такий:

$$\frac{\partial L}{\partial w} = -\frac{2}{n} \sum_i x_i e_i, \qquad \frac{\partial L}{\partial b} = -\frac{2}{n} \sum_i e_i, \qquad e_i = y_i - (w x_i + b)$$

Обидва способи мають дати те саме число. Це і є момент, коли похідна перестає бути магією.

def похідна_чисельно(w, b, крок=1e-6):
    "Наскільки зміниться втрата, якщо зсунути w на крихітний крок праворуч."
    return (втрата(w + крок, b) - втрата(w, b)) / крок


def градієнт(w, b):
    "Обидві похідні одразу: по w і по b. Це і називають градієнтом."
    помилки = y - (w * x + b)
    похідна_по_w = -2 * np.mean(x * помилки)
    похідна_по_b = -2 * np.mean(помилки)
    return похідна_по_w, похідна_по_b


w_проба, b_проба = 1.0, 4.0
чисельна = похідна_чисельно(w_проба, b_проба)
аналітична, _ = градієнт(w_проба, b_проба)

print(f"чисельно  : {чисельна:.6f}")
print(f"формулою  : {аналітична:.6f}")
print(f"різниця   : {abs(чисельна - аналітична):.2e}")

assert np.allclose(чисельна, аналітична, atol=1e-3), "похідні розійшлися!"
print("✅ два незалежні способи дали те саме число")

чисельно  : -12.533991
формулою  : -12.534000
різниця   : 9.18e-06
✅ два незалежні способи дали те саме число

5 · Градієнтний спуск з нуля

Уся суть методу — три рядки всередині циклу: порахуй градієнт, зсунь параметри проти нього, повтори. Важливо: обидва параметри оновлюються зі старих значень, тому спершу рахуємо обидві похідні й лише потім присвоюємо нові w і b.

def градієнтний_спуск(швидкість_навчання, кількість_ітерацій):
    "Повертає знайдені (w, b) та історію втрати на кожній ітерації."
    w, b = 0.0, 0.0
    історія_втрат = []
    for _ in range(кількість_ітерацій):
        похідна_по_w, похідна_по_b = градієнт(w, b)      # спершу обидві похідні
        w = w - швидкість_навчання * похідна_по_w        # і лише потім оновлення
        b = b - швидкість_навчання * похідна_по_b
        # захист від розльоту: далі рахувати немає сенсу, буде переповнення
        if not np.isfinite(w) or not np.isfinite(b) or abs(w) > 1e12:
            історія_втрат.append(np.inf)
            break
        історія_втрат.append(втрата(w, b))
    return w, b, np.array(історія_втрат)


w_спуск, b_спуск, історія = градієнтний_спуск(0.05, 1000)

print(f"після 1000 ітерацій: w = {w_спуск:.6f}   b = {b_спуск:.6f}")
print(f"втрата на старті  : {втрата(0.0, 0.0):.4f}")
print(f"втрата наприкінці : {історія[-1]:.4f}")

після 1000 ітерацій: w = 1.618275   b = 4.236765
втрата на старті  : 83.0359
втрата наприкінці : 6.8029

Подивимось, як падала втрата. Вісь по вертикалі логарифмічна — інакше видно було б лише перший обвал, а вся тонка робота наприкінці злилась би в лінію.

plt.figure(figsize=(7, 3.6))
plt.plot(історія, color="#c2185b")
plt.yscale("log")
plt.xlabel("ітерація")
plt.ylabel("втрата L (лог. шкала)")
plt.title("Як падає втрата під час спуску")
plt.grid(alpha=0.3)
plt.show()

print("втрата на ітераціях 1, 10, 100, 1000:",
      [round(float(історія[i]), 4) for i in [0, 9, 99, 999]])

втрата на ітераціях 1, 10, 100, 1000: [10.5302, 8.9595, 6.8122, 6.8029]

6 · Звірка: спуск, формула і scikit-learn

Тепер найцінніше. Ті самі коефіцієнти знайдемо ще двома способами:

• аналітично — методом найменших квадратів (np.linalg.lstsq розвʼязує систему напряму, без обернення матриці);
• бібліотечно — через LinearRegression.

Якщо всередині бібліотеки немає магії, всі три відповіді мають збігтися.

# матриця ознак: перший стовпець — сама ознака, другий — одиниці для вільного члена
матриця = np.c_[x, np.ones_like(x)]
w_формула, b_формула = np.linalg.lstsq(матриця, y, rcond=None)[0]

модель = LinearRegression().fit(x.reshape(-1, 1), y)
w_бібліотека = модель.coef_[0]
b_бібліотека = модель.intercept_

print(f"{'спосіб':<22}{'w':>12}{'b':>12}{'втрата':>12}")
print("-" * 58)
print(f"{'градієнтний спуск':<22}{w_спуск:>12.6f}{b_спуск:>12.6f}{втрата(w_спуск, b_спуск):>12.6f}")
print(f"{'формула МНК':<22}{w_формула:>12.6f}{b_формула:>12.6f}{втрата(w_формула, b_формула):>12.6f}")
print(f"{'LinearRegression':<22}{w_бібліотека:>12.6f}{b_бібліотека:>12.6f}{втрата(w_бібліотека, b_бібліотека):>12.6f}")

assert np.allclose([w_формула, b_формула], [w_бібліотека, b_бібліотека]), "формула розійшлася з бібліотекою!"
assert np.allclose([w_спуск, b_спуск], [w_формула, b_формула], atol=1e-4), "спуск не дійшов до розвʼязку!"
print("\n✅ усі три способи дали ті самі числа")

спосіб                           w           b      втрата
----------------------------------------------------------
градієнтний спуск         1.618275    4.236765    6.802919
формула МНК               1.618275    4.236765    6.802919
LinearRegression          1.618275    4.236765    6.802919

✅ усі три способи дали ті самі числа

Прочитати результат можна так: кожен рік новизни додає до ціни приблизно w тисяч гривень, а телефон 2018 року випуску коштує в середньому b тисяч.

print(f"кожен рік новизни додає  {w_спуск * 1000:>7.0f} грн")
print(f"телефон 2018 року коштує {b_спуск * 1000:>7.0f} грн")
print(f"типова помилка прогнозу  {np.sqrt(втрата(w_спуск, b_спуск)) * 1000:>7.0f} грн (RMSE)")

кожен рік новизни додає     1618 грн
телефон 2018 року коштує    4237 грн
типова помилка прогнозу     2608 грн (RMSE)

7 · Три швидкості навчання — три долі

Швидкість навчання ми задаємо самі, і це найважливіший важіль. Проженемо три значення на однаковій кількості ітерацій і подивимось, що вийде.

Для квадратичної втрати межу збіжності можна порахувати точно: спуск не розлітається, поки швидкість менша за 2 / c, де c — найбільша власна кривина функції втрат.

# найбільша кривина чаші — найбільше власне число матриці других похідних
матриця_кривин = 2 * np.array([[np.mean(x ** 2), np.mean(x)],
                               [np.mean(x),      1.0]])
найбільша_кривина = np.linalg.eigvalsh(матриця_кривин).max()
межа_збіжності = 2 / найбільша_кривина

print(f"найбільша кривина : {найбільша_кривина:.4f}")
print(f"межа збіжності    : {межа_збіжності:.4f}")

найбільша кривина : 19.7641
межа збіжності    : 0.1012

результати = []
for швидкість in [0.001, 0.05, 0.15]:
    w_знайдене, b_знайдене, історія_швидкості = градієнтний_спуск(швидкість, 300)
    фінальна_втрата = історія_швидкості[-1]
    результати.append({
        "швидкість": швидкість,
        "w": w_знайдене,
        "b": b_знайдене,
        "втрата": фінальна_втрата,
        "ітерацій": len(історія_швидкості),
    })

таблиця = pd.DataFrame(результати)
print(таблиця.to_string(index=False))
print()
print(f"для довідки: правильна відповідь w = {w_формула:.4f}, b = {b_формула:.4f}, "
      f"втрата = {втрата(w_формула, b_формула):.4f}")

 швидкість             w             b   втрата  ітерацій
     0.001  2.435752e+00  1.304038e+00 9.565028       300
     0.050  1.618388e+00  4.236363e+00 6.802919       300
     0.150 -1.401839e+12 -3.941195e+11      inf        40

для довідки: правильна відповідь w = 1.6183, b = 4.2368, втрата = 6.8029

Читаємо таблицю рядок за рядком:

• 0.001 — метод працює правильно, просто повзе: за 300 ітерацій до відповіді так і не дійшов;
• 0.05 — збігся точно в мінімум;
• 0.15 — більша за межу збіжності: кожен крок перестрибує дно й опиняється вище, ніж був. Цикл обірвався достроково, бо числа перестали бути скінченними.

Ось ті самі три долі на одному графіку.

plt.figure(figsize=(7, 3.8))
for швидкість, колір in [(0.001, "#0f766e"), (0.05, "#c2185b"), (0.12, "#c2620f")]:
    _, _, історія_швидкості = градієнтний_спуск(швидкість, 60)
    скінченні = історія_швидкості[np.isfinite(історія_швидкості)]
    plt.plot(скінченні, color=колір, label=f"швидкість {швидкість}")
plt.yscale("log")
plt.xlabel("ітерація")
plt.ylabel("втрата L (лог. шкала)")
plt.title("Замала, добра і завелика швидкість навчання")
plt.legend()
plt.grid(alpha=0.3)
plt.show()

print("бірюзова повзе, рожева падає й лягає на дно, бурштинова летить угору")

бірюзова повзе, рожева падає й лягає на дно, бурштинова летить угору

8 · Стохастичний спуск

Повний градієнт на кожному кроці проходить усі 700 оголошень. Стохастичний бере випадкову партію — одне оголошення або кілька десятків — і одразу робить крок. Напрямок виходить кривий, зате кроків за той самий обчислювальний бюджет набагато більше.

Бюджет міряємо в епохах: одна епоха — це один повний перегляд усіх даних.

def стохастичний_спуск(швидкість_навчання, кількість_епох, розмір_партії, seed=0):
    "Спуск партіями. Повертає (w, b), кількість кроків та історію втрати."
    генератор = np.random.default_rng(seed)
    w, b = 0.0, 0.0
    кроків = 0
    історія_втрат = []
    for _ in range(кількість_епох):
        порядок = генератор.permutation(len(x))       # щоразу новий порядок обʼєктів
        for початок in range(0, len(x), розмір_партії):
            індекси = порядок[початок:початок + розмір_партії]
            x_партії, y_партії = x[індекси], y[індекси]
            помилки = y_партії - (w * x_партії + b)
            # градієнт рахуємо тільки по партії — саме тут і виникає економія
            w = w + швидкість_навчання * 2 * np.mean(x_партії * помилки)
            b = b + швидкість_навчання * 2 * np.mean(помилки)
            кроків += 1
        історія_втрат.append(втрата(w, b))            # втрату міряємо на всіх даних
    return w, b, кроків, np.array(історія_втрат)


w_сгд, b_сгд, кроків_сгд, історія_сгд = стохастичний_спуск(0.005, 5, розмір_партії=1)

print("стохастичний спуск, 5 епох, партія 1")
print(f"кроків зроблено : {кроків_сгд}")
print(f"знайдено        : w = {w_сгд:.4f}, b = {b_сгд:.4f}")
print(f"втрата          : {втрата(w_сгд, b_сгд):.4f}   (мінімум {втрата(w_формула, b_формула):.4f})")

# для порівняння: скільки встиг повний спуск за той самий бюджет
w_за_5, b_за_5, _ = градієнтний_спуск(0.005, 5)
print()
print(f"повний спуск за ті самі 5 епох робить 5 кроків, втрата {втрата(w_за_5, b_за_5):.4f}")

стохастичний спуск, 5 епох, партія 1
кроків зроблено : 3500
знайдено        : w = 2.0815, b = 4.3107
втрата          : 8.9496   (мінімум 6.8029)

повний спуск за ті самі 5 епох робить 5 кроків, втрата 36.1719

Тепер найголовніше порівняння: однаковий бюджет, різна кількість кроків. За 5 епох повний спуск робить 5 кроків, партія по 32 — 110 кроків, партія по 8 — 440, партія по 1 — 3500. Обчислень витрачено стільки ж у всіх чотирьох.

ЕПОХ = 5
ШВИДКІСТЬ = 0.005

порівняння = []
# повний градієнт: одна епоха = один крок
w_повний, b_повний, історія_повного = градієнтний_спуск(ШВИДКІСТЬ, ЕПОХ)
порівняння.append({"режим": "повний (700 на крок)", "кроків": ЕПОХ,
                   "w": w_повний, "b": b_повний, "втрата": втрата(w_повний, b_повний)})

for розмір in [256, 32, 8, 1]:
    w_ч, b_ч, кроків_ч, _ = стохастичний_спуск(ШВИДКІСТЬ, ЕПОХ, розмір_партії=розмір)
    порівняння.append({"режим": f"партія {розмір}", "кроків": кроків_ч,
                       "w": w_ч, "b": b_ч, "втрата": втрата(w_ч, b_ч)})

print(pd.DataFrame(порівняння).to_string(index=False))
print()
print(f"мінімально можлива втрата: {втрата(w_формула, b_формула):.4f}")

               режим  кроків        w        b    втрата
повний (700 на крок)       5 1.041518 0.348820 36.171938
          партія 256      15 2.009923 0.732037 13.638878
           партія 32     110 2.344111 1.708248  8.866733
            партія 8     440 2.015809 3.328319  7.275310
            партія 1    3500 2.081533 4.310730  8.949623

мінімально можлива втрата: 6.8029

Читаємо таблицю згори вниз. Повний спуск за свої 5 кроків майже нічого не встиг. Що менша партія, то більше кроків і то нижча втрата — але лише до певної межі. Партія 1 робить у 700 разів більше кроків за повний спуск і все одно програє партії 8: кожен її крок настільки шумний, що траєкторія перестає прицільно йти на дно й починає тремтіти навколо нього.

Це і є справжня відповідь на питання «яку партію брати»: не найменшу, а таку, де виграш від кількості кроків ще перекриває шкоду від шуму.

Подивись на це збоку — як поводиться втрата за епохами при партії 1 і при партії 32.

plt.figure(figsize=(7, 3.8))
for розмір, колір in [(1, "#c2185b"), (32, "#0f766e")]:
    _, _, _, історія_розміру = стохастичний_спуск(0.005, 30, розмір_партії=розмір)
    plt.plot(range(1, 31), історія_розміру, color=колір, label=f"партія {розмір}")
plt.axhline(втрата(w_формула, b_формула), color="#5b6b7c", linestyle="--", label="мінімум")
plt.xlabel("епоха")
plt.ylabel("втрата L на всіх даних")
plt.title("Партія 1 шумить, партія 32 сідає рівно")
plt.legend()
plt.grid(alpha=0.3)
plt.show()

_, _, _, історія_один = стохастичний_спуск(0.005, 30, розмір_партії=1)
print("партія 1, втрата на епохах 10, 20, 30:",
      [round(float(історія_один[i]), 4) for i in [9, 19, 29]])
print("розмах коливань наприкінці:",
      round(float(історія_один[20:].max() - історія_один[20:].min()), 4))

партія 1, втрата на епохах 10, 20, 30: [7.7201, 7.605, 6.8792]
розмах коливань наприкінці: 1.4394

Рожева крива не сідає в точку: біля мінімуму шум не зникає, і параметри вічно тремтять у невеликій області навколо відповіді. Саме тому на практиці або беруть партію більшу за одиницю, або зменшують швидкість навчання під кінець навчання.

9 · Що ми показали

• похідна, порахована чисельно й за формулою, дає те саме число;
• градієнтний спуск із нуля знаходить ті самі w і b, що й формула найменших квадратів та LinearRegression;
• швидкість навчання має вузьке робоче вікно: замала — повземо, завелика — розліт;
• за однаковий обчислювальний бюджет спуск партіями виграє в повного, але найдрібніша партія не найкраща: занадто шумні кроки з'їдають перевагу.

---

Завдання

🟢 Рівень 1

Додай до градієнтний_спуск критерій зупинки: якщо втрата за ітерацію зменшилась менш ніж на 1e-9, вихід із циклу. Скільки ітерацій знадобиться при швидкості 0.05? А при 0.01? Надрукуй обидва числа.

🟡 Рівень 2

Побудуй карту втрат, як у лекції. Перебери сітку значень w від −1 до 4 і b від −1 до 9 (наприклад, по 120 точок на вісь), порахуй втрату в кожному вузлі й покажи plt.contourf. Поверх карти намалюй траєкторію спуску при швидкості 0.005 за 60 кроків: збери список (w, b) усередині циклу. Опиши словами, чому траєкторія йде не прямо на мінімум.

🔴 Рівень 3

Додай до спуску момент: заведи змінні швидкість_w і швидкість_b, на кожній ітерації онови їх як бета * швидкість + похідна, а параметри зсувай на −швидкість_навчання * швидкість. Порівняй кількість ітерацій до втрати 6.81 для звичайного спуску й для спуску з моментом бета = 0.9 при однаковій швидкості навчання 0.005. Постав assert, що обидва варіанти приходять до тих самих w і b з точністю 1e-3.