2018-08-27

Python: グローバルスコープにあるオブジェクトの del() でインポートしたときの挙動について

今回は Python の __del__() メソッドでちょっと不思議な挙動を目にしてから色々と調べてみた話。具体的には、グローバルスコープにあるオブジェクトの __del__() で別のモジュールをインポートしてるとき、そのオブジェクトがプロセス終了時に破棄されると場合によっては例外になる。ただし、これは Python の仕様かというとかなり微妙で CPython の 3.x 系でしか同じ問題は観測できていない。少なくとも、同じ CPython でも 2.x 系や、同じ Python 3.x 系の実装でも PyPy3 では発生しない。おそらく実装上の都合によるものだと思う。

使った環境は次の通り。

$ cat /etc/lsb-release 
DISTRIB_ID=Ubuntu
DISTRIB_RELEASE=18.04
DISTRIB_CODENAME=bionic
DISTRIB_DESCRIPTION="Ubuntu 18.04.1 LTS"
$ uname -r
4.15.0-20-generic
$ python3 -V
Python 3.6.5
$ python -V
Python 2.7.15rc1

オブジェクトの `del()` メソッドについて

まずは前提となる知識から。 Python のオブジェクトは __del__() という特殊メソッドを定義することで、自身が破棄されるときの挙動をオーバーライドできる。この __del__() メソッドは、デストラクタやファイナライザとも呼ばれる。

例えば、以下のサンプルコードでは Example クラスに __del__() メソッドを定義している。それを main() 関数の中でインスタンス化した後に del 文を使って明示的に破棄している。

#!/usr/bin/env python
# -*- coding: utf-8 -*-


class Example(object):

    def __init__(self):
        """オブジェクトが作られるとき呼び出される"""
        print('born:', id(self))

    def __del__(self):
        """オブジェクトが破棄されるとき呼び出される"""
        print('died:', id(self))


def main():
    print('making')
    o = Example()  # オブジェクトを作る
    print('deleting')
    del o  # オブジェクトを明示的に破棄する
    print('done')


if __name__ == '__main__':
    main()

上記を実行してみよう。

$ python3 explicitdel.py
making
born: 139716100508416
deleting
died: 139716100508416
done

del 文が発行されたタイミングで __del__() メソッドが呼び出されていることが分かる。

また、明示的に del 文を発行しなくても、オブジェクトの寿命と共に呼び出される。次のサンプルコードでは、先ほどとは異なり明示的に del 文を発行していない。ただし、関数スコープの終了と共にオブジェクトは破棄されることが期待できる。

#!/usr/bin/env python
# -*- coding: utf-8 -*-


class Example(object):

    def __init__(self):
        """オブジェクトが作られるとき呼び出される"""
        print('born:', id(self))

    def __del__(self):
        """オブジェクトが破棄されるとき呼び出される"""
        print('died:', id(self))


def main():
    print('making')
    o = Example()
    # 明示的にオブジェクトを破棄しない
    print('done')


if __name__ == '__main__':
    print('start')
    main()
    print('end')

上記を実際に実行してみよう。

$ python3 implicitdel.py 
start
making
born: 140412486589464
done
died: 140412486589464
end

たしかに main() 関数が終了するタイミングで __del__() メソッドが呼び出されていることが分かる。

本題

ここからが本題なんだけど、例えば以下のようなコードがあったとする。特徴としては二つある。まずひとつ目は __del__() メソッドの中で別のモジュールをインポートしているところ。そしてふたつ目が、そのオブジェクトがグローバルスコープにあるところ。

#!/usr/bin/env python
# -*- coding: utf-8 -*-


class Example(object):

    def __init__(self):
        """オブジェクトが作られるとき呼び出される"""
        print('born:', id(self))

    def __del__(self):
        """オブジェクトが破棄されるとき呼び出される"""
        import sys  # __del__() の中で他のモジュールをインポートする
        print('died:', id(self))


# オブジェクトをグローバルスコープに設置する
print('making')
o = Example()
print('done')
# プロセスが終了するタイミングでオブジェクトが破棄される

で、これを 3.x 系の CPython で実行すると、次のような例外になる。

$ python3 gimplicitdel.py 
making
born: 139889227380104
done
Exception ignored in: <bound method Example.__del__ of <__main__.Example object at 0x7f3a7fb4b588>>
Traceback (most recent call last):
  File "globaldel.py", line 13, in __del__
ImportError: sys.meta_path is None, Python is likely shutting down

なんかもう Python が終了しようとしてるからインポートするの無理っすみたいなエラーになってる。

ちなみに、同じ CPython でも 2.x 系であれば例外にならない。

$ python gimplicitdel.py
making
('born:', 140683376878608)
done
('died:', 140683376878608)

回避策 (ワークアラウンド)

この問題が起こる理由は一旦置いといて、とりあえず回避する方法としては以下の三つがある。尚 Python 2.x を使うという選択肢は、もちろんない。

__del__() メソッド内でインポートするのをやめる
オブジェクトを del 文で明示的に破棄する
CPython 以外の実装を使う

ひとつ目の回避策は一見するともっともで、そもそもファイルの先頭以外でのインポートは PEP8 に準拠していない。とはいえ、現実にはそうもいかない場合があって。例えば標準パッケージでもファイルの先頭以外でインポートしている例は見つかる。このようなコードを間接的にでも呼び出すと、同じ問題が起こる。

github.com

上記でインポートしているモジュール dbm は、環境によっては存在しない。そこで、あるときだけ使うためにこのようなコードになっているんだと思う。

ふたつ目の選択肢は、ひとつ目がダメなときは現実的になってくるかもしれない。 Python は atexit というモジュールを使ってインタプリタの終了ハンドラが登録できる。その中で del 文を発行すれば良いと思う。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import atexit


class Example(object):

    def __init__(self):
        """オブジェクトが作られるとき呼び出される"""
        print('born:', id(self))

    def __del__(self):
        """オブジェクトが破棄されるとき呼び出される"""
        import sys  # __del__() の中で他のモジュールをインポートする
        print('died:', id(self))


# オブジェクトをグローバルスコープに設置する
print('making')
o = Example()
print('making done')

def _atexit():
    """オブジェクトを後始末する"""
    print('atexit')
    global o  # グローバルスコープの変数を変更する
    print('deleting')
    del o  # 終了ハンドラの中で明示的にオブジェクトを破棄する
    print('deleting done')


def main():
    # 終了ハンドラに登録する
    atexit.register(_atexit)


if __name__ == '__main__':
    main()

上記を実行してみよう。

$ python3 delatexit.py 
making
born: 140689479625640
making done
atexit
deleting
died: 140689479625640
deleting done

たしかに、今度はエラーにならない。

ただ、ここまで書いておいてなんだけど、そもそも本当に回避する必要はあるのか？という点も議論の余地があるかもしれない。問題は Python のプロセスが終了するタイミングの話なので、放っておいても結局のところオブジェクトとか関係なくメモリは開放される。とはいえ、例外のせいで終了時のステータスコードは非ゼロにセットされるし、対処しないと気持ち悪いのはたしか。

ソースコードから問題について調べる

ここからは、この問題について CPython のソースコードを軽く追ってみた話。

まず、前述した例外を出しているのは以下のようだ。

https://github.com/python/cpython/blob/v3.6.6/Lib/importlib/_bootstrap.py#L873,L876

コメントには PyImport_Cleanup() 関数が呼ばれている最中か、あるいは既に呼ばれたことで起こると書いてある。

ドキュメントを確認すると、この関数は用途が内部利用なもののインポート関連のリソースを開放する目的があるらしい。

モジュールのインポート — Python 3.6.6 ドキュメント

ソースコードでいうと、以下に該当する。

https://github.com/python/cpython/blob/v3.6.6/Python/import.c#L335

