データサイエンティストのための優しいJulia入門（基本文法編）
この記事はJulia Advent Calendar 2018の11日の投稿です。
Juliaと言えば、その高速性やマクロ、多重ディスパッチといった豊富な言語仕様が特徴です。ですが、ここではそのあたりの優位性にあまり突っ込まず、普段PythonやRを使っているデータサイエンティストのためのJulia入門用Jupyter Notebookを整備しました。 実行環境はv1.0.2です。
基本文法を扱っています。
- 変数と型（数値,文字列,ブーリアン,配列,辞書,タプル,集合）
- 制御構文（if,for)
- 関数（function）
Juliaの型システムや多様な構文糖衣の説明も最小限の解説に抑え、データサイエンティストのための実践的な処理（roundとか、absとか、配列の走査とか、for文のenumerateはどう書くか、mapでラムダ関数的に処理する時はどうするかなど）に寄せた内容にしています。
参考にさせて頂いた記事
Julia高速チュートリアル@bicycle1885
St_Hakky’s blog Julia入門
biostatistics Julia
以上です。
JuliaのインストールとIntelliJ IDEAでの開発環境の構築
Juliaでの開発ってJupyter Notebook使いの人が多そうだけど、Pycharmに慣れている人にはIntelliJ IDEAが整っているのでこちらがおすすめです。
Qiitaにも類似の記事（ページ下部に記載）がありますが、メモまでに。
Juliaのinstall
公式からインストーラをダウンロードして任意の場所に配置。
ここではv1.0.2をインストール。
IntelliJ（PyCharm）のinstall
同じく公式からインストーラをダウンロードしてインストール。
Download PyCharm: Python IDE for Professional Developers by JetBrains
ここではCommunity Edition Version 2018.3をインストール。
Juliaプラグインの導入
Preference -> Plugins -> Marketplace -> 検索窓で"Julia"を入力してEnter -> install
PyCharmの再起動を促されるので再起動。
Create New Projectをすると、JuliaのInterpreterを指定できるようになっている。
初回はJuliaのbinのPathを指定してあげる。
New -> Julia File でファイルを作成。
任意のプログラムで実行確認。
以上。
参考
IntelliJ IDEAでJulia用プラグインを入れたときのメモ - Qiita
IntelliJでJuliaを書く(不具合あり) - Qiita
GitHub - ice1000/julia-intellij: Julia Plugin for IntelliJ IDEA ┗ ┛ ┏ ┓ ┗ ┛
JuliaでDeep Learning vol.0
Deep Learningにはすでに各種の実装が用意されています。
ここでは、Juliaでの実装（Juliaで書かれていて、Juliaで記述する実装）である
を使ってみます。
■Mochaとは
MochaのGitHubより。
Mocha is a Deep Learning framework for Julia, inspired by the C++ framework Caffe.
文字通り、Juliaのためのディープラーニングフレームワークであり、また、Caffeに影響を受けて作られたものです。（実際、記述方法などはCaffeと似ています）
■利用できる手法
同じくMochaのGitHubより。
Effcient implementation of general stochastic gradient solvers and common layers in Mocha could be used to train deep / shallow (convolutional) neural networks, with optinal unsupervised pre-training via (stacked) auto-encoders.
現状のディープラーニングは、大きく以下に２分出来るかと思いますが、
Full-Connected系
- いわゆる(?)ディープラーニング
- Pre-Training必須
- 更にRestricted Boltzmann Machine（RBM）系とAutoencoder系がある
Convolutional Neural Network（CNN）系
- 畳み込み+プーリングを多層にしたもの
- Pre-Training不要
上述を読むと、Mochaでは、Autoencoder系とCNN系ができるみたいです。
RBM系はMochaのDocにもなさそうですが、また今度追いたいと思います。
■豊富なチュートリアル
MochaのDocには、チュートリアルが豊富に載っています。
- Training LeNet on MNIST
- Alex's CIFAR-10 tutorial in Mocha
- CNN系のcuda-convnetのモデルをCIFAR-10のデータに適用
- Image Classification with Pre-trained Model
- CNN系のCaffeから引き継いだ学習済みモデルをilsvrc_2012のデータに適用
- Pre-training with Stacked De-noising Auto-encoders
- Autoencoder系のStacked De-noising Auto-encoders(SDA)のモデルをMNISTのデータに適用
■Training LeNet on MNISTの実践
ここでは、チュートリアルの最初に挙げられるTraining LeNet on MNISTに沿ってMochaを使っていきます。
これは0～9の手書き数字文字であるMNISTを、CNN系で分類するタスクです。
□データセットの整備
MochaからマスターごとDownroad.zipして、以下のシェルスクリプトを実行すると、
\Mocha.jl-master\examples\mnist\get-mnist.sh
データの場所を示す「train.txt」「test.txt」の2つのファイルと、 実際のMNISTデータである「train.hdf5」「test.hdf5」が出来て準備完了です。
シェルスクリプトが使えない環境の場合は、シェルスクリプトに書いてある手順を手動を踏んで、データを用意しましょう。
□Mochaのインストール
Pkg.add("Mocha") #初回のみ実行
using Mocha
□ネットワークの各層の定義
ここでは以下の深層ネットワーク構造を構築。
- 畳み込み層１（conv=~）：フィルター＝20, カーネル=5×5
- →プーリング層１（pool=~）：カーネル=2×2, ストライド=2×2
- →畳み込み層２（conv2=~）：フィルター＝50, カーネル=5×5
- →プーリング層２（pool2=~）：カーネル=2×2, ストライド=2×2
上記のフィルター、カーネル、ストライドのパラメータの意味については、この資料のP.26がわかり易いです。
また、層１, 層２間の伝搬関数にはRectified Linear Unitを指定します。
（fc1~のneuron=Neurons.ReLU()の部分）
#まずランダムシードの設定 srand(12345678) #まず学習データを指定。ミニバッチは64に指定 data_layer = HDF5DataLayer(name="train-data", source="data/train.txt", batch_size=64, shuffle=@windows ? false : true) #ネットワークの各層の定義 conv_layer = ConvolutionLayer(name="conv1", n_filter=20, kernel=(5,5), bottoms=[:data], tops=[:conv]) pool_layer = PoolingLayer(name="pool1", kernel=(2,2), stride=(2,2), bottoms=[:conv], tops=[:pool]) conv2_layer = ConvolutionLayer(name="conv2", n_filter=50, kernel=(5,5), bottoms=[:pool], tops=[:conv2]) pool2_layer = PoolingLayer(name="pool2", kernel=(2,2), stride=(2,2), bottoms=[:conv2], tops=[:pool2]) fc1_layer = InnerProductLayer(name="ip1", output_dim=500, neuron=Neurons.ReLU(), bottoms=[:pool2], tops=[:ip1]) fc2_layer = InnerProductLayer(name="ip2", output_dim=10, bottoms=[:ip1], tops=[:ip2]) #学習のための損失関数も定義 loss_layer = SoftmaxLossLayer(name="loss", bottoms=[:ip2,:label])
□CPUバックエンドの設定
ここでGPUBackend()を指定することも可能
backend = CPUBackend() init(backend)
□ネットワークの定義
ネットワーク名、先ほど作った各層、バックエンドを指定してネットワークを定義。
#先ほど作った各層をまとめる
common_layers = [conv_layer, pool_layer, conv2_layer, pool2_layer, fc1_layer, fc2_layer]
#ネットワーク名、各層をまとめたもの、バックエンドを指定
net = Net("MNIST-train", backend, [data_layer, common_layers..., loss_layer])
□SGD(確率的勾配降下法)ソルバーの設定
各パラメータの説明は、Training LeNet on MNISTに記載があります。
max_iterは最大イテレーション数、regu_coefは正則化係数など。
exp_dir = "snapshots"
params = SolverParameters(max_iter=10000, regu_coef=0.0005,
mom_policy=MomPolicy.Fixed(0.9),
lr_policy=LRPolicy.Inv(0.01, 0.0001, 0.75),
load_from=exp_dir)
solver = SGD(params)
□SGDソルバー実行内容、学習状態保存のタイミングを指定
# 1000イテレーションごとに学習の統計情報を返す setup_coffee_lounge(solver, save_into="$exp_dir/statistics.jld", every_n_iter=1000) # 100イテレーションごとに学習のサマリーを返す add_coffee_break(solver, TrainingSummary(), every_n_iter=100) # 5000イテレーションごとにスナップショットを返す add_coffee_break(solver, Snapshot(exp_dir), every_n_iter=5000)
□学習とテスト
テストデータを指定し、精度評価の指標にaccuracyを指定して、やっと学習の実行。
#テストデータを指定。ミニバッチは100に指定
data_layer_test = HDF5DataLayer(name="test-data", source="data/test.txt", batch_size=100)
#テストデータに対するパフォーマンスは、accuracyで測定
acc_layer = AccuracyLayer(name="test-accuracy", bottoms=[:ip2, :label])
#ネットワーク名、各層をまとめたもの、バックエンドを指定
test_net = Net("MNIST-test", backend, [data_layer_test, common_layers..., acc_layer])
# 1000イテレーションごとにテストデータに対する精度を返す
add_coffee_break(solver, ValidationPerformance(test_net), every_n_iter=1000)
# 学習とテストの実行
solve(solver, net)
最大イテレーション数の10000回を行い、精度は以下のようになりました。
なんと98.9%
01-1 03:40:17:INFO:root:## Performance on Validation Set after 10000 iterations 01-1 03:40:17:INFO:root:--------------------------------------------------------- 01-1 03:40:17:INFO:root: Accuracy (avg over 10000) = 98.9100% 01-1 03:40:17:INFO:root:---------------------------------------------------------
□ネットワークとバックエンドを閉じる
destroy(net) destroy(test_net) shutdown(backend)
以上です。他のチュートリアルやパラメータいじりもやっていきたいですね。
■参考文献
Juliaでスパム判定の機械学習分類器を作る
あと１週間でクリスマスイブですね。
Julia Advent Calendar 2014の17日目の記事（@yutajuly）です。
ここでは、スパムデータ分類を例にとり、実データ特有の問題に対応した機械学習モデルの構築をJuliaで実装します。
■分類器構築の流れ
以下のような流れを考え、スパムデータを分類する分類器を構築する。
1. 学習データの取得
2. 前処理の前半戦
2−1. ラベル付け
2−2. 特徴選択・抽出
3. 前処理の後半戦
3−1. データのスケール調整
3−2. 不均衡データ処理
4. 分類器の構築
4−1. パラメータチューニング
4−2. 精度評価
■1.学習データの取得
ここでは、HP研が収集したSpam E-mail Databaseを扱う。Rのkernlabパッケージにあるspamデータです。
- データ数：4601通（spam：1813通, non-spam：2788通）
- 1〜57列目：各メールの特徴量（単語や記号文字の出現頻度、大文字の連なりの長さ etc）
- 58列目：各メールのラベル（spam ,non-spamを表す）
読み込みは以下のように実行。
（Cドライブ直下にspam.csvが置いてある想定）
##データを読み込む
Pkg.add("DataFrames") #初回のみ実行
using DataFrames
spam = readtable("C:spam.csv", header = true)
head(spam) #データ数行確認
■2.前処理の前半戦
よく教科書にあるデータでは、当然のように、ラベルも特徴量も所与。
ここで扱う上述のスパムデータも同様。
しかしながら、機械学習の一番大変なところや本質はここにあるといえる。
2-1. ラベル付け
ラベルとは正例, 負例（spam, nonspam）のこと。もちろんラベルは最初からあるわけではない。
Spam E-mail Databaseでは、HP研がGeorge（メールデータの提供者）にヒアリングをしながら、4601通のメール１つ１つにラベル付けを行ったと想定され、その手間は想像を絶する。
2-2. 特徴選択・抽出
機械学習とパターン認識において最も大切なフェーズ。どんなに優れたなアルゴリズムを使っても、分類に影響する特徴を見ないと分類できるはずがない。
Spam E-mail Databaseでは、本当の最初はメールそのものがあるだけ。メールから色んな単語を取ってきて、spam, non-spamの分類に影響する57種類の特徴量を作ったと思われる。
（参考）みにくいアヒルの子の定理
－何らかの形で特徴に重要性を考えたり、取捨選択しなければ、みにくいアヒルの子と普通のアヒルの子の区別もできない
■3.前処理の後半戦
さて、ここからようやく、ちゃんとJuliaコーディング。
とはいえ、まだまだデータの前処理。
3-1. データのスケール調整
全ての特徴量について、平均0, 分散１に揃える（など）の調整を行う処理。
取りうる値の範囲が特徴量により異なる場合において、範囲が大きい特徴量が分類に対して支配的になりうることを避けるために行う。
##特徴量のスケールを調整する
l = spam[:,58]
f = spam[:,1:57]
spam_scaled = DataFrame() #空き箱を用意
for i in 1:ncol(f)
spam_scaled[i] = (f[:,i]-mean(f[:,i]))/sqrt(var(f[:,i]))
end
spam_scaled[58] = l
3-2. 不均衡データ処理
正例と負例のデータ比に偏りがある場合、アルゴリズムが偏った学習をしてしまう。
これを避けるための処理を不均衡データ処理という。
例えば、10,000件のメールの内、spam10件、non-spam9,990件の場合、すべてnon-spamと判定しても正解率99.9%となるので、アルゴリズムの学習が偏りがちになる。
対処法は大きく以下の２つ。
1. 少ない方を間違えた時のペナルティを、多い方を間違えた時より大きくする
－Weighted SVM など
2. データ数を調整して正例数＝負例数にする
－Over Sampling
－Under Sampling
今回は、Over Samplingで対応。
##オーバーサンプリング srand(123) #乱数シード固定 spam_scaled_p = spam_scaled[spam_scaled[:x58] .== "spam",:] #正例のデータのみ抽出 spam_scaled_f = spam_scaled[spam_scaled[:x58] .== "nonspam",:] #負例のデータのみ抽出 #正例数と負例数の差分を算出 imbalance = nrow(spam_scaled_f) - nrow(spam_scaled_p) #正例データから、先ほどの差分数をサンプリング（IDを決定） samplingID = rand(1:nrow(spam_scaled_p), imbalance) #正例データから、先ほどの差分数をサンプリング（当該IDを抽出） spam_scaled_p_over = spam_scaled_p[samplingID,:] spam_scaled_sampled = rbind(spam_scaled_p, spam_scaled_p_over, spam_scaled_f)
■4. 分類器の構築
4-1. パラメータチューニング
機械学習による分類器には、SVM（Support Vector Machine）を用いる。
パラメータチューニングは、グリッドサーチと交差検証で行う。
まずは準備。パッケージを読みこんで、データをSVMパッケージの形式に合わせる。
##パラメータチューニングの準備
#各種ライブラリの読み込み
Pkg.add("SVM") #初回のみ実行
using SVM
Pkg.add("MLBase") #初回のみ実行
using MLBase
#SVMの入力形式に合わせ、ラベルを"spam", "nonspam"から、1, -1に変更
label = DataFrame(zeros(nrow(spam_scaled_sampled),1))
for i in 1:nrow(spam_scaled_sampled)
if spam_scaled_sampled[i,58] == "spam"
label[i,1] = 1
else
label[i,1] = -1
end
end
feature = spam_scaled_sampled[:,1:57]
では、いよいよ、機械学習。
ここでは以下の２つのパッケージを用いる。
SVM.jl：SVMの中でもPegasosアルゴリズムで双対問題を解く実装を採用しているSVM.jlパッケージ。線形カーネルを用いたオンラインSVMにあたる。ちなみに、JuliaではLIBSVMの実装もある。
MLBase：グリッドサーチと交差検証には、機械学習手続きのフレームワークを提供するMLBaseパッケージを活用。 （14日目のsfchaosさんの記事も、是非ご参照ください）
この２つのパッケージを使って、交差検証とグリッドサーチによるパラメータチューニングを行う。
精度評価の指標は、Accuracyを用いる。
##交差検証によるパラメータのグリッドサーチ
#交差検証は3hold。Kfold関数でデータIDを分割
cv = 3
datanum = nrow(spam_scaled_sampled)
gen = collect(Kfold(datanum, cv))
model_dict = Dict() #交差検証ごとのモデルを格納する空ディクショナリ
score_dict = Dict() #交差検証ごとのスコアを格納する空ディクショナリ
#モデル推定関数の定義
#3holdの交差検証で3種類のデータセットでモデルを構築して格納
function estfun(lambda, k, T)
for k in 1:cv
learnID = gen[k]
feature_learn = transpose(array(feature[learnID,:]))
label_learn = array(label[learnID,:])[:,1]
model_dict[k] = svm(feature_learn, label_learn, lambda = lambda, k = k, T = T)
end
return model_dict
end
#精度評価関数
#3holdの交差検証で3種類のデータセットでAccuracyを算出。3回の平均値で評価
function evalfun(model_dict)
for k in 1:cv
testID = setdiff(1:datanum, gen[k])
feature_test = transpose(array(feature[testID,:]))
label_test = array(label[testID,:])[:,1]
label_pred = predict(model_dict[k], feature_test)
score_dict[k] = correctrate(int(label_test), int(label_pred))
end
return mean(values(score_dict))
end
#グリッドサーチの実行
out = gridtune(estfun, evalfun,
("lambda", [0.0001, 0.001, 0.01, 0.1, 1.0]),
("k", [5, 10, 50, 100, 500, 1000, 5000]),
("T", [100, 500, 1000, 5000, 10000, 50000]);
ord=Forward, # 精度評価関数の評価値は、昇順か、降順かを指定
verbose=true) # 推定毎の結果出力の有無を指定
#ベストモデル、ベストパラメータ、ベストスコアを出力
best_model, best_cfg, best_score = out
4-2. 精度評価
上記の出力の結果を示す。
ベストパラメータ
| lambda | k | T |
|---|---|---|
| 0.01 | 5000 | 50000 |
ベストスコア
| Accuracy |
|---|
| 0.922 |
ちなみに、svm.jlのsrcコードを覗くとわかるが、デフォルトパラメータは、lambda:0.1, k：5, T:100であり、この場合、今回のデータではAccuracyが0.871であった。
パラメータチューニング大事。
なお、MLBaseでは、Accuracy以外にも、混合行列、Precision, Recall, ROC曲線なども精度指標として出力できる。
最後に、ベストパラメータモデルで、混合行列と各種指標を算出する。
##ベストパラメータでの各種評価指標を出力
#ベストパラメータモデルを作成（データは交差検証の1hold目を用いる）
lambda = best_cfg[1]
k = best_cfg[2]
T = best_cfg[3]
#交差検証ごとの混合行列の計算
cmatrix = [0 0; 0 0]
for k in 1:cv
learnID = gen[k]
feature_learn = transpose(array(feature[learnID,:]))
label_learn = array(label[learnID,:])[:,1]
testID = setdiff(1:datanum, gen[k])
feature_test = transpose(array(feature[testID,:]))
label_test = array(label[testID,:])[:,1]
best_model = svm(feature_learn, label_learn, lambda = lambda, k = k, T = T)
pred = int(predict(best_model, feature_test))
#分類の出力が1, 2, …じゃないとダメぽい。-1, 1だったので変換して入力
gt0 = int(label_test)
gt = DataFrame(zeros(length(gt0),1))
for i in 1:length(gt0)
if gt0[i,1] == 1
gt[i,1] = 1
else
gt[i,1] = 2
end
end
pred0 = int(predict(best_model, feature_test))
pred = DataFrame(zeros(length(pred0),1))
for i in 1:length(pred0)
if pred0[i,1] == 1
pred[i,1] = 1
else
pred[i,1] = 2
end
end
cmatrix = cmatrix + confusmat(2, int(array(gt))[:], int(array(pred))[:])
end
#混合行列から、precision, recall, Fvalueを算出
precision = cmatrix[1, 1]/sum(cmatrix[:, 1])
recall = cmatrix[1, 1]/sum(cmatrix[1, :])
Fvalue = 2/(1/precision + 1/recall)
結果は以下。
混合行列
2537 251 182 2606
混合行列から各種指標を算出
| precision | recall | Fvalue |
|---|---|---|
| 0.933 | 0.910 | 0.921 |
どうもありがとうございました。
あと１週間、アドベントカレンダーを楽しみましょう。メリークリスマス！
■参考文献
@chezou, Julia v0.3.0でSVM.jlを使う
@sfchaos, Juliaによる機械学習の予測モデル構築・評価
@sfchaos, 不均衡データのクラス分類
@sleepy_yoshi, SVM実践ガイド (A Practical Guide to Support Vector Classification)
さいごの碧, kernlabパッケージのspamデータ
