このイベントは北陸先端科学技術大学院大学 (JAIST) の私的内部イベントです. そのため，外部の方の参加はご遠慮ください. 多分交通費がもったいないです.

基本的に python で機械学習を行うハンズオンになります. 実践中心で 色々な手法をとりあえず使ってみて, そこから考えるような感じにするので, アルゴリズム的な紹介は結構大雑把に行います(一日で機械学習をきちんと話せって無理でしょ).

• 基本, 時系列データは対象にしません.

事前準備その他

なるべくサンプルを弄ってもらおうと思いますので python3 環境の構築を行っておいてください.

なお, 発表者の環境は以下の通りですので, OS レベルの齟齬に対してはフォローの対象外です.

• python : 3.4.1
• OS : Linux (Fedora 21)

また, 以下のライブラリを使用します.

• ipython
• pandas
• scikit-learn
• matplotlib (ggplot2 になるかも)

mac を使用している場合, homebrew に任せると嫌な思いをするかも(昔,ちょっとね).

• brew は home 以下だから...
• port からの場合, 問題なかったりする時もあります.

windows の人ガンバ. とりしま, pip でなんとかなると思われ.

linux は apt とは dnf とかでなんとかなるよー

• というか, pandas, scikit-learn, matplotlib は C のライブラリ依存なんで pip だけだとミスる場合も
• とりあえず使ってみる場合は apt 等に任せたほうが楽かも
• 結論 Linux が一番楽だと思うの.

最終的にターミナル(コマンドライン) から

$ ipython

とやると python が開いて,

>> import pandas
>> import sklearn
>> import matplotlib

とやってエラーとかワーニングとか吐かれないなら準備終了です.

機械学習とは何か?

機械学習というと, 一般的にはクラスタリングがどうしても目立ってしまいますが, 別段それしかできないわけではありません.

• 回帰 : パラメータの決定 : GMM, GLM
• 次元圧縮 : 特徴量の次元変更 : PCA, NN
• クラスタリング: 分類器 : SVM, K-means, RF

結局のところ, データベースドに何かのパラメータを最適化することがここでいう学習だと 思ってくれればよいです. その結果として, 例えばクラスタリングができたり, 予測モデルを作成できたり, データを低次元化できたりします.

• 具体例に上げた手法は主にこういう目的で使うことが多いよという話で，別に他の目的に使えないわけではありません -- NN とかクラスタリングにも使うし

手法的な大枠 : (入力データから見た場合)

• 教師あり学習: 正解データを元にパラメータの最適化を行う
• 教師なし学習: 入力データのみからパラメータの最適化を行う

当然, 単純な精度では教師あり学習のほうが高い.

教師あり学習に｢驚き｣はない

手法的な大枠 : (アルゴリズムから見た場合)

• 線形モデル
• 非線形モデル

線形性を持っているかが本質だが, 一般的には, 確率分布を前提にするかしないかでもいいかな?

実践編

さて，ではいよいよハンスオンです． ここでは，実際に python をいじってもらって ４種類のアルゴリズムを試してもらいます．

• K-means: 非線形，教師なし
• Random Forest: 非線形, 教師あり
• Gaussian Mixture Models: 線形，教師なし
• Generalized Linear Models (Logistic regression): 線形，教師あり -- ここで紹介したいのは Logistic regression です -- とは言いつつ，これは GLM の一種です

このうち，多分皆様が使いやすいと思われるのは RF, GLM です．

本題以外

まずちょっと本題以外の話をさせてください． つまり，python 周りの準備を話ます．

ココらへんの話は良い感じのオライリー本が出ています． ここではそのうちの本のさわりだけをやるので興味があったら読んでみてください．

• Pythonによるデータ分析入門

ipython (開発環境)

まず，このハンズオンでは ipython というインタラクティブシェルを使用します． デフォルトのインタラクティブシェルよりも明らかに便利なので．

• 公式サイト

私が好きな機能のみを挙げると以下のようなものがあります．