上記の関数が呼ばれるのは、以下の二箇所かな。

https://github.com/python/cpython/blob/v3.6.6/Python/pylifecycle.c#L608

https://github.com/python/cpython/blob/v3.6.6/Python/pylifecycle.c#L881

上記の PyImport_Cleanup() という関数がオブジェクトが破棄されるよりも前に呼び出されているとアウトっぽいことが分かった。

続いては、実際にタイミングを動的解析で調べてみよう。まずは GDB と Python3 のデバッグ用パッケージをインストールする。

$ sudo apt-get install gdb python3-dbg

次に gdb コマンド経由で Python を起動する。

$ gdb --args python3 gimplicitdel.py

お目当ての関数にブレークポイントを打つ。

(gdb) b PyImport_Cleanup
Breakpoint 1 at 0x573b80: file ../Python/import.c, line 336.

プログラムを走らせるとブレークポイントに引っかかった。

(gdb) run
Starting program: /usr/bin/python3 gimplicitdel.py
[Thread debugging using libthread_db enabled]
Using host libthread_db library "/lib/x86_64-linux-gnu/libthread_db.so.1".
making
born: 140737352545616
done

Breakpoint 1, PyImport_Cleanup () at ../Python/import.c:336
336 ../Python/import.c: No such file or directory.

バックトレースはこんな感じ。さっき確認した場所から呼ばれてる。

(gdb) bt
#0  PyImport_Cleanup () at ../Python/import.c:336
#1  0x0000000000426906 in Py_FinalizeEx () at ../Python/pylifecycle.c:608
#2  0x0000000000426b15 in Py_FinalizeEx () at ../Python/pylifecycle.c:740
#3  0x0000000000441c22 in Py_Main (argc=argc@entry=2, argv=argv@entry=0xa8f260)
    at ../Modules/main.c:830
#4  0x0000000000421ff4 in main (argc=2, argv=<optimized out>)
    at ../Programs/python.c:69

プログラムを進めると、今度は前述した例外が上がった。

(gdb) c
Continuing.
Exception ignored in: <bound method Example.__del__ of <__main__.Example object at 0x7ffff7e7b550>>
Traceback (most recent call last):
  File "gimplicitdel.py", line 13, in __del__
ImportError: sys.meta_path is None, Python is likely shutting down
[Inferior 1 (process 11210) exited normally]

たしかにオブジェクトの __del__() メソッドが呼ばれるより前に PyImport_Cleanup() 関数が呼ばれているようだ。

いじょう。

スマートPythonプログラミング: Pythonのより良い書き方を学ぶ

作者: もみじあめ
発売日: 2016/03/12
メディア: Kindle版
この商品を含むブログ (1件) を見る

2018-08-26

Python: scikit-learn のロジスティック回帰を使ってみる

macOS Python scikit-learn 機械学習統計

最近、意外とロジスティック回帰が使われていることに気づいた。もちろん世間にはもっと表現力のある分類器がたくさんあるけど、問題によってどれくらい複雑なモデルが適しているかは異なる。それに、各特徴量がどのように働くか重みから確認したり、単純なモデルなのでスコアをベンチマークとして利用する、といった用途もあるらしい。今回は、そんなロジスティック回帰を scikit-learn の実装で試してみる。

使った環境は次の通り。

$ sw_vers 
ProductName:    Mac OS X
ProductVersion: 10.13.6
BuildVersion:   17G65
$ python -V               
Python 3.6.6
$ pip list --format=columns | grep -i scikit-learn
scikit-learn    0.19.2

下準備

まずは scikit-learn をインストールしておく。

$ pip install scikit-learn

乳がんデータセットをロジスティック回帰で分類してみる

以下にロジスティック回帰を使って乳がんデータセットを分類するサンプルコードを示す。とはいえ scikit-learn は API が統一されているので、分類器がロジスティック回帰になってる以外に特筆すべき点はないかも。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

from sklearn import datasets
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import StratifiedKFold
from sklearn.model_selection import cross_validate


def main():
    # 乳がんデータセットを読み込む
    dataset = datasets.load_breast_cancer()
    X, y = dataset.data, dataset.target

    # ロジスティック回帰
    clf = LogisticRegression()
    # Stratified K-Fold CV で性能を評価する
    skf = StratifiedKFold(shuffle=True)
    scoring = {
        'acc': 'accuracy',
        'auc': 'roc_auc',
    }
    scores = cross_validate(clf, X, y, cv=skf, scoring=scoring)

    print('Accuracy (mean):', scores['test_acc'].mean())
    print('AUC (mean):', scores['test_auc'].mean())


if __name__ == '__main__':
    main()

実行してみよう。だいたい精度 (Accuracy) の平均で 0.947 前後のスコアが得られた。

$ python logistic.py
Accuracy (mean): 0.9472570314675578
AUC (mean): 0.991659762496548

正直そんなに高くないけど、ロジスティック回帰くらい単純なモデルではこれくらいなんだなっていう指標にはなると思う。

ロジスティック回帰の利点

ロジスティック回帰の良いところは、モデルが単純で解釈も容易なところ。例えば、基本的に線形モデルの眷属なので、各特徴量の重み (傾き) が確認できる。

実際に、それを体験してみよう。今回例に挙げた乳がんデータセットは、腫瘍の特徴を記録した 30 次元の教師データだった。

>>> from sklearn import datasets
>>> dataset = datasets.load_breast_cancer()
>>> dataset.feature_names
array(['mean radius', 'mean texture', 'mean perimeter', 'mean area',
       'mean smoothness', 'mean compactness', 'mean concavity',
       'mean concave points', 'mean symmetry', 'mean fractal dimension',
       'radius error', 'texture error', 'perimeter error', 'area error',
       'smoothness error', 'compactness error', 'concavity error',
       'concave points error', 'symmetry error',
       'fractal dimension error', 'worst radius', 'worst texture',
       'worst perimeter', 'worst area', 'worst smoothness',
       'worst compactness', 'worst concavity', 'worst concave points',
       'worst symmetry', 'worst fractal dimension'], dtype='<U23')
>>> len(dataset.feature_names)
30

目的変数は 0 が良性で 1 が悪性となっている。

>>> import numpy as np
>>> np.unique(y)
array([0, 1])

ホールドアウト検証を使ってデータを分割したら、モデルを学習させよう。

>>> from sklearn.model_selection import train_test_split
>>> X, y = dataset.data, dataset.target
>>> X_train, X_test, y_train, y_test = train_test_split(X, y, shuffle=True, random_state=42)
>>> from sklearn.linear_model import LogisticRegression
>>> clf = LogisticRegression()
>>> clf.fit(X_train, y_train)
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,
          intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
          penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
          verbose=0, warm_start=False)

