#背景
日本ディープラーニング協会の資格試験（E資格）受験資格を得るためにラビット・チャレンジの講座プログラムを受けることにした。
当該講座を受講する時、科目ごとにレポートを作成しWebに投稿する必要があるためQiitaで作成することとした。
個人的には各科目の全要点などを記述よりも自分が難しいと思うところまたは理解不足なポイントをまとめる形にしたいと考えている。
ちなみにラビット・チャレンジについては以下から参考できる
#機械学習
機械学習については以下６つの科目でレポートする
１．線形回帰モデル
２．非線形回帰モデル
３．ロジスティクス回帰モデル
４．主成分分析
５．アルゴリズム
６．サポートベクターマシーン
##1. 線形回帰モデル
線形回帰モデルの特徴は以下３つとなる
・回帰問題を解くための機械学習モデルの一つ
・教師あり学習
・入力とm次元パラメータの線形結合を出力するモデル
回帰問題とはある入力から出力を予測する問題でああり、その中に直線で予測できる問題が線形回帰モデルで扱う。つまり入力と出力の関係を1次式で表すことができる。
入力をx，パラメータをw，予測値をŷとし、
x=(x_1,x_2,…,x_m)^T∈ℝ^m \\
w=(w_1,w_2,…,w_m)^T∈R^m
入力、パラメータ、予測値の関係が以下の式で表せる
ŷ =w^Tx+w_0=\sum_{j=1}^{m}w_jx_j+w_0
教師あり学習の場合、訓練データ入力と出力が既知であるため、パラメータ w は最小二乗法により推定できる。
ŵ =arg\min_{w∈ℝ^{m+1}}MSE_{train}
MSE（最小二乗誤差）を最小にするようにMSEをｗに関して微分したものが０となるｗの点を求める。
\frac{∂}{∂w}MSE_{train}=0
計算の結果：
ŵ=(x^{(train)T}x^{train})^{-1}x^{(train)T}y^{train}
したがって、新たな入力値に対して予測値が以下となる
ŷ=xŵ=x(x^{(train)T}x^{train})^{-1}x^{(train)T}y^{train}
###実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
##データの作成
xs = np.linspace(0,1,100)
ys_true = 2*xs + 5
#add noise
ys = ys_true + np.random.normal(scale=.2,size= ys_true.shape)
##学習
def train(xs, ys):
cov = np.cov(xs, ys, ddof=0)
a = cov[0, 1] / cov[0, 0]
b = np.mean(ys) - a*np.mean(xs)
return cov, a, b
cov, a, b = train(xs, ys)
print("cov: {}".format(cov))
print("coef: {}".format(a))
print("intercept: {}".format(b))
##予測
#x_newに対してys_pred = a*x_new + b で予測する
xs_new = np.linspace(0, 1, 100)
ys_pred = a * xs_new + b
plt.scatter(xs, ys, facecolor="none", edgecolor="b", s=50, label="training data")
plt.plot(xs_new, ys_true, label="$2 x + 5$")
plt.plot(xs_new, ys_pred, label="prediction (a={:.2}, b={:.2})".format(a, b))
plt.legend()
plt.show()
#####実行結果：
##2. 非線形回帰モデル
複雑な非線形構造のデータでは線形回帰モデルでは捉えられないため非線形回帰モデルで表現する必要なある。
非線形回帰モデルの基底関数は多項式関数やガウス型基底関数などを用いるこどが多い。
多項式:
ϕ_j = x^j
ガウス型基底:
ϕ_j(x)=exp\left\{
{\frac{(x-u_j)^2}{2h_j}}
\right\}
非線形回帰の場合、入力値に対して非線形であるが、パラメータについては線形のモデルとすることで、wの導出には線形回帰の場合と同様に扱うことができる。
入力値：
x_i=(x_{i1},x_{i2},…,x_{im})^T∈ℝ^m
非線形関数ベクトル：
\phi(x_i)=(\phi_1(x_i),\phi_2(x_i),…,\phi_k(x_i))^T∈ℝ^k
非線形関数の計画行列：
\Phi^{train}=(\phi(x_i),\phi(x_i),…,\phi(x_n))^T∈ℝ^{n \times k}
最尤法による予測値：
ŷ=\Phi(\Phi^{(train)T}\Phi^{train})^{-1}\Phi^{(train)T}y^{train}
非線形回帰モデルにおいて、未学習と過学習を防ぐ、汎化性能を高めるために主に以下３つの手法が挙げれる。
１）訓練データを増やす
２）不要な基底関数を除去する
３）正則法を用いる
###実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
##訓練データの作成
xs = np.linspace(0, 1, 10)
ys_true = np.sin(2 * np.pi * xs)
ys = ys_true + np.random.normal(scale=.25, size=ys_true.shape)
##学習
def polynomial_features(xs, degree=3):
"""多項式特徴ベクトルに変換
X = [[1, x1, x1^2, x1^3],
[1, x2, x2^2, x2^3],
...
[1, xn, xn^2, xn^3]]"""
X = np.ones((len(xs), degree+1))
X_t = X.T #(100, 4)
for i in range(1, degree+1):
X_t[i] = X_t[i-1] * xs
return X_t.T
Phi = polynomial_features(xs)
#print(xs)
#print(Phi)
Phi_inv = np.dot(np.linalg.inv(np.dot(Phi.T, Phi)), Phi.T)
w = np.dot(Phi_inv, ys)
print(w)
#予測
Phi_test = polynomial_features(xs)
ys_pred = np.dot(Phi_test, w)
print(Phi_test)
print(ys_pred)
plt.scatter(xs, ys, facecolor="none", edgecolor="b", s=50, label="training data")
plt.plot(xs, ys_true, label="$\sin(2\pi x)$")
plt.plot(xs, ys_pred, label="prediction")
plt.legend()
plt.show()
#####実行結果：
##3. ロジスティクス回帰モデル
ロジスティクス回帰モデルは主に二値分類に使われるアルゴリズムであり、出力値が２種類しかないような教師あり学習に適用される。
ロジスティクス回帰モデルは次のような式を前提とするモデルである。
P(Y=1|X=x)=\sigma(w_o+\sum_{j=1}^dx_jw_j)
ここでσはシグモイド関数となる。
\sigma=\frac{1}{1+e^{-x}}
シグモイド関数がよく用いられる理由としてはシグモイド関数の出力が０〜１となるため確率として扱えることと、シグモイド関数の微分がそれ自体で表現で最尤法の計算が容易になるからである。
次にY=1とY=0になる確実を求めて表現するためにベルヌーイ分布が用いられる。
P(y)=p^y(1-p)^{1-y}
最尤推定法にしたがって尤度関数が以下の式となる。
\begin{eqnarray}
P(y_1.y_2,…,y_n|w_0,w_1,…w_m)&=&\prod_{i=1}^{n}p_i^{y_i}(1-p_i)^{1-y_i}\\
&=&\prod_{i=1}^{n}\sigma(w^Tx_i)^{y_i}(1-\sigma(w^Tx_i))^{1-y_i}\\
&=&L(w)
\end{eqnarray}
L(w)を最大にするパラメータｗを求めると、
w^{(k+1)}=w^{(k)}+\eta\sum_{i=1}^n(y_i-p_i)x_i
###実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
##訓練データの作成
n_sample = 100
harf_n_sample = 50
var = .2
def gen_data(n_sample, harf_n_sample):
x0 = np.random.normal(size=n_sample).reshape(-1, 2) - 1.
x1 = np.random.normal(size=n_sample).reshape(-1, 2) + 1.
x_train = np.concatenate([x0, x1])
y_train = np.concatenate([np.zeros(harf_n_sample), np.ones(harf_n_sample)]).astype(np.int)
return x_train, y_train
def plt_data(x_train, y_train):
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train, facecolor="none", edgecolor="b", s=50, label="training data")
plt.legend()
#データ作成
x_train, y_train = gen_data(n_sample, harf_n_sample)
#データ表示
#plt_data(x_train, y_train)
##学習
def add_one(x):
return np.concatenate([np.ones(len(x))[:, None], x], axis=1)
def sigmoid(x):
return 1 / (1 + np.exp(-x))
def sgd(X_train, max_iter, eta):
w = np.zeros(X_train.shape[1])
for _ in range(max_iter):
w_prev = np.copy(w)
sigma = sigmoid(np.dot(X_train, w))
grad = np.dot(X_train.T, (sigma - y_train))
w -= eta * grad
if np.allclose(w, w_prev):
return w
return w
X_train = add_one(x_train)
max_iter=100
eta = 0.01
w = sgd(X_train, max_iter, eta)
#print(X_train)
print(w)
#予測
xx0, xx1 = np.meshgrid(np.linspace(-5, 5, 100), np.linspace(-5, 5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T
X_test = add_one(xx)
proba = sigmoid(np.dot(X_test, w))
y_pred = (proba > 0.5).astype(np.int)
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
plt.contourf(xx0, xx1, proba.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
#####実行結果：
##4. 主成分分析
主成分分析は教師なし学習の一種であり次元圧縮の手法となる。特に２次元に圧縮することで多次元のデータを可視化するに役に立つ。
ｍ次元の入力データの次元を圧縮するには以下の線形変換を行う。
s_j=(s_{1j},s_{2j},…,s_{nj})^T=Xa_j
線形変換後の分散は以下の式で表せる。
V_{ar}(s_j)=\frac{1}{n}s_j^Ts_j=a_j^TV_{ar}(X)a_j
ここで解が無限になることを防ぐために、ノルムが1となる制約を入れる。結果は制約付きの最適化問題としてラグランジュ関数を最大にする係数ベクトルを探査することになる。(微分した結果が０になる)
E(a_j)=a_j^TVar(X)a_j-\lambda(a_j^Ta_j-1) \\
\frac{∂E(a_j)}{∂a_j}=0
計算の結果：
Var(X)a_j=\lambda{a_j}
結果として、元のデータの分散共分散行列の固有値、固有ベクトルが解と一致する。
また、累積寄与率とは各主成分の寄与率を大きい順に足し上げていったのもので、そこまでの主成分でデータの持っていた情報量がどのぐらい説明されているかを示す。
###実装演習
####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
#訓練データ生成
n_sample = 100
def gen_data(n_sample):
mean = [0, 0]
cov = [[2, 0.7], [0.7, 1]]
return np.random.multivariate_normal(mean, cov, n_sample)
def plt_data(X):
plt.scatter(X[:, 0], X[:, 1])
plt.xlim(-5, 5)
plt.ylim(-5, 5)
X = gen_data(n_sample)
#print(X)
#plt_data(X)
#学習
n_components=2
def get_moments(X):
mean = X.mean(axis=0)
stan_cov = np.dot((X - mean).T, X - mean) / (len(X) - 1)
return mean, stan_cov
def get_components(eigenvectors, n_components):
# W = eigenvectors[:, -n_components:]
# return W.T[::-1]
W = eigenvectors[:, ::-1][:, :n_components]
return W.T
def plt_result(X, first, second):
plt.scatter(X[:, 0], X[:, 1])
plt.xlim(-5, 5)
plt.ylim(-5, 5)
# 第１主成分
plt.quiver(0, 0, first[0], first[1], width=0.01, scale=6, color='red')
# 第２主成分
plt.quiver(0, 0, second[0], second[1], width=0.01, scale=6, color='green')
#分散共分散行列を標準化
meean, stan_cov = get_moments(X)
#固有値と固有ベクトルを計算
eigenvalues, eigenvectors = np.linalg.eigh(stan_cov)
components = get_components(eigenvectors, n_components)
#plt_result(X, eigenvectors[0, :], eigenvectors[1, :])
#変換（射影）
def transform_by_pca(X, pca):
mean = X.mean(axis=0)
return np.dot(X-mean, components)
Z = transform_by_pca(X, components.T)
#plt.scatter(Z[:, 0], Z[:, 1])
#plt.xlim(-5, 5)
#plt.ylim(-5, 5)
#逆変換
mean = X.mean(axis=0)
X_ = np.dot(Z, components.T) + mean
plt.scatter(X_[:, 0], X_[:, 1])
plt.xlim(-5, 5)
plt.ylim(-5, 5)
#####実行結果：
学習結果
変換結果
逆変換結果
##5. アルゴリズム
###k近傍法
k近傍法は分類問題のための機械学習手法であり、最近傍のデータを個取ってきて、それらがもっとも多く所属するクラスに識別する方法である。kを大きくすると決定境界は滑らかになる。
###k平均法
k平均法もクラスの分類手法の１つとなるが、k近傍法との明確な違いとしては、k平均法はクラスタリング手法という手法であり、教師なし学習の分類問題とみなす。
###k近傍法実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
#訓練データ生成
def gen_data():
x0 = np.random.normal(size=50).reshape(-1, 2) - 1
x1 = np.random.normal(size=50).reshape(-1, 2) + 1.
x_train = np.concatenate([x0, x1])
y_train = np.concatenate([np.zeros(25), np.ones(25)]).astype(np.int)
return x_train, y_train
X_train, ys_train = gen_data()
plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)
#学習
##なし
#予測
def distance(x1, x2):
return np.sum((x1 - x2)**2, axis=1)
def knc_predict(n_neighbors, x_train, y_train, X_test):
y_pred = np.empty(len(X_test), dtype=y_train.dtype)
for i, x in enumerate(X_test):
distances = distance(x, X_train)
nearest_index = distances.argsort()[:n_neighbors]
mode, _ = stats.mode(y_train[nearest_index])
y_pred[i] = mode
return y_pred
def plt_resut(x_train, y_train, y_pred):
xx0, xx1 = np.meshgrid(np.linspace(-5, 5, 100), np.linspace(-5, 5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
plt.contourf(xx0, xx1, y_pred.reshape(100, 100).astype(dtype=np.float), alpha=0.2, levels=np.linspace(0, 1, 3))
n_neighbors = 3
xx0, xx1 = np.meshgrid(np.linspace(-5, 5, 100), np.linspace(-5, 5, 100))
X_test = np.array([xx0, xx1]).reshape(2, -1).T
y_pred = knc_predict(n_neighbors, X_train, ys_train, X_test)
plt_resut(X_train, ys_train, y_pred)
#####実行結果：
###k平均法実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
#訓練データ生成
def gen_data():
x1 = np.random.normal(size=(100, 2)) + np.array([-5, -5])
x2 = np.random.normal(size=(100, 2)) + np.array([5, -5])
x3 = np.random.normal(size=(100, 2)) + np.array([0, 5])
return np.vstack((x1, x2, x3))
##データ作成
X_train = gen_data()
##データ描画
#plt.scatter(X_train[:, 0], X_train[:, 1])
#学習
def distance(x1, x2):
return np.sum((x1 - x2)**2, axis=1)
n_clusters = 3
iter_max = 100
# 各クラスタ中心をランダムに初期化
centers = X_train[np.random.choice(len(X_train), n_clusters, replace=False)]
for _ in range(iter_max):
prev_centers = np.copy(centers)
D = np.zeros((len(X_train), n_clusters))
# 各データ点に対して、各クラスタ中心との距離を計算
for i, x in enumerate(X_train):
D[i] = distance(x, centers)
# 各データ点に、最も距離が近いクラスタを割り当
cluster_index = np.argmin(D, axis=1)
# 各クラスタの中心を計算
for k in range(n_clusters):
index_k = cluster_index == k
centers[k] = np.mean(X_train[index_k], axis=0)
# 収束判定
if np.allclose(prev_centers, centers):
break
#予測
def plt_result(X_train, centers, xx):
# データを可視化
plt.scatter(X_train[:, 0], X_train[:, 1], c=y_pred, cmap='spring')
# 中心を可視化
plt.scatter(centers[:, 0], centers[:, 1], s=200, marker='X', lw=2, c='black', edgecolor="white")
# 領域の可視化
pred = np.empty(len(xx), dtype=int)
for i, x in enumerate(xx):
d = distance(x, centers)
pred[i] = np.argmin(d)
plt.contourf(xx0, xx1, pred.reshape(100, 100), alpha=0.2, cmap='spring')
y_pred = np.empty(len(X_train), dtype=int)
for i, x in enumerate(X_train):
d = distance(x, centers)
y_pred[i] = np.argmin(d)
xx0, xx1 = np.meshgrid(np.linspace(-10, 10, 100), np.linspace(-10, 10, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T
plt_result(X_train, centers, xx)
#####実行結果：
##６. サポートベクターマシーン
サポートベクターマシーンは２クラス分類のための機械学習手法であり、線形判別関数ともっとも近いデータ点との距離（マージン）が最大となる線形判別関数を求める手法となる。
y=w^Tx+b=\sum_{j=1}^mw_jx_j+b
各点と決定境界との距離は以下であり、
\frac{|w^Tx_i+b|}{||w||}=\frac{t_i(w^Tx_i+b)}{||w||}
したがって、マージンとは決定境界と最も距離の近い点との距離なので、下式で表現できる。
\min_i\frac{t_i(w^Tx_i+b)}{||w||}
すべての点の中で決定境界までの距離が一番近いものを最大化したいので目的関数は以下で表すことができる。
\max_{w,b}[\min_i\frac{t_i(w^Tx_i+b)}{||w||}]
また、マージン上の点においては下記の条件が成り立つため、
t_i(w^Tx_i+b)\geqq 1
最終的には以下を求めることになる。
\max_{w,b}\frac{1}{||w||}
###①線形分離可能 実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
#訓練データ生成
def gen_data():
x0 = np.random.normal(size=50).reshape(-1, 2) - 2.
x1 = np.random.normal(size=50).reshape(-1, 2) + 2.
X_train = np.concatenate([x0, x1])
ys_train = np.concatenate([np.zeros(25), np.ones(25)]).astype(np.int)
return X_train, ys_train
X_train, ys_train = gen_data()
#plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)
#学習
t = np.where(ys_train == 1.0, 1.0, -1.0)
n_samples = len(X_train)
# 線形カーネル
K = X_train.dot(X_train.T)
eta1 = 0.01
eta2 = 0.001
n_iter = 500
H = np.outer(t, t) * K
a = np.ones(n_samples)
for _ in range(n_iter):
grad = 1 - H.dot(a)
a += eta1 * grad
a -= eta2 * a.dot(t) * t
a = np.where(a > 0, a, 0)
#予測
index = a > 1e-6
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]
term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()
xx0, xx1 = np.meshgrid(np.linspace(-5, 5, 100), np.linspace(-5, 5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T
X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
y_project[i] += a * sv_t * sv.dot(X_test[i])
y_pred = np.sign(y_project)
# 訓練データを可視化
plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
s=100, facecolors='none', edgecolors='k')
# 領域を可視化
#plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])
# マージンと決定境界を可視化
plt.quiver(0, 0, 0.1, 0.35, width=0.01, scale=1, color='pink')
#####実行結果：
###②線形分離不可能 実装演習
#####コード：
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
#訓練データ生成
factor = .2
n_samples = 50
linspace = np.linspace(0, 2 * np.pi, n_samples // 2 + 1)[:-1]
outer_circ_x = np.cos(linspace)
outer_circ_y = np.sin(linspace)
inner_circ_x = outer_circ_x * factor
inner_circ_y = outer_circ_y * factor
X = np.vstack((np.append(outer_circ_x, inner_circ_x),
np.append(outer_circ_y, inner_circ_y))).T
y = np.hstack([np.zeros(n_samples // 2, dtype=np.intp),
np.ones(n_samples // 2, dtype=np.intp)])
X += np.random.normal(scale=0.15, size=X.shape)
x_train = X
y_train = y
#plt.scatter(x_train[:,0], x_train[:,1], c=y_train)
#学習
def rbf(u, v):
sigma = 0.8
return np.exp(-0.5 * ((u - v)**2).sum() / sigma**2)
X_train = x_train
t = np.where(y_train == 1.0, 1.0, -1.0)
n_samples = len(X_train)
# RBFカーネル
K = np.zeros((n_samples, n_samples))
for i in range(n_samples):
for j in range(n_samples):
K[i, j] = rbf(X_train[i], X_train[j])
eta1 = 0.01
eta2 = 0.001
n_iter = 5000
H = np.outer(t, t) * K
a = np.ones(n_samples)
for _ in range(n_iter):
grad = 1 - H.dot(a)
a += eta1 * grad
a -= eta2 * a.dot(t) * t
a = np.where(a > 0, a, 0)
#予測
index = a > 1e-6
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]
term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()
xx0, xx1 = np.meshgrid(np.linspace(-1.5, 1.5, 100), np.linspace(-1.5, 1.5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T
X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
y_project[i] += a * sv_t * rbf(X_test[i], sv)
y_pred = np.sign(y_project)
# 訓練データを可視化
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
s=100, facecolors='none', edgecolors='k')
# 領域を可視化
plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])
#####実行結果：
#参考書籍
機械学習の勉強では講師先生のおすすめで以下の本を購入した。
基本の応用数学、Python文法及び機械学習アルゴリズムも詳しく書かれている。
個人的にはアルゴリズムの実装も詳しく説明があったので分かりやすかった。