Python 機械学習プログラミング 達人データサイエンティストによる理論と実践 第3章 写経

目次

始めに

この記事では、[第2版]Python 機械学習プログラミング 達人データサイエンティストによる理論と実践

の3章を写経した記録をまとめます。前回はこちらPython 機械学習プログラミング 達人データサイエンティストによる理論と実践 第2章 写経。 実行はJupyter Labにて実行しています。

写経は以下のようにまとめました。

• 写経ではありますが、関数を作成し、コードの再実行をやりやすいようにした部分があります。

• より良いと思われるコードを考えた場合は書き換えてコメントを添えるようにし、変更点をなるべく明示するようにしてあります。

• 個人的に気になった点のメモを残すようにしてあります。同じような疑問を持った方の助けになれば幸いです。

• 以前書いたコードと同じようなコード（例えばグラフの描写等）は効率化のために飛ばしているところもあります。

• 記事内で使用するモジュールなどは一番最初に宣言するようにしてあります。

[第2版]Python 機械学習プログラミング 達人データサイエンティストによる理論と実践を読んでいている際のちょっとした正誤表代わりになればと思います。

3 分類問題ー機械学習ライブラリscikit-learnの活用

使用モジュール

from sklearn import datasets
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import Perceptron
from sklearn.metrics import accuracy_score
from matplotlib.colors import ListedColormap
import matplotlib.pyplot as plt
from IPython.core.pylabtools import figsize
from sklearn.linear_model import LogisticRegression
from sklearn.svm import SVC
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier
from sklearn.neighbors import KNeighborsClassifier

3.2 scikit-learnの活用へのファーストステップ：パーセプトロンのトレーニング

iris = datasets.load_iris()
X = iris.data[:, [2, 3]]
y = iris.target

print('Class labels:', np.unique(y))

Class labels: [0 1 2]

メモ（sklearn.utils.Bunchオブジェクトについて）

irisはklearn.utils.Bunchというクラスのインスタンスになっています。

type(iris)

Output:

    sklearn.utils.Bunch

klearn.utils.Bunchは辞書に似ているオブジェクトで、アトリビュートによってデータを取り出すことができます。以下に簡単な使い方を紹介しておきます。

import sklearn.utils

以下のようにklearn.utils.Bunchインスタンスを作成します。以下では、キーaに1をキーbに2を所持するklearn.utils.Bunchインスタンスを作成します。

bunch = sklearn.utils.Bunch(a=1, b=2)

通常の辞書のように、コンストラクタで指定したキーの文字列によりデータにアクセスることができます。

bunch['a']

Output:

    1

アトリビュートとしてもデータにアクセスすることができます。

bunch.b

Output:

    2

新しいキーをとバリューを追加することもできます。

bunch.c = 3

bunch['c']

Output:

    3

X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.3, random_state=1, stratify=y)

print('Labels counts in y:', np.bincount(y))
print('Labels counts in y_train:', np.bincount(y_train))
print('Labels counts in y_test:', np.bincount(y_test))

Output:

    Labels counts in y: [50 50 50]
    Labels counts in y_train: [35 35 35]
    Labels counts in y_test: [15 15 15]

sc = StandardScaler()
sc.fit(X_train)
X_train_std = sc.transform(X_train)
X_test_std = sc.transform(X_test)

ppn = Perceptron(max_iter=40, eta0=0.1, random_state=1) # n_iterは今後削除される予定であるというwarningが出た。n_iterはmax_iterに引き継がれたよう。
ppn.fit(X_train_std, y_train)

Output:

    Perceptron(alpha=0.0001, class_weight=None, early_stopping=False, eta0=0.1,
          fit_intercept=True, max_iter=40, n_iter=None, n_iter_no_change=5,
          n_jobs=None, penalty=None, random_state=1, shuffle=True, tol=None,
          validation_fraction=0.1, verbose=0, warm_start=False)

y_pred = ppn.predict(X_test_std)
print('Misclassified samples: %d' % (y_test != y_pred).sum())