すると、学習したモデルで切片 (LogisticRegression#intercept_) と重み (LogisticRegression#coef_) が確認できる。重みは各特徴量ごとにある。

>>> clf.intercept_
array([0.40407439])
>>> clf.coef_
array([[ 2.18931343e+00,  1.51512837e-01, -1.57814199e-01,
        -1.03404299e-03, -1.29170075e-01, -4.23805008e-01,
        -6.47620520e-01, -3.37002545e-01, -1.97619418e-01,
        -3.23607668e-02, -6.88409834e-02,  1.48012177e+00,
         4.81243097e-02, -1.05177866e-01, -1.40690243e-02,
        -3.50323361e-02, -7.06715773e-02, -3.93587747e-02,
        -4.81468850e-02, -2.01238862e-03,  1.20675464e+00,
        -3.93262696e-01, -4.96613892e-02, -2.45385329e-02,
        -2.43248181e-01, -1.21314110e+00, -1.60969567e+00,
        -6.01906976e-01, -7.28573372e-01, -1.21974174e-01]])
>>> len(clf.coef_[0])
30

元のデータセットを標準化していないので、重みの値の大小については単純な比較が難しい。ただ、特徴量がプラスであればその特徴量は悪性の方向に、反対にマイナスなら良性の方向に働く。

学習済みモデルの切片と重みから推論内容を確認する

ここからは学習済みモデルの切片と重みから計算した内容が推論と一致することを確認してみる。ロジスティック回帰は線形回帰の式をシグモイド関数で 0 ~ 1 に変換したものになっている。

以下の式は線形回帰の式で、重みが $w$ で切片が $b$ に対応する。ようするに特徴量と重みをかけて切片を足すだけ。

$y = wX + b$

上記をシグモイド関数に放り込むと値が 0 ~ 1 の範囲に収まる。

$z = \frac{1}{1 + e^{-y}}$

これがロジスティック回帰の出力となる。

実際に上記を学習済みモデルで確認してみよう。例えば検証用データの先頭は悪性と判定されている。

>>> clf.predict([X_test[0]])
array([1])
>>> clf.predict_proba([X_test[0]])
array([[0.19165157, 0.80834843]])

悪性の確率は 0.80834843 になる。

正解を確認すると、たしかに悪性のようだ。

>>> y_test[0]
1

それでは学習済みモデルから計算した内容と上記が一致するかを確認してみよう。まずは線形回帰の式を作る。これが上記の $y = wX + b$ に対応する。

>>> import numpy as np
>>> y = np.sum(clf.coef_ * X_test[0]) + clf.intercept_
>>> y
array([1.4393142])

続いてシグモイド関数を定義しておく。

>>> def sigmoid(x):
...     return 1 / (1 + np.exp(-x))
...

あとは、さきほど得られた結果をシグモイド関数に放り込むだけ。

>>> z = sigmoid(y)
>>> z
array([0.80834843])

結果は 0.80834843 となって、見事に先ほど得られた内容と一致している。

ばっちり。

スマートPythonプログラミング: Pythonのより良い書き方を学ぶ

作者: もみじあめ
発売日: 2016/03/12
メディア: Kindle版
この商品を含むブログ (1件) を見る

2018-08-25

Python: scikit-learn の Pipeline 機能をデバッグする

macOS Python pandas scikit-learn 機械学習統計

今回はだいぶ小ネタ。以前にこのブログでも記事にしたことがある scikit-learn の Pipeline 機能について。

blog.amedama.jp

scikit-learn の Pipeline 機能は機械学習に必要となる複数の工程を一つのパイプラインで表現できる。ただ、パイプラインを組んでしまうと途中のフェーズで出力がどうなっているか、とかが確認しにくい問題がある。この問題について調べると以下の StackOverflow が見つかるんだけど、なかなかシンプルな解決方法だった。

stackoverflow.com

先に概要を述べると、特に何もしないフェーズを用意して、そこでデバッグ用の出力をするというもの。

下準備

まずは必要になるパッケージをインストールしておく。

$ pip install pandas scikit-learn scipy numpy

Pipeline に組み込むデバッグ用のオブジェクト

早速だけど以下にサンプルコードを示す。このコードでは Debug という名前でデバッグ用のクラスを用意している。このクラスは scikit-learn の Pipeline に組み込むことができる。実体としては Pipeline が使うメソッドでデバッグ用の出力をするだけの内容になっている。

#!/usr/bin/env python
# -*- coding: utf-8 -*-


from sklearn.pipeline import Pipeline
from sklearn.decomposition import PCA
from sklearn import datasets
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
from sklearn.base import BaseEstimator
from sklearn.base import TransformerMixin

import pandas as pd


class Debug(BaseEstimator, TransformerMixin):
    """Pipelineをデバッグするためのクラス"""

    def transform(self, X):
        # 受け取った X を DataFrame に変換して先頭部分だけを出力する
        print(pd.DataFrame(X).head())
        # データはそのまま横流しする
        return X

    def fit(self, X, y=None, **fit_params):
        # 特に何もしない
        return self


def main():
    # Iris データセットを読み込む
    dataset = datasets.load_iris()
    X, y = dataset.data, dataset.target

    # パイプラインを構成する
    steps = [
        ('pca', PCA()),
        ('debug', Debug()),  # PCA の出力結果を確認する
        ('rf', RandomForestClassifier()),
    ]
    pipeline = Pipeline(steps=steps)

    # 学習
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33, random_state=42)
    pipeline.fit(X_train, y_train)

    # 推論
    y_pred = pipeline.predict(X_test)
    score = accuracy_score(y_test, y_pred)
    print(score)


if __name__ == '__main__':
    main()

上記を実行してみよう。すると、Iris データセットを PCA (主成分分析) した結果の先頭部分が標準出力にプリントされる。二回出力されているのは、モデルの学習 (fit() メソッド) と評価 (predict()) の二回で呼び出されているため。

$ python debug.py                            
          0         1         2         3
0  0.321625 -0.235144  0.057917  0.125637
1  3.355396  0.578542 -0.331641  0.076760
2  0.606069 -0.315582  0.300426  0.187366
3 -2.727847  0.438338  0.013295  0.002542
4  3.455577  0.501194 -0.562918  0.098940
          0         1         2         3
0  0.868341 -0.114257 -0.250542  0.271719
1 -2.233869  0.987378 -0.045914 -0.029639
2  3.746741  0.287862 -0.513685 -0.094163
3  0.760309 -0.111519  0.023542  0.020324
4  1.283430  0.320953 -0.507830 -0.063090
0.96

ばっちり。

スマートPythonプログラミング: Pythonのより良い書き方を学ぶ

作者: もみじあめ
発売日: 2016/03/12
メディア: Kindle版
この商品を含むブログ (1件) を見る

2018-08-25

Python: 層化抽出法を使ったK-分割交差検証 (Stratified K-Fold CV)

macOS Python 機械学習統計 scikit-learn

K-分割交差検証 (K-Fold CV) を用いた機械学習モデルの評価では、元のデータセットを K 個のサブセットに分割する。そして、分割したサブセットの一つを検証用に、残りの K - 1 個を学習用に用いる。

上記の作業で、元のデータセットを K 個のサブセットに分割する工程に着目してみよう。果たして、どのようなルールにもとづいて分割するのが良いのだろうか？このとき、誤ったやり方で分割すると、モデルの学習が上手くいかなかったり、汎化性能を正しく評価できない恐れがある。

今回は、分割方法として層化抽出法を用いたK-分割交差検証 (Stratified K-Fold CV) について書いてみる。この方法を使うと、学習用データと検証用データで目的変数の偏りが少なくなる。実装には scikit-learn の sklearn.model_selection.StratifiedKFold を用いた。

使った環境は次の通り。

$ sw_vers
ProductName:    Mac OS X
ProductVersion: 10.13.6
BuildVersion:   17G65
$ python -V
Python 3.6.6
$ pip list --format=columns | grep -i scikit-learn
scikit-learn 0.19.2

下準備

まずは今回使う Python のパッケージをインストールしておこう。

$ pip install scikit-learn numpy scipy

本当は怖い KFold CV

セクションのタイトルはちょっと煽り気味になっちゃったけど、実際のところ知っていないと怖い。例えば scikit-learn が実装している KFold は、データの分割になかなか大きな落とし穴がある。

次のサンプルコードでは sklearn.model_selection.KFold を使ってデータセットを分割している。問題を単純化するために、データセットには 4 つの要素しか入っていない。そして、目的変数に相当する変数 y には 0 と 1 が 2 つずつ入っている。このデータセットを 2 つに分割 (2-Fold) したとき、どのような結果が得られるだろうか？

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np

from sklearn.model_selection import KFold


