Momentum, RMSProp আর Adam — derive করে আর তুলনা করে

AI-এর ভিত্তি কোর্সের পর্ব ৩৪-এ gradient descent-এর নিয়ম x ← x − η∇f(x), learning rate আর SGD-র ধারণা এসেছিল, আর সেখানে এক লাইনে বলা হয়েছিল যে momentum, RMSProp আর Adam এরই উন্নত রূপ। এই কোর্সের পর্ব ১১-এ Hessian আর পর্ব ১০-এ quadratic form-এর বাটির আকার দেখেছেন। এবার প্রশ্ন: সরু, লম্বা উপত্যকায় সাধারণ gradient descent কেন এত ধীর, আর এই optimizer-গুলো ঠিক কোন কৌশলে, কতটা আর কখন তাকে হারায়? সব কটি optimizer আমরা NumPy-তে নিজে লিখব, তিন-চার লাইনে। চালানো যায় এমন চারটি ব্লক পাতাতেই চলে (ব্রাউজারের Python, NumPy আর matplotlib); প্রথমবার প্যাকেজ নামতে কয়েক সেকেন্ড লাগে।

১ · সমস্যাটি — সরু উপত্যকায় gradient descent-এর এঁকেবেঁকে চলা

এই পর্বের প্রথম পরীক্ষাগার একটি সরু বাটি: f(x) = ½(x₁² + ২৫x₂²). এর Hessian হলো diag(১, ২৫), তাই দুটি eigenvalue: সবচেয়ে ছোট μ = ১ (সমতল দিক, x₁) আর সবচেয়ে বড় L = ২৫ (x₂, খাড়া দিক)। তাদের অনুপাত condition number κ = L/μ = ২৫. পর্ব ১২ (Convexity — কেন কিছু loss সহজ আর কিছু কঠিন)-এ L-কে বলা হচ্ছে smoothness আর μ-কে strong convexity-র মাপ; এখানে quadratic-এ দুটিই সরাসরি eigenvalue.

Gradient descent এই বাটিতে প্রতিটি eigen-দিক আলাদাভাবে চালায়: λ eigenvalue-এর দিকে প্রতি ধাপে দূরত্ব (১ − ηλ) গুণ হয়। খাড়া দিক স্থির রাখতে η < ২/L লাগে; অথচ সমতল দিকের গতি নির্ভর করে ημ-এর ওপর, যা তখন ২/κ-এর চেয়েও ছোট। সেরা স্থির ধাপ η = ২/(L + μ)-এ দুই দিকের গুণক সমান হয়: (κ − ১)/(κ + ১). κ = ২৫-এ তা ২৪/২৬ ≈ ০.৯২৩, অর্থাৎ প্রতি ধাপে ভুল মাত্র ৭.৭% কমে, আর খাড়া দিকে গুণক −০.৯২৩: চিহ্ন প্রতি ধাপে উল্টায়। ওটিই zig-zag. পর্ব ১৩ (Gradient descent-এর পরিবার — batch, SGD, mini-batch আর learning-rate schedule) এই হার, step size-এর সীমা আর feature scaling দিয়ে κ কমানো বিস্তারিত দেখাচ্ছে; এই পর্ব সেখান থেকে শুরু করে প্রশ্ন করে, κ না কমিয়েই কি দ্রুত যাওয়া যায়?

২ · Heavy-ball momentum — সূত্র, derive আর এক ধাপ হাতে

১৯৬৪ সালে সোভিয়েত গণিতবিদ Boris Polyak একটি প্রবন্ধে gradient descent-এ আগের ধাপের একটি অংশ যোগ করার প্রস্তাব দেন; বলটি ভারী হলে জড়তার কারণে যেমন চলে, তাই এর প্রচলিত নাম heavy-ball পদ্ধতি। আধুনিক লেখায় এটি দুই লাইনে:

v_t+1 = β v_t − η ∇f(x_t), x_t+1 = x_t + v_t+1

এখানে v হলো velocity (গতিবেগ, অর্থাৎ শেষ ধাপটি নিজেই), η হলো learning rate, আর β ∈ [০, ১) হলো momentum-এর সহগ, সাধারণত ০.৯. β = ০ হলে এটি সাধারণ gradient descent. v-কে প্রসারিত করলে (v₀ = 0 ধরে) দেখা যায় এটি আসলে অতীতের সব gradient-এর একটি ওজনদার যোগফল:

v_t+1 = −η [ g_t + β g_t−1 + β² g_t−2 + … ]