• 優秀な補完機能
• イントロスペクション
• マジックコマンド
• shell コマンドの使用
• 履歴操作

まあ，ちょっと触って見ましょうか． とりあえず，端末上で以下のコマンドを入力すると ipython が起動します．

$ ipython

• 正しい python のヴァージョンが読み込まれているかのみ注意してください．
• brew とか， apt とかで入れているとデフォルトシェルの python に合わせている場合があります．

で，補完です．補完．起動したIPythonのシェルで import o と入力して キーを押してみましょう． ちゃんと補完できますね．

python のデフォルトシェルは基本補完をしてくれないので， これだけでも便利だということがわかってもらえると思います．

わーい．

• インポートしたモジュールの中まで補完候補となります
• もちろん，自分で宣言した変数名等も補完対象になりますよ

更に、Pythonのオブジェクトやキーワードだけでなく、ファイル名やディレクトリ名も補完対象になります。 例えば filePath = "./ くらいで を押して見たりすればよいかと

• win は知らん．

でね，python って基本，すべてがオブジェクトなので，オブジェクト情報超大事じゃないっすか． これを調べる機能をイントロスペクションと言います．

例えば以下のような変数を作ったとします．

>>> filePath = "./home.txt"
>>> filePath?

するといろんな情報出てきますよね．

わーい

• ?? とするとより詳細な情報が出てきたりします

マジックコマンド とはなんと言うか, ipython 自身を操作するメタコマンドだと思ってくれればいいです. 私の好きな機能のみ列挙しておくと,こんな感じですね.

• %run 外部スクリプトの実行
• %edit エディタを開いてスクリプトを編集
• %paste クリップボードのコードをペーストして実行する
• %save コマンドの実行履歴をファイルに保存
• %hist コマンド履歴を保存
• %reset 変数や名前空間などをクリアする
• %time 実行時間の計測
• %whos 変数のリスト
• %page オブジェクトの中身の確認

多分,今回多用するのは %edit もしくは %run だと思います.

まず, ipython を開いた状態で, %edit と打ってみてください. すると,多分 vi かなんかが開くかと思います(別に好みのエディタに設定可能ですが,ここでは省略). ここに python スクリプトを書いて保存して閉じると, そのスクリプトが実行されます.

例えば, この資料の用に,スクリプトが細切れになっていたり, 各行の先頭に飾りがついている場合にそれを編集するのに便利なほか, ちょっと頭を使って,整理して書きたいとき(関数の作成とか)に便利です.

一方, 例えば tmux 等で画面分割をしてしまい, スクリプトを書きながらテストしていくのが好きな方もいらしゃるでしょう(私はこっち派). その場合便利なのが %run です. これは以下のように使用します.

>>> %run hogehoge.py

こうすると, hogehoge.py が実行され結果が ipython 上に反映されます.

名前空間を確認したい場合は %whos が便利で, 名前空間を綺麗にしたい場合は %reset を使用するがよいです. で, 履歴をファイルに保存したい時には %save hogehoge.py で保存されます.

わーい.

• 他にも色々なコマンドが用意されています.公式サイトとか見るとよいです.

これは多用厳禁ですが, ipython を使用すると python 中で shell コマンドの使用が可能です.

>>> !ls
>>> files = !ls

• いや, perl とか ruby のコマンドも呼び出せるのよ.
• python ってなんなんだろう...

で, 皆さんお気づきのことと思いますが, ipython では In[1]: hogehoge のように実行コマンドに ID が振られるじゃないですか? これ, いつでも呼び出せます.

>>> In[1]
>>> Out[3]

というわけで, ipython を使うと如何に便利になるかがお分かりいただけたと思います.

あ, ちなみに閉じるのは quit と打てばよいです. んで, 閉じたら ipython notebook と入力してみてください. 面白いことがおきますよ.

• でも CUI 好きの私は使わないんだけどね.

pandas (データ処理)

次にデータ操作関連の話をします. 皆さん R 言語って触ったことあります?