def main():
    # 目的変数のつもり
    y = np.array([0, 0, 1, 1])
    # 説明変数のつもり
    X = np.arange(len(y))

    # デフォルトでは先頭からの並びで分割される
    # 目的変数の並びに規則性があると確実に偏りが生じる
    kf = KFold(n_splits=2)
    for train_index, test_index in kf.split(X, y):
        # どう分割されたか確認する
        print('TRAIN:', y[train_index], 'TEST:', y[test_index])


if __name__ == '__main__':
    main()

上記を保存して実行してみよう。学習用とテスト用のサブセットで、目的変数が偏ってしまっている。

$ python kfold.py                                 
TRAIN: [1 1] TEST: [0 0]
TRAIN: [0 0] TEST: [1 1]

仮に、上記のような偏ったデータを機械学習モデルに学習させて評価させた場合を考えてみよう。最初の試行では学習データの目的変数に 1 しかないので 0 のパターンをモデルは覚えることができない。そして、覚えていない 0 だけのデータでモデルが評価されることになる。もちろん、次の試行でも同様に学習データと検証用データが偏ることになる。これでは正しくモデルを学習させて評価することはできない。

どうしてこんなことが起こるかというと、デフォルトで K-Fold はデータの並び順にもとづいて分割するため。例えば、先ほどの例と同じようにデータの並び順に規則性のある Iris データセットを使って検証してみよう。

>>> from sklearn import datasets
>>> dataset = datasets.load_iris()
>>> dataset.target
array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
       0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
       0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
       2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
       2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2])

上記のように Iris データセットはあやめの各品種 (目的変数) ごとに規則性を持ってデータが並んでいる。

このように並び順に規則性を持ったデータセットを scikit-learn の KFold のデフォルトパラメータで分割してみよう。

>>> from sklearn.model_selection import KFold
>>> kf = KFold(n_splits=2)
>>> ite = kf.split(dataset.data, dataset.target)
>>> train_index, test_index = next(ite)
>>> dataset.target[train_index]
array([1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
       2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
       2, 2, 2, 2, 2, 2, 2, 2, 2])
>>> dataset.target[test_index]
array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
       0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
       0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
       1, 1, 1, 1, 1, 1, 1, 1, 1])

データが順番通りに真っ二つにされていることが上記から確認できる。

無作為抽出法を用いたK-分割交差検証

上記のようなデータの偏りを減らす方法として無作為抽出法 (Random Sampling) を使うやり方がある。これは、順番に依存せず無作為にデータを選んでサブセットを作るというもの。

例えば scikit-learn の KFold であれば、オプションに shuffle=True を渡すと無作為抽出になる。次のサンプルコードで試してみよう。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np

from sklearn.model_selection import KFold


def main():
    # 目的変数のつもり
    y = np.array([0, 0, 1, 1])
    # 説明変数のつもり
    X = np.arange(len(y))

    # 無作為抽出法を使って分割する (実行結果は試行によって異なる)
    kf = KFold(n_splits=2, shuffle=True)
    for train_index, test_index in kf.split(X, y):
        # どう分割されたか確認する
        print('TRAIN:', y[train_index], 'TEST:', y[test_index])


if __name__ == '__main__':
    main()

上記を実行してみよう。試行にもよるけど、ちゃんとデータが偏らずに分割されるパターンもある。

$ python rndkfold.py
TRAIN: [0 1] TEST: [0 1]
TRAIN: [0 1] TEST: [0 1]
$ python rndkfold.py
TRAIN: [0 0] TEST: [1 1]
TRAIN: [1 1] TEST: [0 0]

層化抽出法を用いたK-分割交差検証

先ほどの無作為抽出法では試行によってはサブセットに偏りができる場合もあった。もちろん、データセットが大きければ大きいほど大数の法則に従って偏りはできにくくなる。とはいえゼロではないので、そこで登場するのが今回紹介する層化抽出法 (Stratified Sampling) を用いる方法となる。

層化抽出法を使うと、サブセットを作るときに目的変数の比率がなるべく元のままになるように分割できる。次のサンプルコードでは、実装に StratifiedKFold を使うことで層化抽出法を使った分割を実現している。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np

from sklearn.model_selection import StratifiedKFold


def main():
    # 目的変数のつもり
    y = np.array([0, 0, 1, 1])
    # 説明変数のつもり
    X = np.arange(len(y))

    # 層化抽出法を使って分割する
    kf = StratifiedKFold(n_splits=2, shuffle=True)
    for train_index, test_index in kf.split(X, y):
        # どう分割されたか確認する
        print('TRAIN:', y[train_index], 'TEST:', y[test_index])


if __name__ == '__main__':
    main()

上記を実行してみよう。

$ python skfold.py  
TRAIN: [0 1] TEST: [0 1]
TRAIN: [0 1] TEST: [0 1]
$ python skfold.py
TRAIN: [0 1] TEST: [0 1]
TRAIN: [0 1] TEST: [0 1]

何度実行しても偏りができないように分割されることが分かる。

試しに Iris データセットを使ったパターンも確認しておこう。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np

from sklearn.model_selection import StratifiedKFold
from sklearn import datasets


def main():
    # Iris データセットを読み込む
    dataset = datasets.load_iris()
    X, y = dataset.data, dataset.target
    
    # 層化抽出法を使って分割する
    kf = StratifiedKFold(n_splits=2, shuffle=True)
    for train_index, test_index in kf.split(X, y):
        # どう分割されたか確認する
        print('TRAIN:', y[train_index], 'TEST:', y[test_index])


if __name__ == '__main__':
    main()

上記を実行すると分割した結果に偏りがないことが分かる。

$ python skfoldiris.py
TRAIN: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1
 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 2] TEST: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1
 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 2]
TRAIN: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1
 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 2] TEST: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1
 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 2]

めでたしめでたし。

スマートPythonプログラミング: Pythonのより良い書き方を学ぶ

作者: もみじあめ
発売日: 2016/03/12
メディア: Kindle版
この商品を含むブログ (1件) を見る

2018-08-25

ソースコードから Python をインストールするときにビルドされないモジュールを確認する

Linux Ubuntu18.04LTS Unix Python

ソースコードから Python をインストールするとき、環境によってはビルドされないモジュールが出てくる。今回は、どんなモジュールがビルドされなかったかを確認する方法について。先に結論から書くと、ビルドされなかったモジュールがあるときはログにメッセージが残る。

使った環境は次の通り。

$ cat /etc/lsb-release 
DISTRIB_ID=Ubuntu
DISTRIB_RELEASE=18.04
DISTRIB_CODENAME=bionic
DISTRIB_DESCRIPTION="Ubuntu 18.04.1 LTS"
$ uname -r
4.15.0-20-generic

下準備

まずは、あえてビルドされないモジュールが出るように環境を整える。最低限のビルドだけはできるように build-essential パッケージだけインストールしておこう。

$ sudo apt-get -y install build-essential

もちろん使う環境によるけど、これだと結構なモジュールがビルドされないはず。

続いて Python のソースコードを取得する。

$ wget -O - https://www.python.org/ftp/python/3.6.6/Python-3.6.6.tgz | tar zxvf -
$ cd Python-3.6.6

Python をビルドする

あとは ./configure と make を使って Python をビルドする。

$ ./configure
$ make 2>&1 | tee make.log

上記で、ちゃんとビルドのログを残しておくのがポイント。

すると、ビルドしたディレクトリに Python のインタプリタが作られる。

$ file python
python: ELF 64-bit LSB shared object, x86-64, version 1 (SYSV), dynamically linked, interpreter /lib64/ld-linux-x86-64.so.2, for GNU/Linux 3.2.0, BuildID[sha1]=27b252471eb5f2c16941f992d26f2e13ec5e3ada, with debug_info, not stripped

