Parabola JournalMachine Learning

Scikit-learnを使って主成分分析などを中心に遊んだ

Machine Learning Advent Calendar 2014の企画で書いているブログ記事です.最近私生活で使ってたことや,そのときに躓いたことをまとめました.

主成分分析

名前は聞いたことる方が多いと思いますが,できるだけデータの情報を損失することなく新しい軸を作るための手法です.高次元ベクトルデータなどは現実に可視化することは不可能ですが,2次元や3次元ぐらいまでに落としてあげると可視化することが出来ますし,タスクに依っては低次元部分だけのデータで十分なこともあります.
outputよくある例ではこのような二次元データからoutput2このように第一主成分(青破線)と第二主成分(黒破線)を求めます.そうするとデータ間のユークリッド距離は,青線上で測った場合でも概ね元のデータの特徴を保持している,ということが言えます.PCAはscikit-learnに実装されているものを利用するか,自分で固有値分解をすると計算可能です.
x = np.linspace(-3, 3, n)
y = 2.5 * x + np.random.randn(n)

data = np.c_[x, y]
dim  = 1
pca  = sklearn.decomposition.PCA(dim)
result = pca.fit_transform(data)
print data.shape      # (n, 2)
print result.shape    # (n, 1)

完全に余談ですが固有ベクトル計算は縦に出てきます.
cov_data = np.cov(data)
[values, vectors] = scipy.linalg.eigh(cov_data)
sorted_ix = np.argsort(values)[::-1] # 昇順から降順へ(PCAの場合,大きい固有値に対応する固有ベクトルを使う)
values_sorted  = evs[sorted_ix]
vectors_sorted = evc[:,sorted_ix] # ←ココ!!

一週間このミスでハマりました.それはともかく,意外と大事です.データ全部ぶち込んで自動で上手く行けばいいですが,現実にはなかなかそうは上手くいきません.データをじっくり眺めるために敢えて基本的な手法に通してみるのも大事だと思います.

MDS

主成分分析と似たような形でデータを低次元空間で可視化できる手法は他にもたくさんあり,お手軽なものは多次元尺度構成法 – Wikipediaというものがあります(英語版の方が詳しい; Multidimensional scaling – Wikipedia, the free encyclopedia).MDSでは

  • 似ているデータ→近くに置く
  • 似ていないデータ→遠く置く
という凄い単純な原則に則った可視化を行ってくれます.関係パッケージはsklearn.manifoldにあります.同時にsklearn.metrics当たりを使うといい感じに使うことができます.もし距離を計算してある場合には,metrics関係のパッケージは使わなくて良いです.

# 距離行列が計算済み(http://www1.doshisha.ac.jp/~mjin/R/27/27.html)
lbls = ["Hyogo", "Wakayama", "Osaka", "Nara", "Shiga", "Kyoto"]
data = np.matrix([[0, 134,85, 116, 118, 60],
                  [134, 0, 68, 66, 145, 141],
                  [85, 68, 0, 32, 83, 75],
                  [116, 66, 32, 0, 79, 95],
                  [118, 145, 83, 79, 0, 63],
                  [60,141, 75, 95, 63, 0]])
mds = skm.MDS(n_components = 2, \
              dissimilarity='precomputed')
y = mds.fit_transform(data)

こんな図が生成されます.

mds

上のMDSが入っているmanifoldパッケージには他にも面白い手法がたくさん実装されていて,例えば非線形な次元圧縮手法で恐らく一番有名なIsomapが実装されています.いろいろ調べていましたが,結局全部scikit-learnのスタイルで実装されているので,使い方はMDSとだいたい同じです(パラメータを設定,fit, predict/transformをする.終わり).下で使います.

まとめ.scikit-learn内のメソッドはだいたい
  • fit()
  • fit_predict()
  • fit_transform()
といった関数を持っています.ですので何か試したい手法が会った場合,そのクラスコンストラクタに渡す必要があるパラメータをよく読み,fit/fit_predict/fit_transformなどを実行するとだいたい望んだ結果を得る事ができます.またndarrayを操作するとき,スライスや既に実装されているviewの取得方法を上手く使うことで,ムダな計算を省くことができるようです.

オマケ

せっかくなので小ネタとして使ってみます.データとしては今季最高のアニメとごく一部のクラスタで名高い,TVアニメ「結城友奈は勇者である」公式サイトのスペシャルからアイコンダウンロード(ここです→アイコンダウンロード | TVアニメ「結城友奈は勇者である」公式サイト)へ行って取ってきます.全部jpgファイルなので,前処理でpngに変換してあります.

