ラビットチャレンジ stage2 機械学習 その他手法

1.主成分分析

1-1.多変量データの持つ構造をより小数個の指標に圧縮

　　○変量の個数を減らすことに伴う、情報の損失はなるべく小さくしたい
　　○少数変数を利用した分析や可視化(2・3次元の場合)が実現可能
　　
　　　学習データ

$$ \boldsymbol{x}_i = (x _{i1}, x _{i2}, \cdots , x _{im}) \in \mathbb{R}^m
$$

　　　平均（ベクトル）

$$ \bar{x} = \frac{1}{n} \sum_{i=1}^{n} \boldsymbol{x}_i $$

　　　データ行列

$$ \bar{X} = (\boldsymbol{x}_1 - \bar{\boldsymbol{x}}, \cdots ,
\boldsymbol{x}_n - \bar{\boldsymbol{x}} )^T \in \mathbb{R}^{n \times m} $$

　　分散共分散行列(復習用)

$$ \Sigma = Var(\bar{X}) = \frac{1}{n} \bar{X}^T \bar{x} $$

　　　線形変換後のベクトル

$$ s_j = (s _{1j} , \cdots , s _{nj})^T = \bar{X}^Ta_j　　　　a_j \in \mathbb{R}^m

$$
　　

1-2.係数ベクトルが変われば線形変換後の値が変化

　○情報の量を分散の大きさと捉える
　○線形変換後の変数の分散が最大となる射影軸を探索
　
　　線形変換後の分散

  Var(s_j) = \frac{1}{n} s_j^T s_j = \frac{1}{n} (\bar{X} a_j)^T (\bar{X} a_j)
  = \frac{1}{n} a_j^T \bar{X} \bar{X} a_j = a_j^T Var(\bar{X}) a_j

1-3.以下の制約付き最適化問題を解く

　　○ノルムが1となる制約を入れる(制約を入れないと無限に解がある)
　
　　　　目的変数：　$\underset{a \in \mathbb{R}^m}{arg \ max \ }a_j^T Var(\bar{X}) a_j
$

　　　　制約条件： $ a_j^T a_j = 1 $

1-4.制約つき最適化問題の解き方

　　○ラグランジュ関数を最大にする係数ベクトルを探索(微分して0になる点)
　
　　　ラグランジュ関数
$$ E(a_j) = a_j^T Var (\bar{X}) a_j - \lambda (a_j^T a_j - 1) $$

1-5.ラグランジュ関数を微分して最適解を求める

　　○元のデータの分散共分散行列の固有値と固有ベクトルが、上記の制約付き最適化問題の解となる

　　　微分
$$ \frac{\partial E (a_j)}{\partial a_j} = 2Var(\bar{X}) a_j - 2 \lambda
a_j = 0 $$

　　　解
$$ Var(\bar{X}) a_j = \lambda a_j $$

　　○射影先の分散は固有値と一致

$$ Var(s_1) = a_1^T Var (\bar{X}) a_1 = \lambda_1 a_1^T a_1 = \lambda_1 $$

1-6.分散共分散行列は正定値対称行列→固有値は必ず0以上・固有ベクトルは直行

　　・分散共分散行列を計算
　　　　　↓
　　・固有値問題を解く
　　　　　↓
　　・(最大)m個の固有値と固有ベクトルのペアが出現
　
　　○k番目の固有値(昇順)並べ、対応する固有ベクトルを第k主成分と呼ぶ

1-7.寄与率

　　○第1~元次元分の主成分の分散は、元のデータの分散と一致
　　　■2次元のデータを2次元の主成分で表示した時、固有値の和と元のデータの分散が一致
　　　■第k主成分の分散は主成分に対応する固有値
　$$ V_{total} = \sum_{i=1}^{m} \lambda_i $$

　　○寄与率：第k主成分の分散の全分散に対する割合(第k主成分が持つ情報量の割合)
　　○累積寄与率：第1-k主成分まで圧縮した際の情報損失量の割合

$$ c_k = \frac{\lambda_k}{ \sum_{i=1}^{m} \lambda_i} $$

$$ r_k = \frac{\sum_{j=1}^{k} \lambda_j}{\sum_{i=1}^{m} \lambda_i} $$

1-8.ハンズオン

　　●乳がんデータ分析

#必要なデータやライブラリを準備