実行すると、ちゃんと動く。

$ ./python -c "print('Hello, World')"
Hello, World

ただし、上記だとおそらくビルドされないモジュールがある。具体的には、ビルドログを以下の文字列で検索をかければ分かる。

$ grep -A 3 "The necessary bits to build these optional modules were not found:" make.log 
The necessary bits to build these optional modules were not found:
_bz2                  _curses               _curses_panel      
_dbm                  _gdbm                 _lzma              
_sqlite3              _tkinter

色々と出てきた。

例えば _sqlite3 という文字列が見えるので sqlite3 パッケージをインポートしてみよう。

$ ./python -c "import sqlite3"
Traceback (most recent call last):
  File "<string>", line 1, in <module>
  File "/home/vagrant/Python-3.6.6/Lib/sqlite3/__init__.py", line 23, in <module>
    from sqlite3.dbapi2 import *
  File "/home/vagrant/Python-3.6.6/Lib/sqlite3/dbapi2.py", line 27, in <module>
    from _sqlite3 import *
ModuleNotFoundError: No module named '_sqlite3'

想定通りではあるけど ModuleNotFoundError というエラーになってしまった。

もちろん、上記のインポートはディストリビューションにデフォルトで入っている Python ならちゃんと動く。

$ which python3
/usr/bin/python3
$ python3 -c "import sqlite3"

つまりビルドした環境に色々とパッケージが足りていない。

ビルドされないモジュールを無くす

ここからは、先ほどビルドされなかったモジュールが、ちゃんとビルドされるように直していこう。例えば Python のインストールマネージャの pyenv では、各環境用に必要なパッケージを示している。

github.com

ここに記載されているコマンドを、そのまま実行してみよう。

$ sudo apt-get install make build-essential libssl-dev zlib1g-dev libbz2-dev \
  libreadline-dev libsqlite3-dev wget curl llvm libncurses5-dev libncursesw5-dev \
  xz-utils tk-dev libffi-dev liblzma-dev

インストールできたら、再度 Python をビルドする。

$ ./configure && make 2>&1 | tee make.log

もう一度ログを確認すると、ほとんどのモジュールはビルドできたようだ。ただ、ごく一部に関してはまだ残っている。

$ grep -A 3 "The necessary bits to build these optional modules were not found:" make.log 
The necessary bits to build these optional modules were not found:
_dbm                  _gdbm                                    
To find the necessary bits, look in setup.py in detect_modules() for the module's name.

該当するモジュールに必要なパッケージをインストールして再度チャレンジする。

$ sudo apt-get install libgdbm-dev libdb-dev
$ ./configure && make 2>&1 | tee make.log

ログを確認すると、今度は grep に何も引っかからなくなった。

$ grep -A 3 "The necessary bits to build these optional modules were not found:" make.log

先ほどエラーになったインポートなどを試しても、今度は問題なく動いている。

$ ./python -c "import sqlite3"
$ ./python -c "import dbm"
$ ./python -c "from dbm import gnu"

ばっちり。

pyenv を使う場合

ちなみに pyenv を使って Python をインストールする場合にも同じ問題が起こる恐れがある。そのときはインストールするときに -v オプションをつけてビルドログを残しておこう。

$ pyenv install 3.6.6 -v 2>&1 | tee make.log

あとは同じ要領でビルドされなかったモジュールを確認できる。

めでたしめでたし。

スマートPythonプログラミング: Pythonのより良い書き方を学ぶ

作者: もみじあめ
発売日: 2016/03/12
メディア: Kindle版
この商品を含むブログ (1件) を見る

2018-08-20

SSHFS を使ってリモートホストのディレクトリをマウントする

Ubuntu18.04LTS

SSH でログインできるリモートホストとのファイルのやり取りは SCP を使うことが多い。ただ、頻繁にやり取りするときは、それも面倒に感じることがある。ただ、あんまり手間のかかる設定作業はしたくない。そんなときは SSHFS を使うと手軽に楽ができそう。

使った環境は次の通り。

$ cat /etc/lsb-release 
DISTRIB_ID=Ubuntu
DISTRIB_RELEASE=18.04
DISTRIB_CODENAME=bionic
DISTRIB_DESCRIPTION="Ubuntu 18.04.1 LTS"
$ uname -r
4.15.0-20-generic

下準備

まず、ファイルを頻繁にやり取りする先として 192.168.33.10 という IP アドレスのホストがいたとする。

$ ping -c 3 192.168.33.10
PING 192.168.33.10 (192.168.33.10) 56(84) bytes of data.
64 bytes from 192.168.33.10: icmp_seq=1 ttl=64 time=0.761 ms
64 bytes from 192.168.33.10: icmp_seq=2 ttl=64 time=0.473 ms
64 bytes from 192.168.33.10: icmp_seq=3 ttl=64 time=0.646 ms

--- 192.168.33.10 ping statistics ---
3 packets transmitted, 3 received, 0% packet loss, time 2029ms
rtt min/avg/max/mdev = 0.473/0.626/0.761/0.121 ms

このホストには、あらかじめ SSH の公開鍵をインストールしておく。

$ ssh-keygen -t rsa -P '' -f ~/.ssh/id_rsa
Generating public/private rsa key pair.
Your identification has been saved in /home/vagrant/.ssh/id_rsa.
Your public key has been saved in /home/vagrant/.ssh/id_rsa.pub.
The key fingerprint is:
SHA256:Yby6vBiMZBv1q76ncrJYOQAMiVbHOaTEFmcGX09pTXQ vagrant@vagrant
The key's randomart image is:
+---[RSA 2048]----+
|o.o==*.. .=o E   |
|=..oB+..oo ..    |
|o...... =.       |
|.  . . . o       |
|. +   . S        |
| + *   o         |
|  * o o          |
| oo..=..         |
|. .**++.         |
+----[SHA256]-----+
$ ssh-copy-id -i ~/.ssh/id_rsa.pub vagrant@192.168.33.10
/usr/bin/ssh-copy-id: INFO: Source of key(s) to be installed: "/home/vagrant/.ssh/id_rsa.pub"
The authenticity of host '192.168.33.10 (192.168.33.10)' can't be established.
ECDSA key fingerprint is SHA256:8qF9F+rTQA5Mqn+DhSCuAdo6jvL6RrXNBDQAEeuSkRk.
Are you sure you want to continue connecting (yes/no)? yes
/usr/bin/ssh-copy-id: INFO: attempting to log in with the new key(s), to filter out any that are already installed
/usr/bin/ssh-copy-id: INFO: 1 key(s) remain to be installed -- if you are prompted now it is to install the new keys
vagrant@192.168.33.10's password: 

Number of key(s) added: 1

Now try logging into the machine, with:   "ssh 'vagrant@192.168.33.10'"
and check to make sure that only the key(s) you wanted were added.

これで、パスワードを入力しなくてもホストに SSH でログインできるようになった。

$ ssh vagrant@192.168.33.10 "echo Hello, World!"
Hello, World!

続いてローカルのマシンには SSHFS をインストールしておく。

$ sudo apt-get -y install sshfs

これで下準備が整った。

SCP を使ってファイルをやり取りする場合

まず、前提として SCP を使ってファイルをやり取りする場合に必要となる作業を確認しておこう。 scp コマンドを使って、転送元と転送先を指定する。

$ echo "Hello, SCP!" > greet.txt
$ scp greet.txt vagrant@192.168.33.10:~/
greet.txt                                                                                                 100%   12     9.7KB/s   00:00

SSH で確認すると、ちゃんとファイルが転送できている。

$ ssh vagrant@192.168.33.10 "cat greet.txt"
Hello, SCP!

ただ、何度もこれをやるのは結構しんどい。

基本的な使い方

続いては SSHFS を使ってみよう。

