Juliaによる機械学習の予測モデル構築・評価

これは,Julia Advent Calendar 2014 14日目の記事です.MLBaseパッケージを用いて機械学習の予測モデルを構築し,評価する方法について説明します.
以下では,Julia0.3.2,MLBase0.5.1,DecisionTree0.3.4,RDatasets0.1.1を使用しています.

Juliaで使用できる機械学習の手法

Juliaで使用できる機械学習の手法には,以下のようなものがある.

手法 パッケージ
決定木 DecisionTree
ランダムフォレスト DecisionTree, RandomForests(by @bicycle1885さん)
サポートベクタマシン SVM, LIBSVM

他の手法については,Awesome Machine Learningにまとまっている.

ランダムフォレストを試してみる

Juliaでランダムフォレストを実行するためには,DecisionTreeパッケージ,またはRandomForestsパッケージ(by @bicycle1885さん)を用いる.ここでは,DecisionTreeパッケージによる実行方法を示す.
Pima Indianデータセットに対して,ランダムフォレストの予測モデルを構築し,テストデータに対する予測結果を評価している.

using RDatasets
using DecisionTree
srand(123)
# 訓練データの作成
Pima_tr = dataset("MASS", "Pima.tr")
X_tr = array(Pima_tr[:, 1:7])
y_tr = array(Pima_tr[:, 8])
# テストデータの作成
Pima_te = dataset("MASS", "Pima.te")
X_te = array(Pima_te[:, 1:7])
y_te = array(Pima_te[:, 8])

# 予測モデルの構築(特徴量の個数:4, 木の個数:2000, 訓練データの割合:70%)
model = build_forest(y_tr, X_tr, 4, 2000, 0.7)
# テストデータに対する予測
pred = apply_forest(model, array(Pima_te[:, 1:7]))
# 分割表
confusion_matrix(y_te[:, 1], pred[:, 1])

以上を実行すると,下記の結果を得る.

Classes:  {"No","Yes"}
Matrix:   
2x2 Array{Int64,2}:
 195  28
  46  63
Accuracy: 0.7771084337349398
Kappa:    0.4723594347321851

他にも,DecisionTreeパッケージにはクロスバリデーションを実行するnfoldCV_forest関数が用意されている.次の例は,Pima Indianデータセットに対して,10-foldのクロスバリデーションを実行している.

# 10-foldクロスバリデーションの実行(特徴量の個数:4, 木の個数:200, fold数:10, 訓練データの割合:70%)
accuracy = nfoldCV_forest(labels, features, 4, 200, 10, 0.7)

以上のように,DecisionTreeパッケージを使用すると,ランダムフォレストの予測モデルを構築し評価することが可能である.しかし,予測モデルに用いるアルゴリズムが異なると,別のパッケージを使用し異なる関数を用いて上記のタスクを実行する必要がある.アルゴリズムによらず,統一的なフレームワークにより予測モデルを構築し評価することができれば非常に便利である.

MLBaseパッケージによる予測モデル構築・評価の統一的なフレームワーク

MLBaseは,RのcaretやmlrPythonのscikit-learnのように,機械学習の予測モデルを構築・評価する汎用的なパッケージである.このパッケージを使用することにより,クロスバリデーション,ハイパーパラメータのグリッドサーチなどの予測モデルの構築・評価に関わるタスクを統一的なフレームワークで実行することが可能になる.

クロスバリデーション

MLBaseでは,クロスバリデーションはcross_validate関数を用いることにより実行できる.予測モデルを構築する関数,予測モデルを評価する関数の2つを作成して,cross_validate関数に渡す必要がある.
以下の例では,Pima Indianデータセットに対して,10-foldのクロスバリデーションを実行している.

using MLBase
srand(123)
# 訓練データの作成
Pima_tr = dataset("MASS", "Pima.tr")
X_tr = array(Pima_tr[:, 1:7])
y_tr = array(Pima_tr[:, 8])
# テストデータの作成
Pima_te = dataset("MASS", "Pima.te")
X_te = array(Pima_te[:, 1:7])
y_te = array(Pima_te[:, 8])

# 予測モデルを構築する関数
function trainfun(inds)
  # ランダムフォレストによる予測モデルの構築(特徴量の個数:4, 木の個数:2000, 訓練に用いるデータの割合:70%)
  model = build_forest(y_tr[inds, 1], X_tr[inds, :], 4, 2000, 0.7)
  return model