from google.colab import drive
drive.mount('/content/drive')

import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegressionCV
from sklearn.metrics import confusion_matrix
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
%matplotlib inline

cancer_df = pd.read_csv('/content/drive/My Drive/study_ai_ml/data/cancer.csv')
print('cancer df shape: {}'.format(cancer_df.shape))

#実行結果

cancer df shape: (569, 33)

cancer_df

#不要なデータを削除

cancer_df.drop('Unnamed: 32', axis=1, inplace=True)

# 目的変数の抽出
y = cancer_df.diagnosis.apply(lambda d: 1 if d == 'M' else 0)

# 説明変数の抽出
X = cancer_df.loc[:, 'radius_mean':]

# 学習用とテスト用でデータを分離
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)

# 標準化
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

# ロジスティック回帰で学習
logistic = LogisticRegressionCV(cv=10, random_state=0)
logistic.fit(X_train_scaled, y_train)

# 検証
print('Train score: {:.3f}'.format(logistic.score(X_train_scaled, y_train)))
print('Test score: {:.3f}'.format(logistic.score(X_test_scaled, y_test)))
print('Confustion matrix:\n{}'.format(confusion_matrix(y_true=y_test, y_pred=logistic.predict(X_test_scaled))))

#実行結果

Train score: 0.988
Test score: 0.972
Confustion matrix:
[[89  1]
 [ 3 50]]

　　○まずはロジスティック回帰で、学習し検証したところ、testでも97%の精度を出せている。
　　○現時点で、３０個の説明変数を使用して学習しているので、説明変数を減らすことができるかどうか？
　　○主成分分析を行う前に、まずは現状の説明変数に相関があるかを確認する。

#次元圧縮が有効かの確認（相関の確認）

plt.figure(figsize=[20, 20])
seaborn.heatmap(pd.DataFrame(X_train_scaled).corr(), annot=True)

　　○画像の中には相関係数が１のものも幾つか含まれており、相関がある（次元圧縮が有効）であることが分かる。

#寄与率の確認

pca = PCA(n_components=30)
pca.fit(X_train_scaled)
plt.bar([n for n in range(1, len(pca.explained_variance_ratio_)+1)], pca.explained_variance_ratio_)

　　○ここで寄与率を確認したところ、
　　　　第一主成分は約４3％
　　　　第二主成分は約２０％　・・・
　　　　ということが分かる。

　　○実際に第一主成分と、第二主成分の２次元に圧縮し、散布図で確認すると

# PCA
# 次元数2まで圧縮
pca = PCA(n_components=2)
X_train_pca = pca.fit_transform(X_train_scaled)
print('X_train_pca shape: {}'.format(X_train_pca.shape))
# X_train_pca shape: (426, 2)

# 寄与率
print('explained variance ratio: {}'.format(pca.explained_variance_ratio_))
# explained variance ratio: [ 0.43315126  0.19586506]

# 散布図にプロット
temp = pd.DataFrame(X_train_pca)
temp['Outcome'] = y_train.values
b = temp[temp['Outcome'] == 0]
m = temp[temp['Outcome'] == 1]
plt.scatter(x=b[0], y=b[1], marker='o') # 良性は○でマーク
plt.scatter(x=m[0], y=m[1], marker='^') # 悪性は△でマーク
plt.xlabel('PC 1') # 第1主成分をx軸
plt.ylabel('PC 2') # 第2主成分をy軸

2.k近傍法(KNN)

2-1.k近傍法とは

　●分類問題のための機械学習手法
　　○最近傍のデータをk個とってきて、それらが多く所属するクラスに識別

　●kを変化させると結果も変わる

　●kを大きくすると決定境界は滑らかになる

　　識別ルール

  識別クラス = \left\{
  \begin{array}{ll}
   j & \big\{ k \big\} = max \big\{ k_1, \cdots , k_K \big\}  \\
   reject & \big\{ k_i , \cdots , k_j \big\} = max \big\{ k_1, \cdots , k_K 
   \big\}  \\
  \end{array}
  \right.

　　　　※rejectではなく、ランダムに割り当てても良い

　　①学習データ

  \boldsymbol{x}_i = (x _{i1},x _{i2}, \cdots ,x _{im}) \in \mathbb{R}^m

　　②所属するクラス

 \Omega = \big\{ C_1, C_2, \cdots, C_k \big\} 

　　③i番目のデータが所属するクラス

   w_i \in \Omega

　　④入力にもっとも近いk個のデータ集合

   k( \boldsymbol{x} ) = \big\{x _{i1},x _{i2}, \cdots ,x _{ik} \big\}

　　⑤④の中でクラスjに属するデータの数

   k_j

2-2.ハンズオン

#必要なライブラリをインポート

%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の数を５個にしてみると

　　○少しではあるが決定境界線が変化したことも分かる。

3.k-平均法（k-means)