まずはマウントポイントになるディレクトリを用意する。

$ mkdir sshfsmnt

続いて sshfs コマンドを使って、マウントしたいリモートホストのディレクトリと、ローカルホストのマウント先のディレクトリを指定する。

$ sshfs vagrant@192.168.33.10:/home/vagrant/ sshfsmnt

これで df コマンドの結果に SSHFS でマウントした内容が表示されるようになる。

$ df
Filesystem                           1K-blocks     Used Available Use% Mounted on
udev                                    484508        0    484508   0% /dev
tmpfs                                   100916     5440     95476   6% /run
/dev/mapper/vagrant--vg-root          64800356  1849852  59629060   4% /
tmpfs                                   504560        0    504560   0% /dev/shm
tmpfs                                     5120        0      5120   0% /run/lock
tmpfs                                   504560        0    504560   0% /sys/fs/cgroup
vagrant                              244810132 88761060 156049072  37% /vagrant
tmpfs                                   100912        0    100912   0% /run/user/1000
vagrant@192.168.33.10:/home/vagrant/  64800356  1793032  59685880   3% /home/vagrant/sshfsmnt

マウントしたディレクトリの中身を確認すると、先ほど scp コマンドで転送したファイルもちゃんと見えている。

$ ls sshfsmnt/
greet.txt
$ cat sshfsmnt/greet.txt 
Hello, SCP!

試しに、このファイルを書き換えてみよう。

$ echo "Hello, SSHFS!" > sshfsmnt/greet.txt

改めて ssh コマンドを使ってリモートホスト上でファイルの内容を確認する。

$ ssh vagrant@192.168.33.10 "cat greet.txt"
Hello, SSHFS!

ちゃんとリモートホスト上のファイルが書き換わっていることが分かる。

いじょう。

Linuxで動かしながら学ぶTCP/IPネットワーク入門

作者:もみじあめ
発売日: 2020/02/29
メディア: Kindle版

2018-08-18

Python: ベイズ最適化で機械学習モデルのハイパーパラメータを選ぶ

macOS Python 機械学習統計 Matplotlib scikit-learn

機械学習モデルにおいて、人間によるチューニングが必要なパラメータをハイパーパラメータと呼ぶ。ハイパーパラメータをチューニングするやり方は色々とある。例えば、良さそうなパラメータの組み合わせを全て試すグリッドサーチや、無作為に試すランダムサーチなど。

今回は、それとはちょっと違ったベイズ最適化というやり方を試してみる。ベイズ最適化では、過去の試行結果から次に何処を調べれば良いかを確率分布と獲得関数にもとづいて決める。これにより、比較的少ない試行回数でより優れたハイパーパラメータが選べるとされる。

Python でベイズ最適化をするためのパッケージとしては Bayesian Optimization や skopt、GPyOpt などがある。今回は、その中でも Bayesian Optimization を使ってみることにした。

使った環境は次の通り。

$ sw_vers
ProductName:    Mac OS X
ProductVersion: 10.13.6
BuildVersion:   17G65
$ python -V          
Python 3.6.6

下準備

まずは今回使う Python のパッケージを一通りインストールしておく。

$ pip install bayesian-optimization matplotlib tqdm

基本的な考え方

最初に、単純な例を使って基本的な考え方を説明しておきたい。

まず、機械学習モデルにおけるハイパーパラメータの選択は、未知の関数の出力を最大化あるいは最小化する問題と捉えられる。これは、ハイパーパラメータを入力として、それで学習した機械学習モデルを関数、モデルを評価して得られた何らかの性能指標 (出力) を最も良くする (最大化・最小化) ことが目的のため。つまり、未知の関数の出力が最大あるいは最小となる点を見つけることができるなら、それはハイパーパラメータの選択にも適用できる可能性がある。

ごく単純な例として、次のような関数の出力を最大化することを考えてみよう。グラフからは、だいたい x が 2 あたりで y が最大の 1.4 前後となることが分かる。それ以外の場所では、だいたい 0 と 6 あたりに局所解が存在していて、そこでは y が 1 くらいになる。

f:id:momijiame:20180818132951p:plain

ただし、この関数はあくまで本来は未知で、上記のグラフはあらかじめ分からないという前提になる。その状況で、どうやら x が 2 あたりで y が最大になるっぽいぞ、ということを知りたい。

ちなみに、このグラフは次のようなコードで作った。この中にある関数 f() が未知の関数ということになる。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import math

import numpy as np

from matplotlib import pyplot as plt


def f(x):
    """最大値を見つけたい関数"""
    return math.exp(-(x - 2) ** 2) + math.exp(-(x - 6) ** 2 / 10) + 1 / (x ** 2 + 1)


def main():
    # x が -5 ~ +15 の範囲を 0.1 刻みでプロットする
    X = [x for x in np.arange(-5, 15, 0.1)]
    y = [f(x) for x in X]
    plt.plot(X, y)
    plt.xlabel('x')
    plt.ylabel('y')
    plt.grid()
    plt.show()


if __name__ == '__main__':
    main()

上記を適当な名前で保存して実行すると、先ほどのグラフが得られる。

$ python fplot.py

ベイズ最適化で関数の最大値を見つける

先ほどの関数を例にして、実際に Bayesian Optimization を使って関数が最大となる場所を探してみよう。

早速、以下にサンプルコードを示す。大体の使い方はコメントを見れば分かると思う。基本的には、探すパラメータとその範囲、そして初期点と探索するイテレーションの数を指定するだけで良い。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import math

from bayes_opt import BayesianOptimization


def f(x):
    """最大値を見つけたい関数"""
    return math.exp(-(x - 2) ** 2) + math.exp(-(x - 6) ** 2 / 10) + 1 / (x ** 2 + 1)


def main():
    # 探索するパラメータと範囲を決める
    pbounds = {
        'x': (-5, 15),
    }
    # 探索対象の関数と、探索するパラメータと範囲を渡す
    bo = BayesianOptimization(f=f, pbounds=pbounds)
    # 最大化する
    bo.maximize(init_points=3, n_iter=10)
    # 結果を出力する
    print(bo.res['max'])


if __name__ == '__main__':
    main()

実行結果は次の通り。尚、結果は実行する毎に違ったものになる。試行によっては局所解を答えてしまう場合もあるかも。

$ python bo.py             
Initialization
-----------------------------------------
 Step |   Time |      Value |         x | 
    1 | 00m00s |    0.87632 |    4.6393 | 
    2 | 00m00s |    0.12396 |   -2.6651 | 
    3 | 00m00s |    0.75934 |   -0.5850 | 
Bayesian Optimization
-----------------------------------------
 Step |   Time |      Value |         x | 
    4 | 00m00s |    0.00520 |   14.6498 | 
    5 | 00m00s |    0.39154 |    9.1121 | 
    6 | 00m00s |    1.39524 |    2.0885 | 
    7 | 00m00s |    1.17158 |    1.4247 | 
    8 | 00m00s |    0.80735 |    3.1369 | 
    9 | 00m00s |    0.97558 |    6.6884 | 
   10 | 00m00s |    0.04176 |   11.7992 | 
   11 | 00m00s |    0.03847 |   -5.0000 | 
   12 | 00m00s |    1.02148 |    5.7120 | 
   13 | 00m00s |    0.74644 |    7.7734 | 
{'max_val': 1.3952415083439782, 'max_params': {'x': 2.0884969890991476}}

上記の試行では、関数が最大となる場所は x が 2.08 あたりで、値は 1.39 くらいという結果が得られた。先ほどグラフから目視で読み取った内容と整合している。

試行の過程を可視化してみる

先ほどの例では何となく上手くいくことは分かった。ただ、ベイズ最適化が具体的に何処をどう探索しているのかよく分からない。そこで、続いてはその過程を可視化してみることにする。