あれを知っている方はアレの DataFrame と playr の機能が整っていると思えばよいです. まあ, 使ってみるのが速いです.

とりあえずおもむろに以下のコマンドを使用してみましょう.

>>> from pandas import read_csv
>>> df = read_csv('https://dataminingproject.googlecode.com/svn-history/r44/DataMiningApp/datasets/Iris/iris.csv')
>>> df.head()

• ちなみに pandas を使用する際に, どの関数のみを使うのか決めていない場合には以下のように読み込むのが普通です
  • import pandas as pd
  • これはあとで使うので,実行しておいてください.

皆様の予想通り, read_csv は csv ファイルを読み込むことをしています. これはローカルファイルでもいいですし, どっか web 上の URL を直接叩いてもいいです.

• ちなみにエクセル形式にも対応しています(.xlsx に対応している数すくないライブラリ)
• sqlite も OK

これでもう, エクセルとか使わなくていいですね. わーい.

んで, pandas を使用してデータを読み込むと,基本的に DataFrame型として読み込まれます. この型は データを sql 的に(つまり, 真偽式で) 操作をすることに長けています.

例えば, Species が setosa のデータのみを取得するには以下のようにします.

>>> df[df['Species'] == 'setosa']

以上,以下,未満,より大きいはこんな感じです.

>>> df[df['Sepal Length'] >= 6]
>>> df[df['Sepal Length'] > 6]

結構複雑なことも可能です. 例えば, Sepal Length が全体の平均より大きいデータを取得してみましょうか.

>>> df[df['Sepal Length'] > df['Sepal Length'].mean()]

もうちょい直感的に Sepal Length と Sepal Width との差が 4 より大きいものを取得してみます.

>>> df[df['Sepal Length'] - df['Sepal Width'] > 4]

条件式を複数使いたい場合は以下のようにします.

>>> df[(df['Species'] == 'setosa') & (df['Sepal Length'] > 5)]
>>> df[(df['Species'] == 'setosa') | (df['Sepal Length'] > 5)]

なにこれすげーでしょ.

わーい.

さらに, pandas は DataFrame をグループ化して, そのグループごとに関数を適応させることができます.

>>> grouped = df.groupby('Species')
>>> grouped.mean()

ちょっと自分で関数作って見ましょうか. ここではすごい単純に Species ごとにデータがいくつあるのかをカウントします.

>>> def countData(data)
>>>     return len(data)
>>> grouped.agg({'Sepal Length': countData})

ここで, どのカラムに関数を適応させているのかを指定していることに注意してください. 指定しない場合,すべてのカラムに関数が適応されます.

matplotlib と連携すると, データの可視化も楽にしてくれます.

>>> import matplotlib.pyplot as plt
>>> import matplotlib
>>> matplotlib.style.use('ggplot')
>>> df.plot()
>>> plt.show()

ただし, これはカラム別に色分けをしてくるので注意です. 色分けの対象になるカラムを指定したい場合, 先に gropby() を行う必要があります.

>>> grouped["Sepal Length"].plot()

デフォルトでは line をしますが, kind オプションをつけることで任意の図がかけます.

>>> grouped["Sepal Length"].plot(kind='hist', alpha=0.5)

大体,使用したいと思われる図は網羅しているので, 詳しくは以下のページを確認してみてください.

• Plotting

最後にデータの保存や, 他のライブラリが求める形式への型変換に関して説明します.

とりあえず,色々加工したデータを csv に保存するには,以下のようにします.

>>> df.to_csv('hogehoge.csv')

一方で例えばあるカラムのデータをListに変換したい場合は以下のようにします.

>>> df['Sepal Length'].tolist()

では, 複数のカラムを一般的な型に変形するにはどうするのかというと,これは2通りのやり方があります. 連携したいライブラリが numpy 依存の場合 df.as_matrix() がよいでしょう.

>>> df.as_matrix()

一方, そうではない場合, 内包表記でリストを作成するのが楽かと思います.

>>> [df[x].tolist() for x in df.columns]

• ちなみに scikit-learn は DataFrame 型がそのまま使用できます.