3-1.k平均法(k-means)とは

　●教師なし学習
　
　●クラスタリング手法
　
　●与えられたデータをk個のクラスタに分類する
　
　●中心の初期値を変えるとクラスタリング結果も変わりうる
　　○初期値が近い→うまくクラスタリングできない
　　○初期値が離れる→うまくクラスタリングできる
　
　●kの値を変えるとクラスタリング結果も変わる
　
　　①学習データ

   \boldsymbol{x}_i = (x _{i1},x _{i2},\cdots,x _{im}) \in \mathbb{R}^m

　　②各クラスタの中心

   \mu_k = \frac{\sum_{i=1}^{N} q_{ik} \boldsymbol{x}_i}{\sum_{i=1}^{N} 
   q_{ik}}

　　③所属クラスタの識別

  q_{ik} = \left\{
  \begin{array}{ll}
  1 & (where \boldsymbol{x}_i \in M(\mu_k)) \\
  0 & (else)
  \end{array}
  \right.

  \ M = \big\{ \mu_1,\cdots,\mu_K \big\}

　　④各クラスタ内データと中心の距離の総和

  J(q_{ik},\mu_{k}) = \sum_{i=1}^{N} \sum_{j=1}^{K} q_{ik} \  
  \| \boldsymbol{x}_i - \boldsymbol{\mu_k} \|^2

　　⑤中心を変化させたときの④の変化量

   \frac{\partial J (q_{ik},\mu_k)}{\partial \mu_k} = 2 \sum_{i=1}^{N} q_{ik} 
   (\boldsymbol{x}_i - \mu_k) = 0

　

3-2.k平均法(k-means)のアルゴリズム

　　　１）各クラスタ中心の初期値を設定する
　
　　　２）各データ点に対して、各クラスタ中心との距離を計算し、最も距離が近いクラスタを割り当てる
　
　　　３）各クラスタの平均ベクトル（中心）を計算する
　
　　　４）収束するまで2, 3の処理を繰り返す

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn import cluster, preprocessing, datasets

from sklearn.cluster import KMeans

wine = datasets.load_wine()

X = wine.data

X.shape

#実行結果
(178, 13)

y=wine.target

y.shape
(178, )

wine.target_names

#実行結果
array(['class_0', 'class_1', 'class_2'], dtype='<U7')

model = KMeans(n_clusters=3)

labels = model.fit_predict(X)

df = pd.DataFrame({'labels': labels})
type(df)

def species_label(theta):
    if theta == 0:
        return wine.target_names[0]
    if theta == 1:
        return wine.target_names[1]
    if theta == 2:
        return wine.target_names[2]

df['species'] = [species_label(theta) for theta in wine.target]

pd.crosstab(df['labels'], df['species'])

　　　○今回はワインの等級を分類する。
　　　○データは１７８件あり、説明変数が１３個ある。
　　　○k-meansにより３つのクラスタに分けてみるが今回はうまく行っていない。

4.サポートベクトルマシン（SVM)

4-1.SVMとは

　　●２クラス分類問題の代表的な手法の一つ
　　　○未知のデータに対して、高い予測精度を持つ関数を構築できる。
　　　○回帰問題や、教師なし学習にも応用される。

　

4-2.決定関数と分類境界

　　　●一般に2クラス分類問題では、特徴ベクトルxがどちらのクラスに属するか判定するために
　　　　次の決定関数(decision function)と呼ばれる関数f(x)が使われる。
　
$$ f(x) = \boldsymbol{w}^T \boldsymbol{x} + b $$

