PRML演習問題1.1,1.2

poly_curve_fitting.py

import numpy as np
from polys import poly_fit as y  # 式1.1

'''
演習問題 1.1

X_n = [x_n^0, x_n^1, x_n^2, ...]
W   = [w_0, w_1, w_2, ...]

E(W)     = 1/2 * Σ[n]{y(X_n, W) - t_n}^2 = 1/2 * Σ[n]{W * X_n - t_n}^2
dE(W)/dW = Σ[n]{W * X_n - t_n}X_n = 0
上記を変形していくと問題で与えられたΣ[j]{A_ij * w_j} = T_i の線形方程式となり、
A * W = Tの形の連立方程式でWが解ける
'''
def T(xl, tl, n):
    return np.array([(xl**i * tl).sum() for i in range(n)]) # 式1.123

def A_ij(x_n, i, j):
    return (x_n**(i + j)).sum()

def train(xl, tl, n):
    A = np.array([A_ij(xl, i, j) for i in range(n) for j in range(n)]).reshape(n, n) # 式1.123
    return np.dot(np.linalg.inv(A), T(xl, tl, n)) # 式1.122

'''
演習問題 1.2
1.1の正則化項付き

E(W)       = 1/2 * Σ[n]{W * X_n - t_n}^2 + λ/2 * Σ{w_n^2}
dE(W)/dw_i = Σ[n]{W * X_n - t_n} * X_n + λ * w_i = 0
上記を変形していくと、Σ[j]{A_ij * w_j} + λ * w_i = T_i の形になる
'''
def train_reg(xl, tl, n, lp):
    A = np.array([A_ij(xl, i, j) if j != i
                  else A_ij(xl, i, j) + lp
                  for i in range(n) for j in range(n)]).reshape(n, n)

    return np.dot(np.linalg.inv(A), T(xl, tl, n))