যেখানে g_t = ∇f(x_t). পুরোনো gradient-এর ওজন জ্যামিতিক হারে কমে, একে বলে exponentially weighted (চলমান) গড়। এর দুটি সরাসরি ফল। এক: gradient যদি অনেক ধাপ ধরে একই থাকে (g), তাহলে ধাপ দাঁড়ায় η(১ + β + β² + …)g = ηg/(১ − β). β = ০.৯ হলে কার্যকর ধাপ η-এর ১০ গুণ। দুই: gradient যদি প্রতি ধাপে চিহ্ন বদলায় (g, −g, g, …), তাহলে যোগফল η(১ − β + β² − …)g = ηg/(১ + β), অর্থাৎ ধাপ ছোট হয়। একই সূত্র দুই দিকে দুই রকম আচরণ করে: নালা ধরে জোরে, পাড়ের দিকে আস্তে।

নিচের stepper সরু বাটিতে η = ০.০৭ আর β = ০.৫ নিয়ে প্রথম দুটি ধাপ হাতে করে। η = ০.০৭ ইচ্ছা করে বড় নেওয়া, খাড়া দিকের সীমা ২/L = ০.০৮-এর কাছাকাছি, যাতে সাধারণ GD zig-zag করে।

Quadratic-এ গতির হিসাব — কেন (√κ − ১)/(√κ + ১)

এক eigen-দিকে, eigenvalue λ, f = (λ/২)x². দুই লাইনের সূত্র থেকে v বাদ দিলে একটি দুই-ধাপের পুনরাবৃত্তি পাওয়া যায়:

x_t+1 = (১ + β − ηλ) x_t − β x_t−1

এর আচরণ ঠিক করে characteristic সমীকরণ z² − (১ + β − ηλ)z + β = ০-এর দুটি মূল; ভুল প্রতি ধাপে মোটামুটি বড় মূলটির মাপ (|z|) গুণ হয়। মূল দুটি জটিল সংখ্যা হলে তাদের গুণফল β, তাই দুটিরই মাপ ঠিক √β, λ যা-ই হোক। Polyak-এর কৌশল: η আর β এমনভাবে বাছা যাতে μ থেকে L পর্যন্ত সব eigenvalue-তে মূল জটিল (বা ঠিক সীমায়) থাকে। হিসাব করলে আসে η = ৪/(√L + √μ)² আর β = ((√κ − ১)/(√κ + ১))², আর তখন প্রতিটি দিকের হার √β:

একই পুনরাবৃত্তি থেকে স্থিতিশীলতার সীমাও আসে: ০ ≤ β < ১ হলে পদ্ধতি স্থির থাকে যদি ηλ < ২(১ + β). β = ০.৯-এ এটি ηλ < ৩.৮, যেখানে GD-র সীমা ηλ < ২. সতর্কতা: এই সুন্দর হার quadratic-এর জন্য প্রমাণিত। সাধারণ strongly convex function-এ heavy ball একই গতির নিশ্চয়তা দেয় না; Lessard, Recht আর Packard-এর ২০১৬ সালের একটি প্রবন্ধে এমন একটি মসৃণ, strongly convex উদাহরণ দেখানো হয়েছে যেখানে এই tuned heavy ball একটি চক্রে আটকে যায়, অভিসারী হয় না।

Zig-zag কীভাবে থামে, তা কয়েক ধাপ ধরে দেখা যাক। একই η = ০.০৭ আর β = ০.৫, একই শুরু।

৩ · Nesterov — যেখানে যাবেন, সেখানকার ঢাল দেখে পা ফেলা

১৯৮৩ সালে Yurii Nesterov একটি ছোট প্রবন্ধে convex আর মসৃণ function-এর জন্য এমন একটি পদ্ধতি দেন, যার ভুল k ধাপে O(১/k²) হারে কমে, যেখানে সাধারণ gradient descent-এর O(১/k). শুধু gradient ব্যবহার করে এমন পদ্ধতির জন্য এটি তাত্ত্বিকভাবে সেরা হার (Nemirovski আর Yudin-এর lower bound-এর সঙ্গে মেলে)। এর একটি জনপ্রিয় রূপ momentum-এর সঙ্গে প্রায় একই, শুধু gradient নেওয়ার বিন্দুটি আলাদা:

v_t+1 = β v_t − η ∇f(x_t + β v_t), x_t+1 = x_t + v_t+1

x + βv হলো সেই বিন্দু, যেখানে শুধু momentum-ই আপনাকে নিয়ে যাবে। Heavy ball আগে এখানকার ঢাল দেখে, তারপর momentum-এ ভেসে যায়; Nesterov আগে ভেসে গিয়ে সেখানকার ঢাল দেখে সংশোধন করে। রিকশাচালকের উদাহরণে: সামনের মোড় দেখে আগেই ব্রেক কষা, মোড়ে পৌঁছে নয়। Overshoot শুরু হওয়ার আগেই সংশোধন আসে বলে Nesterov সাধারণত কম দোলে।