scikit-learn の使い方

んで, ようやく本題. 機械学習に関する実践です.

みんな, 最初の話, 覚えてます? 本日の本題は機械学習を行うことです.

で python には scikit-learn というライブラリがあって,ここには大体今ある有名ドコロの手法が網羅されています.

• 公式サイト

このライブラリをおすすめする理由は, 複数のアルゴリズムが統一的な記法で試すことができる点. ドキュメントがとてもわかり易い点にあります(まー英語だけど).

とりあえず, 基本的な構文は以下の通りです.

>>> estimator = hogehoge(学習モデルのパラメータ)
>>> estimator.fit(学習データ)       # 学習
>>> estimator.predict(未知のデータ) # 推定

• 最初の hogehoge にアルゴリズムの関数名が入ります.

なお, ここでは先ほど読み込んだ, iris のデータを使用します.

• 以下で具体的な使用例を紹介しますが, チューニング等はざっくりと省略しています.
• 本家のドキュメントに書いてあるので,そちらを確認してください.

教師なし/非線形 K-means

とりあえず,機械学習の Hello world じゃないかな? 実際｢とりあえず,試す｣ときに使われます.

アルゴリズムの説明

実践: K-means

では,まーやってみましょう.

>>> from sklearn.cluster import KMeans
>>> # データの読み込み
>>> data = df[df.columns[0:3]]
>>> target = df['Species']
>>>
>>> # データの可視化
>>> color_wheel = {"setosa": "#0392cf", 
>>>                "versicolor": "#7bc043", 
>>>                "virginica": "#ee4035"}
>>> colors = target.map(lambda x: color_wheel.get(x))
>>> ax = scatter_matrix(data, color=colors, alpha=0.6, diagonal='hist')
>>> plt.show()
>>>
>>> # ここからクラスタリング
>>> estimator = KMeans(n_clusters=3)
>>> estimator.fit(data)
>>> estimator.labels_  # クラスタリング結果

これでクラスタリングができました. 結果を見るとデータは 0,1,2 の三種類に 分類されています.

これは, KMeans(n_clusters = 3) でクラスタ数を指定しているからです. このように, KMeans ではクラスタ数のみは指定してやる必要があります.

このままではどのように分類されたのか見難いので少し可視化してみます. これは先程のスクリプトをちょっと変えれば よいですね.

>>> ans = pd.Series(estimator.labels_)
>>> color_wheel = {0: "#b58900",
>>>                1: "#2aa198",
>>>                2: "#d33682"}
>>> colors = ans.map(lambda x: color_wheel.get(x))
>>> ax = scatter_matrix(data, color=colors, alpha=0.6, diagonal='hist')
>>> plt.show()

結構いい感じに弁別できていることがわかるかと思います.

• ただし, あくまでも K-means は教師なし学習なので, どこにもそれが, 自分たちの分離したいものを分離した結果であるという保証はないので注意

まぁ, 見た目上の一致とかよくわかんないので, きちんと評価用の関数を使用してみましょう. metrics には出来上がった学習結果を評価する関数が色々と入っています.

今回は正解データもあることですし, 単純に一致率を見てみましょう.

>>> from sklearn import metrics
>>> metrics.adjusted_rand_score(target, estimator.labels_)

adjusted とあるように, このアルゴリズムでは偶然の一致に関して重みづけをしています. 完全一致をすることはまずないのですが, 大体 0.7 以上の結果が出たらよい感じだと判断するのが通例(私の感覚だけど)です.

• どこまでの精度を求めるのかは作成するアプリケーションによって異なりますので.

教師あり/非線形: ランダムフォレスト

K-means も 非線形な教師なし学習モデルですが, 単純な精度でいうと, ランダムフォレストモデルが デファクトスタンダードです.

実装の基本は決定木学習です. ただし,いっぱい作成します. 名前に ランダム とあるように, いっぱい作成する決定木はランダム性を持ちます. 何がランダムかといえば, 学習に使用するデータと, どの特徴量を使用するのかがランダムです. それらの木の学習結果を, 大雑把に言えば多数決で選んでいって最終的にまとめていく感じですね.