Output:

    Misclassified samples: 9

print('Accuracy: %.2f' % accuracy_score(y_test, y_pred))

Output:

    Accuracy: 0.80

print('Accuracy: %.2f' % ppn.score(X_test_std, y_test))

Output:

    Accuracy: 0.80

def plot_decision_regions(X, y, classifier, test_idx=None, resolution=0.02):
    markers = ('s', 'x', 'o', '^', 'v')
    colors = ('red', 'blue', 'lightgreen', 'gray', 'cyan')
    cmap = ListedColormap(colors[:len(np.unique(y))])

    x1_min, x1_max = X[:, 0].min() - 1, X[:, 0].max() + 1
    x2_min, x2_max = X[:, 1].min() - 1, X[:, 1].max() + 1
    xx1, xx2 = np.meshgrid(np.arange(x1_min, x1_max, resolution),
                           np.arange(x2_min, x2_max, resolution))
    Z = classifier.predict(np.array([xx1.ravel(), xx2.ravel()]).T)
    Z = Z.reshape(xx1.shape)
    plt.contourf(xx1, xx2, Z, alpha=0.3, cmap=cmap)
    plt.xlim(xx1.min(), xx1.max())
    plt.ylim(xx2.min(), xx2.max())

    for idx, cl in enumerate(np.unique(y)):
        plt.scatter(x=X[y == cl, 0], 
                    y=X[y == cl, 1],
                    alpha=0.8, 
                    c=colors[idx],
                    marker=markers[idx], 
                    label=cl, 
                    edgecolor='black')

    if test_idx:
        X_test, y_test = X[test_idx, :], y[test_idx]

        plt.scatter(X_test[:, 0],
                    X_test[:, 1],
                    c='',
                    edgecolor='black',
                    alpha=1.0,
                    linewidth=1,
                    marker='o',
                    s=100, 
                    label='test set')

X_combined_std = np.vstack((X_train_std, X_test_std))
y_combined = np.hstack((y_train, y_test))

figsize(10, 7)
plt.rcParams['font.size'] = 20