end

# 予測モデルを評価する関数
function evalfun(model, inds)
  # 予測
  pred = apply_forest(model, X_tr[inds, :])
  # 分割表
  conf_mat = confusion_matrix(y_tr[inds, 1], pred).matrix
  # Precisionの算出
  prec = conf_mat[1, 1]/sum(conf_mat[:, 1])
  return prec
end

# サンプル数
const n = size(X_tr, 1)

# クロスバリデーションの実行
scores = cross_validate(
  inds -> trainfun(inds),
  (model, inds) -> evalfun(model, inds),
  n,
  Kfold(n, 10))

以上のコードを実行すると,下記の結果が表示される.

10-element Array{Float64,1}:
 0.785714
 0.705882
 0.8     
 0.615385
 0.692308
 0.733333
 0.642857
 0.857143
 0.846154
 0.947368

cross_validate関数の第4引数に,各foldごとに用いるデータのインデクスを指定する.以上では,10-foldのクロスバリデーションによりデータを分割した.クロスバリデーションのスキームとしては,以下の5つの手法

  • K-foldクロスバリーデション(Kfold関数)
  • 層別K-foldクロスバリデーション(StratifiedKfold関数)
  • 一つ抜きクロスバリデーション(LOOCV関数)
  • ランダムサブサンプリング(RandomSub関数)
  • 層別ランダムサブサンプリング(StratifiedRandomSum関数)

の5つが提供されている.

ハイパーパラメータのグリッドサーチ

MLBaseでハイパーパラメータのグリッドサーチを実行するためには,gridtune関数を使用する.

using MLBase
srand(123)
# 訓練データの作成
Pima_tr = dataset("MASS", "Pima.tr")
X_tr = array(Pima_tr[:, 1:7])
y_tr = array(Pima_tr[:, 8])
# テストデータの作成
Pima_te = dataset("MASS", "Pima.te")
X_te = array(Pima_te[:, 1:7])
y_te = array(Pima_te[:, 8])

# 予測モデルを構築する関数
function estfun(ntree, nfeature)
  # 予測モデルの構築(訓練データの割合:70%)
  model = build_forest(y_tr, X_tr, nfeature, ntree, 0.7)
  return model
end

# 予測精度を評価する関数
function evalfun(model)
  # テストデータに対する予測
  pred = apply_forest(model, X_te)
  # 分割表
  conf_mat = confusion_matrix(y_te[:, 1], pred[:, 1]).matrix
  # 評価指標(Precision)の算出
  prec = conf_mat[1, 1]/sum(conf_mat[:, 1])
  return prec
end

# グリッドサーチによるハイパーパラメータの評価(木の個数:{5,10,15}, 特徴量の個数:{4,5,6})
r = gridtune(estfun, evalfun,
             ("ntree", [5, 10, 15]),
             ("nfeature", [4, 5, 6]);
             verbose=true)

以上の例では,ランダムフォレストの木の個数を5,10,15個,特徴量の個数を4,5,6個として,これらの2つのハイパーパラメータのグリッドサーチを実行している.この結果,下記が表示される.

[ntree=5, nfeature=4] => 0.8034934497816594
[ntree=10, nfeature=4] => 0.8035714285714286
[ntree=15, nfeature=4] => 0.8122270742358079
[ntree=5, nfeature=5] => 0.8301886792452831
[ntree=10, nfeature=5] => 0.8285714285714286
[ntree=15, nfeature=5] => 0.7773279352226721
[ntree=5, nfeature=6] => 0.8109243697478992
[ntree=10, nfeature=6] => 0.8302752293577982
[ntree=15, nfeature=6] => 0.8060344827586207
(Ensemble of Decision Trees
Trees:      10
Avg Leaves: 20.6
Avg Depth:  8.5,(10,6),0.8302752293577982)

この結果,木の個数を10個,特徴量の個数を6個としたときに,Precisionは0.830275と最大値となり,この値にハイパーパラメータを設定するのが良さそうであることが分かる.

ハイパーパラメータのグリッドサーチとクロスバリデーションを組み合わせて実行する関数は提供されていないようだが,以上で説明した事項を組み合わせると実行できる.