Strongly convex আর L-smooth function-এ η = ১/L আর β = (√κ − ১)/(√κ + ১) নিলে Nesterov-এর হার প্রায় (১ − ১/√κ), এবং এটি সব এমন function-এ প্রমাণিত, শুধু quadratic-এ নয়। Quadratic-এ tuned heavy ball তার চেয়ে একটু দ্রুত। নিচের ব্লক তিনটিকে সরু বাটিতে পাশাপাশি চালায়: সেরা স্থির η-এর GD, Polyak-এর tuned heavy ball, আর Nesterov. শুরু (−৯, ১.৫), থামার শর্ত ‖x‖ < ১০^−৬.

import numpy as np
import matplotlib.pyplot as plt

A = np.diag([1.0, 25.0])                  # μ = 1, L = 25, তাই κ = 25
grad = lambda x: A @ x
L, mu, x0 = 25.0, 1.0, np.array([-9.0, 1.5])
r = (np.sqrt(L / mu) - 1) / (np.sqrt(L / mu) + 1)    # (√κ − 1)/(√κ + 1)

def gd(x, v):                             # সাধারণ GD, সেরা স্থির η = 2/(L + μ)
    return x - 2 / (L + mu) * grad(x), v
def heavy_ball(x, v):                     # Polyak: β = r², η = 4/(√L + √μ)²
    v = r**2 * v - 4 / (np.sqrt(L) + np.sqrt(mu))**2 * grad(x)
    return x + v, v
def nesterov(x, v):                       # সামনে তাকিয়ে gradient, β = r, η = 1/L
    v = r * v - grad(x + r * v) / L
    return x + v, v

paths = {}
for name, step in [("GD", gd), ("Heavy ball", heavy_ball), ("Nesterov", nesterov)]:
    x, v, P = x0.copy(), np.zeros(2), [x0]
    for t in range(400):
        x, v = step(x, v)
        P.append(x)
    paths[name] = P = np.array(P)
    d = np.linalg.norm(P, axis=1)
    print(f"{name:10s} steps to |x|<1e-6: {np.argmax(d < 1e-6):3d} | |x| after 10 steps: {d[10]:.1e}")

g1, g2 = np.meshgrid(np.linspace(-10, 10, 200), np.linspace(-3, 3, 120))
plt.contour(g1, g2, 0.5 * (g1**2 + 25 * g2**2), levels=15, colors="0.8")
for name, P in paths.items():
    plt.plot(P[:30, 0], P[:30, 1], ".-", label=name)
plt.title("First 30 steps on f = (x1^2 + 25 x2^2) / 2")
plt.xlabel("x1"); plt.ylabel("x2"); plt.legend(); plt.show()

GD         steps to |x|<1e-6: 201 | |x| after 10 steps: 4.1e+00
Heavy ball steps to |x|<1e-6:  48 | |x| after 10 steps: 8.2e-01
Nesterov   steps to |x|<1e-6:  85 | |x| after 10 steps: 2.9e+00

৪ · প্রতিটি parameter-এর নিজের ধাপ — AdaGrad আর RMSProp

Momentum κ-এর ক্ষতি কমায়, কিন্তু সব parameter-এ একই η দেয়। বাস্তবে parameter-গুলোর gradient-এর মাপ খুব আলাদা হতে পারে: বাজারঘরের ডেটায় “অর্ডারের টাকা” (হাজারের ঘরে) আর “ডেলিভারির দিন” (একের ঘরে) feature-এর ওজনের gradient-এ হাজার গুণ তফাত, আর কোনো বিরল শব্দের feature-এর gradient প্রায় সবসময় শূন্য। আদর্শ সমাধান Newton-এর মতো Hessian দিয়ে ভাগ করা (পর্ব ১৫), কিন্তু তা দামি। সস্তা বিকল্প: প্রতিটি parameter-এর ধাপকে তার নিজের gradient-এর সাধারণ মাপ দিয়ে ভাগ করা। একে একটি কর্ণ-matrix দিয়ে preconditioning বলা যায়।

AdaGrad (Duchi, Hazan আর Singer, Journal of Machine Learning Research, ২০১১) প্রতিটি parameter i-এর জন্য অতীতের সব gradient-এর বর্গ যোগ করে রাখে:

G_t,i = G_t−1,i + g_t,i², x_t+1,i = x_t,i − η g_t,i / (√G_t,i + ε)

ε হলো শূন্য দিয়ে ভাগ ঠেকানোর একটি খুব ছোট সংখ্যা (যেমন ১০^−৮)। যে parameter-এর gradient বড়, তার ধাপ ছোট; যে parameter কদাচিৎ gradient পায় (বিরল feature), তার G ছোট থাকে, তাই এলে বড় ধাপ পায়। Sparse ডেটায় (text, click) এটি খুব কাজের। কিন্তু একটি দুর্বলতা স্পষ্ট: G শুধু বাড়ে। Gradient মোটামুটি স্থির g হলে √G_t ≈ |g|√t, তাই ধাপ ≈ η/√t: সময়ের সঙ্গে শূন্যের দিকে যায়, পথ বাকি থাকলেও। দীর্ঘ, non-convex training-এ এটি সমস্যা।