plot_decision_regions(X=X_combined_std, y=y_combined,
                      classifier=ppn, test_idx=range(105, 150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')

plt.tight_layout()
plt.show()

3.3.3 ADALINE実装をロジスティック回帰のアルゴリズムに変換する

class LogisticRegressionGD(object):
    def __init__(self, eta=0.05, n_iter=100, random_state=1):
        self.eta = eta
        self.n_iter = n_iter
        self.random_state = random_state
        
    def fit(self, X, y):
        rgen = np.random.RandomState(self.random_state)
        self.w_ = rgen.normal(loc=0.0, scale=0.01, size=1 + X.shape[1])
        self.cost_ = []

        for i in range(self.n_iter):
            net_input = self.net_input(X)
            output = self.activation(net_input)
            errors = (y - output)
            self.w_[1:] += self.eta * X.T.dot(errors)
            self.w_[0] += self.eta * errors.sum()
            cost = -y.dot(np.log(output)) - ((1 - y).dot(np.log(1 - output)))
            self.cost_.append(cost)
        return self
    
    def net_input(self, X):
        return np.dot(X, self.w_[1:]) + self.w_[0]

    def activation(self, z):
        return 1. / (1. + np.exp(-np.clip(z, -250, 250)))

    def predict(self, X):
        return np.where(self.net_input(X) >= 0.0, 1, 0)

X_train_01_subset = X_train[(y_train == 0) | (y_train == 1)]
y_train_01_subset = y_train[(y_train == 0) | (y_train == 1)]

lrgd = LogisticRegressionGD(eta=0.05, n_iter=1000, random_state=1)
lrgd.fit(X_train_01_subset,
         y_train_01_subset)

plot_decision_regions(X=X_train_01_subset, 
                      y=y_train_01_subset,
                      classifier=lrgd)

plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')

plt.tight_layout()

3.3.4 scikit-learnを使ったロジスティック回帰モデルのトレーニング

lr = LogisticRegression(C=100.0, random_state=1, solver='liblinear', multi_class='ovr') # solverとmulti_class引数を明示的に指定しないとwarningが出力された。
lr.fit(X_train_std, y_train)

plot_decision_regions(X_combined_std, y_combined,
                      classifier=lr, test_idx=range(105, 150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

3.4 サポートベクトルマシンによる最大マージン分類

3.4.2 スラック変数を使った線形分離不可能なケースへの対処

svm = SVC(kernel='linear', C=1.0, random_state=1)
svm.fit(X_train_std, y_train)

plot_decision_regions(X_combined_std, 
                      y_combined,
                      classifier=svm, 
                      test_idx=range(105, 150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

3.5 カーネルSVMを使った非線形問題の求解

3.5.1 線形分離不可能なデータに対するカーネル手法

np.random.seed(1)
X_xor = np.random.randn(200, 2)
y_xor = np.logical_xor(X_xor[:, 0] > 0,
                       X_xor[:, 1] > 0)
y_xor = np.where(y_xor, 1, -1)

plt.scatter(X_xor[y_xor == 1, 0],
            X_xor[y_xor == 1, 1],
            c='b', marker='x',
            label='1')
plt.scatter(X_xor[y_xor == -1, 0],
            X_xor[y_xor == -1, 1],
            c='r',
            marker='s',
            label='-1')

plt.xlim([-3, 3])
plt.ylim([-3, 3])
plt.legend(loc='best')
plt.tight_layout()
plt.show()

svm = SVC(kernel='rbf', random_state=1, gamma=0.10, C=10.0)
svm.fit(X_xor, y_xor)
plot_decision_regions(X_xor, y_xor,
                      classifier=svm)

plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

svm = SVC(kernel='rbf', random_state=1, gamma=0.2, C=1.0)
svm.fit(X_train_std, y_train)

plot_decision_regions(X_combined_std, y_combined,
                      classifier=svm, test_idx=range(105, 150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

svm = SVC(kernel='rbf', random_state=1, gamma=100.0, C=1.0)
svm.fit(X_train_std, y_train)

plot_decision_regions(X_combined_std, y_combined, 
                      classifier=svm, test_idx=range(105, 150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
plt.tight_layout()
#plt.savefig('images/03_16.png', dpi=300)
plt.show()

3.6 決定木学習

3.6.2 決定木の構築

tree = DecisionTreeClassifier(criterion='gini', 
                              max_depth=4, 
                              random_state=1)
tree.fit(X_train, y_train)

X_combined = np.vstack((X_train, X_test))
y_combined = np.hstack((y_train, y_test))
plot_decision_regions(X_combined, y_combined, 
                      classifier=tree, test_idx=range(105, 150))

plt.xlabel('petal length [cm]')
plt.ylabel('petal width [cm]')
plt.legend(loc='upper left')
plt.tight_layout()

plt.show()

3.6.3 ランダムフォレストを使って複数の決定木を結合する

forest = RandomForestClassifier(criterion='gini',
                                n_estimators=25, 
                                random_state=1,
                                n_jobs=2)
forest.fit(X_train, y_train)

plot_decision_regions(X_combined, y_combined, 
                      classifier=forest, test_idx=range(105, 150))

plt.xlabel('petal length [cm]')
plt.ylabel('petal width [cm]')
plt.legend(loc='upper left')
plt.tight_layout()

plt.show()

3.7 k近傍法：怠惰学習アルゴリズム

knn = KNeighborsClassifier(n_neighbors=5, 
                           p=2, 
                           metric='minkowski')
knn.fit(X_train_std, y_train)

plot_decision_regions(X_combined_std, y_combined, 
                      classifier=knn, test_idx=range(105, 150))

plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
plt.tight_layout()

plt.show()

今回の写経は以上です。 ここまで読んでいただきありがとうございました。