https://upload.wikimedia.org/wikipedia/ja/thumb/9/91/Id3_result.gif/200px-Id3_result.gif

実践: ランダムフォレスト

ランダムフォレストの場合, 教師情報を与えるので,必然的にいくつに分けるのかは明示されます. また, 森の中にいくつ木があるのか指定してやる必要があります.

>>> from sklearn.ensemble import RandomForestClassifier
>>> # データの作成
>>> # 今度は教師ありなので, テストデータと評価データを分けてみる
>>> train = df.sample(100)
>>> train_data = train[train.columns[0:4]]
>>> train_target = train['Species']
>>> 
>>> test = df.ix[set(df.index) - set(train.index)]
>>> test_data = test[test.columns[0:4]]
>>> test_target = test['Species']
>>> 
>>> estimator = RandomForestClassifier(n_estimators=100)
>>> model = estimator.fit(train_data, train_target)
>>> test_ans = model.predict(test_data)
>>> print(test_ans)

これで学習結果ができましたね. 先ほどと同じく評価をしてみましょう.

>>> metrics.adjusted_rand_score(test_target, test_ans)

k-means と比べて精度がよくなっています.

わーい.

教師なし学習/線形: GMM

ここからは線形モデルの話です. 線形モデルの難しい部分はデータに関していろいろ仮定を置く必要がある部分です.

GMM はある特徴量はいくつかのクラスタの混合で生成されていると仮定します. その上で,各クラスタはすべてガウス分布にしたがって生成されると仮定します.

GMM のイメージ図としてこんなものを考えてみます. なんでもいいですが, とりあえず人間の身長とかにしてみましょう. ここで, 入力データとして与えられるものは身長の値だけです. これをヒストグラムにしてみると緑のような感じになったと考えてください.

学習器に与えられる情報はこれだけです. でも, 人間の目からすると,この図はなんとなく,2つの山に見えますよね. このなんとなく2つの山がどこにあるのかを推定するのがこの手法です.

ここで,人間には｢なんとなく｣ですが, 機械には何かの基準が必要です. そこで, ｢山とはガウス分布である｣という前提をおいてやるわけです. すると, 機械は入力データのどこにガウス分布が存在していて(平均), どういう形(分散)なのかを考えればよいことになります.

で, その考えた結果が, 曲線で示した部分です. 後は, この山が何なのかを解釈する必要があります(男女差なのかもしれませんし, 大人と子供の差なのかもしれません). これは人間の仕事です.

実践: GMM

>>> from sklearn.mixture import GMM
>>> 
>>> train_data_1D = train_data["Sepal Length"]
>>> estimator = GMM(n_components=3)
>>> estimator.fit(train_data_1D.reshape(-1,1))
>>> 
>>> scale = np.arange(train_data_1D.min(), train_data_1D.max(), 0.001)
>>> for v in zip(estimator.covars_, estimator.means_, estimator.weights_):
>>>     y = (1./np.sqrt(2*np.pi*v[0])) * np.exp(-(scale - v[1])**2/2/v[0]) * v[2]
>>>     plt.plot(x, y)
>>> plt.show()

今回はわかりやすさのため一次元にしましたが, 当然多次元可能です.

>>> estimator2 = GMM(n_components=3)
>>> estimator2.fit(train_data)
>>> test_ans = estimator2.predict(test_data)
>>> metrics.adjusted_rand_score(test_target, test_ans)

ここで, 何気なくクラスタ数を教えていることに注意してください. この手法でもクラスタ数は教えてやる必要があります.

しかし, ここで線形モデルのいいところが, 分布を仮定しているからこそ評価できる項目があるということです. 正規分布の場合, モデルの性能を測る尺度として AIC, BIC というものがあります. これらの基準は, モデル自体の精度から, そのモデルの複雑さを引いたものとして計算されます.

• なお, 値が低いモデルが良いモデルとされます.