RMSProp এই সমস্যার একটি সরল সমাধান: যোগফলের জায়গায় g²-এর exponentially weighted গড়, যাতে পুরোনো ইতিহাস ভুলে যাওয়া যায়:

E_t,i = ρ E_t−1,i + (১ − ρ) g_t,i², x_t+1,i = x_t,i − η g_t,i / (√E_t,i + ε)

ρ সাধারণত ০.৯. নাম এসেছে root mean square থেকে: √E হলো সাম্প্রতিক gradient-এর RMS. Gradient স্থির থাকলে ধাপ ≈ η, কমে না। এর ইতিহাস একটু অদ্ভুত: RMSProp কোনো প্রকাশিত প্রবন্ধে আসেনি; সাধারণত সূত্র হিসেবে উল্লেখ করা হয় Geoffrey Hinton-এর ২০১২ সালের Coursera কোর্স “Neural Networks for Machine Learning”-এর একটি lecture-এর slide (lecture 6e, অর্থাৎ ষষ্ঠ সপ্তাহের পঞ্চম অংশ)। ফলে বিভিন্ন library-তে ε কোথায় বসে (√-এর ভেতরে না বাইরে) বা momentum যোগ হয় কি না, তাতে ছোটখাটো তফাত আছে।

৫ · Adam — momentum আর RMSProp একসঙ্গে, আর bias correction derive

Diederik Kingma আর Jimmy Ba ২০১৪ সালের ডিসেম্বরে arXiv-এ Adam (adaptive moment estimation) প্রকাশ করেন, আর ২০১৫ সালের ICLR-এ উপস্থাপন করেন। ধারণা: gradient-এর প্রথম moment (গড়, অর্থাৎ momentum-এর দিক) আর দ্বিতীয় moment (g²-এর গড়, অর্থাৎ RMSProp-এর মাপ) দুটিরই exponentially weighted অনুমান রাখা:

m_t = β₁m_t−1 + (১ − β₁)g_t, v_t = β₂v_t−1 + (১ − β₂)g_t²

m̂_t = m_t/(১ − β₁^t), v̂_t = v_t/(১ − β₂^t), x_t = x_t−1 − α m̂_t/(√v̂_t + ε)

প্রবন্ধের প্রস্তাবিত default: α (learning rate) = ০.০০১, β₁ = ০.৯, β₂ = ০.৯৯৯, ε = 10⁻⁸. Exponentially weighted গড়ের “স্মৃতির দৈর্ঘ্য” মোটামুটি ১/(১ − β): m মনে রাখে প্রায় শেষ ১০টি gradient, v প্রায় শেষ ১০০০টি।

Bias correction কোথা থেকে আসে

m₀ = 0 দিয়ে শুরু করলে সূত্রটি খুলে লেখা যায়: m_t = (১ − β₁) ∑_i=1..t β₁^t−i g_i. ধরা যাক gradient-এর প্রত্যাশিত মান প্রতি ধাপে একই, E[g]. তাহলে

E[m_t] = E[g] · (১ − β₁) ∑_i=1..t β₁^t−i = E[g] · (১ − β₁^t)

শেষ ধাপে জ্যামিতিক ধারার যোগফল ব্যবহার হয়েছে। অর্থাৎ শুরুর দিকে m শূন্যের দিকে ঝুঁকে থাকে: t = ১-এ সত্যিকারের মানের মাত্র ১০%. (১ − β₁^t) দিয়ে ভাগ করলে এই ঝোঁক ঠিক বাতিল হয়। v-এর জন্য একই যুক্তিতে (১ − β₂^t), যা t = ১-এ মাত্র ০.০০১: v নিজে সত্যিকারের মানের হাজার ভাগের এক ভাগ। সূক্ষ্ম বিষয়টি দুটির অনুপাত। Correction ছাড়া ধাপ হবে α·m/√v, যার প্রত্যাশিত মাপ α·(১ − β₁^t)/√(১ − β₂^t); t = ১-এ তা α·০.১/√০.০০১ ≈ ৩.১৬α. অর্থাৎ correction না থাকলে প্রথম ধাপগুলো ছোট নয়, বড় হয়, কারণ v-এর ঝোঁক m-এর চেয়ে অনেক বেশি। Kingma আর Ba-র প্রবন্ধেও এই কথা আছে: β₂ এক-এর কাছে হলে correction ছাড়া শুরুর ধাপ অনেক বড় হয়।