wget ICON_URL_BASE{001..124}.jpg
for v in `ls`; do;
    convert $v ${v%.*}.png;
done;

numpy/scipyから画像を読み込む方法などについてはドキュメントがある(6. 画像の取扱い — 2012.ProgramingLanguage 0.1 documentation),それを利用してみようと思います.

yuyuyu1


Parabola Journal

IBIS2014(11/17-19@名古屋大学)に参加しました

名古屋大学で17日から19日まで開催されたIBIS2014に参加しました.IBISワークショップ本体は,2011年2012年に引き続き3回目の参加です.毎年来年こそは発表しよう(ポスターセッション)と思うのですが,毎年失敗しています.

企画セッション1: 離散アルゴリズムの機械学習応用

  • モンテカルロ木探索≒コンピュータ囲碁で最初に使われ始めたゲーム木探索方法
  • 電聖戦 http://entcog.c.ooco.jp/entcog/densei/
  • 木探索にプレイアウトをくっつけていい感じに勝ちそうな手を探索する
  • ヒューリスティックをくっつけて定石を回避
  • 囲碁だけじゃなくて木探索を扱う問題や,様々なゲーム木に応用可能になった
  • 離散構造と離散分布
  • SAT/MCと確率計算の対応関係から本質的に難しい
  • 離散における列挙≒連続における積分
  • 指数個存在する真理値表のパラメータをどう効率化するか
  • ベイジアンネットワーク≒条件付き独立+CPT
  • 扱いやすくモデルを作る(BDD,ZDD,KC, etc.)
  • 初期のLifted Inference≒命題変数を述語(F.O.)を利用してまとめて効率化
  • 全体的に難しすぎるので簡単に計算可能っぽい部分だけを抽出する
  • Exchangeable Variable Model(PDF): http://homes.cs.washington.edu/~pedrod/papers/mlc14.pdf
  • 大規模グラフ解析のための乱択スケッチ技法
  • 秋葉さんのスライド(slideshare): http://www.slideshare.net/iwiwi/ss-41752585
  • ネオコグニトロン
  • 参考: http://dbnst.nii.ac.jp/pro/detail/498
  • S-cell/C-cell(当時の脳科学的知見より)
  • Sub-sampling(Max-pooring)との関係
  • add-it-self; 仮想ペナルティノードを入れる
  • あまりまじめにメモしてなかった
  • Deep learning: scaling and applications
  • 何もメモが残っていない


企画セッション2: 学習理論

  • Wasserstein幾何とφ-正規分布族
  • iPadで手書き講演(pure math)
  • 講演と大体同内容の講演者の梗概PDF: http://www.kurims.kyoto-u.ac.jp/~kyodo/kokyuroku/contents/pdf/1916-11.pdf
  • Wasserstein distance = earth mover’s distance
  • Current and Future Trends in Computer Vision
  • メモがない


企画セッション3: ビッグデータ利用の社会的側面

  • 産業領域におけるデータ活用への期待と現状
  • ビッグデータ→高解像度リアルタイム化
  • お見合いおばさん=データインタプリタ
  • 痩せるフォーク
  • ゲノムプライバシの保護と個別化医療への展開
  • 低コスト化されたゲノム解析/SNP解析を基とした様々なサービスに関する諸問題
  • 通常のデータ:識別子と属性
  • ゲノム:識別子と個人識別が可能な属性←むずい
  • 永続的に再利用可能な個人識別子をどう活用し,どう安全に補完したら良いのか
  • →新しい安全性定義・秘密計算・統計的評価・情報開示理論
  • 差分プライバシーは現実的アプリケーションでは強すぎる
  • MLを用いた攻撃,評価


とても楽しかったです,ありがとうございました.またそのうち真面目に内容を補完します.やはりポスターセッションがすごく活発でいいですよね.

PythonMachine Learning

円の方程式/人工データ作成/スペクトラムクラスタリング

よくある例題データみたいなのを使ってscikit-learn(sklearn)を使ってスペクトラムクラスタリングを(英語版Wikipedia)を試そうと思っていろいろやっていたのをメモしておく.

円の方程式

凄い当たり前だけど,半径rの円上の点を生成するには三角関数(cos/sin)を使えばいい.媒介変数とし角度θ(0≦θ≦2π)として,θの基準となる軸からθだけ回転した点の(x,y)座標は (r \cos(\theta), r \sin(\theta))