そのため いくつかクラスタ数の異なるモデルを作成して,そのうちで最もよいモデルを選ぶという手法があります.

>>> estimators = [GMM(n_components = x).fit(train_data) for x in range(1,20)]
>>> aics = pd.Series([est.aic(train_data) for est in estimators])
>>> estimators = [GMM(n_components = x).fit(train_data) for x in range(1,20)]
>>> aics = [est.aic(train_data) for est in estimators]
>>> gmmBest = estimators[np.argmin(aics)]

• ただ,よほどデータが綺麗でないと, 結構細かくなってしまいます.

最新の手法として DPGMM という拡張がなされています. これはチャイニーズレストランプロセスというアルゴリズムで,クラスタ数の推定を行います. これまた, ざっくりとした説明ですが, 簡単に言えば, 学習時に新しいデータが来ると そのデータは今まで見たことあるデータか否かを判断します. んで, 今まで見たことないデータであるのならクラスタ数を増やすという処理を繰り返します. すると,論理上は無限個のクラスタ数が作成可能なクラスタ数の推定が行えます.

>>> from sklearn.mixture import DPGMM
>>> estimator = DPGMM()
>>> gmmDP = estimator.fit(train_data)
>>> len(gmmDP.means_[0])

ある程度はクラスタ数の推定もできましたね. わーい

教師あり/線形: ロジスティクス回帰

• 学習データは二項分布に従う
• 二項目が教師データ

要は真偽の分類を行います. この時に,それがどの程度なのかという感覚を推定するモデルです.

実践: ロジスティクス回帰

例によってまずは一次元でやってみます.

>>> from sklearn.linear_model import LogisticRegression
>>> train_target_1D = train["Species"]
>>> train_target.ix[~(train_target_1D == "setosa")] = 0
>>> train_target.ix[train_target_1D == "setosa"] = 1
>>> logistic = LogisticRegression(C=100)
>>> logistic.fit(train_data_1D.reshape(-1,1), train_target_1D)
>>> 
>>> xvals = np.arange(train_data_1D.min(), train_data_1D.max(),1)
>>> predictions = logistic.predict_proba(train_data2=xvals[:,np.newaxis])
>>> probs = [y for [x, y] in predictions]
>>> plt.plot(xvals, probs)
>>> plt.show()

続いて二次元.

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from sklearn import linear_model, datasets
>>> 
>>> # import some data to play with
>>> iris = datasets.load_iris()
>>> X = iris.data[:, :2]  # we only take the first two features.
>>> Y = iris.target
>>> 
>>> h = .02  # step size in the mesh
>>> 
>>> logreg = linear_model.LogisticRegression(C=1e5)
>>> 
>>> # we create an instance of Neighbours Classifier and fit the data.
>>> logreg.fit(X, Y)
>>> 
>>> # Plot the decision boundary. For that, we will assign a color to each
>>> # point in the mesh [x_min, m_max]x[y_min, y_max].
>>> x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
>>> y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
>>> xx, yy = np.meshgrid(np.arange(x_min, x_max, h), np.arange(y_min, y_max, h))
>>> Z = logreg.predict(np.c_[xx.ravel(), yy.ravel()])
>>> 
>>> # Put the result into a color plot
>>> Z = Z.reshape(xx.shape)
>>> plt.figure(1, figsize=(4, 3))
>>> plt.pcolormesh(xx, yy, Z, cmap=plt.cm.Paired)
>>> 
>>> # Plot also the training points
>>> plt.scatter(X[:, 0], X[:, 1], c=Y, edgecolors='k', cmap=plt.cm.Paired)
>>> plt.xlabel('Sepal length')
>>> plt.ylabel('Sepal width')
>>> 
>>> plt.xlim(xx.min(), xx.max())
>>> plt.ylim(yy.min(), yy.max())
>>> plt.xticks(())
>>> plt.yticks(())
>>> 
>>> plt.show()

• ごめん.いい加減眠い
• ちょっとさぼった

以上.