তাই Adam-এর α-কে “প্রতি ধাপে parameter মোটামুটি কতটা সরবে” হিসেবে পড়া যায়। নিচের ব্লক একটি parameter-এর noisy gradient (গড় ০.৫) নিয়ে correction-সহ আর ছাড়া ধাপের মাপ α-এর গুণিতকে ছাপে, ১০০০ ধাপ পর্যন্ত।

import numpy as np

rng = np.random.default_rng(42)
g = 0.5 + 0.1 * rng.normal(size=1000)      # একটি parameter-এর gradient: গড় 0.5, সামান্য noise
b1, b2, lr, eps = 0.9, 0.999, 0.001, 1e-8  # Adam-এর প্রচলিত default
m = v = 0.0
raw_steps = []
print("   t      g       m   m_hat         v   v_hat   raw  corr")
for t in range(1, 1001):
    m = b1 * m + (1 - b1) * g[t - 1]
    v = b2 * v + (1 - b2) * g[t - 1] ** 2
    m_hat, v_hat = m / (1 - b1**t), v / (1 - b2**t)
    raw = m / (np.sqrt(v) + eps)               # bias correction ছাড়া (lr-এর গুণিতক)
    corr = m_hat / (np.sqrt(v_hat) + eps)      # bias correction-সহ
    raw_steps.append(raw)
    if t in (1, 2, 3, 10, 100, 1000):
        print(f"{t:4d} {g[t-1]:6.3f} {m:7.4f} {m_hat:7.4f} {v:9.6f} {v_hat:7.4f} {raw:5.2f} {corr:5.2f}")
k = int(np.argmax(raw_steps))
print(f"largest raw step: {raw_steps[k]:.2f} x lr at t = {k + 1};  corrected step stays near 1 x lr")

   t      g       m   m_hat         v   v_hat   raw  corr
   1  0.530  0.0530  0.5305  0.000281  0.2814  3.16  1.00
   2  0.396  0.0873  0.4597  0.000438  0.2191  4.17  0.98
   3  0.575  0.1361  0.5023  0.000768  0.2563  4.91  0.99
  10  0.415  0.3006  0.4615  0.002243  0.2254  6.35  0.97
 100  0.409  0.4556  0.4556  0.023859  0.2506  2.95  0.91
1000  0.582  0.5381  0.5381  0.162590  0.2571  1.33  1.06
largest raw step: 6.61 x lr at t = 14;  corrected step stays near 1 x lr

৬ · Rosenbrock-এ AdaGrad, RMSProp আর Adam — নিজে লিখে তুলনা

সরু বাটি সহজ, কারণ সেটি quadratic. Rosenbrock function f(x, y) = (১ − x)² + ১০০(y − x²)² (Howard Rosenbrock, ১৯৬০) একটি বাঁকা, কলার মতো উপত্যকা: তলায় নামা সহজ, কিন্তু উপত্যকার বাঁক ধরে সর্বনিম্ন বিন্দু (১, ১)-এ পৌঁছানো কঠিন, আর বিভিন্ন জায়গায় বাঁকের মাপ আলাদা, তাই কোনো একটি স্থির κ নেই। Optimizer পরীক্ষার এটি একটি পুরোনো মানদণ্ড; পর্ব ১৫-এ Newton আর L-BFGS-ও এখানেই পরীক্ষা হবে।

নিচের ব্লকে প্রতিটি optimizer একটি function, যা gradient, একটি ছোট state dictionary আর ধাপ-সংখ্যা t নিয়ে শুধু “কত সরতে হবে” ফেরত দেয়। চারটি learning rate-এ (০.০০১ থেকে ১) প্রতিটি ২০০০ ধাপ চলে, শুরু (−১.৫, ২). ছকে থাকে (১, ১) থেকে শেষ দূরত্ব, আর ছবিতে প্রতিটি optimizer-এর নিজের সেরা lr-এর পথ। একটি optimizer-কে শুধু একটি lr-এ চালিয়ে বিচার করা অন্যায্য, তাই পরিসর।

import numpy as np
import matplotlib.pyplot as plt

f = lambda p: (1 - p[0])**2 + 100 * (p[1] - p[0]**2)**2       # Rosenbrock, সর্বনিম্ন (1, 1)
grad = lambda p: np.array([-2 * (1 - p[0]) - 400 * p[0] * (p[1] - p[0]**2), 200 * (p[1] - p[0]**2)])
def gd(g, s, t, lr): return lr * g
def adagrad(g, s, t, lr):                                      # সব অতীত g²-এর যোগফল
    s["G"] = s.get("G", 0) + g**2;              return lr * g / (np.sqrt(s["G"]) + 1e-8)
def rmsprop(g, s, t, lr):                                      # g²-এর চলমান গড়
    s["E"] = 0.9 * s.get("E", 0) + 0.1 * g**2;  return lr * g / (np.sqrt(s["E"]) + 1e-8)
