Blog

#!/usr/bin/env python3
"""경사하강법으로 최소제곱과 같은 해에 도달하는지 확인한다."""
import numpy as np
from data import load

ext, size, entropy, ratio, tr, te = load()
X = np.column_stack([entropy[tr], np.ones(len(tr))])
y = ratio[tr]
n = len(y)

def mse(t):
    return float(np.mean((X @ t - y) ** 2))

def grad(t):                      # dL/dtheta = 2/n * X^T (X theta - y)
    return 2.0 / n * (X.T @ (X @ t - y))

theta = np.zeros(2)
lr = 0.02
print(f"{'step':>6}{'MSE':>12}{'w':>10}{'b':>10}{'|grad|':>12}")
for step in range(4001):
    if step in (0, 10, 100, 500, 1000, 2000, 4000):
        g = np.linalg.norm(grad(theta))
        print(f"{step:>6}{mse(theta):>12.6f}{theta[0]:>10.4f}{theta[1]:>10.4f}{g:>12.2e}")
    theta -= lr * grad(theta)

closed = np.linalg.lstsq(X, y, rcond=None)[0]
print(f"\n닫힌 해      : w={closed[0]:+.6f} b={closed[1]:+.6f}  MSE={mse(closed):.8f}")
print(f"경사하강 4000: w={theta[0]:+.6f} b={theta[1]:+.6f}  MSE={mse(theta):.8f}")

#!/usr/bin/env python3
"""학습률만 바꿔 가며 1000스텝 뒤 상태를 비교하고, 이론적 상한과 맞춰본다."""
import numpy as np
from data import load

ext, size, entropy, ratio, tr, te = load()
X = np.column_stack([entropy[tr], np.ones(len(tr))])
y = ratio[tr]
n = len(y)
best = np.linalg.lstsq(X, y, rcond=None)[0]
best_mse = float(np.mean((X @ best - y) ** 2))

# 손실 L(t)=|Xt-y|^2/n 의 헤세 행렬은 (2/n) X^T X 로 상수다.
# 경사하강이 발산하지 않을 조건은 lr < 2 / (헤세의 최대 고윳값).
H = 2.0 / n * (X.T @ X)
eig = np.linalg.eigvalsh(H)
lr_max = 2.0 / eig.max()
print(f"헤세 고윳값     : {eig.min():.4f} ~ {eig.max():.4f}")
print(f"이론적 학습률 상한: {lr_max:.6f}")
print(f"최소 MSE(닫힌 해) : {best_mse:.8f}\n")

def train(lr, steps=1000):
    t = np.zeros(2)
    with np.errstate(over="ignore", invalid="ignore"):
        for s in range(steps):
            t = t - lr * (2.0 / n * (X.T @ (X @ t - y)))
            if not np.all(np.isfinite(t)) or np.abs(t).max() > 1e12:
                return None, s
    return t, steps

print(f"{'학습률':>9}{'1000스텝 MSE':>15}{'w':>9}{'b':>9}  상태")
print("-" * 55)
for lr in (0.0001, 0.001, 0.01, 0.02, 0.028, 0.030, 0.04):
    t, s = train(lr)
    if t is None:
        print(f"{lr:>9}{'발산':>15}{'-':>9}{'-':>9}  {s}스텝 만에 폭발")
        continue
    m = float(np.mean((X @ t - y) ** 2))
    gap = m - best_mse
    state = "수렴" if gap < 1e-7 else ("거의 수렴" if gap < 1e-4 else "느림")
    print(f"{lr:>9}{m:>15.8f}{t[0]:>9.4f}{t[1]:>9.4f}  {state}")

#!/usr/bin/env python3
"""특성 스케일이 경사하강에 미치는 영향을 조건수와 함께 본다."""
import numpy as np
from data import load

ext, size, entropy, ratio, tr, te = load()
y = ratio[tr]
n = len(y)

raw  = np.column_stack([entropy[tr], size[tr], np.ones(n)])
mu, sd = raw[:, :2].mean(0), raw[:, :2].std(0)
std  = np.column_stack([(raw[:, :2] - mu) / sd, np.ones(n)])

def analyze(name, X):
    H = 2.0 / n * (X.T @ X)
    eig = np.linalg.eigvalsh(H)
    lr = 0.9 * 2.0 / eig.max()                 # 상한의 90%로 고정
    best = float(np.mean((X @ np.linalg.lstsq(X, y, rcond=None)[0] - y) ** 2))
    t = np.zeros(X.shape[1])
    hit = None
    with np.errstate(over="ignore", invalid="ignore"):
        for s in range(1, 200001):
            t = t - lr * (2.0 / n * (X.T @ (X @ t - y)))
            if hit is None and np.isfinite(t).all():
                if float(np.mean((X @ t - y) ** 2)) < best * 1.001:
                    hit = s
                    break
    print(f"[{name}]")
    print(f"    조건수 cond(X^T X) = {np.linalg.cond(X.T @ X):.3e}   쓸 수 있는 학습률 = {lr:.3e}")
    print(f"    최소 MSE의 100.1%에 도달한 스텝: {hit if hit else '200000 안에 도달 못함'}")

analyze("원본 스케일 (entropy, size)", raw)
analyze("표준화 (평균 0 / 표준편차 1)", std)

#!/usr/bin/env python3
"""같은 에폭 예산에서 배치 크기가 수렴 속도와 흔들림에 미치는 영향을 본다."""
import numpy as np
from data import load

ext, size, entropy, ratio, tr, te = load()
y = ratio[tr]
n = len(y)
raw = np.column_stack([entropy[tr], size[tr]])
X = np.column_stack([(raw - raw.mean(0)) / raw.std(0), np.ones(n)])
best = float(np.mean((X @ np.linalg.lstsq(X, y, rcond=None)[0] - y) ** 2))

def mse(t):
    return float(np.mean((X @ t - y) ** 2))

def run(batch, lr, epochs=30, seed=0):
    rng = np.random.default_rng(seed)
    t = np.zeros(X.shape[1])
    hist = []
    for _ in range(epochs):
        order = rng.permutation(n)
        for i in range(0, n, batch):
            j = order[i:i + batch]
            t = t - lr * (2.0 / len(j) * (X[j].T @ (X[j] @ t - y[j])))
            hist.append(mse(t))
    return t, np.array(hist)

print(f"최소 MSE = {best:.8f}, 에폭 예산 30회 (학습률은 모두 0.05로 고정)\n")
print(f"{'배치':>8}{'에폭당 갱신':>11}{'최종 MSE':>13}{'최소값 대비':>12}{'마지막 10회 표준편차':>20}")
print("-" * 66)
for batch in (n, 64, 16, 4, 1):
    t, hist = run(batch, 0.05)
    label = f"{batch}(전체)" if batch == n else str(batch)
    print(f"{label:>8}{int(np.ceil(n / batch)):>11}{hist[-1]:>13.8f}"
          f"{hist[-1] / best:>12.4f}{hist[-10:].std():>20.2e}")