　　ここでwは、特徴ベクトルxと同じ次元の数値ベクトルで、 $ w^T x $ は次のように計算する

   \boldsymbol{w}^T \boldsymbol{x} = (w_1, \cdots, w_n) 
   \begin{pmatrix}
     x_1 \\
     \vdots \\ 
     x_n
   \end{pmatrix}
   = w_1 x_1 + w_2 x_2 + \cdots + w_n x_n
   = \sum_{i=1}^{n} w_i x_i

　　●$ w^T $はwの転置ベクトルを表しており、$ w^T x $はいわゆるwとxの内積を表している。

　　●wとxは共にn次元のベクトルだが、ベクトルxの入力に対して、
　　　スカラーの値を返す関数になっている。

4-3.決定関数f(x)の２クラス分類方法

　
　　●ある入力データxに対して、決定関数f（x）を計算し、その符号により２クラスを分類する。

  y = sgn \ f(x) =
  \left\{
  \begin{array}{ll}
    +1 &  ( f(x) > 0 ) \\
    -1 &  ( f(x) < 0 )
  \end{array}
  \right.

　　●変数yはラベルで、sgn(・)は符号関数と呼ばれる関数で、
　　　　引数が正の場合は＋１、負の場合には−１を返す関数である。

　　例えば、特徴ベクトルxが２次元の場合、決定関数f(x)は図のように平面を表す方程式となる。

　　●上の図で、f(x) > 0 の場合（青の領域） では y = 1
　　　　f(x) < 0 の場合（赤の領域) では　y = -1 となる。

　　●図の中にある f(x) = 0 の線は２つのクラスを分ける境界線となる。
　　　　一般的にこのような境界線を　分類境界　と呼ぶ。

4-3.線形サポートベクトル分類(ハードマージン)

　

　　●式で登場したパラメータ w や b を適切に決定するためには、
　　　既にラベルの値がわかっている特徴ベクトルが複数必要になる。
　　　このような特徴ラベルとラベルのセットを訓練データと呼ぶ。

　　●一般に訓練データを分離できる分類境界は複数存在する。
　　　SVMでは、それぞれのクラスのデータが分類境界から
　　　なるべく離れるようにして「最適な」分類境界を決定する。

　　
　　●分類境界を挟んで２つのクラスがどのくらい離れているかをマージンと呼ぶ。

　　●SVMではなるべくこのマージンが大きくなる分類境界を探す
　　　大きいマージンを持つ分類境界を探すことをマージン最大化と呼ぶ。

　　●類境界f(x) = 0と分類境界から最も近くにあるデータ $ x_i $ との距離は、
　　　次のように表現する

  \underset{w,b}{min} \ \frac{1}{2} \  ||w||^2 \ , \  y_i[w^T x + b] \geq
  1　　(i=1, \cdots, n)

4-4.線形サポートベクトル分類(ソフトマージン)

　　●ハードマージンでは訓練データを完璧に分類できるf(x)があるという仮定を置き分類していくが、
　　　　ソフトマージンでは制約条件を緩和し、完璧に分類できないデータでも分類できるようにする。

　●制約条件を緩和したことによって、発生した誤分類されたデータに対する誤差を表す変数
　　　スラック変数を用いて、ハードマージンと同様に、マージンを最大化しつつ、誤差を小さくしていく。

　●式で表すと以下のようになる。

   \underset{w,b}{min} \ \Big[ \frac{1}{2} \  ||w||^2 + C \sum_{i=1}^{n} \xi_i 
   \Big] \ , \ y_i[w^T x + b] \geq 1 - \xi_i \ , \  \xi_i \geq　　
   (i=1, \cdots, n)

　　●係数Cは正則化係数と呼ばれる、正の定数である。
　　　Cが大きいほどハードマージンの場合に近づき、
　　　C→∞ではiが有限の値を持つと目的関数が発散してしまうため、必ずi= 0となる。

4-5.ハンズオン

#必要なライブラリのインポート

%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)

　　●特徴空間上で線形なモデル 𝑦(𝑥)=𝑤𝜙(𝑥)+𝑏 を用い、その正負によって2値分類を行う

　　●最適化問題を解いていくが、ラグランジュ乗数方を用いる
　　　最適化問題にはラグランジュ乗数 𝑎(≥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')

　　●新しいデータ点 𝑥 に対しては、

$$ y(x)=w \phi(X) + b = \sum_{i=1}^n a_i t_i k(x,x_i) + b $$

　　の正負によって分類する。