次のサンプルコードでは、ベイズ最適化の過程をグラフとしてプロットしている。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import math

import numpy as np

from bayes_opt import BayesianOptimization

from matplotlib import pyplot as plt


def f(x):
    """最大値を見つけたい関数"""
    return math.exp(-(x - 2) ** 2) + math.exp(-(x - 6) ** 2 / 10) + 1 / (x ** 2 + 1)


def plot_bo(bo):
    # プロット範囲 (決め打ち)
    X = [x for x in np.arange(-5, 15, 0.1)]

    # 真の関数
    y = [f(x) for x in X]
    plt.plot(X, y, label='true')

    # サンプル点
    xs = [p['x'] for p in bo.res['all']['params']]
    ys = bo.res['all']['values']
    plt.scatter(xs, ys, c='green', s=20, zorder=10, label='sample')

    # 予測結果
    mean, sigma = bo.gp.predict(np.array(X).reshape(-1, 1), return_std=True)
    plt.plot(X, mean, label='pred')  # 推定した関数
    plt.fill_between(X, mean + sigma, mean - sigma, alpha=0.1)  # 標準偏差

    # 最大値
    max_x = bo.res['max']['max_params']['x']
    max_y = bo.res['max']['max_val']
    plt.scatter([max_x], [max_y], c='red', s=50, zorder=10, label='pred_max')


def main():
    # 探索するパラメータと範囲を決める
    pbounds = {
        'x': (-5, 15),
    }

    # 探索対象の関数と、探索するパラメータと範囲を渡す
    bo = BayesianOptimization(f=f, pbounds=pbounds)
    # 最大化する
    bo.maximize(init_points=3, n_iter=10)

    # 結果をグラフに描画する
    plot_bo(bo)

    # グラフを表示する
    plt.legend()
    plt.grid()
    plt.show()


if __name__ == '__main__':
    main()

適当な名前をつけて実行してみよう。

$ python boplot.py 
Initialization
-----------------------------------------
 Step |   Time |      Value |         x | 
    1 | 00m00s |    0.06077 |   -3.9330 | 
    2 | 00m00s |    0.74088 |    7.7945 | 
    3 | 00m00s |    0.12176 |   10.6689 | 
Bayesian Optimization
-----------------------------------------
 Step |   Time |      Value |         x | 
    4 | 00m00s |    1.12790 |    2.6188 | 
    5 | 00m00s |    0.87412 |    4.6285 | 
    6 | 00m00s |    1.04842 |    0.0374 | 
    7 | 00m01s |    1.14774 |    1.3868 | 
    8 | 00m01s |    1.39986 |    1.9524 | 
    9 | 00m00s |    0.00473 |   15.0000 | 
   10 | 00m00s |    0.26949 |   -1.6584 | 
   11 | 00m00s |    1.02154 |    6.1962 | 
   12 | 00m00s |    0.01546 |   12.8301 | 
   13 | 00m00s |    0.39699 |    9.0894 |

すると、例えば以下のようなグラフが得られる。ここで true となっているのが真の関数で pred がベイズ最適化で推定した関数となる。推定するのに用いた点は sample で、見つけた最大値と思われる場所は pred_max で図示している。網掛けになっているのは、その周辺を調べていないことから、まだ推定にバラつきが大きいことを示している。

f:id:momijiame:20180818135603p:plain

上記を見ると、いくつかの点を調べることで的確に真の関数に近いものを推定し、最大となる箇所を見つけていることが分かる。

獲得関数 (Acquisition Function)

ベイズ最適化では、次に何処を探すのかを確率分布にもとづいて獲得関数が決める。使う獲得関数によって、まだ調べていない場所の探索に重きを置くのか、それとも着実な改善が見込める場所を探すのか、といった味付けが変わる。

Bayesian Optimization で使える獲得関数は以下の三つで、デフォルトは "ucb" になっている。

Upper Confidence Bound (ucb)
Probability of Improvement (poi)
Expected Improvement (ei)

ちなみに Probability of Improvement は局所解に陥りやすいため、残り二つのどちらかを使うのが良さそう。

以下のサンプルコードでは、それぞれの獲得関数がどのように試行するのかをグラフに可視化してみた。獲得関数のパラメータも二種類ずつ設定している。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import math

import numpy as np

from bayes_opt import BayesianOptimization

from matplotlib import pyplot as plt

from tqdm import tqdm


def f(x):
    """最大値を見つけたい関数"""
    return math.exp(-(x - 2) ** 2) + math.exp(-(x - 6) ** 2 / 10) + 1 / (x ** 2 + 1)


def plot_bo(bo):
    # プロット範囲 (決め打ち)
    X = [x for x in np.arange(-5, 15, 0.1)]

    # 真の関数
    y = [f(x) for x in X]
    plt.plot(X, y, label='true')

    # サンプル点
    xs = [p['x'] for p in bo.res['all']['params']]
    ys = bo.res['all']['values']
    plt.scatter(xs, ys, c='green', s=20, zorder=10, label='sample')

    # 予測結果
    mean, sigma = bo.gp.predict(np.array(X).reshape(-1, 1), return_std=True)
    plt.plot(X, mean, label='pred')  # 推定した関数
    plt.fill_between(X, mean + sigma, mean - sigma, alpha=0.1)  # 標準偏差

    # 最大値
    max_x = bo.res['max']['max_params']['x']
    max_y = bo.res['max']['max_val']
    plt.scatter([max_x], [max_y], c='red', s=50, zorder=10, label='pred_max')


def main():
    # 探索するパラメータと範囲を決める
    pbounds = {
        'x': (-5, 15),
    }

    acq_params = [
        ('ucb_kappa_1', {
            'acq': 'ucb',
            'kappa': 1,
        }),
        ('ucb_kappa_10', {
            'acq': 'ucb',
            'kappa': 10,
        }),
        ('ei_xi_1e-4', {
            'acq': 'ei',
            'xi': 1e-4,
        }),
        ('ei_xi_1e-1', {
            'acq': 'ei',
            'xi': 1e-1,
        }),
        ('poi_xi_1e-4', {
            'acq': 'poi',
            'xi': 1e-4,
        }),
        ('poi_xi_1e-1', {
            'acq': 'poi',
            'xi': 1e-1,
        }),
    ]

    plt.figure(figsize=(8, 12))

    for i, (name, acq_param) in tqdm(list(enumerate(acq_params))):
        # 探索対象の関数と、探索するパラメータと範囲を渡す
        bo = BayesianOptimization(f=f, pbounds=pbounds, verbose=0)
        # 最大化する
        bo.maximize(init_points=3, n_iter=10, **acq_param)

        # 結果をグラフに描画する
        plt.subplot(3, 2, i + 1)
        plt.title(name)
        plot_bo(bo)

    # グラフを表示する
    plt.show()


if __name__ == '__main__':
    main()

パラメータは、まだちゃんと理解できていないものの、次のように解釈できるっぽい。

kappa
- 獲得関数が "ucb" のときに有効なパラメータ
- 大きくするほど「探索」に重きを置く
xi
- 獲得関数が "poi" と "ei" で有効なパラメータ
- 大きくするほど既知の最大値から改善が見込める幅が大きな場所を探す

上記を適当な名前で保存して実行しよう。終わるのには結構時間がかかるし、いくらか警告が出るかもしれない。

$ python boacq.py
...

すると、次のようなグラフが得られる。

f:id:momijiame:20180818201854p:plain

獲得関数が "poi" でパラメータの "xi" が 1e-4 のパターンでは局所解に陥っていることが分かる。

機械学習モデルにベイズ最適化を適用する