def adam(g, s, t, lr, b1=0.9, b2=0.999):                       # momentum + RMSProp + bias correction
    s["m"] = b1 * s.get("m", 0) + (1 - b1) * g
    s["v"] = b2 * s.get("v", 0) + (1 - b2) * g**2
    return lr * s["m"] / (1 - b1**t) / (np.sqrt(s["v"] / (1 - b2**t)) + 1e-8)

lrs, best = [0.001, 0.01, 0.1, 1.0], {}
print("distance to (1, 1) after 2000 steps, start (-1.5, 2)\nlr      " + "".join(f"{lr:>9}" for lr in lrs))
for name, opt in [("GD", gd), ("AdaGrad", adagrad), ("RMSProp", rmsprop), ("Adam", adam)]:
    row = ""
    for lr in lrs:
        p, s, P = np.array([-1.5, 2.0]), {}, []
        for t in range(1, 2001):
            p = p - opt(grad(p), s, t, lr); P.append(p)
            if f(p) > 1e6: break                                  # ছিটকে গেছে, থামুন
        d = np.linalg.norm(p - 1) if f(p) < 1e6 else np.inf
        row += f"{d:9.1e}" if d < np.inf else " diverged"
        if d < best.get(name, (np.inf,))[0]: best[name] = (d, np.array(P), lr)
    print(f"{name:8s}" + row)

u, w = np.meshgrid(np.linspace(-2, 2, 200), np.linspace(-1, 3, 200))
plt.contour(u, w, np.log10(f((u, w)) + 1e-3), levels=20, colors="0.8")
for name, (d, P, lr) in best.items():
    plt.plot(P[:, 0], P[:, 1], lw=1.2, label=f"{name} (lr={lr})")
plt.plot(1, 1, "k*", ms=12); plt.legend(); plt.xlabel("x"); plt.ylabel("y")
plt.title("Rosenbrock: each optimizer at its best lr in the grid"); plt.show()

distance to (1, 1) after 2000 steps, start (-1.5, 2)
lr          0.001     0.01      0.1      1.0
GD        8.2e-01 diverged diverged diverged
AdaGrad   2.7e+00  2.5e+00  1.0e+00  5.7e-01
RMSProp   1.6e+00  1.2e-01  1.2e+00  4.9e-01
Adam      2.4e+00  1.9e-01  1.7e-11  1.3e+00

৭ · L2 না weight decay — AdamW আর একটি ন্যায্য তুলনা

AI-এর ভিত্তি কোর্সের পর্ব ৬১-এ L2 regularization এসেছে: loss-এ (λ/২)‖w‖² যোগ, ফলে gradient-এ λw যোগ হয়। সাধারণ gradient descent-এ তখন হালনাগাদ w ← w − η(∇L + λw) = (১ − ηλ)w − η∇L: প্রতি ধাপে ওজন একটি স্থির অনুপাতে শূন্যের দিকে ছোট হয়। এই “ছোট করা”-কেই বলে weight decay, আর SGD-তে দুটি কথা একই জিনিস (পর্ব ৬৮-এ সেখানে এক লাইনে এর ইঙ্গিত ছিল)।

Adam-এ তারা আলাদা হয়ে যায়। Coupled L2 (Adam + L2): g = ∇L + λw-কে Adam-এ দিলে λw-ও m আর v-এর ভেতরে যায়, আর শেষে √v̂ দিয়ে ভাগ হয়। যে ওজনের gradient ঐতিহাসিকভাবে বড় (v̂ বড়), তার শাস্তি কার্যত ছোট; যার gradient ছোট, তার বড়। Decoupled (AdamW): Ilya Loshchilov আর Frank Hutter ২০১৭ সালে arXiv-এ (পরে ICLR ২০১৯-এ) প্রস্তাব দেন, শাস্তিকে Adam-এর হিসাব থেকে বের করে সরাসরি ওজনে বসাতে:

w ← w − η ( m̂/(√v̂ + ε) + λw )

এখানে m̂ আর v̂ শুধু ∇L থেকে আসে। এই AdamW প্রতিটি ওজনকে একই অনুপাতে (ηλ) ছোট করে। তাঁদের image classification পরীক্ষায় AdamW-এর generalization ভালো ছিল; আজ বহু transformer training-এ এটি প্রচলিত পছন্দ। তবে এটি deep network-এর অভিজ্ঞতা; প্রতিটি সমস্যায় স্বয়ংক্রিয় নিশ্চয়তা নয়।