人工データの作成

適当に小さい円(r=1.5)と大きい円(r=3.0)に関して人工データを作成する.どうやって二次元ベクトルを作るかよく分からなかったのでnumpy.r_[X,Y]というベクトルをくっつけてくれるっぽい感じのメソッドをそのまま使った.媒介変数θは0から2πの間を適当に分割して作った.こういう昨日はMatlabと一緒でnumpy.linspaceがやってくれる.2つ作った店の集合(small/large)を縦にくっつけるのにnumpy.r_とは反対でnumpy.c_を使った.ついでに(標準)正規分布の乱数を適当に付けた(numpy.random.randn)
import numpy as np
import scipy as sp
import matplotlib.pyplot as plt
np.random.seed(0)

N = 500
theta = np.linspace(0, 2 * np.pi, N / 2)
r1, r2 = 1.5, 3.0
small = np.c_[r1 * np.cos(theta), r1 * np.sin(theta)]
large = np.c_[r2 * np.cos(theta), r2 * np.sin(theta)]
dataX = np.r_[small, large] + np.random.randn(N, 2) * 0.25

# plot
plt.figure()
plt.plot(dataX[:,0], dataX[:,1], 'bo')
plt.savefig("sample.png")

プロットするとよくある図みたいにこうなる.
sample

スペクトラムクラスタリング

スペクトラムクラスタリングに関してはリンク先が解説してくれている(スペクトラルクラスタリングの話 – おいしいお米の話).numpy/scipyだけではなく,ここからはscikit-learnを使う.スペクトラムクラスタリングに相当する機能はsklearn.clusterに納められているクラスタリング関係のパッケージ群の中に

  • SpectralClustering (クラス)
  • spectral_clustering (関数)
の2つが何故か用意されていて,関数で直接計算する場合とクラスにいろいろとパラメータを設定してfit→predictするという,いつものsklearn方式が使える,と思いきやfit_predictという関数に統合されているという謎な感じになっている.スペクトラムクラスタリングに必要な行列W(affinity matrix)をどう計算するかという問題があるけど,データをnumpy.ndarrayの形(例えば人工データのdataX)で入れると,クラスSpectralClustering自体が計算してくれる.パラメータを設定すると,事前に自分が計算した行列を使うことも可能.以下のように使う.

# spectral clustering by scikit-learn
from sklearn.cluster import SpectralClustering, spectral_clustering
sc = SpectralClustering(n_clusters = 2, eigen_solver='arpack', affinity = 'rbf', gamma = 10.0)
assigned_labels = sc.fit_predict(dataX)