だいぶ回り道したけど、ここからが今回の本題になる。続いては、ベイズ最適化を使ってサポートベクターマシンのハイパーパラメータを選んでみよう。 RBF カーネルのサポートベクターマシンには "C" と "gamma" という二つのハイパーパラメータがある。これをベイズ最適化を使って選ぶことにする。

以下にサンプルコードを示す。基本的にやっていることは先ほどと変わらない。まず、返り値を最大化したい関数を f() として用意している。この関数は、渡されたハイパーパラメータを使ってサポートベクターマシンを学習する。そして、学習したモデルを 3-Fold CV を使って精度 (Accuracy) で評価して、各結果の平均を返すことになる。あとは、その関数の精度がなるべく高くなるようにベイズ最適化でハイパーパラメータを選ぶ。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import functools

from sklearn.svm import SVC
from sklearn.model_selection import KFold
from sklearn.model_selection import cross_val_score
from sklearn import datasets


import numpy as np

from bayes_opt import BayesianOptimization

from matplotlib import pyplot as plt
from matplotlib import cm


def f(X, y, **params):
    """最大化したい関数 (モデルを交差検証して得たスコアを返す)"""
    svm = SVC(kernel='rbf', **params)
    kf = KFold(n_splits=4, shuffle=True, random_state=42)
    scores = cross_val_score(svm, X=X, y=y, cv=kf)
    return scores.mean()


def main():
    dataset = datasets.load_iris()
    X, y = dataset.data, dataset.target

    # 最大化したい関数にデータセットを部分適用する
    pf = functools.partial(f, X, y)

    pbounds = {
        'C': (1e+0, 1e+2),
        'gamma': (1e-2, 1e+1),
    }
    bo = BayesianOptimization(f=pf, pbounds=pbounds)
    # 最大化する
    bo.maximize(init_points=3, n_iter=20)
    # 結果を出力する
    print(bo.res['max'])

    # 試行をプロットする
    xs = [p['C'] for p in bo.res['all']['params']]
    ys = [p['gamma'] for p in bo.res['all']['params']]
    zs = np.array(bo.res['all']['values'])
    s_zs = (zs - zs.min()) / (zs.max() - zs.min())  # 0 ~ 1 の範囲で標準化する
    sc = plt.scatter(xs, ys, c=s_zs, s=20, zorder=10, cmap=cm.cool)
    plt.colorbar(sc)

    plt.xlabel('C')
    plt.xscale('log')
    plt.ylabel('gamma')
    plt.yscale('log')
    plt.grid()
    plt.show()


if __name__ == '__main__':
    main()

上記に適当な名前をつけて実行してみよう。

$ python bosvm.py
...(省略)...
{'max_val': 0.9797297297297298, 'max_params': {'C': 31.622294778743424, 'gamma': 0.011103966665153005}}

上記から、最も汎化性能が高くなるのは C が 31.6 前後で gamma が 0.01 前後のとき、ということが分かった。

同時に、次のようなグラフが得られる。これは、各ハイパーパラメータの値ごとに汎化性能を対数スケールの散布図で表したもの。ハイパーパラメータによる精度の違いは元の数値ではぶっちゃけ微々たるものなので、数値は 0 ~ 1 に正規化してある。また、汎化性能の高低は色で表現している。

f:id:momijiame:20180818232846p:plain

参考までに、グリッドサーチで調べた結果も以下に示す。

f:id:momijiame:20180819173530p:plain

上記のグラフを描くのに使ったソースコードは次の通り。

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np

from sklearn import datasets
from sklearn.svm import SVC
from sklearn.model_selection import GridSearchCV

from matplotlib import pyplot as plt
from matplotlib import cm


def main():
    # アイリスデータセットを読み込む
    dataset = datasets.load_iris()
    X, y = dataset.data, dataset.target

    # 優れたハイパーパラメータを見つけたいモデル
    clf = SVC(kernel='rbf')

    # 試行するパラメータを羅列する
    params = {
        'C': [1e+0, 1e+1, 1e+2],
        'gamma': [1e-2, 1e-1, 1e+0, 1e+1],
    }
    
    # グリッドサーチで優れたハイパーパラメータを探す
    grid_search = GridSearchCV(clf,
                               param_grid=params,
                               cv=3)
    grid_search.fit(X, y)

    # スコアとパラメータをグラフにプロットする
    score_and_params = np.array([(params['C'], params['gamma'], mean) for params, mean, _ in grid_search.grid_scores_])
    scores = score_and_params[:, 2]
    standard_scores = (scores - scores.min()) / (scores.max() - scores.min())
    sc = plt.scatter(score_and_params[:, 0], score_and_params[:, 1], c=standard_scores, s=20, zorder=10, cmap=cm.cool)
    plt.colorbar(sc)

    plt.xlabel('C')
    plt.xscale('log')
    plt.ylabel('gamma')
    plt.yscale('log')
    plt.grid()
    plt.show()


if __name__ == '__main__':
    main()

まとめ

ベイズ最適化を使うと未知の関数の出力が最大あるいは最小になる場所を見つけることができる
機械学習モデルにおけるハイパーパラメータの選択は未知の関数の出力を最大化あるいは最小化する問題と捉えられる
つまり、機械学習モデルのハイパーパラメータ選びにベイズ最適化を適用できる場合がある

スマートPythonプログラミング: Pythonのより良い書き方を学ぶ

作者: もみじあめ
発売日: 2016/03/12
メディア: Kindle版
この商品を含むブログ (1件) を見る

CUBE SUGAR CONTAINER

技術系のこと書きます。

Python: グローバルスコープにあるオブジェクトの del() でインポートしたときの挙動について

オブジェクトの `del()` メソッドについて

本題

回避策 (ワークアラウンド)

ソースコードから問題について調べる

Python: scikit-learn のロジスティック回帰を使ってみる

下準備

乳がんデータセットをロジスティック回帰で分類してみる

ロジスティック回帰の利点

学習済みモデルの切片と重みから推論内容を確認する

Python: scikit-learn の Pipeline 機能をデバッグする

下準備

Pipeline に組み込むデバッグ用のオブジェクト

Python: 層化抽出法を使ったK-分割交差検証 (Stratified K-Fold CV)

下準備

本当は怖い KFold CV

無作為抽出法を用いたK-分割交差検証

層化抽出法を用いたK-分割交差検証

ソースコードから Python をインストールするときにビルドされないモジュールを確認する

下準備

Python をビルドする

ビルドされないモジュールを無くす

pyenv を使う場合

SSHFS を使ってリモートホストのディレクトリをマウントする

下準備

SCP を使ってファイルをやり取りする場合

基本的な使い方

Python: ベイズ最適化で機械学習モデルのハイパーパラメータを選ぶ

下準備

基本的な考え方

ベイズ最適化で関数の最大値を見つける

試行の過程を可視化してみる

獲得関数 (Acquisition Function)

機械学習モデルにベイズ最適化を適用する

まとめ

オブジェクトの __del__() メソッドについて

本題

回避策 (ワークアラウンド)

ソースコードから問題について調べる

下準備

乳がんデータセットをロジスティック回帰で分類してみる

ロジスティック回帰の利点

学習済みモデルの切片と重みから推論内容を確認する

下準備

Pipeline に組み込むデバッグ用のオブジェクト

下準備

本当は怖い KFold CV

無作為抽出法を用いたK-分割交差検証

層化抽出法を用いたK-分割交差検証

下準備

Python をビルドする

ビルドされないモジュールを無くす

pyenv を使う場合

下準備

SCP を使ってファイルをやり取りする場合

基本的な使い方

下準備

基本的な考え方

ベイズ最適化で関数の最大値を見つける

試行の過程を可視化してみる

獲得関数 (Acquisition Function)

機械学習モデルにベイズ最適化を適用する

まとめ

オブジェクトの `del()` メソッドについて