একটি গুরুত্বপূর্ণ ফল: AdamW আর loss + (λ/২)‖w‖²-এর সর্বনিম্ন বিন্দুতে যায় না, কারণ তার স্থির বিন্দুর শর্ত আলাদা। নিচের ব্লক এটি সংখ্যায় দেখায়। কাল্পনিক বাজারঘরের ৩০০টি অর্ডারের ফেরত (logistic regression, চারটি feature, ইচ্ছা করে scale খুব আলাদা: ১, ১, ১০, ০.১), objective J(w) = গড় log loss + (০.০৫/২)‖w‖². Newton-এর ১০ ধাপে সঠিক সর্বনিম্ন বের করে (পর্ব ১১-এর মতো), তারপর Adam + L2 আর AdamW ৩০০০ ধাপ। দ্বিতীয় অংশে GD, momentum (β = ০.৯) আর Adam-কে চারটি learning rate আর পাঁচটি seed-এ ৩০০ ধাপ চালিয়ে J − J*-এর গড়।

import numpy as np
sig = lambda z: 0.5 * (1 + np.tanh(z / 2))              # sigmoid, overflow ছাড়া
def data(seed, n=300):                                  # কাল্পনিক বাজারঘর: অর্ডার ফেরত এল কি না
    rng = np.random.default_rng(seed)
    X = rng.normal(size=(n, 4)) * [1, 1, 10, 0.1]       # feature-গুলোর scale খুব আলাদা
    y = (rng.random(n) < sig(X @ [1.0, -2.0, 0.2, 5.0])).astype(float)
    return X, y
lam, lrs = 0.05, [0.001, 0.01, 0.1, 1.0]
grad = lambda w, X, y: X.T @ (sig(X @ w) - y) / len(y)
J = lambda w, X, y: np.mean(np.logaddexp(0, X @ w) - y * (X @ w)) + lam / 2 * w @ w
def exact(X, y, l=lam, w=np.zeros(4)):                  # Newton-এর ১০ ধাপ: আসল সর্বনিম্ন
    for _ in range(10):
        p = sig(X @ w)
        H = X.T @ (X * (p * (1 - p))[:, None]) / len(y) + l * np.eye(4)
        w = w - np.linalg.solve(H, grad(w, X, y) + l * w)
    return w
def run(opt, X, y, lr, steps=300, decoupled=False):
    w, m, v = np.zeros(4), np.zeros(4), np.zeros(4)
    for t in range(1, steps + 1):
        g = grad(w, X, y) + (0 if decoupled else lam * w)        # coupled: L2 যায় gradient-এ
        if opt == "GD": w = w - lr * g
        elif opt == "Momentum": m = 0.9 * m + g; w = w - lr * m
        else:
            m, v = 0.9 * m + 0.1 * g, 0.999 * v + 0.001 * g**2
            step = m / (1 - 0.9**t) / (np.sqrt(v / (1 - 0.999**t)) + 1e-8)
            w = w - lr * (step + (lam * w if decoupled else 0))   # decoupled: আলাদা করে ছোট করা
    return w
X, y = data(0)
for name, w in [("exact, lam=0.05", exact(X, y)), ("exact, lam=0", exact(X, y, 0.0)),
                ("Adam + L2", run("Adam", X, y, 0.01, 3000)), ("AdamW", run("Adam", X, y, 0.01, 3000, True))]:
    print(f"{name:15s} w = {np.round(w, 3)}  J = {J(w, X, y):.5f}")
print("\nJ - J* after 300 steps, mean of 5 seeds\nlr       " + "".join(f"{lr:>9}" for lr in lrs))
for opt in ["GD", "Momentum", "Adam"]:
    gap = lambda lr: np.mean([J(run(opt, *d, lr), *d) - J(exact(*d), *d) for d in map(data, range(5))])
    print(f"{opt:9s}" + "".join(f"{gap(lr):9.1e}" for lr in lrs))

exact, lam=0.05 w = [ 0.688 -1.086  0.106  0.083]  J = 0.46774
exact, lam=0    w = [ 1.181 -1.885  0.141  3.246]  J = 0.78156
Adam + L2       w = [ 0.688 -1.086  0.106  0.083]  J = 0.46774
AdamW           w = [ 1.169 -1.861  0.14   3.181]  J = 0.76830

J - J* after 300 steps, mean of 5 seeds
lr           0.001     0.01      0.1      1.0
GD         1.2e-01  4.0e-02  2.2e-05  3.8e+00
Momentum   4.0e-02  1.6e-05  5.5e-15  6.0e+00
Adam       6.9e-02  6.1e-06  2.3e-15  4.4e-14

৮ · কোনটি কখন — সৎ তুলনা আর জানা ব্যর্থতা

উপরের run-গুলো দুটি ছোট function আর একটি convex model, পূর্ণ gradient-সহ; এখান থেকে সাধারণ সিদ্ধান্ত টানা ভুল হবে। বাস্তব training-এ gradient noisy (mini-batch থেকে আসে, অর্থাৎ SGD-র পরিবেশ), loss non-convex, আর মাপকাঠি training loss নয়, নতুন ডেটায় ফল। কয়েকটি সতর্ক পর্যবেক্ষণ:

Adam সবসময় অভিসারী নয়। Sashank Reddi, Satyen Kale আর Sanjiv Kumar-এর ২০১৮ সালের ICLR প্রবন্ধ “On the Convergence of Adam and Beyond” দেখায় যে Adam-এর মূল প্রবন্ধের অভিসারের প্রমাণে একটি ফাঁক ছিল, আর একটি সরল এক-মাত্রিক online convex উদাহরণ দেয়, যেখানে বেশিরভাগ ধাপে ছোট gradient এক দিকে আর মাঝে মাঝে একটি বড় gradient অন্য দিকে আসে। v-এর স্মৃতি সীমিত বলে বড় gradient-এর তথ্য দ্রুত মুছে যায়, আর Adam ভুল প্রান্তে গিয়ে স্থির হয়। তাঁদের প্রস্তাবিত AMSGrad v̂-কে কখনো কমতে দেয় না (অতীতের সর্বোচ্চ রাখে)। বাস্তবে এই ব্যর্থতা কত ঘন ঘন ঘটে, তা নিয়ে মতভেদ আছে; PyTorch-এর Adam-এ amsgrad=True একটি বিকল্প হিসেবে আছে।

Generalization নিয়ে বিতর্ক। Ashia Wilson ও সহকর্মীদের ২০১৭ সালের NeurIPS প্রবন্ধ “The Marginal Value of Adaptive Gradient Methods in Machine Learning” কয়েকটি পরীক্ষায় দেখায় যে ভালোভাবে tune করা SGD + momentum-এর test ফল adaptive পদ্ধতির চেয়ে ভালো বা সমান ছিল। আবার পরে অনেক বড় benchmark গবেষণা (যেমন Schmidt, Schneider আর Hennig-এর ২০২১ সালের “Descending through a Crowded Valley”) মোটামুটি এই সিদ্ধান্তে আসে যে কোনো একটি optimizer সব কাজে সেরা নয়, আর tuning-এর বাজেট প্রায়ই optimizer বাছাইয়ের চেয়ে বেশি গুরুত্বপূর্ণ। অভিজ্ঞতার কথা হিসেবে বলা যায়: transformer-জাতীয় model-এ AdamW প্রচলিত, কিছু image model-এ SGD + momentum আজও জনপ্রিয়; কিন্তু এগুলো প্রচলন, প্রমাণিত নিয়ম নয়।

scikit-learn-এ এই optimizer-গুলো সরাসরি দেখা যায় MLPClassifier আর MLPRegressor-এ: default solver="adam", সঙ্গে beta_1=0.9, beta_2=0.999, epsilon=1e-8, আর learning_rate_init=0.001, অর্থাৎ Kingma ও Ba-র প্রস্তাবিত মানগুলোই। কিন্তু alpha সেখানে loss-এ যোগ হওয়া L2 শাস্তি, AdamW-এর decoupled decay নয়। আর SGDClassifier-এ কোনো momentum নেই; সেখানে learning-rate schedule আছে (পর্ব ১৩)। Deep learning framework-এর optimizer আসবে Deep Learning কোর্সে।

৯ · AI-এর যুগে — AI-এর লেখা optimizer কোড আর দাবি যাচাই

“Adam দিয়ে training loop লিখে দাও” বললে AI assistant কয়েক সেকেন্ডে পরিষ্কার কোড দেয়। ঝুঁকি থাকে খুঁটিনাটিতে, যা চালিয়ে না দেখলে ধরা পড়ে না: t শূন্য থেকে শুরু করলে bias correction-এ শূন্য দিয়ে ভাগ, weight decay ভুল জায়গায়, বা একটি run দেখে “Adam সেরা” সিদ্ধান্ত।

AI-এর দেওয়া optimizer কোড হাতে পেলে — যাচাইয়ের চেকলিস্ট

• প্রথম দুই-তিনটি ধাপ হাতে (বা একটি ছোট stepper-এর মতো) হিসাব করে কোডের সঙ্গে মিলিয়েছি
• t = ১ থেকে শুরু হচ্ছে, আর correction-সহ প্রথম ধাপ প্রায় lr কি না দেখেছি
• Weight decay coupled না decoupled, তা নিজের library-র version-এর docs থেকে নিশ্চিত হয়েছি
• একটি সঠিক উত্তর জানা সমস্যায় (quadratic বা Newton দিয়ে সমাধান করা logistic) optimizer যাচাই করেছি
• তুলনা করেছি একাধিক lr আর seed-এ, প্রতিটির নিজের সেরা lr দিয়ে, validation ডেটায়
• Test set তালা দেওয়া রেখেছি; lr বা optimizer test ফল দেখে বাছিনি
• কোনো আসল গ্রাহকের ব্যক্তিগত ডেটা AI-তে পেস্ট করিনি