assigned_labels(SpectralClustering#fit_predictの返り値)にはfitしたデータ(dataX)の各データ点がどのクラスタ番号に割り振られたか,というラベル情報が返ってくるので,それとnumpy.whereを使ってデータを分けることが可能.

plt.figure()
c1 = dataX[np.where(assigned_labels == 1)]
c0 = dataX[np.where(assigned_labels == 0)]
plt.plot(c1[:,0], c1[:,1], 'bo')
plt.plot(c0[:,0], c0[:,1], 'rx')
plt.savefig("clustered.png")

プロットするとよくある感じの結果が得られた.

clustered

Parabola Journal

イラストで学ぶ機械学習のMATLABコードをMacのOctaveで実行する

 

こちらにも情報がある.(春山 征吾のくけー : 「イラストで学ぶ機械学習」のMATLABコードを全部試した – livedoor Blog(ブログ)).ただ第1刷でのミスは僕が買った第3刷ではだいたい修正済みの様子.全部は見てないからわからないけれど.

しかし
まったく最小二乗法は最高だぜ!な「イラストで学ぶ機械学習」を読み終えた。 – EchizenBlog-Zwei
に書いてあるとおりイラストがまったく何のために挿入されているか分からないという稀有な本であった.

 
MATLABコードを実行するために,MacOS用のOctaveを使う.インストール方法はここを見るとMacportsとか使ってもいけるようだけど,僕は普通にパッケージ版3.8.0がインストールされていた.

Octave for MacOS X – Octave

流れ

  • CUI用コンソールを起動する(Octave-cli.app)
  • pwd + cdでスクリプト置き場に移動する
  • ファイル名を打ち込む(source.mってファイルに書いたならsourceとうつ)

最初のソースコードの例.setenv打ち込まないとグラフ描画が出来ない.たぶん環境変数に打ち込むのが良いんだろうな.

setenv("GNUTERM", "qt");

n = 50;
N = 1000;
x = linspace(-3, 3, n)';
X = linspace(-3, 3, N)';
pix = pi*x;
y   = sin(pix)./(pix)+0.1*x+0.05*randn(n, 1);

p(:, 1) = ones(n, 1);
P(:, 1) = ones(N, 1);

for j=1:15
  p(:,2*j) = sin(j/2*x); p(:,2*j+1) = cos(j/2*x);
  P(:,2*j) = sin(j/2*X); P(:,2*j+1) = cos(j/2*X);
end

t = p\y; F = P*t;

figure(1); 
clf;
hold on;

axis([-2.8 2.8 -0.5 1.2]);
plot(X, F, 'g-');
plot(x, y, 'bo');



実行結果の例

sample

Parabola JournalプログラミングMachine Learning

Subgroup DiscoveryとPython Orangeで遊んでいた

 
Machine Learning Advent Calendar 2013の一環で書いてます.動かしてみるつもりだったんですが,諸般の事情で動かなかったものもあったりして・・・.機械学習関係で比較的メジャーなタスクは分類と回帰だと思うのですが,RのEcdatパッケージに含まれているHousing data(data(Housing)で使う)を線形回帰してみます.

入力
output

出力(回帰直線を重ねたもの)
output_line

わーい嬉しい.まぁ細かいところを見てみたらこのデータに線形回帰を突っ込むのがそもそもどうなのか,という問題はありますが,別に2次の項を入れた曲線で近似したところでどうこうなる問題でもないでしょう.そこで,結局横軸(x: lotsize)と縦軸(y: price)がダメなんでしょ?という話を考えたい.いえ,別に積極的に考えたくはないのですが.あと問題にもよりますが.

なんとなく線形回帰ダメっぽいという話からすると,このグラフ(軸のとり方)がそもそもあかんのや,という話を考えるしかないでしょう.別にどこかの特徴空間に写像してもらって,うまく回帰するならそれはそれでもちろんいいですけどね,そのあたりはよく分かりません.Housingデータはもともと,price, lotsize, bedrooms, bathrms, stories(謎), driveway, … と他にも属性ありますので,そのあたりもコミコミで計算したらもっといい感じの回帰式にはなると思います(何がいい回帰式なのか,はとりあえず考えない).

  1. とりあえずデータがあります(離散値,整数値,カテゴリ値,二値,etc.で構成される)
  2. なんか処理したい(←今ココ
という状況を思い浮かべていただきたいと思います.こういう状況になるとさすがにパッケージにぶち込んで\(^o^)/というわけにはいかないので,データの整形・特徴抽出・正規化などなどをやる必要があるかもしれません.結局データ事態が言うほどキレイでない,というのが問題なんですね.データの生成がそもそもi.i.d.でないとか,いろんなデータの種類が混じっているとか,要因はいろいろでしょうが,なかなか難しい.

全部ごちゃまぜのデータを扱うのももちろんそれはそれで良いのですが,出来ればなんかキレイなデータにしてから処理した方が,扱いやすくなります.例えば先にクラスタリング(何でもいいので,レコード間の距離か類似度を決めて,例えば階層クラスタリングにぶち込んでグループ分け,みたいなことする)などの手法をかけたデータの前処理を職人芸でやって,さらに機械学習の諸手法で大域モデルを学習すれば,そのグループの中での精度は高くなるでしょう,仮にただの線形回帰だとしても.

ここでようやく表題なのですが,Subgroup Discovery(日本語訳不明)というのはデータマイニングの諸手法の一つで,簡単に言えば,与えられたデータの中で特徴的な部分集合を発見してあげる手法の総称です.(教師あり学習+教師なし学習)/2みたいなポジションです.なんでかというと,目的変数(target; 決定木とかでも決めているあれ,いわゆる目標変数)を決めた上で,マイニングとルール生成と評価を繰り返していくからでしょうね.もちろんこういう手法だけでなく,RやPythonの可視化(Numpy/Scipy->Matplotlib), Wekaほげほげを利用して可視化しながらデータをよーく眺めるというのが,まず最初にやる作業となるとは思いますが,前処理にデータマイニング処理使ってもいいじゃないか.そこで今回はPythonのOrangeという環境でそれをやろうとしました,という話です.データマイナーを名乗るための道のりは長く険しいですね.

パッケージとしてはよく分かりませんが,Python Orange用(http://kt.ijs.si/petra_kralj/SubgroupDiscovery/),Rapid Miner用(http://kt.ijs.si/petra_kralj/SubgroupDiscovery/rm.html),R用(詳細未確認)(http://www.subgroup.cc/)をとりあえず見つけました.

特徴的な部分集合とはどういうことか

線形回帰の例を続けます.全体の回帰式Y_{\mathrm{Price},i} = A + B X_{\mathrm{lotsize}, i} + \varepsilon_iが得られているとします.次にやることは,ある部分集合D (つまり,subgroup)を検索します.この部分集合Dにおいて回帰式を学習すると,例えば別の回帰式Y_{D,\mathrm{Price},i} = A_D + B_D X_{\mathrm{lotsize}, i} + \varepsilon_iが得られるとします.もし,2つの例えば回帰直線の傾きが統計的に有意に異なるなら(つまり帰無仮説H_0\: B_D = B)をチェックする.もし有意に異なるなら,今見ている部分集合Dを,何らかの特徴によって全体と異なる部分だ,と見なす(´・ω・`)

部分集合の作り方: 簡単なイメージ

例えば離散値属性がある場合(性別: 男/女,年齢層:10代/20代/30代/…),これをそのまま条件としたサブグループというのを作ることができる.
  • サブグループ1: 男性
  • サブグループ2: 10代男性
  • サブグループ3: 10代男性既婚
ある条件が有意に異なるサブグループを発見した場合,さらに条件を加えて条件1∧条件2,条件1∧条件2∧条件3,…と頑張って条件を難しくしていく.この探索空間はデータベースの次元(属性数,か)が上がると当然のように手に負えなくなるので,近似的な探索(例えばグループを表現するルールの精度やカバー率を定義して,その値が大きいものから決められた個数だけ発展的に検索するなど.要するにビーム探索)をすることも.というか,基本的に全列挙は諦めていると思う.

上のHousingデータの属性pref_area(Yes/No二値)を例えば対象とすると,Yes(好ましい)とNo(好ましくない)という+/-がデータ中の各レコードに付いていると考えることができる.次にやるのは上で考えた空間の探索で,簡単な規則から難しい規則を生成して,空間を走査していく.

ある条件(ルール)が得られた時,ルールで表現されるレコードのうち,pref_areaが,例えばYesが75%以上であれば,そのルールを特徴的な部分を表現していると見なす.アイテム集合発見における出現割合(# of occurrences / The entire size of DB)みたいなもの.もちろんそういう評価の仕方だけでなく,いろんな評価があってもいいのだけど,空間内のルールと評価関数においてdownward-closure propertyが成立していた方が,探索は簡単になる.

必要なもの

真面目に考えると,Subgroup Discoveryを適用するには次のことを決める必要がある.
  • 何を対象として集合を分けるのか
  • 決定木作るときの対象属性と同じ.Subgroup Discoveryは基本的に1つの対象変数のみを扱う.2つ以上を同時に扱うような拡張は,一部でExceptional Model Miningと呼ばれる.
  • そのようにしてサブグループを表現するのか
  • 命題論理式,述語論理式,数値区間,その他いろいろを使う.定式化的には探索空間を効率的に検索するための精密化演算子が上手く定義出来,下で述べる基準について単調性が保証される方が望ましい.
  • サブグループ選択のための基準
  • ルールを評価するために使う.カバー率や頻度,エントロピーなど.
  • 探索をどのように効率化するか
  • いくら枝刈りを頑張っても全列挙はだいたいしんどいので近似探索手法も同時に使う.

Orangeで動かしてみる

Subgroup Discovery用のWidgetをインストールしたOrangeにデータをぶち込む
Orange_1

見たい場合はData Tableに通す
Orange_2

インストールしたSubgroup Discovery Widget群を使って
Orange_3

Subgroup Discoveryを動かすと
Orange_4

楽しい✌(‘ω’✌ )三✌(‘ω’)✌三( ✌’ω’)✌ここでは相関規則系のSubgroup Discovery手法を利用しているので,出てきているサブグループの記述もそういう形です.(条件→嬉しい✌(‘ω’✌ )三✌(‘ω’)✌三( ✌’ω’)✌)とりあえず家賃400(単位はポンド・・・か?)あたりでいいらしいです.
Orange_5

データマイナーへの道のりは長く険しいですねぇ(二度目).機械学習手法を用いるとき,自分のデータがどういうデータなのかに注意をはらって,楽しい学習ライフを送りましょうということで,ひとつ.

参考文献
  • F. Herrera et al. An overview on subgroup discovery: foundations and applications

参考リンク